Remove charge groups from domdec and localtop
[alexxy/gromacs.git] / src / gromacs / domdec / domdec.cpp
index d46ee6ce00a1ed41325f6f217dbce733bb558239..4ba6cc34846d2cf8ea4dd16d9ce07fe844799759 100644 (file)
@@ -320,8 +320,6 @@ void dd_move_x(gmx_domdec_t             *dd,
 
     comm = dd->comm;
 
-    const gmx::RangePartitioning &atomGrouping = dd->atomGrouping();
-
     nzone   = 1;
     nat_tot = comm->atomRanges.numHomeAtoms();
     for (int d = 0; d < dd->ndim; d++)
@@ -340,46 +338,37 @@ void dd_move_x(gmx_domdec_t             *dd,
             int                        n          = 0;
             if (!bPBC)
             {
-                for (int g : ind.index)
+                for (int j : ind.index)
                 {
-                    for (int j : atomGrouping.block(g))
-                    {
-                        sendBuffer[n] = x[j];
-                        n++;
-                    }
+                    sendBuffer[n] = x[j];
+                    n++;
                 }
             }
             else if (!bScrew)
             {
-                for (int g : ind.index)
+                for (int j : ind.index)
                 {
-                    for (int j : atomGrouping.block(g))
+                    /* We need to shift the coordinates */
+                    for (int d = 0; d < DIM; d++)
                     {
-                        /* We need to shift the coordinates */
-                        for (int d = 0; d < DIM; d++)
-                        {
-                            sendBuffer[n][d] = x[j][d] + shift[d];
-                        }
-                        n++;
+                        sendBuffer[n][d] = x[j][d] + shift[d];
                     }
+                    n++;
                 }
             }
             else
             {
-                for (int g : ind.index)
+                for (int j : ind.index)
                 {
-                    for (int j : atomGrouping.block(g))
-                    {
-                        /* Shift x */
-                        sendBuffer[n][XX] = x[j][XX] + shift[XX];
-                        /* Rotate y and z.
-                         * This operation requires a special shift force
-                         * treatment, which is performed in calc_vir.
-                         */
-                        sendBuffer[n][YY] = box[YY][YY] - x[j][YY];
-                        sendBuffer[n][ZZ] = box[ZZ][ZZ] - x[j][ZZ];
-                        n++;
-                    }
+                    /* Shift x */
+                    sendBuffer[n][XX] = x[j][XX] + shift[XX];
+                    /* Rotate y and z.
+                     * This operation requires a special shift force
+                     * treatment, which is performed in calc_vir.
+                     */
+                    sendBuffer[n][YY] = box[YY][YY] - x[j][YY];
+                    sendBuffer[n][ZZ] = box[ZZ][ZZ] - x[j][ZZ];
+                    n++;
                 }
             }
 
@@ -433,8 +422,6 @@ void dd_move_f(gmx_domdec_t             *dd,
 
     comm = dd->comm;
 
-    const gmx::RangePartitioning &atomGrouping = dd->atomGrouping();
-
     nzone   = comm->zones.n/2;
     nat_tot = comm->atomRanges.end(DDAtomRanges::Type::Zones);
     for (int d = dd->ndim-1; d >= 0; d--)
@@ -487,16 +474,13 @@ void dd_move_f(gmx_domdec_t             *dd,
             int n = 0;
             if (!bShiftForcesNeedPbc)
             {
-                for (int g : ind.index)
+                for (int j : ind.index)
                 {
-                    for (int j : atomGrouping.block(g))
+                    for (int d = 0; d < DIM; d++)
                     {
-                        for (int d = 0; d < DIM; d++)
-                        {
-                            f[j][d] += receiveBuffer[n][d];
-                        }
-                        n++;
+                        f[j][d] += receiveBuffer[n][d];
                     }
+                    n++;
                 }
             }
             else if (!bScrew)
@@ -504,43 +488,37 @@ void dd_move_f(gmx_domdec_t             *dd,
                 /* fshift should always be defined if this function is
                  * called when bShiftForcesNeedPbc is true */
                 assert(nullptr != fshift);
-                for (int g : ind.index)
+                for (int j : ind.index)
                 {
-                    for (int j : atomGrouping.block(g))
+                    for (int d = 0; d < DIM; d++)
                     {
-                        for (int d = 0; d < DIM; d++)
-                        {
-                            f[j][d] += receiveBuffer[n][d];
-                        }
-                        /* Add this force to the shift force */
-                        for (int d = 0; d < DIM; d++)
-                        {
-                            fshift[is][d] += receiveBuffer[n][d];
-                        }
-                        n++;
+                        f[j][d] += receiveBuffer[n][d];
+                    }
+                    /* Add this force to the shift force */
+                    for (int d = 0; d < DIM; d++)
+                    {
+                        fshift[is][d] += receiveBuffer[n][d];
                     }
+                    n++;
                 }
             }
             else
             {
-                for (int g : ind.index)
+                for (int j : ind.index)
                 {
-                    for (int j : atomGrouping.block(g))
+                    /* Rotate the force */
+                    f[j][XX] += receiveBuffer[n][XX];
+                    f[j][YY] -= receiveBuffer[n][YY];
+                    f[j][ZZ] -= receiveBuffer[n][ZZ];
+                    if (fshift)
                     {
-                        /* Rotate the force */
-                        f[j][XX] += receiveBuffer[n][XX];
-                        f[j][YY] -= receiveBuffer[n][YY];
-                        f[j][ZZ] -= receiveBuffer[n][ZZ];
-                        if (fshift)
+                        /* Add this force to the shift force */
+                        for (int d = 0; d < DIM; d++)
                         {
-                            /* Add this force to the shift force */
-                            for (int d = 0; d < DIM; d++)
-                            {
-                                fshift[is][d] += receiveBuffer[n][d];
-                            }
+                            fshift[is][d] += receiveBuffer[n][d];
                         }
-                        n++;
                     }
+                    n++;
                 }
             }
         }
@@ -571,8 +549,6 @@ void dd_atom_spread_real(gmx_domdec_t *dd, real v[])
 
     comm = dd->comm;
 
-    const gmx::RangePartitioning &atomGrouping = dd->atomGrouping();
-
     nzone   = 1;
     nat_tot = comm->atomRanges.numHomeAtoms();
     for (int d = 0; d < dd->ndim; d++)
@@ -587,12 +563,9 @@ void dd_atom_spread_real(gmx_domdec_t *dd, real v[])
             DDBufferAccess<gmx::RVec> sendBufferAccess(comm->rvecBuffer, ind.nsend[nzone + 1]);
             gmx::ArrayRef<real>       sendBuffer = realArrayRefFromRvecArrayRef(sendBufferAccess.buffer);
             int                       n          = 0;
-            for (int g : ind.index)
+            for (int j : ind.index)
             {
-                for (int j : atomGrouping.block(g))
-                {
-                    sendBuffer[n++] = v[j];
-                }
+                sendBuffer[n++] = v[j];
             }
 
             DDBufferAccess<gmx::RVec> receiveBufferAccess(comm->rvecBuffer2, cd->receiveInPlace ? 0 : ind.nrecv[nzone + 1]);
@@ -634,8 +607,6 @@ void dd_atom_sum_real(gmx_domdec_t *dd, real v[])
 
     comm = dd->comm;
 
-    const gmx::RangePartitioning &atomGrouping = dd->atomGrouping();
-
     nzone   = comm->zones.n/2;
     nat_tot = comm->atomRanges.end(DDAtomRanges::Type::Zones);
     for (int d = dd->ndim-1; d >= 0; d--)
@@ -677,13 +648,10 @@ void dd_atom_sum_real(gmx_domdec_t *dd, real v[])
                        sendBuffer, receiveBuffer);
             /* Add the received forces */
             int n = 0;
-            for (int g : ind.index)
+            for (int j : ind.index)
             {
-                for (int j : atomGrouping.block(g))
-                {
-                    v[j] += receiveBuffer[n];
-                    n++;
-                }
+                v[j] += receiveBuffer[n];
+                n++;
             }
         }
         nzone /= 2;
@@ -699,7 +667,7 @@ real dd_cutoff_multibody(const gmx_domdec_t *dd)
     comm = dd->comm;
 
     r = -1;
-    if (comm->bInterCGBondeds)
+    if (comm->haveInterDomainMultiBodyBondeds)
     {
         if (comm->cutoff_mbody > 0)
         {
@@ -2187,8 +2155,6 @@ static void set_dd_limits_and_grid(const gmx::MDLogger &mdlog,
     /* Allocate the charge group/atom sorting struct */
     comm->sort = std::make_unique<gmx_domdec_sort_t>();
 
-    comm->bCGs = (ncg_mtop(mtop) < mtop->natoms);
-
     /* We need to decide on update groups early, as this affects communication distances */
     comm->useUpdateGroups = false;
     if (ir->cutoff_scheme == ecutsVERLET)
@@ -2197,16 +2163,10 @@ static void set_dd_limits_and_grid(const gmx::MDLogger &mdlog,
         setupUpdateGroups(mdlog, *mtop, *ir, cutoffMargin, cr->nnodes, comm);
     }
 
-    comm->bInterCGBondeds = ((ncg_mtop(mtop) > gmx_mtop_num_molecules(*mtop)) ||
-                             mtop->bIntermolecularInteractions);
-    if (comm->bInterCGBondeds)
-    {
-        comm->bInterCGMultiBody = (multi_body_bondeds_count(mtop) > 0);
-    }
-    else
-    {
-        comm->bInterCGMultiBody = FALSE;
-    }
+    // TODO: Check whether all bondeds are within update groups
+    comm->haveInterDomainBondeds          = (mtop->natoms > gmx_mtop_num_molecules(*mtop) ||
+                                             mtop->bIntermolecularInteractions);
+    comm->haveInterDomainMultiBodyBondeds = (multi_body_bondeds_count(mtop) > 0);
 
     if (comm->useUpdateGroups)
     {
@@ -2265,7 +2225,7 @@ static void set_dd_limits_and_grid(const gmx::MDLogger &mdlog,
      *       Note that we would need to improve the pairlist buffer case.
      */
 
-    if (comm->bInterCGBondeds)
+    if (comm->haveInterDomainBondeds)
     {
         if (options.minimumCommunicationRange > 0)
         {
@@ -2401,7 +2361,7 @@ static void set_dd_limits_and_grid(const gmx::MDLogger &mdlog,
                            !isDlbDisabled(comm),
                            options.dlbScaling,
                            comm->cellsize_limit, comm->cutoff,
-                           comm->bInterCGBondeds);
+                           comm->haveInterDomainBondeds);
 
         if (dd->nc[XX] == 0)
         {
@@ -2500,7 +2460,7 @@ static void set_dd_limits_and_grid(const gmx::MDLogger &mdlog,
         comm->slb_frac[ZZ] = get_slb_frac(mdlog, "z", dd->nc[ZZ], options.cellSizeZ);
     }
 
-    if (comm->bInterCGBondeds && comm->cutoff_mbody == 0)
+    if (comm->haveInterDomainBondeds && comm->cutoff_mbody == 0)
     {
         if (comm->bBondComm || !isDlbDisabled(comm))
         {
@@ -2669,16 +2629,12 @@ static void writeSettings(gmx::TextWriter       *log,
     const bool haveInterDomainVsites =
         (countInterUpdategroupVsites(*mtop, comm->updateGroupingPerMoleculetype) != 0);
 
-    if (comm->bInterCGBondeds ||
+    if (comm->haveInterDomainBondeds ||
         haveInterDomainVsites ||
         dd->splitConstraints || dd->splitSettles)
     {
         std::string decompUnits;
-        if (comm->bCGs)
-        {
-            decompUnits = "charge groups";
-        }
-        else if (comm->useUpdateGroups)
+        if (comm->useUpdateGroups)
         {
             decompUnits = "atom groups";
         }
@@ -2707,7 +2663,7 @@ static void writeSettings(gmx::TextWriter       *log,
             }
         }
 
-        if (comm->bInterCGBondeds)
+        if (comm->haveInterDomainBondeds)
         {
             log->writeLineFormatted("%40s  %-7s %6.3f nm",
                                     "two-body bonded interactions", "(-rdd)",
@@ -2863,7 +2819,7 @@ gmx_bool dd_bonded_molpbc(const gmx_domdec_t *dd, int ePBC)
      * or we use domain decomposition for each periodic dimension,
      * we do not need to take pbc into account for the bonded interactions.
      */
-    return (ePBC != epbcNONE && dd->comm->bInterCGBondeds &&
+    return (ePBC != epbcNONE && dd->comm->haveInterDomainBondeds &&
             !(dd->nc[XX] > 1 &&
               dd->nc[YY] > 1 &&
               (dd->nc[ZZ] > 1 || ePBC == epbcXY)));