Merge release-2019 branch into master
[alexxy/gromacs.git] / src / gromacs / pulling / pull.cpp
index 1f9722efc9dc5953a8c00d9665e634316e51bdf8..2b6017d56efc99399746ccf8d5847970a76738a5 100644 (file)
@@ -65,6 +65,7 @@
 #include "gromacs/mdtypes/inputrec.h"
 #include "gromacs/mdtypes/md_enums.h"
 #include "gromacs/mdtypes/mdatom.h"
+#include "gromacs/mdtypes/state.h"
 #include "gromacs/pbcutil/pbc.h"
 #include "gromacs/topology/mtop_lookup.h"
 #include "gromacs/topology/topology.h"
@@ -1691,6 +1692,21 @@ void dd_make_local_pull_groups(const t_commrec *cr, struct pull_t *pull)
 
             comm->bParticipate = bWillParticipate;
             comm->nparticipate = count[0];
+
+            /* When we use the previous COM for PBC, we need to broadcast
+             * the previous COM to ranks that have joined the communicator.
+             */
+            for (pull_group_work_t &group : pull->group)
+            {
+                if (group.epgrppbc == epgrppbcPREVSTEPCOM)
+                {
+                    GMX_ASSERT(comm->bParticipate || !MASTER(cr),
+                               "The master rank has to participate, as it should pass an up to date prev. COM "
+                               "to bcast here as well as to e.g. checkpointing");
+
+                    gmx_bcast(sizeof(dvec), group.x_prev_step, cr);
+                }
+            }
         }
     }
 
@@ -1877,7 +1893,7 @@ init_pull(FILE *fplog, const pull_params_t *pull_params, const t_inputrec *ir,
         /* Set up the global to local atom mapping for PBC atoms */
         for (pull_group_work_t &group : pull->group)
         {
-            if (group.epgrppbc == epgrppbcREFAT)
+            if (group.epgrppbc == epgrppbcREFAT || group.epgrppbc == epgrppbcPREVSTEPCOM)
             {
                 /* pbcAtomSet consists of a single atom */
                 group.pbcAtomSet = gmx::compat::make_unique<gmx::LocalAtomSet>(atomSets->add({&group.params.pbcatom, &group.params.pbcatom + 1}));
@@ -1893,7 +1909,6 @@ init_pull(FILE *fplog, const pull_params_t *pull_params, const t_inputrec *ir,
     pull->bFOutAverage = pull_params->bFOutAverage;
 
     GMX_RELEASE_ASSERT(pull->group[0].params.nat == 0, "pull group 0 is an absolute reference group and should not contain atoms");
-    pull->group[0].x_prev_step[XX] = NAN;
 
     pull->numCoordinatesWithExternalPotential = 0;
 
@@ -2160,8 +2175,6 @@ init_pull(FILE *fplog, const pull_params_t *pull_params, const t_inputrec *ir,
             init_pull_group_index(fplog, cr, g, pgrp,
                                   bConstraint, pulldim_con,
                                   mtop, ir, lambda);
-
-            pgrp->x_prev_step[XX] = NAN;
         }
         else
         {
@@ -2226,7 +2239,7 @@ init_pull(FILE *fplog, const pull_params_t *pull_params, const t_inputrec *ir,
     }
 
     comm->pbcAtomBuffer.resize(pull->group.size());
-    comm->comBuffer.resize(pull->group.size()*DIM);
+    comm->comBuffer.resize(pull->group.size()*c_comBufferStride);
     if (pull->bCylinder)
     {
         comm->cylinderBuffer.resize(pull->coord.size()*c_cylinderBufferStride);
@@ -2253,6 +2266,35 @@ static void destroy_pull(struct pull_t *pull)
     delete pull;
 }
 
+void preparePrevStepPullCom(const t_inputrec *ir, const t_mdatoms *md, t_state *state, const t_state *state_global, const t_commrec *cr, bool startingFromCheckpoint)
+{
+    if (!ir->pull || !ir->pull->bSetPbcRefToPrevStepCOM)
+    {
+        return;
+    }
+    allocStatePrevStepPullCom(state, ir->pull_work);
+    if (startingFromCheckpoint)
+    {
+        if (MASTER(cr))
+        {
+            state->pull_com_prev_step = state_global->pull_com_prev_step;
+        }
+        if (PAR(cr))
+        {
+            /* Only the master rank has the checkpointed COM from the previous step */
+            gmx_bcast(sizeof(double) * state->pull_com_prev_step.size(), &state->pull_com_prev_step[0], cr);
+        }
+        setPrevStepPullComFromState(ir->pull_work, state);
+    }
+    else
+    {
+        t_pbc pbc;
+        set_pbc(&pbc, ir->ePBC, state->box);
+        initPullComFromPrevStep(cr, ir->pull_work, md, &pbc, state->x.rvec_array());
+        updatePrevStepPullCom(ir->pull_work, state);
+    }
+}
+
 void finish_pull(struct pull_t *pull)
 {
     check_external_potential_registration(pull);