From: Berk Hess <hess@kth.se>
Date: Tue, 23 Jun 2020 19:53:11 +0000 (+0200)
Subject: Change grompp disres cost to linear
X-Git-Url: http://biod.pnpi.spb.ru/gitweb/?a=commitdiff_plain;h=5c4a4519fb8694d1651b1dbb9b9ffe58420498b5;p=alexxy%2Fgromacs.git

Change grompp disres cost to linear

The cost of adding distance restraints (actually their parameters)
was quadratic in the number of restraints and has now been changed
to linear.

Fixes #3457
---

diff --git a/docs/release-notes/2021/major/performance.rst b/docs/release-notes/2021/major/performance.rst
index 0c895e7637..0931f3843f 100644
--- a/docs/release-notes/2021/major/performance.rst
+++ b/docs/release-notes/2021/major/performance.rst
@@ -12,3 +12,11 @@ Extend supported use-cases for GPU version of update and constraints
 
 GPU version of update and constraints can now be used for FEP, except mass and constraints
 free-energy perturbation.
+       
+Reduce time spent in grompp with large numbers of distance restraints
+"""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
+
+The time `gmx grompp` spent processing distance restraint has been
+changed from quadratic in the number of restraints to linear.
+       
+:issue:`3457`
diff --git a/src/gromacs/gmxpreprocess/convparm.cpp b/src/gromacs/gmxpreprocess/convparm.cpp
index c9efb733cb..225086c52f 100644
--- a/src/gromacs/gmxpreprocess/convparm.cpp
+++ b/src/gromacs/gmxpreprocess/convparm.cpp
@@ -451,7 +451,6 @@ static int enter_params(gmx_ffparams_t*           ffparams,
                         bool                      bAppend)
 {
     t_iparams newparam;
-    int       type;
     int       rc;
 
     if ((rc = assign_param(ftype, &newparam, forceparams, comb, reppow)) < 0)
@@ -462,21 +461,33 @@ static int enter_params(gmx_ffparams_t*           ffparams,
 
     if (!bAppend)
     {
-        for (type = start; (type < ffparams->numTypes()); type++)
+        if (ftype != F_DISRES)
         {
-            if (ffparams->functype[type] == ftype)
+            for (int type = start; type < ffparams->numTypes(); type++)
             {
-                if (memcmp(&newparam, &ffparams->iparams[type], static_cast<size_t>(sizeof(newparam))) == 0)
+                // Note that the first condition is always met by starting the loop at start
+                if (ffparams->functype[type] == ftype
+                    && memcmp(&newparam, &ffparams->iparams[type], static_cast<size_t>(sizeof(newparam))) == 0)
                 {
                     return type;
                 }
             }
         }
+        else
+        {
+            // Distance restraints should have unique labels and pairs with the same label
+            // should be consecutive, so we here we only need to check the last type in the list.
+            // This changes the complexity from quadratic to linear in the number of restraints.
+            const int type = ffparams->numTypes() - 1;
+            if (type >= 0 && ffparams->functype[type] == ftype
+                && memcmp(&newparam, &ffparams->iparams[type], static_cast<size_t>(sizeof(newparam))) == 0)
+            {
+                return type;
+            }
+        }
     }
-    else
-    {
-        type = ffparams->numTypes();
-    }
+
+    const int type = ffparams->numTypes();
 
     ffparams->iparams.push_back(newparam);
     ffparams->functype.push_back(ftype);