Fix simulated annealing and add test

Refs #2871

Change-Id: Ibde5226a664d3d4d5558477fdf9b81e1219a5295

diff --git a/src/gromacs/gmxpreprocess/readir.cpp b/src/gromacs/gmxpreprocess/readir.cpp
index 126c8f735233b6de3b62a0aa577c6411a11f5a8f..6fd7e16227252611473fe4fb1f146471bb47ba53 100644 (file)
--- a/src/gromacs/gmxpreprocess/readir.cpp
+++ b/src/gromacs/gmxpreprocess/readir.cpp
@@ -3324,9 +3324,9 @@ void do_index(const char* mdparin, const char *ndx,
         simulatedAnnealingGroupNames.resize(0);
     }
     if (!simulatedAnnealingGroupNames.empty() &&
-        simulatedAnnealingGroupNames.size() != size_t(nr))
+        gmx::ssize(simulatedAnnealingGroupNames) != nr)
     {
-        gmx_fatal(FARGS, "Not enough annealing values: %zu (for %d groups)\n",
+        gmx_fatal(FARGS, "Wrong number of annealing values: %zu (for %d groups)\n",
                   simulatedAnnealingGroupNames.size(), nr);
     }
     else
@@ -3372,7 +3372,8 @@ void do_index(const char* mdparin, const char *ndx,
                               simulatedAnnealingPoints.size(), simulatedAnnealingGroupNames.size());
                 }
                 convertInts(wi, simulatedAnnealingPoints, "annealing points", ir->opts.anneal_npoints);
