Fix perturbation to Verlet scheme when not perturbing
authorMark Abraham <mark.j.abraham@gmail.com>
Sun, 16 Mar 2014 12:05:28 +0000 (13:05 +0100)
committerGerrit Code Review <gerrit@gerrit.gromacs.org>
Mon, 17 Mar 2014 15:43:27 +0000 (16:43 +0100)
A .top that has B-state parameters for non-bonded interactions
should run the same with free-energy = no as a .top that has no
B-state parameters. Renamed variable to better reflect its use,
and used it to correctly trigger the FEP code path.

Change-Id: I47656b29ef4ca2462af0fb587e0424d0eab3d5c4

src/gromacs/mdlib/forcerec.c

index 95fbb19bab2d167d4f3c6400ff650d1dc907bae6..3f4c0d8ddc20b4558ef5830bd95f2b1f2f2ce64c 100644 (file)
@@ -606,7 +606,7 @@ static cginfo_mb_t *init_cginfo_mb(FILE *fplog, const gmx_mtop_t *mtop,
     int                  *a_con;
     int                   ftype;
     int                   ia;
-    gmx_bool              bId, *bExcl, bExclIntraAll, bExclInter, bHaveVDW, bHaveQ, bFEP;
+    gmx_bool              bId, *bExcl, bExclIntraAll, bExclInter, bHaveVDW, bHaveQ, bHavePerturbedAtoms;
 
     ncg_tot = ncg_mtop(mtop);
     snew(cginfo_mb, mtop->nmolblock);
@@ -720,11 +720,11 @@ static cginfo_mb_t *init_cginfo_mb(FILE *fplog, const gmx_mtop_t *mtop,
                 /* bExclIntraAll: all intra cg interactions excluded
                  * bExclInter:    any inter cg interactions excluded
                  */
-                bExclIntraAll = TRUE;
-                bExclInter    = FALSE;
-                bHaveVDW      = FALSE;
-                bHaveQ        = FALSE;
-                bFEP          = FALSE;
+                bExclIntraAll       = TRUE;
+                bExclInter          = FALSE;
+                bHaveVDW            = FALSE;
+                bHaveQ              = FALSE;
+                bHavePerturbedAtoms = FALSE;
                 for (ai = a0; ai < a1; ai++)
                 {
                     /* Check VDW and electrostatic interactions */
@@ -733,7 +733,7 @@ static cginfo_mb_t *init_cginfo_mb(FILE *fplog, const gmx_mtop_t *mtop,
                     bHaveQ  = bHaveQ    || (molt->atoms.atom[ai].q != 0 ||
                                             molt->atoms.atom[ai].qB != 0);
 
-                    bFEP    = bFEP || (PERTURBED(molt->atoms.atom[ai]) != 0);
+                    bHavePerturbedAtoms = bHavePerturbedAtoms || (PERTURBED(molt->atoms.atom[ai]) != 0);
 
                     /* Clear the exclusion list for atom ai */
                     for (aj = a0; aj < a1; aj++)
@@ -795,7 +795,7 @@ static cginfo_mb_t *init_cginfo_mb(FILE *fplog, const gmx_mtop_t *mtop,
                 {
                     SET_CGINFO_HAS_Q(cginfo[cgm+cg]);
                 }
-                if (bFEP)
+                if (bHavePerturbedAtoms && fr->efep != efepNO)
                 {
                     SET_CGINFO_FEP(cginfo[cgm+cg]);
                     *bFEP_NonBonded = TRUE;