Generic tpr handling for electric field

Remove special code to store inputrec values to tpr files from
IInputRecExtension, and instead store them as a generic KeyValueTree.

Currently, there is no support for changing the storage format, other
than by adding new options. We may not need this immediately, and we
think it should not be any extra work to introduce it later (rather
than right now). With the present approach, tpr backwards
compatibility is essentially the same as mdp backwards compatibility.

There are now two paths to setting up modules for use, depending
whether the data is coming from mdp parsing in e.g. grompp or tests,
or from tpr parsing in e.g. mdrun, dump or check (but not tpbconv).

Change-Id: I0a0f552041b2a7401bc82972c9648d1f3a64af34

diff --git a/src/gromacs/applied-forces/electricfield.cpp b/src/gromacs/applied-forces/electricfield.cpp
index d5a5f533fc05bf5f81508e0740d2a08cc3b54b46..964ad3411a78184eb6ec23aaf0e3714fab553d9d 100644 (file)
--- a/src/gromacs/applied-forces/electricfield.cpp
+++ b/src/gromacs/applied-forces/electricfield.cpp
@@ -47,9 +47,6 @@
 
 #include "gromacs/commandline/filenm.h"
 #include "gromacs/fileio/gmxfio.h"
-#include "gromacs/fileio/gmxfio-xdr.h"
-#include "gromacs/fileio/readinp.h"
-#include "gromacs/fileio/warninp.h"
 #include "gromacs/fileio/xvgr.h"
 #include "gromacs/gmxlib/network.h"
 #include "gromacs/math/units.h"
@@ -165,7 +162,6 @@ class ElectricField : public IInputRecExtension, public IForceProvider
         ElectricField() : fpField_(nullptr) {}
 
         // From IInputRecExtension
-        virtual void doTpxIO(t_fileio *fio, bool bRead);
         virtual void initMdpTransform(IKeyValueTreeTransformRules *transform);
         virtual void initMdpOptions(IOptionsContainerWithSections *options);
         virtual void broadCast(const t_commrec *cr);
@@ -247,54 +243,6 @@ class ElectricField : public IInputRecExtension, public IForceProvider
         FILE             *fpField_;
 };
 
