Decouple update-group aspects of vsites and constraints

Update groups should not be managed from within domdec module, so
shift the responsiblity to runner temporarily to simplify future
changes.

This made clear that some internal fields of domdec module were
inaccurately named. Update groups preclude the possibility of
constraints or vsites being split across domains, but when update
groups are not in use such splits are merely possible, not
assured. The code doesn't change, because we always have to account
for the possiblity that some are split.

Refs #4004

diff --git a/src/gromacs/domdec/domdec.cpp b/src/gromacs/domdec/domdec.cpp
index e6e3df58e390d4f74105f7b8d876339b7f4c6dc9..d684d8057e27002fe6b638bb87eddb725e42b0c8 100644 (file)
--- a/src/gromacs/domdec/domdec.cpp
+++ b/src/gromacs/domdec/domdec.cpp
@@ -982,9 +982,9 @@ int dd_pme_maxshift_y(const gmx_domdec_t& dd)
     }
 }
 
-bool ddHaveSplitConstraints(const gmx_domdec_t& dd)
+bool ddMayHaveSplitConstraints(const gmx_domdec_t& dd)
 {
-    return dd.comm->systemInfo.haveSplitConstraints;
+    return dd.comm->systemInfo.mayHaveSplitConstraints;
 }
 
 bool ddUsesUpdateGroups(const gmx_domdec_t& dd)
