src/gromacs/modularsimulator/velocityscalingtemperaturecoupling.h

   1 /*
   2  * This file is part of the GROMACS molecular simulation package.
   3  *
   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.
   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 /*! \internal \file
  36  * \brief Declares a velocity-scaling temperature coupling element for
  37  * the modular simulator
  38  *
  39  * \author Pascal Merz <pascal.merz@me.com>
  40  * \ingroup module_modularsimulator
  41  *
  42  * This header is only used within the modular simulator module
  43  */
  44
  45 #ifndef GMX_MODULARSIMULATOR_VELOCITYSCALINGTEMPERATURECOUPLING_H
  46 #define GMX_MODULARSIMULATOR_VELOCITYSCALINGTEMPERATURECOUPLING_H
  47
  48 #include "gromacs/utility/arrayref.h"
  49
  50 #include "energydata.h"
  51 #include "modularsimulatorinterfaces.h"
  52 #include "propagator.h"
  53
  54 struct t_commrec;
  55
  56 namespace gmx
  57 {
  58 class ITemperatureCouplingImpl;
  59 class LegacySimulatorData;
  60 struct TemperatureCouplingData;
  61
  62 //! Enum describing whether the thermostat is using full or half step kinetic energy
  63 enum class UseFullStepKE
  64 {
  65     Yes,
  66     No,
  67     Count
  68 };
  69
  70 //! Enum describing whether the thermostat is reporting conserved energy from the previous step
  71 enum class ReportPreviousStepConservedEnergy
  72 {
  73     Yes,
  74     No,
  75     Count
  76 };
  77
  78 //! Typedef to match current use of ints as types.
  79 using TemperatureCouplingType = int;
  80
  81 /*! \internal
  82  * \ingroup module_modularsimulator
  83  * \brief Element implementing the a velocity-scaling thermostat
  84  *
  85  * This element takes a callback to the propagator and updates the velocity
  86  * scaling factor according to the internal temperature coupling implementation.
  87  *
  88  * Note that the concrete implementation is handled by the concrete
  89  * implementations of the ITemperatureCouplingImpl interface, while the element
  90  * handles the scheduling and interfacing with other elements.
  91  */
  92 class VelocityScalingTemperatureCoupling final : public ISimulatorElement, public ICheckpointHelperClient
  93 {
  94 public:
  95     //! Constructor
  96     VelocityScalingTemperatureCoupling(int                               nstcouple,
  97                                        int                               offset,
  98                                        UseFullStepKE                     useFullStepKE,
  99                                        ReportPreviousStepConservedEnergy reportPreviousConservedEnergy,
 100                                        int64_t                           seed,
 101                                        int                               numTemperatureGroups,
 102                                        double                            couplingTimeStep,
 103                                        const real*                       referenceTemperature,
 104                                        const real*                       couplingTime,
 105                                        const real*                       numDegreesOfFreedom,
 106                                        EnergyData*                       energyData,
 107                                        TemperatureCouplingType           couplingType);
 108
 109     /*! \brief Register run function for step / time
 110      *
 111      * \param step                 The step number
 112      * \param time                 The time
 113      * \param registerRunFunction  Function allowing to register a run function
 114      */
 115     void scheduleTask(Step step, Time time, const RegisterRunFunction& registerRunFunction) override;
 116
 117     //! Sanity check at setup time
 118     void elementSetup() override;
 119     //! No element teardown needed
 120     void elementTeardown() override {}
 121
 122     //! Contribution to the conserved energy (called by energy data)
 123     [[nodiscard]] real conservedEnergyContribution() const;
 124
 125     //! Connect this to propagator
 126     void connectWithPropagator(const PropagatorThermostatConnection& connectionData);
 127
 128     //! ICheckpointHelperClient write checkpoint implementation
 129     void saveCheckpointState(std::optional<WriteCheckpointData> checkpointData, const t_commrec* cr) override;
 130     //! ICheckpointHelperClient read checkpoint implementation
 131     void restoreCheckpointState(std::optional<ReadCheckpointData> checkpointData, const t_commrec* cr) override;
 132     //! ICheckpointHelperClient key implementation
 133     const std::string& clientID() override;
 134
 135     /*! \brief Factory method implementation
 136      *
 137      * \param legacySimulatorData  Pointer allowing access to simulator level data
 138      * \param builderHelper  ModularSimulatorAlgorithmBuilder helper object
 139      * \param statePropagatorData  Pointer to the \c StatePropagatorData object
 140      * \param energyData  Pointer to the \c EnergyData object
 141      * \param freeEnergyPerturbationData  Pointer to the \c FreeEnergyPerturbationData object
 142      * \param globalCommunicationHelper  Pointer to the \c GlobalCommunicationHelper object
 143      * \param offset  The step offset at which the thermostat is applied
 144      * \param useFullStepKE  Whether full step or half step KE is used
 145      * \param reportPreviousStepConservedEnergy  Report the previous or the current step conserved energy
 146      *
 147      * \return  Pointer to the element to be added. Element needs to have been stored using \c storeElement
 148      */
 149     static ISimulatorElement*
 150     getElementPointerImpl(LegacySimulatorData*                    legacySimulatorData,
 151                           ModularSimulatorAlgorithmBuilderHelper* builderHelper,
 152                           StatePropagatorData*                    statePropagatorData,
 153                           EnergyData*                             energyData,
 154                           FreeEnergyPerturbationData*             freeEnergyPerturbationData,
 155                           GlobalCommunicationHelper*              globalCommunicationHelper,
 156                           int                                     offset,
 157                           UseFullStepKE                           useFullStepKE,
 158                           ReportPreviousStepConservedEnergy reportPreviousStepConservedEnergy);
 159
 160 private:
 161     //! The frequency at which the thermostat is applied
 162     const int nstcouple_;
 163     //! If != 0, offset the step at which the thermostat is applied
 164     const int offset_;
 165     //! Whether we're using full step kinetic energy
 166     const UseFullStepKE useFullStepKE_;
 167     //! Whether we are reporting the conserved energy from the previous step
 168     const ReportPreviousStepConservedEnergy reportPreviousConservedEnergy_;
 169
 170     //! The number of temperature groups
 171     const int numTemperatureGroups_;
 172     //! The coupling time step - simulation time step x nstcouple_
 173     const double couplingTimeStep_;
 174     //! Coupling temperature per group
 175     const std::vector<real> referenceTemperature_;
 176     //! Coupling time per group
 177     const std::vector<real> couplingTime_;
 178     //! Number of degrees of freedom per group
 179     const std::vector<real> numDegreesOfFreedom_;
 180     //! Work exerted by thermostat per group
 181     std::vector<double> temperatureCouplingIntegral_;
 182     //! Work exerted by thermostat per group (backup from previous step)
 183     std::vector<double> temperatureCouplingIntegralPreviousStep_;
 184
 185     // TODO: Clarify relationship to data objects and find a more robust alternative to raw pointers (#3583)
 186     //! Pointer to the energy data (for ekindata)
 187     EnergyData* energyData_;
 188
 189     //! Callback to let propagator know that we updated lambda
 190     PropagatorCallback propagatorCallback_;
 191
 192     //! Set new lambda value (at T-coupling steps)
 193     void setLambda(Step step);
 194
 195     //! The temperature coupling implementation
 196     std::unique_ptr<ITemperatureCouplingImpl> temperatureCouplingImpl_;
 197
 198     //! CheckpointHelper identifier
 199     const std::string identifier_ = "VelocityScalingTemperatureCoupling";
 200     //! Helper function to read from / write to CheckpointData
 201     template<CheckpointDataOperation operation>
 202     void doCheckpointData(CheckpointData<operation>* checkpointData);
 203 };
 204
 205 } // namespace gmx
 206
 207 #endif // GMX_MODULARSIMULATOR_VELOCITYSCALINGTEMPERATURECOUPLING_H