Merge branch 'rotation-4-5' into rotation

diff --git a/src/gmxlib/checkpoint.c b/src/gmxlib/checkpoint.c
index c31137433a262beef2f9b0c62ec2206ee72fed9b..2218171eb63c82191ad532b4c281288c7806b801 100644 (file)
--- a/src/gmxlib/checkpoint.c
+++ b/src/gmxlib/checkpoint.c
@@ -1745,8 +1745,16 @@ static void read_checkpoint(const char *fn,FILE **pfplog,
                 if (_locking(fileno(gmx_fio_getfp(chksum_file)), _LK_NBLCK, LONG_MAX)==-1)
 #endif
                 {
-                    gmx_fatal(FARGS,"Failed to lock: %s. Already running "
-                        "simulation?", outputfiles[i].filename);
+                    if (errno!=EACCES && errno!=EAGAIN)
+                    {
+                        gmx_fatal(FARGS,"Failed to lock: %s. %s.",
+                                  outputfiles[i].filename, strerror(errno));
+                    }
+                    else 
+                    {
+                        gmx_fatal(FARGS,"Failed to lock: %s. Already running "
+                                  "simulation?", outputfiles[i].filename);
+                    }
                 }
             }
             

diff --git a/src/gmxlib/futil.c b/src/gmxlib/futil.c
index f06e09892c0c043c648a333adfcd82d10cc11cd7..1c58f6d73e17a80610bc5c2afa3fe370bd230eef 100644 (file)
--- a/src/gmxlib/futil.c
+++ b/src/gmxlib/futil.c
@@ -454,7 +454,7 @@ FILE *ffopen(const char *file,const char *mode)
 
     bRead= (mode[0]=='r'&&mode[1]!='+');
     strcpy(buf,file);