@@ -2051,14 +2051,14 @@ static DDSystemInfo getSystemInfo(const gmx::MDLogger&           mdlog,
 
     if (systemInfo.useUpdateGroups)
     {
-        systemInfo.haveSplitConstraints = false;
-        systemInfo.haveSplitSettles     = false;
+        systemInfo.mayHaveSplitConstraints = false;
+        systemInfo.mayHaveSplitSettles     = false;
     }
     else
     {
-        systemInfo.haveSplitConstraints = (gmx_mtop_ftype_count(mtop, F_CONSTR) > 0
-                                           || gmx_mtop_ftype_count(mtop, F_CONSTRNC) > 0);
-        systemInfo.haveSplitSettles     = (gmx_mtop_ftype_count(mtop, F_SETTLE) > 0);
+        systemInfo.mayHaveSplitConstraints = (gmx_mtop_ftype_count(mtop, F_CONSTR) > 0
+                                              || gmx_mtop_ftype_count(mtop, F_CONSTRNC) > 0);
+        systemInfo.mayHaveSplitSettles     = (gmx_mtop_ftype_count(mtop, F_SETTLE) > 0);
     }
 
     if (ir.rlist == 0)
@@ -2180,7 +2180,7 @@ static DDSystemInfo getSystemInfo(const gmx::MDLogger&           mdlog,
     }
 
     systemInfo.constraintCommunicationRange = 0;
-    if (systemInfo.haveSplitConstraints && options.constraintCommunicationRange <= 0)
+    if (systemInfo.mayHaveSplitConstraints && options.constraintCommunicationRange <= 0)
     {
         /* There is a cell size limit due to the constraints (P-LINCS) */
         systemInfo.constraintCommunicationRange = gmx::constr_r_max(mdlog, &mtop, &ir);
@@ -2225,7 +2225,7 @@ static void checkDDGridSetup(const DDGridSetup&   ddGridSetup,
 {
     if (options.numCells[XX] <= 0 && (ddGridSetup.numDomains[XX] == 0))
     {
-        const bool  bC = (systemInfo.haveSplitConstraints
+        const bool  bC = (systemInfo.mayHaveSplitConstraints
                          && systemInfo.constraintCommunicationRange > systemInfo.minCutoffForMultiBody);
         std::string message =
                 gmx::formatString("Change the number of ranks or mdrun option %s%s%s",
@@ -2577,7 +2577,7 @@ static void writeSettings(gmx::TextWriter*   log,
             (countInterUpdategroupVsites(mtop, comm->systemInfo.updateGroupingsPerMoleculeType) != 0);
 
     if (comm->systemInfo.haveInterDomainBondeds || haveInterDomainVsites
-        || comm->systemInfo.haveSplitConstraints || comm->systemInfo.haveSplitSettles)
+        || comm->systemInfo.mayHaveSplitConstraints || comm->systemInfo.mayHaveSplitSettles)
     {
         std::string decompUnits;
         if (comm->systemInfo.useUpdateGroups)
@@ -2631,7 +2631,7 @@ static void writeSettings(gmx::TextWriter*   log,
         {
             log->writeLineFormatted("%40s  %-7s %6.3f nm", "virtual site constructions", "(-rcon)", limit);
         }
-        if (comm->systemInfo.haveSplitConstraints || comm->systemInfo.haveSplitSettles)
+        if (comm->systemInfo.mayHaveSplitConstraints || comm->systemInfo.mayHaveSplitSettles)
         {
             std::string separation =
                     gmx::formatString("atoms separated by up to %d constraints", 1 + ir.nProjOrder);

diff --git a/src/gromacs/domdec/domdec.h b/src/gromacs/domdec/domdec.h
index c866788c32bde9a4bd18db52375498bdb82ecd9b..b969eab4767e17e93bdca3307e7f62206fd18d36 100644 (file)
--- a/src/gromacs/domdec/domdec.h
+++ b/src/gromacs/domdec/domdec.h
@@ -151,8 +151,8 @@ int dd_pme_maxshift_x(const gmx_domdec_t& dd);
 /*! \brief Returns the maximum shift for coordinate communication in PME, dim y */
 int dd_pme_maxshift_y(const gmx_domdec_t& dd);
 
-/*! \brief Return whether constraints, not including settles, cross domain boundaries */
-bool ddHaveSplitConstraints(const gmx_domdec_t& dd);
+/*! \brief Return whether constraints, not including settles, may cross domain boundaries */
+bool ddMayHaveSplitConstraints(const gmx_domdec_t& dd);
 
 /*! \brief Return whether update groups are used */
 bool ddUsesUpdateGroups(const gmx_domdec_t& dd);

diff --git a/src/gromacs/domdec/domdec_constraints.cpp b/src/gromacs/domdec/domdec_constraints.cpp
index 1ba64f5c9f6f1b0451af6a561ef10f11ca33a058..cd6fec86c7c0f867b5029c052ead5db9b378423d 100644 (file)
--- a/src/gromacs/domdec/domdec_constraints.cpp
+++ b/src/gromacs/domdec/domdec_constraints.cpp
@@ -407,7 +407,7 @@ int dd_make_local_constraints(gmx_domdec_t*                  dd,
     // This code should not be called unless this condition is true,
     // because that's the only time init_domdec_constraints is
     // called...
-    GMX_RELEASE_ASSERT(dd->comm->systemInfo.haveSplitConstraints || dd->comm->systemInfo.haveSplitSettles,
+    GMX_RELEASE_ASSERT(dd->comm->systemInfo.mayHaveSplitConstraints || dd->comm->systemInfo.mayHaveSplitSettles,
                        "dd_make_local_constraints called when there are no local constraints");
     // ... and init_domdec_constraints always sets
     // dd->constraint_comm...
@@ -439,7 +439,7 @@ int dd_make_local_constraints(gmx_domdec_t*                  dd,
 
     gmx::ArrayRef<const std::vector<int>> at2settle_mt;
     /* When settle works inside charge groups, we assigned them already */
-    if (dd->comm->systemInfo.haveSplitSettles)
+    if (dd->comm->systemInfo.mayHaveSplitSettles)
     {
         // TODO Perhaps gmx_domdec_constraints_t should keep a valid constr?
         GMX_RELEASE_ASSERT(constr != nullptr, "Must have valid constraints object");

diff --git a/src/gromacs/domdec/domdec_internal.h b/src/gromacs/domdec/domdec_internal.h
index 019fb362073da25b9e7eda427d3cfad1343bcd97..ddfa6b847ffa2b027bcf03c0207639774599efb1 100644 (file)
--- a/src/gromacs/domdec/domdec_internal.h
+++ b/src/gromacs/domdec/domdec_internal.h
@@ -442,9 +442,9 @@ struct DDSystemInfo
     real cellsizeLimit = 0;
 
     //! Can atoms connected by constraints be assigned to different domains?
-    bool haveSplitConstraints = false;
+    bool mayHaveSplitConstraints = false;
     //! Can atoms connected by settles be assigned to different domains?
-    bool haveSplitSettles = false;
+    bool mayHaveSplitSettles = false;
     //! Estimated communication range needed for constraints
     real constraintCommunicationRange = 0;
 

diff --git a/src/gromacs/domdec/domdec_topology.cpp b/src/gromacs/domdec/domdec_topology.cpp
index 9c85eb1a07b97991a93b016744649cda9e1c383a..5bacabd878df49d2723d1e25a6e98a1599272214 100644 (file)
--- a/src/gromacs/domdec/domdec_topology.cpp
+++ b/src/gromacs/domdec/domdec_topology.cpp
@@ -800,8 +800,8 @@ void dd_make_reverse_top(FILE*                           fplog,
                                || ddBondedChecking == DDBondedChecking::All,
                        "Invalid enum value for mdrun -ddcheck");
     const ReverseTopOptions rtOptions(ddBondedChecking,
-                                      !dd->comm->systemInfo.haveSplitConstraints,
-                                      !dd->comm->systemInfo.haveSplitSettles);
+                                      !dd->comm->systemInfo.mayHaveSplitConstraints,
+                                      !dd->comm->systemInfo.mayHaveSplitSettles);
 
     dd->reverse_top = std::make_unique<gmx_reverse_top_t>(
             mtop, inputrec.efep != FreeEnergyPerturbationType::No, rtOptions);
@@ -832,7 +832,7 @@ void dd_make_reverse_top(FILE*                           fplog,
         init_domdec_vsites(dd, numInterUpdategroupVirtualSites);
     }
 
-    if (dd->comm->systemInfo.haveSplitConstraints || dd->comm->systemInfo.haveSplitSettles)
+    if (dd->comm->systemInfo.mayHaveSplitConstraints || dd->comm->systemInfo.mayHaveSplitSettles)
     {
         init_domdec_constraints(dd, mtop);
     }

diff --git a/src/gromacs/domdec/partition.cpp b/src/gromacs/domdec/partition.cpp
index 20513be8a4ac6e06338cb9057b0456e5fc4aef96..6bd203df67f4adadb3428f6a02c3bfc0d78278dd 100644 (file)
--- a/src/gromacs/domdec/partition.cpp
+++ b/src/gromacs/domdec/partition.cpp
@@ -3222,7 +3222,7 @@ void dd_partition_system(FILE*                     fplog,
                 }
                 break;
             case DDAtomRanges::Type::Constraints:
-                if (dd->comm->systemInfo.haveSplitConstraints || dd->comm->systemInfo.haveSplitSettles)
+                if (dd->comm->systemInfo.mayHaveSplitConstraints || dd->comm->systemInfo.mayHaveSplitSettles)
                 {
                     /* Only for inter-cg constraints we need special code */
                     n = dd_make_local_constraints(dd,

diff --git a/src/gromacs/mdlib/constr.cpp b/src/gromacs/mdlib/constr.cpp
index ecaaa2c85f5a0b1e92e2c8f8bd819e71afc5407d..ef7be1b06d0a5b873afb67f63a54d8b713b13a1d 100644 (file)
--- a/src/gromacs/mdlib/constr.cpp
+++ b/src/gromacs/mdlib/constr.cpp
@@ -104,6 +104,7 @@ public:
          pull_t*               pull_work,
          FILE*                 log_p,
          const t_commrec*      cr_p,
+         bool                  useUpdateGroups,
          const gmx_multisim_t* ms,
          t_nrnb*               nrnb,
          gmx_wallcycle*        wcycle_p,
@@ -1080,13 +1081,14 @@ Constraints::Constraints(const gmx_mtop_t&     mtop,
                          pull_t*               pull_work,
                          FILE*                 log,
                          const t_commrec*      cr,
+                         const bool            useUpdateGroups,
                          const gmx_multisim_t* ms,
                          t_nrnb*               nrnb,
                          gmx_wallcycle*        wcycle,
                          bool                  pbcHandlingRequired,
                          int                   numConstraints,
                          int                   numSettles) :
-    impl_(new Impl(mtop, ir, pull_work, log, cr, ms, nrnb, wcycle, pbcHandlingRequired, numConstraints, numSettles))
+    impl_(new Impl(mtop, ir, pull_work, log, cr, useUpdateGroups, ms, nrnb, wcycle, pbcHandlingRequired, numConstraints, numSettles))
 {
 }
 
@@ -1095,6 +1097,7 @@ Constraints::Impl::Impl(const gmx_mtop_t&     mtop_p,
                         pull_t*               pull_work,
                         FILE*                 log_p,
                         const t_commrec*      cr_p,
+                        const bool            mayHaveSplitConstraints,
                         const gmx_multisim_t* ms_p,
                         t_nrnb*               nrnb_p,
                         gmx_wallcycle*        wcycle_p,
@@ -1151,23 +1154,20 @@ Constraints::Impl::Impl(const gmx_mtop_t&     mtop_p,
             }
         }
 
+        // When there are multiple PP domains and update groups are not in use,
+        // the constraints might be split across them.
         if (ir.eConstrAlg == ConstraintAlgorithm::Lincs)
         {
-            lincsd = init_lincs(log,
-                                mtop,
-                                nflexcon,
-                                at2con_mt,
-                                DOMAINDECOMP(cr) && ddHaveSplitConstraints(*cr->dd),
-                                ir.nLincsIter,
-                                ir.nProjOrder);
+            lincsd = init_lincs(
+                    log, mtop, nflexcon, at2con_mt, mayHaveSplitConstraints, ir.nLincsIter, ir.nProjOrder);
         }
 
         if (ir.eConstrAlg == ConstraintAlgorithm::Shake)
         {
-            if (DOMAINDECOMP(cr) && ddHaveSplitConstraints(*cr->dd))
+            if (mayHaveSplitConstraints)
             {
                 gmx_fatal(FARGS,
-                          "SHAKE is not supported with domain decomposition and constraint that "
+                          "SHAKE is not supported with domain decomposition and constraints that "
                           "cross domain boundaries, use LINCS");
             }
             if (nflexcon)

diff --git a/src/gromacs/mdlib/constr.h b/src/gromacs/mdlib/constr.h
index 03eb4772ac0b2f581117f0a01b3a0d9d0e24999e..f5bc84ea2ce427aa9f4f51ef734cf0e69d77ac3d 100644 (file)
--- a/src/gromacs/mdlib/constr.h
+++ b/src/gromacs/mdlib/constr.h
@@ -107,6 +107,7 @@ private:
                 pull_t*               pull_work,
                 FILE*                 log,
                 const t_commrec*      cr,
+                bool                  haveHaveSplitConstraints,
                 const gmx_multisim_t* ms,
                 t_nrnb*               nrnb,
                 gmx_wallcycle*        wcycle,

diff --git a/src/gromacs/mdlib/tests/energyoutput.cpp b/src/gromacs/mdlib/tests/energyoutput.cpp
index c60df687e74c9b26003ff86c116bd9834954ae6f..51a06b9f769d93837b2399ba68fab5371452be6f 100644 (file)
--- a/src/gromacs/mdlib/tests/energyoutput.cpp
+++ b/src/gromacs/mdlib/tests/energyoutput.cpp
@@ -371,7 +371,7 @@ public:
         // TODO This object will always return zero as RMSD value.
         //      It is more relevant to have non-zero value for testing.
         constraints_ = makeConstraints(
-                mtop_, inputrec_, nullptr, false, nullptr, &cr_, nullptr, nullptr, nullptr, false);
+                mtop_, inputrec_, nullptr, false, nullptr, &cr_, false, nullptr, nullptr, nullptr, false);
     }
 
     /*! \brief Helper function to generate synthetic data to output

diff --git a/src/gromacs/mdlib/vsite.cpp b/src/gromacs/mdlib/vsite.cpp
index 8621154d7f7d60b154ce52f9d20549516b251a16..18aa3fd88e6fe447221cc6700b8b1de26c7d4b5b 100644 (file)
--- a/src/gromacs/mdlib/vsite.cpp
+++ b/src/gromacs/mdlib/vsite.cpp
@@ -232,7 +232,10 @@ class VirtualSitesHandler::Impl
 {
 public:
     //! Constructor, domdec should be nullptr without DD
-    Impl(const gmx_mtop_t& mtop, gmx_domdec_t* domdec, PbcType pbcType);
+    Impl(const gmx_mtop_t&                 mtop,
+         gmx_domdec_t*                     domdec,
+         PbcType                           pbcType,
+         ArrayRef<const RangePartitioning> updateGroupingPerMoleculeType);
 
     //! Returns the number of virtual sites acting over multiple update groups
     int numInterUpdategroupVirtualSites() const { return numInterUpdategroupVirtualSites_; }
@@ -2574,7 +2577,8 @@ int countInterUpdategroupVsites(const gmx_mtop_t&                           mtop
 
 std::unique_ptr<VirtualSitesHandler> makeVirtualSitesHandler(const gmx_mtop_t& mtop,
                                                              const t_commrec*  cr,
-                                                             PbcType           pbcType)
+                                                             PbcType           pbcType,
+                                                             ArrayRef<const RangePartitioning> updateGroupingPerMoleculeType)
 {
     GMX_RELEASE_ASSERT(cr != nullptr, "We need a valid commrec");
 
@@ -2603,7 +2607,7 @@ std::unique_ptr<VirtualSitesHandler> makeVirtualSitesHandler(const gmx_mtop_t& m
         return vsite;
     }
 
-    return std::make_unique<VirtualSitesHandler>(mtop, cr->dd, pbcType);
+    return std::make_unique<VirtualSitesHandler>(mtop, cr->dd, pbcType, updateGroupingPerMoleculeType);
 }
 
 ThreadingInfo::ThreadingInfo() : numThreads_(gmx_omp_nthreads_get(ModuleMultiThread::VirtualSite))
@@ -2632,27 +2636,21 @@ ThreadingInfo::ThreadingInfo() : numThreads_(gmx_omp_nthreads_get(ModuleMultiThr
     }
 }
 
-//! Returns the number of inter update-group vsites
-static int getNumInterUpdategroupVsites(const gmx_mtop_t& mtop, const gmx_domdec_t* domdec)
-{
-    gmx::ArrayRef<const gmx::RangePartitioning> updateGroupingsPerMoleculeType;
-    if (domdec)
-    {
-        updateGroupingsPerMoleculeType = getUpdateGroupingsPerMoleculeType(*domdec);
-    }
-
-    return countInterUpdategroupVsites(mtop, updateGroupingsPerMoleculeType);
-}
-
-VirtualSitesHandler::Impl::Impl(const gmx_mtop_t& mtop, gmx_domdec_t* domdec, const PbcType pbcType) :
-    numInterUpdategroupVirtualSites_(getNumInterUpdategroupVsites(mtop, domdec)),
+VirtualSitesHandler::Impl::Impl(const gmx_mtop_t&                       mtop,
+                                gmx_domdec_t*                           domdec,
+                                const PbcType                           pbcType,
+                                const ArrayRef<const RangePartitioning> updateGroupingPerMoleculeType) :
+    numInterUpdategroupVirtualSites_(countInterUpdategroupVsites(mtop, updateGroupingPerMoleculeType)),
     domainInfo_({ pbcType, pbcType != PbcType::No && numInterUpdategroupVirtualSites_ > 0, domdec }),
     iparams_(mtop.ffparams.iparams)
 {
 }
 
-VirtualSitesHandler::VirtualSitesHandler(const gmx_mtop_t& mtop, gmx_domdec_t* domdec, const PbcType pbcType) :
-    impl_(new Impl(mtop, domdec, pbcType))
+VirtualSitesHandler::VirtualSitesHandler(const gmx_mtop_t&                       mtop,
+                                         gmx_domdec_t*                           domdec,
+                                         const PbcType                           pbcType,
+                                         const ArrayRef<const RangePartitioning> updateGroupingPerMoleculeType) :
+    impl_(new Impl(mtop, domdec, pbcType, updateGroupingPerMoleculeType))
 {
 }
 

diff --git a/src/gromacs/mdlib/vsite.h b/src/gromacs/mdlib/vsite.h
index f6a4807b7fc2412b14a108f424c691ab82d384d2..a0783f875c33dcf6fb8dfc0edc5396dc53be22f8 100644 (file)
--- a/src/gromacs/mdlib/vsite.h
+++ b/src/gromacs/mdlib/vsite.h
@@ -96,7 +96,10 @@ class VirtualSitesHandler
 {
 public:
     //! Constructor, used only be the makeVirtualSitesHandler() factory function
-    VirtualSitesHandler(const gmx_mtop_t& mtop, gmx_domdec_t* domdec, PbcType pbcType);
+    VirtualSitesHandler(const gmx_mtop_t&                 mtop,
+                        gmx_domdec_t*                     domdec,
+                        PbcType                           pbcType,
+                        ArrayRef<const RangePartitioning> updateGroupingPerMoleculeType);
 
     ~VirtualSitesHandler();
 
@@ -187,14 +190,18 @@ int countInterUpdategroupVsites(const gmx_mtop_t&                 mtop,
 
 /*! \brief Create the virtual site handler
  *
- * \param[in] mtop      The global topology
- * \param[in] cr        The communication record
- * \param[in] pbcType   The type of PBC
+ * \param[in] mtop                           The global topology
+ * \param[in] cr                             The communication record
+ * \param[in] pbcType                        The type of PBC
+ * \param[in] updateGroupingPerMoleculeType  Update grouping per molecule type, pass
+ *                                           empty when not using update groups
  * \returns A valid vsite handler object or nullptr when there are no virtual sites
  */
-std::unique_ptr<VirtualSitesHandler> makeVirtualSitesHandler(const gmx_mtop_t& mtop,
-                                                             const t_commrec*  cr,
-                                                             PbcType           pbcType);
+std::unique_ptr<VirtualSitesHandler>
+makeVirtualSitesHandler(const gmx_mtop_t&                 mtop,
+                        const t_commrec*                  cr,
+                        PbcType                           pbcType,
+                        ArrayRef<const RangePartitioning> updateGroupingPerMoleculeType);
 
 } // namespace gmx
 

diff --git a/src/gromacs/mdrun/runner.cpp b/src/gromacs/mdrun/runner.cpp
index e7874c8d0c0f5ce275f171de08ae184beec75ecd..e4c8c1ea787da1e3af5a07fe0aa8e904ce6d37c2 100644 (file)
--- a/src/gromacs/mdrun/runner.cpp
+++ b/src/gromacs/mdrun/runner.cpp
@@ -1725,7 +1725,12 @@ int Mdrunner::mdrunner()
         }
 
         /* Initialize the virtual site communication */
-        vsite = makeVirtualSitesHandler(mtop, cr, fr->pbcType);
+        gmx::ArrayRef<const gmx::RangePartitioning> updateGroupingsPerMoleculeType;
+        if (DOMAINDECOMP(cr))
+        {
+            updateGroupingsPerMoleculeType = getUpdateGroupingsPerMoleculeType(*cr->dd);
+        }
+        vsite = makeVirtualSitesHandler(mtop, cr, fr->pbcType, updateGroupingsPerMoleculeType);
 
         calc_shifts(box, fr->shift_vec);
 
@@ -1922,8 +1927,17 @@ int Mdrunner::mdrunner()
         }
 
         /* Let makeConstraints know whether we have essential dynamics constraints. */
-        auto constr = makeConstraints(
-                mtop, *inputrec, pull_work, doEssentialDynamics, fplog, cr, ms, &nrnb, wcycle.get(), fr->bMolPBC);
+        auto constr = makeConstraints(mtop,
+                                      *inputrec,
+                                      pull_work,
+                                      doEssentialDynamics,
+                                      fplog,
+                                      cr,
+                                      DOMAINDECOMP(cr) && ddMayHaveSplitConstraints(*cr->dd),
+                                      ms,
+                                      &nrnb,
+                                      wcycle.get(),
+                                      fr->bMolPBC);
 
         /* Energy terms and groups */
         gmx_enerdata_t enerd(mtop.groups.groups[SimulationAtomGroupType::EnergyOutput].size(),