-void ElectricField::doTpxIO(t_fileio *fio, bool bRead)
-{
-    // The content of the tpr file for this feature has
-    // been the same since gromacs 4.0 that was used for
-    // developing.
-    for (int j = 0; (j < DIM); j++)
-    {
-        int n = 0, nt = 0;
-        if (!bRead)
-        {
-            n = 1;
-            if (omega(j) != 0 || sigma(j) != 0 || t0(j) != 0)
-            {
-                nt = 1;
-            }
-        }
-        gmx_fio_do_int(fio, n);
-        gmx_fio_do_int(fio, nt);
-        std::vector<real> aa, phi, at, phit;
-        if (!bRead)
-        {
-            aa.push_back(a(j));
-            phi.push_back(t0(j));
-            at.push_back(omega(j));
-            phit.push_back(sigma(j));
-        }
-        else
-        {
-            aa.resize(n+1);
-            phi.resize(nt+1);
-            at.resize(nt+1);
-            phit.resize(nt+1);
-        }
-        gmx_fio_ndo_real(fio, aa.data(),  n);
-        gmx_fio_ndo_real(fio, phi.data(), n);
-        gmx_fio_ndo_real(fio, at.data(),  nt);
-        gmx_fio_ndo_real(fio, phit.data(), nt);
-        if (bRead && n > 0)
-        {
-            setFieldTerm(j, aa[0], at[0], phi[0], phit[0]);
-            if (n > 1 || nt > 1)
-            {
-                gmx_fatal(FARGS, "Can not handle tpr files with more than one electric field term per direction.");
-            }
-        }
-    }
-}
-
 //! Converts static parameters from mdp format to E0.
 real convertStaticParameters(const std::string &value)
 {

diff --git a/src/gromacs/applied-forces/tests/electricfield.cpp b/src/gromacs/applied-forces/tests/electricfield.cpp
index 3e0d0ccff6e4baec5884b6cf314cc34344b62699..2bbb7ea0fb28b600907e1c82c60851584043dc58 100644 (file)
--- a/src/gromacs/applied-forces/tests/electricfield.cpp
+++ b/src/gromacs/applied-forces/tests/electricfield.cpp
@@ -100,7 +100,7 @@ class ElectricFieldTest : public ::testing::Test
                 .keyMatchType("/", gmx::StringCompareType::CaseAndDashInsensitive);
             module.initMdpTransform(transform.rules());
             auto result = transform.transform(mdpValues.build(), nullptr);
-            module.assignOptionsToModules(result.object(), nullptr);
+            module.assignOptionsToModulesFromMdp(result.object(), nullptr);
 
             t_mdatoms        md;
             PaddedRVecVector f = { { 0, 0, 0 } };

diff --git a/src/gromacs/fileio/tpxio.cpp b/src/gromacs/fileio/tpxio.cpp
index cdd4a22389ba916c820c762a5303e4f3aee98909..a4ef8718bd4297a9b2e91e43bb68c3f988b0d09d 100644 (file)
--- a/src/gromacs/fileio/tpxio.cpp
+++ b/src/gromacs/fileio/tpxio.cpp
@@ -69,6 +69,8 @@
 #include "gromacs/utility/fatalerror.h"
 #include "gromacs/utility/futil.h"
 #include "gromacs/utility/gmxassert.h"
+#include "gromacs/utility/keyvaluetreebuilder.h"
+#include "gromacs/utility/keyvaluetreeserializer.h"
 #include "gromacs/utility/smalloc.h"
 #include "gromacs/utility/snprintf.h"
 #include "gromacs/utility/txtdump.h"
@@ -112,6 +114,7 @@ enum tpxv {
     tpxv_RemoveTwinRange,                                    /**< removed support for twin-range interactions */
     tpxv_ReplacePullPrintCOM12,                              /**< Replaced print-com-1, 2 with pull-print-com */
     tpxv_PullExternalPotential,                              /**< Added pull type external potential */
+    tpxv_GenericParamsForElectricField,                      /**< Introduced KeyValueTree and moved electric field parameters */
     tpxv_Count                                               /**< the total number of tpxv versions */
 };
 
@@ -919,6 +922,40 @@ static void do_swapcoords_tpx(t_fileio *fio, t_swapcoords *swap, gmx_bool bRead,
 
 }
 
+static void do_legacy_efield(t_fileio *fio, gmx::KeyValueTreeObjectBuilder *root)
+{
+    const char *const dimName[] = { "x", "y", "z" };
+
+    auto              appliedForcesObj = root->addObject("applied-forces");
+    auto              efieldObj        = appliedForcesObj.addObject("electric-field");
+    // The content of the tpr file for this feature has
+    // been the same since gromacs 4.0 that was used for
+    // developing.
+    for (int j = 0; j < DIM; ++j)
+    {
+        int n, nt;
+        gmx_fio_do_int(fio, n);
+        gmx_fio_do_int(fio, nt);
+        std::vector<real> aa(n+1), phi(nt+1), at(nt+1), phit(nt+1);
+        gmx_fio_ndo_real(fio, aa.data(),  n);
+        gmx_fio_ndo_real(fio, phi.data(), n);
+        gmx_fio_ndo_real(fio, at.data(),  nt);
+        gmx_fio_ndo_real(fio, phit.data(), nt);
+        if (n > 0)
+        {
+            if (n > 1 || nt > 1)
+            {
+                gmx_fatal(FARGS, "Can not handle tpr files with more than one electric field term per direction.");
+            }
+            auto dimObj = efieldObj.addObject(dimName[j]);
+            dimObj.addValue<real>("E0", aa[0]);
+            dimObj.addValue<real>("omega", at[0]);
+            dimObj.addValue<real>("t0", phi[0]);
+            dimObj.addValue<real>("sigma", phit[0]);
+        }
+    }
+}
+
 
 static void do_inputrec(t_fileio *fio, t_inputrec *ir, gmx_bool bRead,
                         int file_version, real *fudgeQQ)
@@ -939,6 +976,9 @@ static void do_inputrec(t_fileio *fio, t_inputrec *ir, gmx_bool bRead,
         return;
     }
 
+    gmx::KeyValueTreeBuilder       paramsBuilder;
+    gmx::KeyValueTreeObjectBuilder paramsObj = paramsBuilder.rootObject();
+
     /* Basic inputrec stuff */
     gmx_fio_do_int(fio, ir->eI);
     if (file_version >= 62)
@@ -1626,7 +1666,10 @@ static void do_inputrec(t_fileio *fio, t_inputrec *ir, gmx_bool bRead,
         ir->wall_ewald_zfac  = 3;
     }
     /* Cosine stuff for electric fields */
-    ir->efield->doTpxIO(fio, bRead);
+    if (file_version < tpxv_GenericParamsForElectricField)
+    {
+        do_legacy_efield(fio, &paramsObj);
+    }
 
     /* Swap ions */
     if (file_version >= tpxv_ComputationalElectrophysiology)
@@ -1681,6 +1724,26 @@ static void do_inputrec(t_fileio *fio, t_inputrec *ir, gmx_bool bRead,
         }
         /* end of QMMM stuff */
     }
