#include "gromacs/mdtypes/awh_params.h"
#include "gromacs/utility/stringutil.h"
+#include "gromacs/applied_forces/awh/tests/awh_setup.h"
#include "testutils/refdata.h"
#include "testutils/testasserts.h"
{
//! The number of lambda states to use in the tests.
-const int numLambdaStates = 16;
+static constexpr int c_numLambdaStates = 16;
-/*! \internal \brief
- * Struct that gathers all input for setting up and using a Bias
- */
-struct AwhFepLambdaStateTestParameters
-{
- AwhFepLambdaStateTestParameters() = default;
- //! Move constructor
- AwhFepLambdaStateTestParameters(AwhFepLambdaStateTestParameters&& o) noexcept :
- beta(o.beta),
- awhDimParams(o.awhDimParams),
- awhBiasParams(o.awhBiasParams),
- awhParams(o.awhParams),
- dimParams(std::move(o.dimParams))
- {
- awhBiasParams.dimParams = &awhDimParams;
- awhParams.awhBiasParams = &awhBiasParams;
- }
- double beta; //!< 1/(kB*T)
-
- AwhDimParams awhDimParams; //!< Dimension parameters pointed to by \p awhBiasParams
- AwhBiasParams awhBiasParams; //!< Bias parameters pointed to by \[ awhParams
- AwhParams awhParams; //!< AWH parameters, this is the struct to actually use
-
- std::vector<DimParams> dimParams; //!< Dimension parameters for setting up Bias
-};
-
-//! Helper function to set up the C-style AWH parameters for the test
-static AwhFepLambdaStateTestParameters getAwhFepLambdaTestParameters(int eawhgrowth, int eawhpotential)
-{
- AwhFepLambdaStateTestParameters params;
-
- params.beta = 0.4;
-
- AwhDimParams& awhDimParams = params.awhDimParams;
-
- awhDimParams.period = 0;
- // Correction for removal of GaussianGeometryFactor/2 in histogram size
- awhDimParams.diffusion = 1e-4 / (0.12927243028700 * 2);
- awhDimParams.origin = 0;
- awhDimParams.end = numLambdaStates - 1;
- awhDimParams.coordValueInit = awhDimParams.origin;
- awhDimParams.coverDiameter = 0;
- awhDimParams.eCoordProvider = eawhcoordproviderFREE_ENERGY_LAMBDA;
-
- AwhBiasParams& awhBiasParams = params.awhBiasParams;
-
- awhBiasParams.ndim = 1;
- awhBiasParams.dimParams = &awhDimParams;
- awhBiasParams.eTarget = eawhtargetCONSTANT;
- awhBiasParams.targetBetaScaling = 0;
- awhBiasParams.targetCutoff = 0;
- awhBiasParams.eGrowth = eawhgrowth;
- awhBiasParams.bUserData = FALSE;
- awhBiasParams.errorInitial = 1.0 / params.beta;
- awhBiasParams.shareGroup = 0;
- awhBiasParams.equilibrateHistogram = FALSE;
-
- int64_t seed = 93471803;
-
- params.dimParams.push_back(DimParams::fepLambdaDimParams(numLambdaStates, params.beta));
-
- AwhParams& awhParams = params.awhParams;
-
- awhParams.numBias = 1;
- awhParams.awhBiasParams = &awhBiasParams;
- awhParams.seed = seed;
- awhParams.nstOut = 0;
- awhParams.nstSampleCoord = 1;
- awhParams.numSamplesUpdateFreeEnergy = 10;
- awhParams.ePotential = eawhpotential;
- awhParams.shareBiasMultisim = FALSE;
-
- return params;
-}
//! Convenience typedef: growth type enum, potential type enum, disable update skips
-typedef std::tuple<int, int, BiasParams::DisableUpdateSkips> BiasTestParameters;
+typedef std::tuple<AwhHistogramGrowthType, AwhPotentialType, BiasParams::DisableUpdateSkips> BiasTestParameters;
/*! \brief Test fixture for testing Bias updates
*/
class BiasFepLambdaStateTest : public ::testing::TestWithParam<BiasTestParameters>
{
+private:
+ //! Storage for test parameters.
+ std::unique_ptr<AwhTestParameters> params_;
+
public:
//! Random seed for AWH MC sampling
int64_t seed_;
* and disableUpdateSkips do not affect the point state.
* But the reference data will also ensure this.
*/
- int eawhgrowth;
- int eawhpotential;
+ AwhHistogramGrowthType eawhgrowth;
+ AwhPotentialType eawhpotential;
BiasParams::DisableUpdateSkips disableUpdateSkips;
std::tie(eawhgrowth, eawhpotential, disableUpdateSkips) = GetParam();
/* Set up a basic AWH setup with a single, 1D bias with parameters
* such that we can measure the effects of different parameters.
*/
- const AwhFepLambdaStateTestParameters params =
- getAwhFepLambdaTestParameters(eawhgrowth, eawhpotential);
-
- seed_ = params.awhParams.seed;
+ constexpr AwhCoordinateProviderType coordinateProvider = AwhCoordinateProviderType::FreeEnergyLambda;
+ constexpr int coordIndex = 0;
+ constexpr double origin = 0;
+ constexpr double end = c_numLambdaStates - 1;
+ constexpr double period = 0;
+ // Correction for removal of GaussianGeometryFactor/2 in histogram size
+ constexpr double diffusion = 1e-4 / (0.12927243028700 * 2);
+ const auto awhDimBuffer =
+ awhDimParamSerialized(coordinateProvider, coordIndex, origin, end, period, diffusion);
+ const auto awhDimArrayRef = gmx::arrayRefFromArray(&awhDimBuffer, 1);
+ params_ = std::make_unique<AwhTestParameters>(getAwhTestParameters(
+ eawhgrowth, eawhpotential, awhDimArrayRef, false, 0.4, true, 1.0, c_numLambdaStates));
+
+ seed_ = params_->awhParams.seed();
double mdTimeStep = 0.1;
- bias_ = std::make_unique<Bias>(-1, params.awhParams, params.awhBiasParams, params.dimParams,
- params.beta, mdTimeStep, 1, "", Bias::ThisRankWillDoIO::No,
+ bias_ = std::make_unique<Bias>(-1,
+ params_->awhParams,
+ params_->awhParams.awhBiasParams()[0],
+ params_->dimParams,
+ params_->beta,
+ mdTimeStep,
+ 1,
+ "",
+ Bias::ThisRankWillDoIO::No,
disableUpdateSkips);
}
};
int nSteps = 501;
/* Some energies to use as base values (to which some noise is added later on). */
- std::vector<double> neighborLambdaEnergies(numLambdaStates);
- std::vector<double> neighborLambdaDhdl(numLambdaStates);
+ std::vector<double> neighborLambdaEnergies(c_numLambdaStates);
+ std::vector<double> neighborLambdaDhdl(c_numLambdaStates);
const double magnitude = 12.0;
- for (int i = 0; i < numLambdaStates; i++)
+ for (int i = 0; i < c_numLambdaStates; i++)
{
neighborLambdaEnergies[i] = magnitude * std::sin(i * 0.1);
neighborLambdaDhdl[i] = magnitude * std::cos(i * 0.1);
int umbrellaGridpointIndex = bias.state().coordState().umbrellaGridpoint();
awh_dvec coordValue = { bias.getGridCoordValue(umbrellaGridpointIndex)[0], 0, 0, 0 };
double potential = 0;
- gmx::ArrayRef<const double> biasForce = bias.calcForceAndUpdateBias(
- coordValue, neighborLambdaEnergies, neighborLambdaDhdl, &potential, &potentialJump,
- nullptr, nullptr, step * mdTimeStep, step, seed_, nullptr);
+ gmx::ArrayRef<const double> biasForce = bias.calcForceAndUpdateBias(coordValue,
+ neighborLambdaEnergies,
+ neighborLambdaDhdl,
+ &potential,
+ &potentialJump,
+ nullptr,
+ nullptr,
+ step * mdTimeStep,
+ step,
+ seed_,
+ nullptr);
force.push_back(biasForce[0]);
pot.push_back(potential);
}
std::vector<double> pointBias, logPmfsum;
- for (auto& point : bias.state().points())
+ for (const auto& point : bias.state().points())
{
pointBias.push_back(point.bias());
logPmfsum.push_back(point.logPmfSum());
*/
INSTANTIATE_TEST_CASE_P(WithParameters,
BiasFepLambdaStateTest,
- ::testing::Combine(::testing::Values(eawhgrowthLINEAR, eawhgrowthEXP_LINEAR),
- ::testing::Values(eawhpotentialUMBRELLA),
+ ::testing::Combine(::testing::Values(AwhHistogramGrowthType::Linear,
+ AwhHistogramGrowthType::ExponentialLinear),
+ ::testing::Values(AwhPotentialType::Umbrella),
::testing::Values(BiasParams::DisableUpdateSkips::yes,
BiasParams::DisableUpdateSkips::no)));
// Test that we detect coverings and exit the initial stage at the correct step
TEST(BiasFepLambdaStateTest, DetectsCovering)
{
- const AwhFepLambdaStateTestParameters params =
- getAwhFepLambdaTestParameters(eawhgrowthEXP_LINEAR, eawhpotentialUMBRELLA);
+ constexpr AwhCoordinateProviderType coordinateProvider = AwhCoordinateProviderType::FreeEnergyLambda;
+ constexpr int coordIndex = 0;
+ constexpr double origin = 0;
+ constexpr double end = c_numLambdaStates - 1;
+ constexpr double period = 0;
+ constexpr double diffusion = 1e-4 / (0.12927243028700 * 2);
+ auto awhDimBuffer =
+ awhDimParamSerialized(coordinateProvider, coordIndex, origin, end, period, diffusion);
+ auto awhDimArrayRef = gmx::arrayRefFromArray(&awhDimBuffer, 1);
+ const AwhTestParameters params(getAwhTestParameters(AwhHistogramGrowthType::ExponentialLinear,
+ AwhPotentialType::Umbrella,
+ awhDimArrayRef,
+ false,
+ 0.4,
+ true,
+ 1.0,
+ c_numLambdaStates));
const double mdTimeStep = 0.1;
- Bias bias(-1, params.awhParams, params.awhBiasParams, params.dimParams, params.beta, mdTimeStep,
- 1, "", Bias::ThisRankWillDoIO::No);
+ Bias bias(-1,
+ params.awhParams,
+ params.awhParams.awhBiasParams()[0],
+ params.dimParams,
+ params.beta,
+ mdTimeStep,
+ 1,
+ "",
+ Bias::ThisRankWillDoIO::No);
const int64_t exitStepRef = 320;
bool inInitialStage = bias.state().inInitialStage();
/* Some energies to use as base values (to which some noise is added later on). */
- std::vector<double> neighborLambdaEnergies(numLambdaStates);
- std::vector<double> neighborLambdaDhdl(numLambdaStates);
+ std::vector<double> neighborLambdaEnergies(c_numLambdaStates);
+ std::vector<double> neighborLambdaDhdl(c_numLambdaStates);
const double magnitude = 12.0;
- for (int i = 0; i < numLambdaStates; i++)
+ for (int i = 0; i < c_numLambdaStates; i++)
{
neighborLambdaEnergies[i] = magnitude * std::sin(i * 0.1);
neighborLambdaDhdl[i] = magnitude * std::cos(i * 0.1);
double potential = 0;
double potentialJump = 0;
- bias.calcForceAndUpdateBias(coordValue, neighborLambdaEnergies, neighborLambdaDhdl,
- &potential, &potentialJump, nullptr, nullptr, step, step,
- params.awhParams.seed, nullptr);
+ bias.calcForceAndUpdateBias(coordValue,
+ neighborLambdaEnergies,
+ neighborLambdaDhdl,
+ &potential,
+ &potentialJump,
+ nullptr,
+ nullptr,
+ step,
+ step,
+ params.awhParams.seed(),
+ nullptr);
inInitialStage = bias.state().inInitialStage();
if (!inInitialStage)
biod.pnpi.spb.ru / alexxy / gromacs.git / blobdiff