Fixed the execution order in OpenMM so that the generated trajectory frames correspon...
authorRossen Apostolov <rossen@cbr.su.se>
Wed, 18 Aug 2010 14:53:46 +0000 (16:53 +0200)
committerRossen Apostolov <rossen@cbr.su.se>
Wed, 18 Aug 2010 14:53:46 +0000 (16:53 +0200)
src/kernel/md_openmm.c
src/kernel/openmm_wrapper.cpp
src/kernel/openmm_wrapper.h

index 2e7c0938c65d4743227ff7d0ab4c2d7358c56bbd..2aefbbb927eb9f75146322bcfdde6d89f9e23e46 100644 (file)
@@ -515,26 +515,15 @@ double do_md_openmm(FILE *fplog,t_commrec *cr,int nfile,const t_filenm fnm[],
         };
         do_ene = (do_per_step(step,ir->nstenergy) || bLastStep);
 
-        if (mdof_flags != 0)
+        if (mdof_flags != 0 || do_ene || do_log)
         {
+            wallcycle_start(wcycle,ewcTRAJ);
+            bF = (mdof_flags & MDOF_F);
             bX = (mdof_flags & (MDOF_X | MDOF_XTC | MDOF_CPT));
             bV = (mdof_flags & (MDOF_V | MDOF_CPT));
 
-            wallcycle_start(wcycle,ewcTRAJ);
-            openmm_copy_state(openmmData, state, &t, f, enerd, bX, bV, 0, 0);
-            wallcycle_stop(wcycle,ewcTRAJ);
-        }
-
-        openmm_take_one_step(openmmData);
+            openmm_copy_state(openmmData, state, &t, f, enerd, bX, bV, bF, do_ene);
 