+
+    if (file_version >= tpxv_GenericParamsForElectricField)
+    {
+        gmx::FileIOXdrSerializer serializer(fio);
+        if (bRead)
+        {
+            paramsObj.mergeObject(
+                    gmx::deserializeKeyValueTree(&serializer));
+        }
+        else
+        {
+            GMX_RELEASE_ASSERT(ir->params != nullptr,
+                               "Parameters should be present when writing inputrec");
+            gmx::serializeKeyValueTree(*ir->params, &serializer);
+        }
+    }
+    if (bRead)
+    {
+        ir->params = new gmx::KeyValueTreeObject(paramsBuilder.build());
+    }
 }
 
 

diff --git a/src/gromacs/gmxpreprocess/readir.cpp b/src/gromacs/gmxpreprocess/readir.cpp
index bad314a1a870d54dd24fab2978dcbe8a4e496a04..fc848135b1086f8e5bf057e3158a38af5fbfb536 100644 (file)
--- a/src/gromacs/gmxpreprocess/readir.cpp
+++ b/src/gromacs/gmxpreprocess/readir.cpp
@@ -2251,7 +2251,7 @@ void get_ir(const char *mdparin, const char *mdparout,
         auto                         result
             = transform.transform(convertedValues, &errorHandler);
         errorHandler.setBackMapping(result.backMapping());
-        mdModules->assignOptionsToModules(result.object(), &errorHandler);
+        mdModules->assignOptionsToModulesFromMdp(result.object(), &errorHandler);
     }
 
     /* Ion/water position swapping ("computational electrophysiology") */

diff --git a/src/gromacs/mdlib/broadcaststructs.cpp b/src/gromacs/mdlib/broadcaststructs.cpp
index fcebb03e97f1ebe171ff73ac87d756149d66e8db..35f0e9e6da4edab3cbe17b27805cd0901119d41e 100644 (file)
--- a/src/gromacs/mdlib/broadcaststructs.cpp
+++ b/src/gromacs/mdlib/broadcaststructs.cpp
@@ -668,6 +668,7 @@ static void bc_inputrec(const t_commrec *cr, t_inputrec *inputrec)
     gmx::IInputRecExtension *eptr = inputrec->efield;
     block_bc(cr, *inputrec);
     inputrec->efield = eptr;
+    inputrec->params = nullptr;
 
     bc_grpopts(cr, &(inputrec->opts));
 

