From 9453ea8cfaffca76d999a3bb9659a4916519f012 Mon Sep 17 00:00:00 2001
From: Paul Bauer <paul.bauer.q@gmail.com>
Date: Thu, 28 Feb 2019 10:39:57 +0100
Subject: [PATCH] Fix segmentation fault with simulated annealing

Fix contributed by Daniel Fragiadakis.

Fixes #2871

Change-Id: I237d49f9996eb2ae869d59e4c39a6c2d33c9ece7
---
 docs/release-notes/2019/2019.2.rst   |  8 ++++++++
 src/gromacs/gmxpreprocess/readir.cpp | 21 +++++++++++++++++----
 src/gromacs/mdtypes/inputrec.cpp     |  7 ++++++-
 3 files changed, 31 insertions(+), 5 deletions(-)

diff --git a/docs/release-notes/2019/2019.2.rst b/docs/release-notes/2019/2019.2.rst
index 257a08c716..b0b7d9acfa 100644
--- a/docs/release-notes/2019/2019.2.rst
+++ b/docs/release-notes/2019/2019.2.rst
@@ -19,6 +19,14 @@ Fixes where mdrun could behave incorrectly
 Fixes for ``gmx`` tools
 ^^^^^^^^^^^^^^^^^^^^^^^
 
+Fix segmentation fault when preparing simulated annealing inputs
+""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
+
+grompp was unable to prepare tpr files for inputs containing simulated annealing
+procedures. The code has been fixed to allow the generation of those files again.
+
+:issue:`2871`
+
 Fixes that affect portability
 ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
 
diff --git a/src/gromacs/gmxpreprocess/readir.cpp b/src/gromacs/gmxpreprocess/readir.cpp
index 7c0a259fc9..6b82ca99b2 100644
--- a/src/gromacs/gmxpreprocess/readir.cpp
+++ b/src/gromacs/gmxpreprocess/readir.cpp
@@ -3328,7 +3328,7 @@ void do_index(const char* mdparin, const char *ndx,
     if (!simulatedAnnealingGroupNames.empty() &&
         simulatedAnnealingGroupNames.size() != size_t(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
@@ -3374,7 +3374,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++)
+                int k = 0;
+                for (i = 0; i < nr; i++)
                 {
                     if (ir->opts.anneal_npoints[i] == 1)
                     {
@@ -3386,6 +3387,7 @@ void do_index(const char* mdparin, const char *ndx,
                 }
 
                 auto simulatedAnnealingTimes = gmx::splitString(is->anneal_time);
+
                 if (simulatedAnnealingTimes.size() != size_t(k))
                 {
                     gmx_fatal(FARGS, "Found %zu annealing-time values, wanted %d\n",
@@ -3398,12 +3400,21 @@ void do_index(const char* mdparin, const char *ndx,
                               simulatedAnnealingTemperatures.size(), k);
                 }
 
-                convertReals(wi, simulatedAnnealingTimes, "anneal-time", ir->opts.anneal_time[i]);
-                convertReals(wi, simulatedAnnealingTemperatures, "anneal-temp", ir->opts.anneal_temp[i]);
+                real *allSimulatedAnnealingTimes;
+                real *allSimulatedAnnealingTemperatures;
+                snew(allSimulatedAnnealingTimes, k);
+                snew(allSimulatedAnnealingTemperatures, k);
+                convertReals(wi, simulatedAnnealingTimes, "anneal-time", allSimulatedAnnealingTimes);
+                convertReals(wi, simulatedAnnealingTemperatures, "anneal-temp", allSimulatedAnnealingTemperatures);
+                /* TODO the code here could use some more general clean up and restructuring.
+                 * This will also help in understanding the use of variables in the loop.
+                 */
                 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))
@@ -3427,6 +3438,8 @@ void do_index(const char* mdparin, const char *ndx,
                         k++;
                     }
                 }
+                sfree(allSimulatedAnnealingTimes);
+                sfree(allSimulatedAnnealingTemperatures);
                 /* Print out some summary information, to make sure we got it right */
                 for (i = 0, k = 0; i < nr; i++)
                 {
diff --git a/src/gromacs/mdtypes/inputrec.cpp b/src/gromacs/mdtypes/inputrec.cpp
index 6edc255c4f..1afe02a840 100644
--- a/src/gromacs/mdtypes/inputrec.cpp
+++ b/src/gromacs/mdtypes/inputrec.cpp
@@ -3,7 +3,7 @@
  *
  * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
  * Copyright (c) 2001-2010, The GROMACS development team.
- * Copyright (c) 2012,2014,2015,2016,2017,2018, by the GROMACS development team, led by
+ * Copyright (c) 2012,2014,2015,2016,2017,2018,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.
@@ -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);
-- 
2.11.4.GIT