Merge branch release-2021
[alexxy/gromacs.git] / src / gromacs / gmxlib / nonbonded / nb_free_energy.cpp
index 5c62efc57d4cd660dc1293eab9e59117579c209e..3c9b933758a5b4d5f77f0812f8cf01f867e69113 100644 (file)
@@ -375,6 +375,10 @@ static void nb_free_energy_kernel(const t_nblist* gmx_restrict nlist,
     real* gmx_restrict f      = &(forceWithShiftForces->force()[0][0]);
     real* gmx_restrict fshift = &(forceWithShiftForces->shiftForces()[0][0]);
 
+    const real rlistSquared = gmx::square(fr->rlist);
+
+    int numExcludedPairsBeyondRlist = 0;
+
     for (int n = 0; n < nri; n++)
     {
         int npair_within_cutoff = 0;
@@ -431,6 +435,11 @@ static void nb_free_energy_kernel(const t_nblist* gmx_restrict nlist,
             }
             npair_within_cutoff++;
 
+            if (rSq > rlistSquared)
+            {
+                numExcludedPairsBeyondRlist++;
+            }
+
             if (rSq > 0)
             {
                 /* Note that unlike in the nbnxn kernels, we do not need
@@ -722,7 +731,7 @@ static void nb_free_energy_kernel(const t_nblist* gmx_restrict nlist,
                 }
             }
 
-            if (vdwInteractionTypeIsEwald && r < rVdw)
+            if (vdwInteractionTypeIsEwald && (r < rVdw || !bPairIncluded))
             {
                 /* See comment in the preamble. When using LJ-Ewald interactions
                  * (unless we use a switch modifier) we subtract the reciprocal-space
@@ -839,6 +848,20 @@ static void nb_free_energy_kernel(const t_nblist* gmx_restrict nlist,
      */
 #pragma omp atomic
     inc_nrnb(nrnb, eNR_NBKERNEL_FREE_ENERGY, nlist->nri * 12 + nlist->jindex[nri] * 150);
+
+    if (numExcludedPairsBeyondRlist > 0)
+    {
+        gmx_fatal(FARGS,
+                  "There are %d perturbed non-bonded pair interactions beyond the pair-list cutoff "
+                  "of %g nm, which is not supported. This can happen because the system is "
+                  "unstable or because intra-molecular interactions at long distances are "
+                  "excluded. If the "
+                  "latter is the case, you can try to increase nstlist or rlist to avoid this."
+                  "The error is likely triggered by the use of couple-intramol=no "
+                  "and the maximal distance in the decoupled molecule exceeding rlist.",
+                  numExcludedPairsBeyondRlist,
+                  fr->rlist);
+    }
 }
 
 typedef void (*KernelFunction)(const t_nblist* gmx_restrict nlist,