diff --git a/src/gromacs/mdrunutility/mdmodules.cpp b/src/gromacs/mdrunutility/mdmodules.cpp
index 01ad57f58e0d9c730d20752033c7555f7f4f46bf..c7bcd9a5be4d41c784affcf2f92c93418eca6b0e 100644 (file)
--- a/src/gromacs/mdrunutility/mdmodules.cpp
+++ b/src/gromacs/mdrunutility/mdmodules.cpp
@@ -78,6 +78,14 @@ class MDModules::Impl
             }
         }
 
+        void makeModuleOptions(Options *options)
+        {
+            // Create a section for applied-forces modules
+            auto appliedForcesOptions = options->addSection(OptionSection("applied-forces"));
+            field_->initMdpOptions(&appliedForcesOptions);
+            // In future, other sections would also go here.
+        }
+
         IInputRecExtensionPtr  field_;
         t_inputrec            *ir_;
 };
@@ -108,16 +116,32 @@ void MDModules::initMdpTransform(IKeyValueTreeTransformRules *rules)
     impl_->field_->initMdpTransform(rules);
 }
 
-void MDModules::assignOptionsToModules(const KeyValueTreeObject  &optionValues,
-                                       IKeyValueTreeErrorHandler *errorHandler)
+void MDModules::assignOptionsToModulesFromMdp(const KeyValueTreeObject  &mdpOptionValues,
+                                              IKeyValueTreeErrorHandler *errorHandler)
 {
-    Options options;
-    // Create a section for applied-forces modules
-    auto    appliedForcesOptions = options.addSection(OptionSection("applied-forces"));
-    impl_->field_->initMdpOptions(&appliedForcesOptions);
-    // In future, other sections would also go here.
+    Options moduleOptions;
+    impl_->makeModuleOptions(&moduleOptions);
+
+    KeyValueTreeObject keyValueParameters(mdpOptionValues);
+    impl_->ir_->params = new KeyValueTreeObject(adjustKeyValueTreeFromOptions(keyValueParameters, moduleOptions));
+    // The actual output is in the data fields of the modules that
+    // were set up in the module options.
+    assignOptionsFromKeyValueTree(&moduleOptions, *impl_->ir_->params, errorHandler);
+}
 
-    assignOptionsFromKeyValueTree(&options, optionValues, errorHandler);
+void MDModules::assignOptionsToModulesFromTpr()
+{
+    Options moduleOptions;
+    impl_->makeModuleOptions(&moduleOptions);
+
+    // Note that impl_->ir_->params was set up during tpr reading, so
+    // all we need to do here is integrate that with the module
+    // options, which e.g. might have changed between versions.
+    // The actual output is in the data fields of the modules that
+    // were set up in the module options.
+    //
+    // TODO error handling
+    assignOptionsFromKeyValueTree(&moduleOptions, *impl_->ir_->params, nullptr);
 }
 
 } // namespace gmx

diff --git a/src/gromacs/mdrunutility/mdmodules.h b/src/gromacs/mdrunutility/mdmodules.h
index 4b5a6a19269416876dbaf3e3eb3d770ecc5b8ade..916caa59d3b95c653d6ed8e1872fe685f1a3e4e3 100644 (file)
--- a/src/gromacs/mdrunutility/mdmodules.h
+++ b/src/gromacs/mdrunutility/mdmodules.h
@@ -82,9 +82,11 @@ class IKeyValueTreeTransformRules;
  * loops over a container of these interfaces, and/or groups of them (e.g.
  * applied forces), instead of the current single pointer.
  *
- * The assignOptionsToModules() method of this class also takes
- * responsibility for wiring up the options (and their defaults) for
- * each
+ * The assignOptionsToModules() and
+ * assignOptionsToModulesFromInputrec() methods of this class also
+ * take responsibility for wiring up the options (and their defaults)
+ * for each module, respectively for mdp- and tpr-style input of those
+ * options.
  *
  * \inlibraryapi
  * \ingroup module_mdrunutility
