Rename DD constraint booleans
authorBerk Hess <hess@kth.se>
Tue, 18 Sep 2018 12:46:48 +0000 (14:46 +0200)
committerMark Abraham <mark.j.abraham@gmail.com>
Sat, 22 Sep 2018 12:29:39 +0000 (14:29 +0200)
This is preparation for adding update groups

Change-Id: I135abe196569c6e1a7dde535c8e616e1a479bcac

src/gromacs/domdec/domdec.cpp
src/gromacs/domdec/domdec_constraints.cpp
src/gromacs/domdec/domdec_struct.h
src/gromacs/domdec/domdec_topology.cpp
src/gromacs/domdec/partition.cpp
src/gromacs/mdlib/constr.cpp

index 099b2c07b6b2d8e34eca3a87567b1d265ea297a8..5f5aef8697c453a513cf9d1e76293aea298ce128 100644 (file)
@@ -2111,8 +2111,8 @@ static void set_dd_limits_and_grid(const gmx::MDLogger &mdlog,
         comm->bInterCGMultiBody = FALSE;
     }
 
-    dd->bInterCGcons    = gmx::inter_charge_group_constraints(*mtop);
-    dd->bInterCGsettles = gmx::inter_charge_group_settles(*mtop);
+    dd->splitConstraints = gmx::inter_charge_group_constraints(*mtop);
+    dd->splitSettles     = gmx::inter_charge_group_settles(*mtop);
 
     if (ir->rlist == 0)
     {
@@ -2229,7 +2229,7 @@ static void set_dd_limits_and_grid(const gmx::MDLogger &mdlog,
     }
 
     real rconstr = 0;
-    if (dd->bInterCGcons && options.constraintCommunicationRange <= 0)
+    if (dd->splitConstraints && options.constraintCommunicationRange <= 0)
     {
         /* There is a cell size limit due to the constraints (P-LINCS) */
         rconstr = gmx::constr_r_max(mdlog, mtop, ir);
@@ -2243,7 +2243,7 @@ static void set_dd_limits_and_grid(const gmx::MDLogger &mdlog,
     }
     else if (options.constraintCommunicationRange > 0)
     {
-        /* Here we do not check for dd->bInterCGcons,
+        /* Here we do not check for dd->splitConstraints.
          * because one can also set a cell size limit for virtual sites only
          * and at this point we don't know yet if there are intercg v-sites.
          */
@@ -2299,7 +2299,7 @@ static void set_dd_limits_and_grid(const gmx::MDLogger &mdlog,
         if (dd->nc[XX] == 0)
         {
             char     buf[STRLEN];
-            gmx_bool bC = (dd->bInterCGcons && rconstr > r_bonded_limit);
+            gmx_bool bC = (dd->splitConstraints && rconstr > r_bonded_limit);
             sprintf(buf, "Change the number of ranks or mdrun option %s%s%s",
                     !bC ? "-rdd" : "-rcon",
                     comm->dlbState != DlbState::offUser ? " or -dds" : "",
@@ -2560,7 +2560,7 @@ static void writeSettings(gmx::TextWriter       *log,
 
     if (comm->bInterCGBondeds ||
         bInterCGVsites ||
-        dd->bInterCGcons || dd->bInterCGsettles)
+        dd->splitConstraints || dd->splitSettles)
     {
         log->writeLine("The maximum allowed distance for charge groups involved in interactions is:");
         log->writeLineFormatted("%40s  %-7s %6.3f nm", "non-bonded interactions", "", comm->cutoff);
@@ -2596,7 +2596,7 @@ static void writeSettings(gmx::TextWriter       *log,
             log->writeLineFormatted("%40s  %-7s %6.3f nm",
                                     "virtual site constructions", "(-rcon)", limit);
         }
-        if (dd->bInterCGcons || dd->bInterCGsettles)
+        if (dd->splitConstraints || dd->splitSettles)
         {
             std::string separation = gmx::formatString("atoms separated by up to %d constraints",
                                                        1+ir->nProjOrder);
index 4e2d634dbf024ae8335fda6d6d063de40d0e98f5..f9c8364336d56f8a0b9c6d684a47c0868cb94c8f 100644 (file)
@@ -428,12 +428,12 @@ int dd_make_local_constraints(gmx_domdec_t *dd, int at_start,
     // 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->bInterCGcons || dd->bInterCGsettles, "dd_make_local_constraints called when there are no local constraints");
+    GMX_RELEASE_ASSERT(dd->splitConstraints || dd->splitSettles, "dd_make_local_constraints called when there are no local constraints");
     // ... and init_domdec_constraints always sets
     // dd->constraint_comm...
     GMX_RELEASE_ASSERT(dd->constraint_comm, "Invalid use of dd_make_local_constraints before construction of constraint_comm");
     // ... which static analysis needs to be reassured about, because
-    // otherwise, when dd->bInterCGsettles is
+    // otherwise, when dd->splitSettles is
     // true. dd->constraint_comm is unilaterally dereferenced before
     // the call to atoms_to_settles.
 
@@ -461,7 +461,7 @@ int dd_make_local_constraints(gmx_domdec_t *dd, int at_start,
 
     gmx::ArrayRef < const std::vector < int>> at2settle_mt;
     /* When settle works inside charge groups, we assigned them already */
-    if (dd->bInterCGsettles)
+    if (dd->splitSettles)
     {
         // TODO Perhaps gmx_domdec_constraints_t should keep a valid constr?
         GMX_RELEASE_ASSERT(constr != nullptr, "Must have valid constraints object");
index 226beb930dddd8cd541189dbc84ea849ae1dcf38..85425ff3c626c42a4914806504ef59fd99258396 100644 (file)
@@ -169,9 +169,10 @@ struct gmx_domdec_t { //NOLINT(clang-analyzer-optin.performance.Padding)
     /* Only available on the master node */
     std::unique_ptr<AtomDistribution> ma;
 
-    /* Are there inter charge group constraints */
-    gmx_bool bInterCGcons;
-    gmx_bool bInterCGsettles;
+    /* Can atoms connected by constraints be assigned to different domains? */
+    bool splitConstraints;
+    /* Can atoms connected by settles be assigned to different domains? */
+    bool splitSettles;
 
     /* Global atom number to interaction list */
     gmx_reverse_top_t  *reverse_top;
index 04966c882f546e1f225805e7500370faca91f799..3e98a16f5261f07113f5cd7641126e5d6806de6c 100644 (file)
@@ -770,7 +770,7 @@ void dd_make_reverse_top(FILE *fplog,
     dd->reverse_top  = new gmx_reverse_top_t;
     *dd->reverse_top =
         make_reverse_top(mtop, ir->efep != efepNO, vsitePbcPerMoltype,
-                         !dd->bInterCGcons, !dd->bInterCGsettles,
+                         !dd->splitConstraints, !dd->splitSettles,
                          bBCheck, &dd->nbonded_global);
 
     gmx_reverse_top_t *rt = dd->reverse_top;
@@ -822,7 +822,7 @@ void dd_make_reverse_top(FILE *fplog,
         init_domdec_vsites(dd, vsite->n_intercg_vsite);
     }
 
-    if (dd->bInterCGcons || dd->bInterCGsettles)
+    if (dd->splitConstraints || dd->splitSettles)
     {
         init_domdec_constraints(dd, mtop);
     }
index f1ddb85993d19d10a1f759ce3cadd6c8a48075e1..31662ca546bb22dac4d083939e7d04d0068dd181 100644 (file)
@@ -3476,7 +3476,7 @@ void dd_partition_system(FILE                *fplog,
                 }
                 break;
             case DDAtomRanges::Type::Constraints:
-                if (dd->bInterCGcons || dd->bInterCGsettles)
+                if (dd->splitConstraints || dd->splitSettles)
                 {
                     /* Only for inter-cg constraints we need special code */
                     n = dd_make_local_constraints(dd, n, top_global, fr->cginfo,
index ca77a6ad273fa1f70b3a6b265bbb0f24d32122ae..f1ec807362582d86c24f4d07275cf225af602575 100644 (file)
@@ -1045,13 +1045,13 @@ Constraints::Impl::Impl(const gmx_mtop_t     &mtop_p,
         {
             lincsd = init_lincs(log, mtop,
                                 nflexcon, at2con_mt,
-                                DOMAINDECOMP(cr) && cr->dd->bInterCGcons,
+                                DOMAINDECOMP(cr) && cr->dd->splitConstraints,
                                 ir.nLincsIter, ir.nProjOrder);
         }
 
         if (ir.eConstrAlg == econtSHAKE)
         {
-            if (DOMAINDECOMP(cr) && cr->dd->bInterCGcons)
+            if (DOMAINDECOMP(cr) && cr->dd->splitConstraints)
             {
                 gmx_fatal(FARGS, "SHAKE is not supported with domain decomposition and constraint that cross charge group boundaries, use LINCS");
             }