-    if (gmx_fexist(buf) || !bRead) {
+    if (!bRead || gmx_fexist(buf)) {
         if ((ff=fopen(buf,mode))==NULL)
             gmx_file(buf);
         where();

diff --git a/src/gmxlib/main.c b/src/gmxlib/main.c
index b50ded75b68478fffa99c586a93830d0a08e6ec9..aeda0e0a795081d4d37aa32a236e46c9db611bbe 100644 (file)
--- a/src/gmxlib/main.c
+++ b/src/gmxlib/main.c
@@ -238,17 +238,13 @@ void gmx_log_open(const char *lognm,const t_commrec *cr,gmx_bool bMasterOnly,
   
     debug_gmx();
 
-    if (!bMasterOnly)
+    if (!bMasterOnly && !MASTER(cr))
     {
         /* Since log always ends with '.log' let's use this info */
         par_fn(tmpnm,efLOG,cr,FALSE,!bMasterOnly,buf,255);
         fp = gmx_fio_fopen(buf, bAppend ? "a+" : "w+" );
     }
-#ifdef GMX_FAHCORE
     else if (!bAppend)
-#else
-    else
-#endif
     {
         fp = gmx_fio_fopen(tmpnm, bAppend ? "a+" : "w+" );
     }

diff --git a/src/kernel/grompp.c b/src/kernel/grompp.c
index ee14674d5c8dcd66988efadd7bf35c851cd5a1c4..17772bc9a9b9474871b45791c10a424724f2bfda 100644 (file)
--- a/src/kernel/grompp.c
+++ b/src/kernel/grompp.c
@@ -1529,7 +1529,11 @@ int main (int argc, char *argv[])
     set_pull_init(ir,sys,state.x,state.box,oenv,opts->pull_start);
   
   if (ir->bRot)
-      set_reference_positions(ir->rot,sys,state.x,state.box,opt2fn("-ref",NFILE,fnm),wi);
+  {
+      set_reference_positions(ir->rot,sys,state.x,state.box,
+                              opt2fn("-ref",NFILE,fnm),opt2bSet("-ref",NFILE,fnm),
+                              wi);
+  }
 
   /*  reset_multinr(sys); */
   

diff --git a/src/kernel/mdrun.c b/src/kernel/mdrun.c
index 5cc5a964e305243dc0816aaab2527e8959830156..ed6db41ded10ef3fd1efc8b05b102283a7f6655e 100644 (file)
--- a/src/kernel/mdrun.c
+++ b/src/kernel/mdrun.c
@@ -566,6 +566,11 @@ int main(int argc,char *argv[])
   sim_part_fn = sim_part;
   if (opt2bSet("-cpi",NFILE,fnm))
   {
+      if (bSepPot && bAppendFiles)
+      {
+          gmx_fatal(FARGS,"Output file appending is not supported with -seppot");
+      }
+
       bAppendFiles =
                 read_checkpoint_simulation_part(opt2fn_master("-cpi", NFILE,
                                                               fnm,cr),

diff --git a/src/kernel/readir.h b/src/kernel/readir.h
index a4730f9224f9b16f7833d1b3b983fd30ed6a0795..3f99910756986cca0381b202dd243f1c84a47303 100644 (file)
--- a/src/kernel/readir.h
+++ b/src/kernel/readir.h
@@ -144,6 +144,6 @@ extern void make_rotation_groups(t_rot *rot,char **rotgnames,
 /* Process the rotation parameters after reading the index groups */
 
 extern void set_reference_positions(t_rot *rot, gmx_mtop_t *mtop, rvec *x, matrix box,
-        const char *fn, warninp_t wi);
+        const char *fn, gmx_bool bSet, warninp_t wi);
 
 #endif /* _readir_h */

diff --git a/src/kernel/readrot.c b/src/kernel/readrot.c
index f2873a8d998cb28588499570de726627b1a2e609..dfde11268839e46370471877c76d8442192bb3b1 100644 (file)
--- a/src/kernel/readrot.c
+++ b/src/kernel/readrot.c
@@ -42,6 +42,8 @@
 #include "trnio.h"
 #include "txtdump.h"
 
+static char *RotStr = {"Enforced rotation:"};
+
 
 static char s_vec[STRLEN];
 
@@ -113,8 +115,8 @@ extern char **read_rotparams(int *ninp_p,t_inpfile **inp_p,t_rot *rot,
             sprintf(warn_buf, "rot_vec%d = 0", g);
             warning_error(wi, warn_buf);
         }
-        fprintf(stderr, "Enforced rotation: Group %d (%s) normalized rot. vector: %f %f %f\n", 
-                g, erotg_names[rotg->eType], vec[0], vec[1], vec[2]);
+        fprintf(stderr, "%s Group %d (%s) normalized rot. vector: %f %f %f\n",
+                RotStr, g, erotg_names[rotg->eType], vec[0], vec[1], vec[2]);
         for(m=0; m<DIM; m++)
             rotg->vec[m] = vec[m];
         
@@ -179,7 +181,7 @@ extern char **read_rotparams(int *ninp_p,t_inpfile **inp_p,t_rot *rot,
 
 
 /* Check whether the box is unchanged */
-static void check_box(matrix f_box, matrix box, char fn[],warninp_t wi)
+static void check_box(matrix f_box, matrix box, char fn[], warninp_t wi)
 {
     int i,ii;
     gmx_bool bSame=TRUE;
@@ -192,7 +194,8 @@ static void check_box(matrix f_box, matrix box, char fn[],warninp_t wi)
                 bSame = FALSE;
     if (!bSame)
     {
-        sprintf(warn_buf, "Enforced rotation: Box size in reference file %s differs from actual box size!", fn);
+        sprintf(warn_buf, "%s Box size in reference file %s differs from actual box size!",
+                RotStr, fn);
         warning(wi, warn_buf);
         pr_rvecs(stderr,0,"Your box is:",box  ,3);
         pr_rvecs(stderr,0,"Box in file:",f_box,3);
@@ -201,8 +204,9 @@ static void check_box(matrix f_box, matrix box, char fn[],warninp_t wi)
 
 
 /* Extract the reference positions for the rotation group(s) */
-extern void set_reference_positions(t_rot *rot, gmx_mtop_t *mtop, rvec *x, matrix box,
-        const char *fn,warninp_t wi)
+extern void set_reference_positions(
+        t_rot *rot, gmx_mtop_t *mtop, rvec *x, matrix box,
+        const char *fn, gmx_bool bSet, warninp_t wi)
 {
     int g,i,ii;
     t_rotgrp *rotg;
@@ -218,16 +222,28 @@ extern void set_reference_positions(t_rot *rot, gmx_mtop_t *mtop, rvec *x, matri
     strcpy(extension,extpos+1);
     *extpos = '\0';
 
+
     for (g=0; g<rot->ngrp; g++)
      {
          rotg = &rot->grp[g];
-         fprintf(stderr, "Enforced rotation: group %d has %d reference positions.",g,rotg->nat);
+         fprintf(stderr, "%s group %d has %d reference positions.\n",RotStr,g,rotg->nat);
          snew(rotg->x_ref, rotg->nat);
          
+         /* Construct the name for the file containing the reference positions for this group: */
          sprintf(reffile, "%s.%d.%s", base,g,extension);
+
+         /* If the base filename for the reference position files was explicitly set by
+          * the user, we issue a fatal error if the group file can not be found */
+         if (bSet && !gmx_fexist(reffile))
+         {
+             gmx_fatal(FARGS, "%s The file containing the reference positions was not found.\n"
+                              "Expected the file '%s' for group %d.\n",
+                              RotStr, reffile, g);
+         }
+
          if (gmx_fexist(reffile))
          {
-             fprintf(stderr, " Reading them from %s.\n", reffile);
+             fprintf(stderr, "  Reading them from %s.\n", reffile);
              read_trnheader(reffile, &header);
              if (rotg->nat != header.natoms)
                  gmx_fatal(FARGS,"Number of atoms in file %s (%d) does not match the number of atoms in rotation group (%d)!\n",

diff --git a/src/mdlib/pull_rotation.c b/src/mdlib/pull_rotation.c
index 00dda8707ecdec44174d7a1e8a5ea7cad1828713..243a00f1ce92e58d12f2e3bbb4a521e9c85d6a44 100644 (file)
--- a/src/mdlib/pull_rotation.c
+++ b/src/mdlib/pull_rotation.c
@@ -1155,9 +1155,8 @@ static double opt_angle_analytic(
 }
 
 
-/* Determine actual angle of this slab by RMSD fit */
+/* Determine actual angle of this slab by RMSD fit to the reference */
 /* Not parallelized, call this routine only on the master */
-/* TODO: make this routine mass-weighted */
 static void flex_fit_angle(
         int  g,
         t_rotgrp *rotg,
@@ -1176,15 +1175,20 @@ static void flex_fit_angle(
     rvec        coord;
     real        scal;
     gmx_enfrotgrp_t erg;     /* Pointer to enforced rotation group data */
+    real        OOm_av;      /* 1/average_mass of a rotation group atom */
+    real        m_rel;       /* Relative mass of a rotation group atom  */
 
     
     erg=rotg->enfrotgrp;
 
+    /* Average mass of a rotation group atom: */
+    OOm_av = erg->invmass*rotg->nat;
+
     /**********************************/
     /* First collect the data we need */
     /**********************************/
 
-    /* Loop over slabs */
+    /* Collect the data for the individual slabs */
     for (n = erg->slab_first; n <= erg->slab_last; n++)
     {
         islab = n - erg->slab_first; /* slab index */
@@ -1199,7 +1203,7 @@ static void flex_fit_angle(
             copy_rvec(erg->xc[l], curr_x);
             
             /* The (unrotated) reference position of this atom is copied to ref_x.
-             * Beware, the xc coords have been sorted in do_flex! */
+             * Beware, the xc coords have been sorted in do_flexible */
             copy_rvec(erg->xc_ref_sorted[l], ref_x);
             
             /* Save data for doing angular RMSD fit later */
@@ -1207,17 +1211,23 @@ static void flex_fit_angle(
             copy_rvec(curr_x, sd->x[ind]);
             /* Save the corresponding reference position */
             copy_rvec(ref_x , sd->ref[ind]);
-            /* Save the weight for this atom in this slab */
-            sd->weight[ind] = gaussian_weight(curr_x, rotg, n);
             
+            /* Maybe also mass-weighting was requested. If yes, additionally
+             * multiply the weights with the relative mass of the atom. If not,
+             * multiply with unity. */
+            m_rel = erg->mc_sorted[l]*OOm_av;
+
+            /* Save the weight for this atom in this slab */
+            sd->weight[ind] = gaussian_weight(curr_x, rotg, n) * m_rel;
+
             /* Next atom in this slab */
             ind++;
         }
     }
 
-    /* Get the geometrical center of the whole rotation group: */
-    get_center(erg->xc, NULL, rotg->nat, act_center);
-
+    /* Get the center of the whole rotation group. Note, again, the erg->xc have
+     * been sorted in do_flexible */
+    get_center(erg->xc, erg->mc_sorted, rotg->nat, act_center);
 
     /******************************/
     /* Now do the fit calculation */
@@ -1230,7 +1240,7 @@ static void flex_fit_angle(
          * prior to performing the fit */
         for (i=0; i<rotg->nat; i++)
         {
-            /* First put geometrical center of positions into origin */
+            /* First put the center of the positions into the origin */
             rvec_sub(erg->xc[i], act_center, coord);
             /* Determine the scaling factor for the length: */
             scal = erg->xc_ref_length[erg->xc_sortind[i]] / norm(coord);
@@ -1245,7 +1255,8 @@ static void flex_fit_angle(
     }
     /* Note that from the point of view of the current positions, the reference has rotated backwards,
      * but we want to output the angle relative to the fixed reference, therefore the minus sign. */
-    fitangle = -opt_angle_analytic(erg->xc_ref_sorted, fitcoords, NULL, rotg->nat, erg->xc_ref_center, act_center, rotg->vec);    
+    fitangle = -opt_angle_analytic(erg->xc_ref_sorted, fitcoords, erg->mc_sorted,
+                                   rotg->nat, erg->xc_ref_center, act_center, rotg->vec);
     fprintf(fp, "%12.3e%6d%12.3f%12.3lf", t, g, degangle, fitangle);