@@ -117,14 +119,23 @@ class MDModules
          */
         void initMdpTransform(IKeyValueTreeTransformRules *rules);
 
-        /*! \brief Use \c optionValues to set the options for each module.
+        /*! \brief Use \c mdpOptionValues to set the options (e.g.read
+         * from mdp input) for each module.
          *
-         * \param[in] optionValues Contains keys and values from user
+         * \param[in] mdpOptionValues Contains keys and values from user
          *     input (and defaults) to configure modules that have
          *     registered options with those keys.
          * \param[out] errorHandler  Called to report errors. */
-        void assignOptionsToModules(const KeyValueTreeObject  &optionValues,
-                                    IKeyValueTreeErrorHandler *errorHandler);
+        void assignOptionsToModulesFromMdp(const KeyValueTreeObject  &mdpOptionValues,
+                                           IKeyValueTreeErrorHandler *errorHandler);
+
+        /*! \brief
+         * Initializes modules based on inputrec values read from tpr file.
+         *
+         * This needs to be called after read_tpx_state() if the modules need
+         * to be accessed.
+         */
+        void assignOptionsToModulesFromTpr();
 
     private:
         class Impl;

diff --git a/src/gromacs/mdtypes/inputrec.cpp b/src/gromacs/mdtypes/inputrec.cpp
index a10bc47d68ab95c539de1958b8a38c4d541f2066..c73c2fa255bb2494759988b83033e126b2ce01d5 100644 (file)
--- 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, by the GROMACS development team, led by
+ * Copyright (c) 2012,2014,2015,2016,2017, 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.
@@ -51,6 +51,7 @@
 #include "gromacs/utility/compare.h"
 #include "gromacs/utility/cstringutil.h"
 #include "gromacs/utility/fatalerror.h"
+#include "gromacs/utility/keyvaluetree.h"
 #include "gromacs/utility/smalloc.h"
 #include "gromacs/utility/snprintf.h"
 #include "gromacs/utility/stringutil.h"
@@ -307,6 +308,7 @@ void done_inputrec(t_inputrec *ir)
         done_pull_params(ir->pull);
         sfree(ir->pull);
     }
+    delete ir->params;
 }
 
 static void pr_qm_opts(FILE *fp, int indent, const char *title, const t_grpopts *opts)

diff --git a/src/gromacs/mdtypes/inputrec.h b/src/gromacs/mdtypes/inputrec.h
index b68f83ffe56e57a3e6f0bdfc55f0e6ba1783d5d5..3a43a39eecc244addfa2bb7360ab77d3898aabe5 100644 (file)
--- a/src/gromacs/mdtypes/inputrec.h
+++ b/src/gromacs/mdtypes/inputrec.h
@@ -67,6 +67,7 @@ namespace gmx
 
 class IKeyValueTreeTransformRules;
 class IOptionsContainerWithSections;
+class KeyValueTreeObject;
 
 /*! \libinternal \brief
  * Inputrec extension interface for a mdrun module.
@@ -91,15 +92,6 @@ class IInputRecExtension
     public:
         virtual ~IInputRecExtension() {}
 
-        /*! \brief Read or write tpr file
-         *
-         * Read or write the necessary data from a tpr file. The routine is responsible
-         * for consistency, such that all data belonging to this module is read or written.
-         * \param[inout] fio Gromacs file descriptor
-         * \param[in]    bRead boolean determines whether we are reading or writing
-         */
-        virtual void doTpxIO(t_fileio *fio, bool bRead) = 0;
-
         /*! \brief
          * Initializes a transform from mdp values to sectioned options.
          *
@@ -492,8 +484,10 @@ struct t_inputrec
     real                     scalefactor;   /* factor for scaling the MM charges in QM calc.*/
 
     /* Fields for removed features go here (better caching) */
