src/gromacs/applied_forces/awh/tests/biasstate.cpp

   1 /*
   2  * This file is part of the GROMACS molecular simulation package.
   3  *
   4  * Copyright (c) 2017,2018,2019,2020, by the GROMACS development team, led by
   5  * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
   6  * and including many others, as listed in the AUTHORS file in the
   7  * top-level source directory and at http://www.gromacs.org.
   8  *
   9  * GROMACS is free software; you can redistribute it and/or
  10  * modify it under the terms of the GNU Lesser General Public License
  11  * as published by the Free Software Foundation; either version 2.1
  12  * of the License, or (at your option) any later version.
  13  *
  14  * GROMACS is distributed in the hope that it will be useful,
  15  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  16  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  17  * Lesser General Public License for more details.
  18  *
  19  * You should have received a copy of the GNU Lesser General Public
  20  * License along with GROMACS; if not, see
  21  * http://www.gnu.org/licenses, or write to the Free Software Foundation,
  22  * Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA.
  23  *
  24  * If you want to redistribute modifications to GROMACS, please
  25  * consider that scientific software is very special. Version
  26  * control is crucial - bugs must be traceable. We will be happy to
  27  * consider code for inclusion in the official distribution, but
  28  * derived work must not be called official GROMACS. Details are found
  29  * in the README & COPYING files - if they are missing, get the
  30  * official version at http://www.gromacs.org.
  31  *
  32  * To help us fund GROMACS development, we humbly ask that you cite
  33  * the research papers on the package. Check out http://www.gromacs.org.
  34  */
  35 #include "gmxpre.h"
  36
  37 #include "gromacs/applied_forces/awh/biasstate.h"
  38
  39 #include <cmath>
  40
  41 #include <memory>
  42 #include <vector>
  43
  44 #include <gmock/gmock.h>
  45 #include <gtest/gtest.h>
  46
  47 #include "gromacs/applied_forces/awh/biasgrid.h"
  48 #include "gromacs/applied_forces/awh/pointstate.h"
  49 #include "gromacs/math/functions.h"
  50 #include "gromacs/mdtypes/awh_params.h"
  51 #include "gromacs/utility/arrayref.h"
  52 #include "gromacs/utility/smalloc.h"
  53
  54 #include "testutils/testasserts.h"
  55 #include "testutils/testfilemanager.h"
  56
  57 namespace gmx
  58 {
  59
  60 namespace test
  61 {
  62
  63 /*! \internal \brief
  64  * Struct that gathers all input for setting up and using a Bias
  65  */
  66 struct AwhTestParameters
  67 {
  68     double beta; //!< 1/(kB*T)
  69
  70     AwhDimParams  awhDimParams[2]; //!< Dimension parameters pointed to by \p awhBiasParams
  71     AwhBiasParams awhBiasParams;   //!< Bias parameters pointed to by \[ awhParams
  72     AwhParams     awhParams;       //!< AWH parameters, this is the struct to actually use
  73 };
  74
  75 //! Helper function to set up the C-style AWH parameters for the test
  76 static AwhTestParameters getAwhTestParameters()
  77 {
  78     AwhTestParameters params;
  79
  80     params.beta = 1.0;
  81
  82     AwhParams& awhParams = params.awhParams;
  83     snew(params.awhParams.awhBiasParams, 1);
  84     AwhBiasParams& awhBiasParams = params.awhParams.awhBiasParams[0];
  85     snew(awhBiasParams.dimParams, 2);
  86
  87     AwhDimParams& awhDimParams0 = awhBiasParams.dimParams[0];
  88
  89     awhDimParams0.period         = 0;
  90     awhDimParams0.diffusion      = 0.1;
  91     awhDimParams0.origin         = 0.5;
  92     awhDimParams0.end            = 1.5;
  93     awhDimParams0.coordValueInit = awhDimParams0.origin;
  94     awhDimParams0.coverDiameter  = 0;
  95     awhDimParams0.eCoordProvider = eawhcoordproviderPULL;
  96
  97     AwhDimParams& awhDimParams1 = awhBiasParams.dimParams[1];
  98
  99     awhDimParams1.period         = 0;
 100     awhDimParams1.diffusion      = 0.1;
 101     awhDimParams1.origin         = 0.8;
 102     awhDimParams1.end            = 1.3;
 103     awhDimParams1.coordValueInit = awhDimParams1.origin;
 104     awhDimParams1.coverDiameter  = 0;
 105     awhDimParams1.eCoordProvider = eawhcoordproviderPULL;
 106
 107     awhBiasParams.ndim                 = 2;
 108     awhBiasParams.eTarget              = eawhtargetCONSTANT;
 109     awhBiasParams.targetBetaScaling    = 0;
 110     awhBiasParams.targetCutoff         = 0;
 111     awhBiasParams.eGrowth              = eawhgrowthLINEAR;
 112     awhBiasParams.bUserData            = TRUE;
 113     awhBiasParams.errorInitial         = 0.5;
 114     awhBiasParams.shareGroup           = 0;
 115     awhBiasParams.equilibrateHistogram = FALSE;
 116
 117     awhParams.numBias                    = 1;
 118     awhParams.seed                       = 93471803;
 119     awhParams.nstOut                     = 0;
 120     awhParams.nstSampleCoord             = 1;
 121     awhParams.numSamplesUpdateFreeEnergy = 10;
 122     awhParams.ePotential                 = eawhpotentialCONVOLVED;
 123     awhParams.shareBiasMultisim          = FALSE;
 124
 125     return params;
 126 }
 127
 128 /*! \brief Test fixture for testing Bias updates
 129  */
 130 class BiasStateTest : public ::testing::TestWithParam<const char*>
 131 {
 132 public:
 133     std::unique_ptr<BiasState> biasState_; //!< The bias state
 134
 135     BiasStateTest()
 136     {
 137         AwhTestParameters      params        = getAwhTestParameters();
 138         const AwhParams&       awhParams     = params.awhParams;
 139         const AwhBiasParams&   awhBiasParams = awhParams.awhBiasParams[0];
 140         std::vector<DimParams> dimParams;
 141         dimParams.push_back(DimParams::pullDimParams(1.0, 15.0, params.beta));
 142         dimParams.push_back(DimParams::pullDimParams(1.0, 15.0, params.beta));
 143         BiasGrid   grid(dimParams, awhBiasParams.dimParams);
 144         BiasParams biasParams(
 145                 awhParams, awhBiasParams, dimParams, 1.0, 1.0, BiasParams::DisableUpdateSkips::no, 1, grid.axis(), 0);
 146         biasState_ = std::make_unique<BiasState>(awhBiasParams, 1.0, dimParams, grid);
 147
 148         // Here we initialize the grid point state using the input file
 149         std::string filename = gmx::test::TestFileManager::getInputFilePath(GetParam());
 150         biasState_->initGridPointState(
 151                 awhBiasParams, dimParams, grid, biasParams, filename, params.awhParams.numBias);
 152
 153         sfree(params.awhParams.awhBiasParams[0].dimParams);
 154         sfree(params.awhParams.awhBiasParams);
 155     }
 156 };
 157
 158 TEST_P(BiasStateTest, InitializesFromFile)
 159 {
 160     gmx::ArrayRef<const PointState> points = biasState_->points();
 161
 162     /* Compute the mean square deviation from the expected values in the file.
 163      * The PMF values are spaced by 0.5 per points and logPmfsum has opposite sign.
 164      * The target is (index + 1)/120.
 165      */
 166     double msdPmf = 0;
 167     for (index i = 0; i < points.ssize(); i++)
 168     {
 169         msdPmf += gmx::square(points[i].logPmfSum() - points[0].logPmfSum() + 0.5 * i) / points.size();
 170         EXPECT_DOUBLE_EQ(points[i].target(), (i + 1) / 120.0);
 171     }
 172
 173     EXPECT_NEAR(0.0, msdPmf, 1e-31);
 174 }
 175
 176 // Test that Bias initialization open and reads the correct initialization
 177 // files and the correct PMF and target distribution is set.
 178 INSTANTIATE_TEST_CASE_P(WithParameters,
 179                         BiasStateTest,
 180                         ::testing::Values("pmf_target_format0.xvg", "pmf_target_format1.xvg"));
 181
 182 } // namespace test
 183 } // namespace gmx