Move densityfitting options and output into seperate file
authorChristian Blau <cblau@gwdg.de>
Thu, 8 Aug 2019 09:15:41 +0000 (11:15 +0200)
committerChristian Blau <cblau@gwdg.de>
Thu, 15 Aug 2019 11:56:22 +0000 (13:56 +0200)
Exposes density fitting options and output interface for testing.

Added test for density fitting options class.

refs: #2282

Change-Id: I0107217b5af86abf1ed298dd77553c4c24b1d3f2

src/gromacs/applied_forces/CMakeLists.txt
src/gromacs/applied_forces/densityfitting.cpp
src/gromacs/applied_forces/densityfittingforceprovider.h
src/gromacs/applied_forces/densityfittingoptions.cpp [new file with mode: 0644]
src/gromacs/applied_forces/densityfittingoptions.h [new file with mode: 0644]
src/gromacs/applied_forces/densityfittingoutputprovider.cpp [new file with mode: 0644]
src/gromacs/applied_forces/densityfittingoutputprovider.h [new file with mode: 0644]
src/gromacs/applied_forces/densityfittingparameters.h
src/gromacs/applied_forces/tests/CMakeLists.txt
src/gromacs/applied_forces/tests/densityfitting.cpp
src/gromacs/applied_forces/tests/densityfittingoptions.cpp [new file with mode: 0644]

index c2454d9c207e7a4cdca97f9d80cad6ab2c876d88..d5a61a6e07a67f062fd3a7a5a1b4d7298e4a11e6 100644 (file)
@@ -36,6 +36,8 @@ gmx_add_libgromacs_sources(
     densityfitting.cpp
     densityfittingamplitudelookup.cpp
     densityfittingforceprovider.cpp
+    densityfittingoptions.cpp
+    densityfittingoutputprovider.cpp
     electricfield.cpp
     )
 
index b18486d714f8275fa32e36bfc00d2cdc21ae2ed0..be7dbbc1d4973a8bc7043e76ee4975c5f2493b32 100644 (file)
 
 #include "densityfitting.h"
 
-#include "gromacs/mdtypes/iforceprovider.h"
 #include "gromacs/mdtypes/imdmodule.h"
-#include "gromacs/mdtypes/imdoutputprovider.h"
-#include "gromacs/mdtypes/imdpoptionprovider.h"
-#include "gromacs/options/basicoptions.h"
-#include "gromacs/options/optionsection.h"
-#include "gromacs/utility/keyvaluetreebuilder.h"
-#include "gromacs/utility/keyvaluetreetransform.h"
-#include "gromacs/utility/stringutil.h"
 
 #include "densityfittingforceprovider.h"
