Fixed pull coord init with rate!=0, t_start!=0
authorBerk Hess <hess@kth.se>
Tue, 2 Dec 2014 10:14:22 +0000 (11:14 +0100)
committerGerrit Code Review <gerrit@gerrit.gromacs.org>
Mon, 8 Dec 2014 21:12:02 +0000 (22:12 +0100)
When the start time of a simulation is non-zero a pull setup with
coordinates with rate!=0 and pull-start, grompp would generate
incorrect pull init reference values.
Also fixed the layout of the grompp pull coord information table.

Change-Id: I992a40b098853aba0922ec1d596430eb229d64a8

src/gromacs/gmxpreprocess/readpull.c

index 4ae7480f0f32d3bb46c33ee5747e14bcbc81b1a1..26c71986049d080054e1f75f8e86ca10d65ce7a0 100644 (file)
@@ -395,10 +395,10 @@ void set_pull_init(t_inputrec *ir, gmx_mtop_t *mtop, rvec *x, matrix box, real l
     t_pull_group *pgrp0, *pgrp1;
     t_pbc         pbc;
     int           c, m;
-    double        t_start, tinvrate;
-    real          init;
+    double        t_start;
+    real          init = 0;
     dvec          dr;
-    double        dev;
+    double        dev, value;
 
     init_pull(NULL, ir, 0, NULL, mtop, NULL, oenv, lambda, FALSE, 0);
     md = init_mdatoms(NULL, mtop, ir->efep);
@@ -415,7 +415,7 @@ void set_pull_init(t_inputrec *ir, gmx_mtop_t *mtop, rvec *x, matrix box, real l
 
     pull_calc_coms(NULL, pull, md, &pbc, t_start, x, NULL);
 
-    fprintf(stderr, "Pull group  natoms  pbc atom  distance at start     reference at t=0\n");
+    fprintf(stderr, "Pull group  natoms  pbc atom  distance at start  reference at t=0\n");
     for (c = 0; c < pull->ncoord; c++)
     {
         pcrd  = &pull->coord[c];
@@ -427,28 +427,21 @@ void set_pull_init(t_inputrec *ir, gmx_mtop_t *mtop, rvec *x, matrix box, real l
         fprintf(stderr, "%8d  %8d  %8d ",
                 pcrd->group[1], pgrp1->nat, pgrp1->pbcatom+1);
 
-        init       = pcrd->init;
-        pcrd->init = 0;
-
-        if (pcrd->rate == 0)
-        {
-            tinvrate = 0;
-        }
-        else
+        if (bStart)
         {
-            tinvrate = t_start/pcrd->rate;
+            init       = pcrd->init;
+            pcrd->init = 0;
         }
-        get_pull_coord_distance(pull, c, &pbc, 0, dr, &dev);
-        fprintf(stderr, " %6.3f ", dev);
+
+        get_pull_coord_distance(pull, c, &pbc, t_start, dr, &dev);
+        /* Calculate the value of the coordinate at t_start */
+        value = pcrd->init + t_start*pcrd->rate + dev;
+        fprintf(stderr, " %10.3f nm", value);
 
         if (bStart)
         {
-            pcrd->init = init + dev - tinvrate;
-        }
-        else
-        {
-            pcrd->init = init;
+            pcrd->init = value + init;
         }
-        fprintf(stderr, " %6.3f\n", pcrd->init);
+        fprintf(stderr, "     %10.3f nm\n", pcrd->init);
     }
 }