2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 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.
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.
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.
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.
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.
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.
36 * \brief Defines the free energy perturbation element for the modular simulator
38 * \author Pascal Merz <pascal.merz@me.com>
39 * \ingroup module_modularsimulator
44 #include "freeenergyperturbationdata.h"
46 #include "gromacs/domdec/domdec_network.h"
47 #include "gromacs/mdlib/freeenergyparameters.h"
48 #include "gromacs/mdlib/md_support.h"
49 #include "gromacs/mdlib/mdatoms.h"
50 #include "gromacs/mdtypes/checkpointdata.h"
51 #include "gromacs/mdtypes/commrec.h"
52 #include "gromacs/mdtypes/inputrec.h"
53 #include "gromacs/mdtypes/mdatom.h"
54 #include "gromacs/mdtypes/state.h"
56 #include "modularsimulator.h"
57 #include "simulatoralgorithm.h"
62 FreeEnergyPerturbationData::FreeEnergyPerturbationData(FILE* fplog
, const t_inputrec
* inputrec
, MDAtoms
* mdAtoms
) :
63 element_(std::make_unique
<Element
>(this, inputrec
->fepvals
->delta_lambda
)),
71 // The legacy implementation only filled the lambda vector in state_global, which is only
72 // available on master. We have the lambda vector available everywhere, so we pass a `true`
73 // for isMaster on all ranks. See #3647.
74 initialize_lambdas(fplog_
, *inputrec_
, true, ¤tFEPState_
, lambda_
);
77 void FreeEnergyPerturbationData::Element::scheduleTask(Step step
,
79 const RegisterRunFunction
& registerRunFunction
)
83 registerRunFunction([this, step
]() { freeEnergyPerturbationData_
->updateLambdas(step
); });
87 void FreeEnergyPerturbationData::updateLambdas(Step step
)
89 // at beginning of step (if lambdas change...)
90 lambda_
= currentLambdas(step
, *(inputrec_
->fepvals
), currentFEPState_
);
94 ArrayRef
<real
> FreeEnergyPerturbationData::lambdaView()
99 ArrayRef
<const real
> FreeEnergyPerturbationData::constLambdaView()
104 int FreeEnergyPerturbationData::currentFEPState()
106 return currentFEPState_
;
109 void FreeEnergyPerturbationData::updateMDAtoms()
111 update_mdatoms(mdAtoms_
->mdatoms(), lambda_
[efptMASS
]);
117 * \brief Enum describing the contents FreeEnergyPerturbationData::Element writes to modular checkpoint
119 * When changing the checkpoint content, add a new element just above Count, and adjust the
120 * checkpoint functionality.
122 enum class CheckpointVersion
124 Base
, //!< First version of modular checkpointing
125 Count
//!< Number of entries. Add new versions right above this!
127 constexpr auto c_currentVersion
= CheckpointVersion(int(CheckpointVersion::Count
) - 1);
130 template<CheckpointDataOperation operation
>
131 void FreeEnergyPerturbationData::Element::doCheckpointData(CheckpointData
<operation
>* checkpointData
)
133 checkpointVersion(checkpointData
, "FreeEnergyPerturbationData version", c_currentVersion
);
135 checkpointData
->scalar("current FEP state", &freeEnergyPerturbationData_
->currentFEPState_
);
136 checkpointData
->arrayRef("lambda vector",
137 makeCheckpointArrayRef
<operation
>(freeEnergyPerturbationData_
->lambda_
));
140 void FreeEnergyPerturbationData::Element::saveCheckpointState(std::optional
<WriteCheckpointData
> checkpointData
,
145 doCheckpointData
<CheckpointDataOperation::Write
>(&checkpointData
.value());
149 void FreeEnergyPerturbationData::Element::restoreCheckpointState(std::optional
<ReadCheckpointData
> checkpointData
,
154 doCheckpointData
<CheckpointDataOperation::Read
>(&checkpointData
.value());
156 if (DOMAINDECOMP(cr
))
158 dd_bcast(cr
->dd
, sizeof(int), &freeEnergyPerturbationData_
->currentFEPState_
);
159 dd_bcast(cr
->dd
, freeEnergyPerturbationData_
->lambda_
.size() * sizeof(real
),
160 freeEnergyPerturbationData_
->lambda_
.data());
164 const std::string
& FreeEnergyPerturbationData::Element::clientID()
169 FreeEnergyPerturbationData::Element::Element(FreeEnergyPerturbationData
* freeEnergyPerturbationElement
,
170 double deltaLambda
) :
171 freeEnergyPerturbationData_(freeEnergyPerturbationElement
),
172 lambdasChange_(deltaLambda
!= 0)
176 void FreeEnergyPerturbationData::Element::elementSetup()
178 freeEnergyPerturbationData_
->updateMDAtoms();
181 FreeEnergyPerturbationData::Element
* FreeEnergyPerturbationData::element()
183 return element_
.get();
186 ISimulatorElement
* FreeEnergyPerturbationData::Element::getElementPointerImpl(
187 LegacySimulatorData gmx_unused
* legacySimulatorData
,
188 ModularSimulatorAlgorithmBuilderHelper gmx_unused
* builderHelper
,
189 StatePropagatorData gmx_unused
* statePropagatorData
,
190 EnergyData gmx_unused
* energyData
,
191 FreeEnergyPerturbationData
* freeEnergyPerturbationData
,
192 GlobalCommunicationHelper gmx_unused
* globalCommunicationHelper
)
194 return freeEnergyPerturbationData
->element();