From ff41a4a4eccc017a8e408d5ed7efeacc9b23cd6a Mon Sep 17 00:00:00 2001 From: Kevin Boyd Date: Fri, 3 Jul 2020 13:16:11 +0000 Subject: [PATCH] Fix segfault in gmx genrestr. A gmx_mtop_t going out of scope freed memory in a related structure. Extended the life of the topology structure, and fixed memory leaks that popped up in the ASAN build. --- docs/release-notes/2020/2020.3.rst | 8 ++ src/gromacs/gmxpreprocess/genrestr.cpp | 32 ++++---- src/gromacs/gmxpreprocess/tests/CMakeLists.txt | 3 +- src/gromacs/gmxpreprocess/tests/genrestr.cpp | 86 ++++++++++++++++++++++ .../GenRestrTest_SimpleRestraintsGenerated.xml | 24 ++++++ 5 files changed, 138 insertions(+), 15 deletions(-) create mode 100644 src/gromacs/gmxpreprocess/tests/genrestr.cpp create mode 100644 src/gromacs/gmxpreprocess/tests/refdata/GenRestrTest_SimpleRestraintsGenerated.xml diff --git a/docs/release-notes/2020/2020.3.rst b/docs/release-notes/2020/2020.3.rst index 6ccebd353c..050153c2da 100644 --- a/docs/release-notes/2020/2020.3.rst +++ b/docs/release-notes/2020/2020.3.rst @@ -60,6 +60,14 @@ Users can configure the default path for dssp using GMX_DSSP_PROGRAM_PATH. :issue:`3520` +Avoid segmentation fault in gmx genrestr +"""""""""""""""""""""""""""""""""""""""" + +The tool could fail when running simple inputs due to memory access errors +caused by accessing free'd memory. + +:issue:`3582` + Fixes that affect portability ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ diff --git a/src/gromacs/gmxpreprocess/genrestr.cpp b/src/gromacs/gmxpreprocess/genrestr.cpp index 137bef2ab2..74515a7638 100644 --- a/src/gromacs/gmxpreprocess/genrestr.cpp +++ b/src/gromacs/gmxpreprocess/genrestr.cpp @@ -3,7 +3,7 @@ * * Copyright (c) 1991-2000, University of Groningen, The Netherlands. * Copyright (c) 2001-2004, The GROMACS development team. - * Copyright (c) 2013,2014,2015,2017,2018,2019, by the GROMACS development team, led by + * Copyright (c) 2013,2014,2015,2017,2018,2019,2020, 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. @@ -142,9 +142,7 @@ int gmx_genrestr(int argc, char* argv[]) FILE* out; int igrp; real d, dd, lo, hi; - int* ind_grp; const char * xfn, *nfn; - char* gn_grp; matrix box; gmx_bool bFreeze; rvec dx, *x = nullptr, *v = nullptr; @@ -159,6 +157,7 @@ int gmx_genrestr(int argc, char* argv[]) { return 0; } + output_env_done(oenv); bFreeze = opt2bSet("-of", NFILE, fnm) || opt2parg_bSet("-freeze", asize(pa), pa); bDisre = bDisre || opt2parg_bSet("-disre_dist", npargs, pa); @@ -181,10 +180,13 @@ int gmx_genrestr(int argc, char* argv[]) const char* title = ""; bool haveTopology = false; + gmx_mtop_t mtop; + int* indexGroups = nullptr; + char* indexGroupNames = nullptr; + if (xfn != nullptr) { fprintf(stderr, "\nReading structure file\n"); - gmx_mtop_t mtop; readConfAndTopology(xfn, &haveTopology, &mtop, nullptr, &x, &v, box); title = *mtop.name; atoms = gmx_mtop_global_atoms(&mtop); @@ -216,18 +218,18 @@ int gmx_genrestr(int argc, char* argv[]) else if ((bDisre || bConstr) && x) { printf("Select group to generate %s matrix from\n", bConstr ? "constraint" : "distance restraint"); - get_index(&atoms, nfn, 1, &igrp, &ind_grp, &gn_grp); + get_index(&atoms, nfn, 1, &igrp, &indexGroups, &indexGroupNames); out = ftp2FILE(efITP, NFILE, fnm, "w"); if (bConstr) { - fprintf(out, "; constraints for %s of %s\n\n", gn_grp, title); + fprintf(out, "; constraints for %s of %s\n\n", indexGroupNames, title); fprintf(out, "[ constraints ]\n"); fprintf(out, ";%4s %5s %1s %10s\n", "i", "j", "tp", "dist"); } else { - fprintf(out, "; distance restraints for %s of %s\n\n", gn_grp, title); + fprintf(out, "; distance restraints for %s of %s\n\n", indexGroupNames, title); fprintf(out, "[ distance_restraints ]\n"); fprintf(out, ";%4s %5s %1s %5s %10s %10s %10s %10s %10s\n", "i", "j", "?", "label", "funct", "lo", "up1", "up2", "weight"); @@ -236,11 +238,11 @@ int gmx_genrestr(int argc, char* argv[]) { for (j = i + 1; j < igrp; j++, k++) { - rvec_sub(x[ind_grp[i]], x[ind_grp[j]], dx); + rvec_sub(x[indexGroups[i]], x[indexGroups[j]], dx); d = norm(dx); if (bConstr) { - fprintf(out, "%5d %5d %1d %10g\n", ind_grp[i] + 1, ind_grp[j] + 1, 2, d); + fprintf(out, "%5d %5d %1d %10g\n", indexGroups[i] + 1, indexGroups[j] + 1, 2, d); } else { @@ -256,8 +258,8 @@ int gmx_genrestr(int argc, char* argv[]) } lo = std::max(0.0_real, d - dd); hi = d + dd; - fprintf(out, "%5d %5d %1d %5d %10d %10g %10g %10g %10g\n", ind_grp[i] + 1, - ind_grp[j] + 1, 1, k, 1, lo, hi, hi + disre_up2, 1.0); + fprintf(out, "%5d %5d %1d %5d %10d %10g %10g %10g %10g\n", indexGroups[i] + 1, + indexGroups[j] + 1, 1, k, 1, lo, hi, hi + disre_up2, 1.0); } } } @@ -267,15 +269,15 @@ int gmx_genrestr(int argc, char* argv[]) else { printf("Select group to position restrain\n"); - get_index(&atoms, nfn, 1, &igrp, &ind_grp, &gn_grp); + get_index(&atoms, nfn, 1, &igrp, &indexGroups, &indexGroupNames); out = ftp2FILE(efITP, NFILE, fnm, "w"); - fprintf(out, "; position restraints for %s of %s\n\n", gn_grp, title); + fprintf(out, "; position restraints for %s of %s\n\n", indexGroupNames, title); fprintf(out, "[ position_restraints ]\n"); fprintf(out, ";%3s %5s %9s %10s %10s\n", "i", "funct", "fcx", "fcy", "fcz"); for (i = 0; i < igrp; i++) { - fprintf(out, "%4d %4d %10g %10g %10g\n", ind_grp[i] + 1, 1, fc[XX], fc[YY], fc[ZZ]); + fprintf(out, "%4d %4d %10g %10g %10g\n", indexGroups[i] + 1, 1, fc[XX], fc[YY], fc[ZZ]); } gmx_ffclose(out); } @@ -285,6 +287,8 @@ int gmx_genrestr(int argc, char* argv[]) sfree(v); done_atom(&atoms); } + sfree(indexGroupNames); + sfree(indexGroups); return 0; } diff --git a/src/gromacs/gmxpreprocess/tests/CMakeLists.txt b/src/gromacs/gmxpreprocess/tests/CMakeLists.txt index c164166a0b..3be1db5528 100644 --- a/src/gromacs/gmxpreprocess/tests/CMakeLists.txt +++ b/src/gromacs/gmxpreprocess/tests/CMakeLists.txt @@ -1,7 +1,7 @@ # # This file is part of the GROMACS molecular simulation package. # -# Copyright (c) 2014,2015,2017,2018,2019, by the GROMACS development team, led by +# Copyright (c) 2014,2015,2017,2018,2019,2020, 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. @@ -36,6 +36,7 @@ gmx_add_unit_test(GmxPreprocessTests gmxpreprocess-test editconf.cpp genconf.cpp genion.cpp + genrestr.cpp gpp_atomtype.cpp gpp_bond_atomtype.cpp insert_molecules.cpp diff --git a/src/gromacs/gmxpreprocess/tests/genrestr.cpp b/src/gromacs/gmxpreprocess/tests/genrestr.cpp new file mode 100644 index 0000000000..ad8daf19dd --- /dev/null +++ b/src/gromacs/gmxpreprocess/tests/genrestr.cpp @@ -0,0 +1,86 @@ +/* + * This file is part of the GROMACS molecular simulation package. + * + * Copyright (c) 2020, 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. + * + * GROMACS is free software; you can redistribute it and/or + * modify it under the terms of the GNU Lesser General Public License + * as published by the Free Software Foundation; either version 2.1 + * of the License, or (at your option) any later version. + * + * GROMACS is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU + * Lesser General Public License for more details. + * + * You should have received a copy of the GNU Lesser General Public + * License along with GROMACS; if not, see + * http://www.gnu.org/licenses, or write to the Free Software Foundation, + * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. + * + * If you want to redistribute modifications to GROMACS, please + * consider that scientific software is very special. Version + * control is crucial - bugs must be traceable. We will be happy to + * consider code for inclusion in the official distribution, but + * derived work must not be called official GROMACS. Details are found + * in the README & COPYING files - if they are missing, get the + * official version at http://www.gromacs.org. + * + * To help us fund GROMACS development, we humbly ask that you cite + * the research papers on the package. Check out http://www.gromacs.org. + */ +/*! \internal \file + * \brief + * Tests for gmx genrestr. + * + * \author Kevin Boyd + */ + +#include "gmxpre.h" + +#include "gromacs/gmxpreprocess/genrestr.h" + +#include "testutils/cmdlinetest.h" +#include "testutils/filematchers.h" +#include "testutils/refdata.h" +#include "testutils/stdiohelper.h" +#include "testutils/testfilemanager.h" +#include "testutils/textblockmatchers.h" + +namespace gmx +{ +namespace test +{ + +class GenRestrTest : public CommandLineTestBase +{ +public: + void runTest(const std::string& interactiveCommandLineInput) + { + StdioTestHelper stdIoHelper(&fileManager()); + stdIoHelper.redirectStringToStdin(interactiveCommandLineInput.c_str()); + + CommandLine& cmdline = commandLine(); + // Provide the name of the module to call + std::string module[] = { "genrestr" }; + cmdline.merge(CommandLine(module)); + + ASSERT_EQ(0, gmx_genrestr(cmdline.argc(), cmdline.argv())); + checkOutputFiles(); + } +}; + +TEST_F(GenRestrTest, SimpleRestraintsGenerated) +{ + setInputFile("-f", "lysozyme.pdb"); + ExactTextMatch settings; + setOutputFile("-o", "restraints.itp", TextFileMatch(settings)); + // Select c-alphas from default index options. + std::string selection = "3"; + runTest(selection); +} +} // namespace test +} // namespace gmx diff --git a/src/gromacs/gmxpreprocess/tests/refdata/GenRestrTest_SimpleRestraintsGenerated.xml b/src/gromacs/gmxpreprocess/tests/refdata/GenRestrTest_SimpleRestraintsGenerated.xml new file mode 100644 index 0000000000..ea751d45d0 --- /dev/null +++ b/src/gromacs/gmxpreprocess/tests/refdata/GenRestrTest_SimpleRestraintsGenerated.xml @@ -0,0 +1,24 @@ + + + + + + + + + -- 2.11.4.GIT