-#include "densityfittingparameters.h"
+#include "densityfittingoptions.h"
+#include "densityfittingoutputprovider.h"
 
 namespace gmx
 {
 
+class IMdpOptionProvider;
+class DensityFittingForceProvider;
+
 namespace
 {
 
-/*! \internal
- * \brief Input data storage for density fitting
- */
-class DensityFittingOptions : public IMdpOptionProvider
-{
-    public:
-        DensityFittingOptions()
-        { }
-
-        //! From IMdpOptionProvider
-        void initMdpTransform(IKeyValueTreeTransformRules * /*transform*/ ) override
-        {}
-
-        /*! \brief
-         * Build mdp parameters for density fitting to be output after pre-processing.
-         * \param[in, out] builder the builder for the mdp options output KV-tree.
-         * \note This should be symmetrical to option initialization without
-         *       employing manual prefixing with the section name string once
-         *       the legacy code blocking this design is removed.
-         */
-        void buildMdpOutput(KeyValueTreeObjectBuilder *builder) const override
-        {
-            builder->addValue<bool>(DensityFittingModuleInfo::name_ + "-" + c_activeTag_,
-                                    active_);
-        }
-
-        /*! \brief
-         * Connect option name and data.
-         */
-        void initMdpOptions(IOptionsContainerWithSections *options) override
-        {
-            auto section = options->addSection(OptionSection(DensityFittingModuleInfo::name_.c_str()));
-            section.addOption(
-                    BooleanOption(c_activeTag_.c_str()).store(&active_));
-        }
-
-        //! Report if this set of options is active
-        bool active() const
-        {
-            return active_;
-        }
-
-        //! Process input options to parameters, including input file reading.
-        DensityFittingParameters buildParameters()
-        {
-            // read map to mrc-header and data vector
-            // convert mrcheader and data vector to grid data
-            return {};
-        }
-
-    private:
-        //! The name of the density-fitting module
-
-        const std::string     c_activeTag_             = "active";
-        bool                  active_                  = false;
-};
-
 /*! \internal
  * \brief Density fitting
  *
  * Class that implements the density fitting forces and potential
  * \note the virial calculation is not yet implemented
  */
-class DensityFitting final : public IMDModule,
-                             public IMDOutputProvider
+class DensityFitting final : public IMDModule
 {
     public:
         DensityFitting() = default;
@@ -146,19 +84,13 @@ class DensityFitting final : public IMDModule,
         }
 
         //! This MDModule provides its own output
-        IMDOutputProvider *outputProvider() override { return this; }
-
-        //! Initialize output
-        void initOutput(FILE * /*fplog*/, int /*nfile*/, const t_filenm /*fnm*/[],
-                        bool /*bAppendFiles*/, const gmx_output_env_t * /*oenv*/) override
-        {}
-
-        //! Finalizes output from a simulation run.
-        void finishOutput() override {}
+        IMDOutputProvider *outputProvider() override { return &densityFittingOutputProvider_; }
 
     private:
+        //! The output provider
+        DensityFittingOutputProvider                 densityFittingOutputProvider_;
         //! The options provided for density fitting
-        DensityFittingOptions densityFittingOptions_;
+        DensityFittingOptions                        densityFittingOptions_;
         //! Object that evaluates the forces
         std::unique_ptr<DensityFittingForceProvider> forceProvider_;
 };
index a5aab572a5b6e84027fe4803a24091bbda4006d7..0e170d34ffac352627c5faa1dc9990bb81ff04e3 100644 (file)
 
 namespace gmx
 {
-class DensityFittingParameters;
+struct DensityFittingParameters;
 
 class DensityFittingForceProvider final : public IForceProvider
 {
     public:
-        DensityFittingForceProvider(const DensityFittingParameters &paramters);
+        DensityFittingForceProvider(const DensityFittingParameters &parameters);
         void calculateForces(const ForceProviderInput &forceProviderInput,
                              ForceProviderOutput      *forceProviderOutput) override;
 
diff --git a/src/gromacs/applied_forces/densityfittingoptions.cpp b/src/gromacs/applied_forces/densityfittingoptions.cpp
new file mode 100644 (file)
index 0000000..e0ce323
--- /dev/null
@@ -0,0 +1,132 @@
+/*
+ * This file is part of the GROMACS molecular simulation package.
+ *
+ * Copyright (c) 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.
+ *
+ * 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
+ * Implements force provider for density fitting
+ *
+ * \author Christian Blau <blau@kth.se>
+ * \ingroup module_applied_forces
+ */
+#include "gmxpre.h"
+
+#include "densityfittingoptions.h"
+
+#include "gromacs/applied_forces/densityfitting.h"
+#include "gromacs/options/basicoptions.h"
+#include "gromacs/options/optionsection.h"
+#include "gromacs/utility/keyvaluetreebuilder.h"
+#include "gromacs/utility/keyvaluetreetransform.h"
+#include "gromacs/utility/strconvert.h"
+
+namespace gmx
+{
+
+namespace
+{
+
+/*! \brief Helper to declare mdp transform rules.
+ *
+ * Enforces uniform mdp options that are always prepended with the correct
+ * string for the densityfitting mdp options.
+ *
+ * \tparam ToType type to be transformed to
+ * \tparam TransformWithFunctionType type of transformation function to be used
+ *
+ * \param[in] rules KVT transformation rules
+ * \param[in] transformationFunction the function to transform the flat kvt tree
+ * \param[in] optionTag string tag that describes the mdp option, appended to the
+ *                      default string for the density guided simulation
+ */
+template <class ToType, class TransformWithFunctionType>
+void densityfittingMdpTransformFromString(IKeyValueTreeTransformRules * rules,
+                                          TransformWithFunctionType     transformationFunction,
+                                          const std::string            &optionTag)
+{
+    rules->addRule()
+        .from<std::string>("/" + DensityFittingModuleInfo::name_ + "-" + optionTag)
+        .to<ToType>("/" + DensityFittingModuleInfo::name_ +"/" + optionTag)
+        .transformWith(transformationFunction);
+}
+/*! \brief Helper to declare mdp output.
+ *
+ * Enforces uniform mdp options output sting that are always prepended with the
+ * correct string for the densityfitting mdp options and is consistent with the
+ * options name and transformation.
+ *
+ * \tparam OptionType the type of the mdp option
+ * \param[in] builder the KVT builder to generate the output
+ * \param[in] option the mdp option
+ * \param[in] optionTag string tag that describes the mdp option, appended to the
+ *                      default string for the density guided simulation
+ */
+template <class OptionType>
+void addDensityFittingMdpOutputValue(KeyValueTreeObjectBuilder *builder,
+                                     const OptionType          &option,
+                                     const std::string         &optionTag)
+{
+    builder->addValue<OptionType>(DensityFittingModuleInfo::name_ + "-" + optionTag,
+                                  option);
+}
+
+}   // namespace
+
+void DensityFittingOptions::initMdpTransform(IKeyValueTreeTransformRules * rules)
+{
+    densityfittingMdpTransformFromString<bool>(rules, &fromStdString<bool>, c_activeTag_);
+}
+
+void DensityFittingOptions::buildMdpOutput(KeyValueTreeObjectBuilder *builder) const
+{
+    addDensityFittingMdpOutputValue(builder, parameters_.active_, c_activeTag_);
+}
+
+void DensityFittingOptions::initMdpOptions(IOptionsContainerWithSections *options)
+{
+    auto section = options->addSection(OptionSection(DensityFittingModuleInfo::name_.c_str()));
+
+    section.addOption(BooleanOption(c_activeTag_.c_str()).store(&parameters_.active_));
+}
+
+bool DensityFittingOptions::active() const
+{
+    return parameters_.active_;
+}
+
+const DensityFittingParameters &DensityFittingOptions::buildParameters()
+{
+    return parameters_;
+}
+
+} // namespace gmx
diff --git a/src/gromacs/applied_forces/densityfittingoptions.h b/src/gromacs/applied_forces/densityfittingoptions.h
new file mode 100644 (file)
index 0000000..523c41f
--- /dev/null
@@ -0,0 +1,90 @@
+/*
+ * This file is part of the GROMACS molecular simulation package.
+ *
+ * Copyright (c) 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.
+ *
+ * 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
+ * Declares options for density fitting
+ *
+ * \author Christian Blau <blau@kth.se>
+ * \ingroup module_applied_forces
+ */
+#ifndef GMX_APPLIED_FORCES_DENSITYFITTINGOPTIONS_H
+#define GMX_APPLIED_FORCES_DENSITYFITTINGOPTIONS_H
+
+#include <string>
+
+#include "gromacs/mdtypes/imdpoptionprovider.h"
+
+#include "densityfittingparameters.h"
+
+namespace gmx
+{
+
+/*! \internal
+ * \brief Input data storage for density fitting
+ */
+class DensityFittingOptions final : public IMdpOptionProvider
+{
+    public:
+        //! From IMdpOptionProvider
+        void initMdpTransform(IKeyValueTreeTransformRules * rules) override;
+
+        /*! \brief
+         * Build mdp parameters for density fitting to be output after pre-processing.
+         * \param[in, out] builder the builder for the mdp options output KV-tree.
+         * \note This should be symmetrical to option initialization without
+         *       employing manual prefixing with the section name string once
+         *       the legacy code blocking this design is removed.
+         */
+        void buildMdpOutput(KeyValueTreeObjectBuilder *builder) const override;
+
+        /*! \brief
+         * Connect option name and data.
+         */
+        void initMdpOptions(IOptionsContainerWithSections *options) override;
+
+        //! Report if this set of options is active
+        bool active() const;
+
+        //! Process input options to parameters, including input file reading.
+        const DensityFittingParameters &buildParameters();
+
+    private:
+        const std::string        c_activeTag_ = "active";
+        DensityFittingParameters parameters_;
+};
+
+} // namespace gmx
+
+#endif
diff --git a/src/gromacs/applied_forces/densityfittingoutputprovider.cpp b/src/gromacs/applied_forces/densityfittingoutputprovider.cpp
new file mode 100644 (file)
index 0000000..c32f49a
--- /dev/null
@@ -0,0 +1,56 @@
+/*
+ * This file is part of the GROMACS molecular simulation package.
+ *
+ * Copyright (c) 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.
+ *
+ * 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
+ * Declares data structure for density fitting output
+ *
+ * \author Christian Blau <blau@kth.se>
+ * \ingroup module_applied_forces
+ */
+
+#include "gmxpre.h"
+
+#include "densityfittingoutputprovider.h"
+namespace gmx
+{
+
+void DensityFittingOutputProvider::initOutput(FILE * /*fplog*/, int /*nfile*/,
+                                              const t_filenm /*fnm*/[], bool /*bAppendFiles*/, const gmx_output_env_t * /*oenv*/)
+{}
+
+void DensityFittingOutputProvider::finishOutput()
+{}
+
+} // namespace gmx
diff --git a/src/gromacs/applied_forces/densityfittingoutputprovider.h b/src/gromacs/applied_forces/densityfittingoutputprovider.h
new file mode 100644 (file)
index 0000000..6bbc3c6
--- /dev/null
@@ -0,0 +1,65 @@
+/*
+ * This file is part of the GROMACS molecular simulation package.
+ *
+ * Copyright (c) 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.
+ *
+ * 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
+ * Declares output to file for density fitting
+ *
+ * \author Christian Blau <blau@kth.se>
+ * \ingroup module_applied_forces
+ */
+#ifndef GMX_APPLIED_FORCES_DENSITYFITTINGOUTPUTPROVIDER_H
+#define GMX_APPLIED_FORCES_DENSITYFITTINGOUTPUTPROVIDER_H
+
+#include "gromacs/mdtypes/imdoutputprovider.h"
+
+namespace gmx
+{
+
+/*! \internal
+ * \brief Handle file output for density guided simulations.
+ */
+class DensityFittingOutputProvider final : public IMDOutputProvider
+{
+    public:
+        //! Initialize output
+        void initOutput(FILE * /*fplog*/, int /*nfile*/, const t_filenm /*fnm*/[],
+                        bool /*bAppendFiles*/, const gmx_output_env_t * /*oenv*/) override;
+        //! Finalizes output from a simulation run.
+        void finishOutput() override;
+};
+
+} // namespace gmx
+
+#endif
index 06af6320e518a939fbb5b917fc0bed021745cd57..e2df3712916255df9b9177bf054c4a6527a5ade1 100644 (file)
 namespace gmx
 {
 
-class DensityFittingParameters
+/*! \internal
+ * \brief Holding all directly user-provided parameters for density fitting.
+ *
+ * Also used for setting all default parameters.
+ */
+struct DensityFittingParameters
 {
+    //! Indicate if density fitting is active
+    bool active_ = false;
 };
 
 }      // namespace gmx
index ff580f70cafdea28a5f06c6377eed1969a891151..d9cd94c6da3c8c78bd1beeadc3ceb2b96b49c9d2 100644 (file)
@@ -33,7 +33,7 @@
 # the research papers on the package. Check out http://www.gromacs.org.
 
 gmx_add_unit_test(AppliedForcesUnitTest applied_forces-test
-                  densityfitting.cpp
+                  densityfittingoptions.cpp
                   densityfittingamplitudelookup.cpp
                   electricfield.cpp
                 )
index 0791ca8aad6005e25942d109f4f9c09510d82c45..ddf950237273d4d4bf5762e54fa63eed487c72e9 100644 (file)
@@ -71,27 +71,27 @@ namespace gmx
 namespace
 {
 
-TEST(DensityFittingTest, Options)
+TEST(DensityFittingTest, ForceOnSingleOption)
 {
-    auto densityFittingModule(gmx::DensityFittingModuleInfo::create());
+    auto densityFittingModule(DensityFittingModuleInfo::create());
 
     // Prepare MDP inputs
-    gmx::KeyValueTreeBuilder mdpValueBuilder;
-    mdpValueBuilder.rootObject().addValue("density-guided-simulation-active", "yes");
+    KeyValueTreeBuilder      mdpValueBuilder;
+    mdpValueBuilder.rootObject().addValue("density-guided-simulation-active", std::string("yes"));
     KeyValueTreeObject       densityFittingMdpValues = mdpValueBuilder.build();
 
     // set up options
-    gmx::Options densityFittingModuleOptions;
+    Options densityFittingModuleOptions;
     densityFittingModule->mdpOptionProvider()->initMdpOptions(&densityFittingModuleOptions);
 
     // Add rules to transform mdp inputs to densityFittingModule data
-    gmx::KeyValueTreeTransformer transform;
-    transform.rules()->addRule().keyMatchType("/", gmx::StringCompareType::CaseAndDashInsensitive);
+    KeyValueTreeTransformer transform;
+    transform.rules()->addRule().keyMatchType("/", StringCompareType::CaseAndDashInsensitive);
     densityFittingModule->mdpOptionProvider()->initMdpTransform(transform.rules());
 
     // Execute the transform on the mdpValues
     auto transformedMdpValues = transform.transform(densityFittingMdpValues, nullptr);
-    gmx::assignOptionsFromKeyValueTree(&densityFittingModuleOptions, transformedMdpValues.object(), nullptr);
+    assignOptionsFromKeyValueTree(&densityFittingModuleOptions, transformedMdpValues.object(), nullptr);
 
     // Build the force provider, once all input data is gathered
     ForceProviders densityFittingForces;
@@ -100,18 +100,18 @@ TEST(DensityFittingTest, Options)
     // Build a minimal simulation system.
     t_mdatoms               mdAtoms;
     mdAtoms.homenr = 1;
-    PaddedVector<gmx::RVec> x             = {{0, 0, 0}};
+    PaddedVector<RVec>      x             = {{0, 0, 0}};
     matrix                  simulationBox = {{0, 0, 0}, {0, 0, 0}, {0, 0, 0}};
     const double            t             = 0.0;
     t_commrec              *cr            = init_commrec();
-    gmx::ForceProviderInput forceProviderInput(x, mdAtoms, t, simulationBox, *cr);
+    ForceProviderInput      forceProviderInput(x, mdAtoms, t, simulationBox, *cr);
 
     // The forces that the force-provider is to update
-    PaddedVector<gmx::RVec>  f = {{0, 0, 0}};
-    gmx::ForceWithVirial     forceWithVirial(f, false);
+    PaddedVector<RVec>       f = {{0, 0, 0}};
+    ForceWithVirial          forceWithVirial(f, false);
 
     gmx_enerdata_t           energyData(1, 0);
-    gmx::ForceProviderOutput forceProviderOutput(&forceWithVirial, &energyData);
+    ForceProviderOutput      forceProviderOutput(&forceWithVirial, &energyData);
 
     // update the forces
     densityFittingForces.calculateForces(forceProviderInput, &forceProviderOutput);
diff --git a/src/gromacs/applied_forces/tests/densityfittingoptions.cpp b/src/gromacs/applied_forces/tests/densityfittingoptions.cpp
new file mode 100644 (file)
index 0000000..7e04304
--- /dev/null
@@ -0,0 +1,123 @@
+/*
+ * This file is part of the GROMACS molecular simulation package.
+ *
+ * Copyright (c) 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.
+ *
+ * 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 density fitting module options.
+ *
+ * \author Christian Blau <blau@kth.se>
+ * \ingroup module_applied_forces
+ */
+#include "gmxpre.h"
+
+#include "gromacs/applied_forces/densityfittingoptions.h"
+
+#include <gtest/gtest.h>
+
+#include "gromacs/gmxpreprocess/keyvaluetreemdpwriter.h"
+#include "gromacs/options/options.h"
+#include "gromacs/options/treesupport.h"
+#include "gromacs/utility/keyvaluetreebuilder.h"
+#include "gromacs/utility/keyvaluetreetransform.h"
+#include "gromacs/utility/stringcompare.h"
+#include "gromacs/utility/stringstream.h"
+#include "gromacs/utility/textwriter.h"
+
+#include "testutils/testasserts.h"
+#include "testutils/testmatchers.h"
+
+
+namespace gmx
+{
+
+namespace
+{
+
+TEST(DensityFittingOptionsTest, DefaultParameters)
+{
+    DensityFittingOptions densityFittingOptions;
+    EXPECT_FALSE(densityFittingOptions.buildParameters().active_);
+}
+
+TEST(DensityFittingOptionsTest, OptionSetsActive)
+{
+    DensityFittingOptions densityFittingOptions;
+    EXPECT_FALSE(densityFittingOptions.buildParameters().active_);
+
+    // Prepare MDP inputs
+    KeyValueTreeBuilder mdpValueBuilder;
+    mdpValueBuilder.rootObject().addValue("density-guided-simulation-active",
+                                          std::string("yes"));
+    KeyValueTreeObject  densityFittingMdpValues = mdpValueBuilder.build();
+
+    // set up options
+    Options densityFittingModuleOptions;
+    densityFittingOptions.initMdpOptions(&densityFittingModuleOptions);
+
+    // Add rules to transform mdp inputs to densityFittingModule data
+    KeyValueTreeTransformer transform;
+    transform.rules()->addRule().keyMatchType("/", StringCompareType::CaseAndDashInsensitive);
+
+    densityFittingOptions.initMdpTransform(transform.rules());
+
+    // Execute the transform on the mdpValues
+    auto transformedMdpValues = transform.transform(densityFittingMdpValues, nullptr);
+    assignOptionsFromKeyValueTree(&densityFittingModuleOptions, transformedMdpValues.object(), nullptr);
+
+    EXPECT_TRUE(densityFittingOptions.buildParameters().active_);
+}
+
+TEST(DensityFittingOptionsTest, OutputDefaultValues)
+{
+    // Transform module data into a flat key-value tree for output.
+
+    StringOutputStream        stream;
+    KeyValueTreeBuilder       builder;
+    KeyValueTreeObjectBuilder builderObject = builder.rootObject();
+
+    DensityFittingOptions     densityFittingOptions;
+    densityFittingOptions.buildMdpOutput(&builderObject);
+
+    {
+        TextWriter writer(&stream);
+        writeKeyValueTreeAsMdp(&writer, builder.build());
+    }
+    stream.close();
+
+    EXPECT_EQ(stream.toString(), std::string("density-guided-simulation-active = false\n"));
+}
+
+} // namespace
+
+} // namespace gmx