-        if (mdof_flags != 0 || do_ene || do_log)
-        {
-            wallcycle_start(wcycle,ewcTRAJ);
-            bF = (mdof_flags & MDOF_F);
-            if (bF || do_ene || do_log)
-            {
-                openmm_copy_state(openmmData, state, &t, f, enerd, 0, 0, bF, do_ene);
-            }
             upd_mdebin(mdebin, NULL,TRUE,
                        t,mdatoms->tmass,enerd,state,lastbox,
                        shake_vir,force_vir,total_vir,pres,
@@ -656,6 +645,8 @@ double do_md_openmm(FILE *fplog,t_commrec *cr,int nfile,const t_filenm fnm[],
         bStartingFromCpt = FALSE;
         step++;
         step_rel++;
+
+        openmm_take_one_step(openmmData);
     }
     /* End of main MD loop */
     debug_gmx();
index 2fd6de3f58862942abf306f38eb7f49d5573e980..5a2aa22e0839d3bbaad290c9791c2808a35dcb1f 100644 (file)
@@ -857,7 +857,7 @@ void* openmm_init(FILE *fplog, const char *platformOptStr,
         if (ir->nstcomm > 0)
             sys->addForce(new CMMotionRemover(ir->nstcomm));
 
-        // Set bonded force field terms.
+        /* Set bonded force field terms. */
         const int* bondAtoms = (int*) idef.il[F_BONDS].iatoms;
         HarmonicBondForce* bondForce = new HarmonicBondForce();
         sys->addForce(bondForce);
@@ -870,7 +870,8 @@ void* openmm_init(FILE *fplog, const char *platformOptStr,
             bondForce->addBond(atom1, atom2,
                                idef.iparams[type].harmonic.rA, idef.iparams[type].harmonic.krA);
         }
-        // Urey-Bradley includes both the angle and bond potential for 1-3 interactions
+
+        /* Urey-Bradley includes both the angle and bond potential for 1-3 interactions */
         const int* ubAtoms = (int*) idef.il[F_UREY_BRADLEY].iatoms;
         HarmonicBondForce* ubBondForce = new HarmonicBondForce();
         HarmonicAngleForce* ubAngleForce = new HarmonicAngleForce();
@@ -888,6 +889,8 @@ void* openmm_init(FILE *fplog, const char *platformOptStr,
             ubAngleForce->addAngle(atom1, atom2, atom3, 
                     idef.iparams[type].u_b.theta*M_PI/180.0, idef.iparams[type].u_b.ktheta);
         }
+
+               /* Set the angle force field terms */
         const int* angleAtoms = (int*) idef.il[F_ANGLES].iatoms;
         HarmonicAngleForce* angleForce = new HarmonicAngleForce();
         sys->addForce(angleForce);
@@ -901,6 +904,8 @@ void* openmm_init(FILE *fplog, const char *platformOptStr,
             angleForce->addAngle(atom1, atom2, atom3, 
                     idef.iparams[type].harmonic.rA*M_PI/180.0, idef.iparams[type].harmonic.krA);
         }
+
+               /* Set proper dihedral terms */
         const int* periodicAtoms = (int*) idef.il[F_PDIHS].iatoms;
         PeriodicTorsionForce* periodicForce = new PeriodicTorsionForce();
         sys->addForce(periodicForce);
@@ -917,6 +922,8 @@ void* openmm_init(FILE *fplog, const char *platformOptStr,
                                       idef.iparams[type].pdihs.phiA*M_PI/180.0, 
                                       idef.iparams[type].pdihs.cpA);
         }
+
+               /* Set improper dihedral terms that are represented by a periodic function (as in AMBER FF) */
         const int* periodicImproperAtoms = (int*) idef.il[F_PIDIHS].iatoms;
         PeriodicTorsionForce* periodicImproperForce = new PeriodicTorsionForce();
         sys->addForce(periodicImproperForce);
@@ -933,6 +940,8 @@ void* openmm_init(FILE *fplog, const char *platformOptStr,
                                       idef.iparams[type].pdihs.phiA*M_PI/180.0,
                                       idef.iparams[type].pdihs.cpA);
         }
+
+        /* Ryckaert-Bellemans dihedrals */
         const int* rbAtoms = (int*) idef.il[F_RBDIHS].iatoms;
         RBTorsionForce* rbForce = new RBTorsionForce();
         sys->addForce(rbForce);
@@ -950,6 +959,7 @@ void* openmm_init(FILE *fplog, const char *platformOptStr,
                                 idef.iparams[type].rbdihs.rbcA[4], idef.iparams[type].rbdihs.rbcA[5]);
         }
 
+               /* Set improper dihedral terms (as in CHARMM FF) */
         const int* improperDihAtoms = (int*) idef.il[F_IDIHS].iatoms;
                CustomTorsionForce* improperDihForce = new CustomTorsionForce("2.0*k*asin(sin((theta-theta0)/2))^2");
         sys->addForce(improperDihForce);
@@ -970,7 +980,7 @@ void* openmm_init(FILE *fplog, const char *platformOptStr,
                                 improperDihParameters);
         }
 
-        // Set nonbonded parameters and masses.
+        /* Set nonbonded parameters and masses. */
         int ntypes = fr->ntype;
         int* types = mdatoms->typeA;
         real* nbfp = fr->nbfp;
index 359e5df2d4f0b895b2783186ddaa90c82d741366..16dd5856b8b4ebb3427d77c6cfdf5140d7f89d8a 100644 (file)
@@ -49,6 +49,8 @@ void* openmm_init(FILE *fplog, const char *platformOptStr,
 
 void openmm_take_one_step(void* data);
 
+void openmm_take_steps(void* data, int nsteps);
+
 void openmm_copy_state(void *data,
                         t_state *state, double *time,
                         rvec f[], gmx_enerdata_t *enerd,
@@ -65,6 +67,8 @@ void* openmm_init(FILE *fplog, const char *platformOptStr,
 
 void openmm_take_one_step(void* data){}
 
+void openmm_take_steps(void* data, int nsteps){}
+
 void openmm_copy_state(void *data,
                         t_state *state, double *time,
                         rvec f[], gmx_enerdata_t *enerd,