-    gmx_bool        bAdress;       // Whether AdResS is enabled - always false if a valid .tpr was read
-    gmx_bool        useTwinRange;  // Whether twin-range scheme is active - always false if a valid .tpr was read
+    gmx_bool                 bAdress;      // Whether AdResS is enabled - always false if a valid .tpr was read
+    gmx_bool                 useTwinRange; // Whether twin-range scheme is active - always false if a valid .tpr was read
+
+    gmx::KeyValueTreeObject *params;
 };
 
 int ir_optimal_nstcalcenergy(const t_inputrec *ir);

diff --git a/src/gromacs/tools/check.cpp b/src/gromacs/tools/check.cpp
index 5b916196ea07ebba5147dcf08cd72b1fa647f4a8..fc61d8fd0e302d3393413914139093a57923769d 100644 (file)
--- a/src/gromacs/tools/check.cpp
+++ b/src/gromacs/tools/check.cpp
@@ -103,6 +103,7 @@ static void comp_tpx(const char *fn1, const char *fn2,
     for (i = 0; i < (fn2 ? 2 : 1); i++)
     {
         read_tpx_state(ff[i], ir[i], &state[i], &(mtop[i]));
+        mdModules[i].assignOptionsToModulesFromTpr();
     }
     if (fn2)
     {

diff --git a/src/gromacs/tools/dump.cpp b/src/gromacs/tools/dump.cpp
index 6b9e31c30d4e32470b278d3c5ccaa21236f651a6..4cc89d9bf9b270964968c49debb3c482c88e0c51 100644 (file)
--- a/src/gromacs/tools/dump.cpp
+++ b/src/gromacs/tools/dump.cpp
@@ -95,6 +95,10 @@ static void list_tpx(const char *fn,
                    ir,
                    &state,
                    tpx.bTop ? &mtop : nullptr);
+    if (tpx.bIr)
+    {
+        mdModules.assignOptionsToModulesFromTpr();
+    }
 
     if (mdpfn && tpx.bIr)
     {

diff --git a/src/gromacs/utility/keyvaluetreebuilder.h b/src/gromacs/utility/keyvaluetreebuilder.h
index f5856df371ddb35a4b5d521e9aedf12d49f79a57..d8ce0d9edd2e7fe6b032ac31119ebeeb06ed02b5 100644 (file)
--- a/src/gromacs/utility/keyvaluetreebuilder.h
+++ b/src/gromacs/utility/keyvaluetreebuilder.h
@@ -1,7 +1,7 @@
 /*
  * This file is part of the GROMACS molecular simulation package.
  *
- * Copyright (c) 2016, by the GROMACS development team, led by
+ * Copyright (c) 2016,2017, 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.
@@ -215,7 +215,10 @@ class KeyValueTreeObjectBuilder
         }
         void mergeObject(KeyValueTreeValue &&value)
         {
-            KeyValueTreeObject &obj = value.asObject();
+            mergeObject(std::move(value.asObject()));
+        }
+        void mergeObject(KeyValueTreeObject &&obj)
+        {
             for (auto &prop : obj.valueMap_)
             {
                 addRawValue(prop.first, std::move(prop.second));

diff --git a/src/programs/mdrun/runner.cpp b/src/programs/mdrun/runner.cpp
index 63327cbd434942f18959a08936335142daa4740f..d17c5f9f4a16e4a661bff574e1ba9aba605bad16 100644 (file)
--- a/src/programs/mdrun/runner.cpp
+++ b/src/programs/mdrun/runner.cpp
@@ -779,6 +779,7 @@ int mdrunner(gmx_hw_opt_t *hw_opt,
     {
         /* Read (nearly) all data required for the simulation */
         read_tpx_state(ftp2fn(efTPR, nfile, fnm), inputrec, state, mtop);
+        mdModules.assignOptionsToModulesFromTpr();
 
         if (inputrec->cutoff_scheme == ecutsVERLET)
         {