-                for (k = 0, i = 0; i < nr; i++)
+                size_t numSimulatedAnnealingFields = 0;
+                for (i = 0; i < nr; i++)
                 {
                     if (ir->opts.anneal_npoints[i] == 1)
                     {
@@ -3380,28 +3381,33 @@ void do_index(const char* mdparin, const char *ndx,
                     }
                     snew(ir->opts.anneal_time[i], ir->opts.anneal_npoints[i]);
                     snew(ir->opts.anneal_temp[i], ir->opts.anneal_npoints[i]);
-                    k += ir->opts.anneal_npoints[i];
+                    numSimulatedAnnealingFields += ir->opts.anneal_npoints[i];
                 }
 
                 auto simulatedAnnealingTimes = gmx::splitString(is->anneal_time);
-                if (simulatedAnnealingTimes.size() != size_t(k))
+
+                if (simulatedAnnealingTimes.size() != numSimulatedAnnealingFields)
                 {
-                    gmx_fatal(FARGS, "Found %zu annealing-time values, wanted %d\n",
-                              simulatedAnnealingTimes.size(), k);
+                    gmx_fatal(FARGS, "Found %zu annealing-time values, wanted %zu\n",
+                              simulatedAnnealingTimes.size(), numSimulatedAnnealingFields);
                 }
                 auto simulatedAnnealingTemperatures = gmx::splitString(is->anneal_temp);
-                if (simulatedAnnealingTemperatures.size() != size_t(k))
+                if (simulatedAnnealingTemperatures.size() != numSimulatedAnnealingFields)
                 {
-                    gmx_fatal(FARGS, "Found %zu annealing-temp values, wanted %d\n",
-                              simulatedAnnealingTemperatures.size(), k);
+                    gmx_fatal(FARGS, "Found %zu annealing-temp values, wanted %zu\n",
+                              simulatedAnnealingTemperatures.size(), numSimulatedAnnealingFields);
                 }
 
-                convertReals(wi, simulatedAnnealingTimes, "anneal-time", ir->opts.anneal_time[i]);
-                convertReals(wi, simulatedAnnealingTemperatures, "anneal-temp", ir->opts.anneal_temp[i]);
+                std::vector<real> allSimulatedAnnealingTimes(numSimulatedAnnealingFields);
+                std::vector<real> allSimulatedAnnealingTemperatures(numSimulatedAnnealingFields);
+                convertReals(wi, simulatedAnnealingTimes, "anneal-time", allSimulatedAnnealingTimes.data());
+                convertReals(wi, simulatedAnnealingTemperatures, "anneal-temp", allSimulatedAnnealingTemperatures.data());
                 for (i = 0, k = 0; i < nr; i++)
                 {
                     for (j = 0; j < ir->opts.anneal_npoints[i]; j++)
                     {
+                        ir->opts.anneal_time[i][j] = allSimulatedAnnealingTimes[k];
+                        ir->opts.anneal_temp[i][j] = allSimulatedAnnealingTemperatures[k];
                         if (j == 0)
                         {
                             if (ir->opts.anneal_time[i][0] > (ir->init_t+GMX_REAL_EPS))
@@ -3426,7 +3432,7 @@ void do_index(const char* mdparin, const char *ndx,
                     }
                 }
                 /* Print out some summary information, to make sure we got it right */
-                for (i = 0, k = 0; i < nr; i++)
+                for (i = 0; i < nr; i++)
                 {
                     if (ir->opts.annealing[i] != eannNO)
                     {

diff --git a/src/gromacs/mdtypes/inputrec.cpp b/src/gromacs/mdtypes/inputrec.cpp
index 46d58cf5021d698259f747247bd26d079a4dc9e6..5f7f16aa5f21437e41be6ad9fda12a95c3a3b9d2 100644 (file)
--- a/src/gromacs/mdtypes/inputrec.cpp
+++ b/src/gromacs/mdtypes/inputrec.cpp
@@ -307,6 +307,11 @@ void done_inputrec(t_inputrec *ir)
 {
     sfree(ir->opts.nrdf);
     sfree(ir->opts.ref_t);
+    for (int i = 0; i < ir->opts.ngtc; i++)
+    {
+        sfree(ir->opts.anneal_time[i]);
+        sfree(ir->opts.anneal_temp[i]);
+    }
     sfree(ir->opts.annealing);
     sfree(ir->opts.anneal_npoints);
     sfree(ir->opts.anneal_time);

diff --git a/src/programs/mdrun/tests/grompp.cpp b/src/programs/mdrun/tests/grompp.cpp
index afd93175ba235b13b07ed664f38554d5fc3e505e..c81c1ee8c68386430c7cdb037948769d8a398574 100644 (file)
--- a/src/programs/mdrun/tests/grompp.cpp
+++ b/src/programs/mdrun/tests/grompp.cpp
@@ -1,7 +1,7 @@
 /*
  * This file is part of the GROMACS molecular simulation package.
  *
- * Copyright (c) 2015, by the GROMACS development team, led by
+ * Copyright (c) 2015,2019, by the GROMACS development team, led by
  * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
  * and including many others, as listed in the AUTHORS file in the
  * top-level source directory and at http://www.gromacs.org.
@@ -75,4 +75,29 @@ TEST_F(GromppTest, EmptyMdpFileWorks)
     runTest();
 }
 
+/* Test for making sure grompp can handle simulated annealing data */
+TEST_F(GromppTest, SimulatedAnnealingWorks)
+{
+    runner_.useStringAsMdpFile("annealing = periodic\n"
+                               "annealing-npoints = 4\n"
+                               "annealing-time = 0 2 4 6\n"
+                               "annealing-temp = 298 320 320 298\n"
+                               );
+    runTest();
+}
+
+TEST_F(GromppTest, SimulatedAnnealingWorksWithMultipleGroups)
+{
+    runner_.useStringAsMdpFile("tc-grps = Methanol SOL\n"
+                               "tau-t = 0.1 0.1\n"
+                               "ref_t = 298 298\n"
+                               "annealing = single periodic\n"
+                               "annealing-npoints = 3 4\n"
+                               "annealing-time = 0 3 6 0 2 4 6\n"
+                               "annealing-temp = 298 280 270 298 320 320 298\n"
+                               );
+    runTest();
+}
+
+
 } // namespace