Pass ref_t to initialize_lambdas
authorJoe Jordan <ejjordan12@gmail.com>
Sun, 11 Apr 2021 00:48:56 +0000 (00:48 +0000)
committerAndrey Alekseenko <al42and@gmail.com>
Sun, 11 Apr 2021 00:48:56 +0000 (00:48 +0000)
src/gromacs/mdrun/md.cpp
src/gromacs/mdrun/mimic.cpp
src/gromacs/mdrun/minimize.cpp
src/gromacs/mdrun/rerun.cpp
src/gromacs/mdtypes/state.cpp
src/gromacs/mdtypes/state.h
src/gromacs/modularsimulator/freeenergyperturbationdata.cpp
src/gromacs/modularsimulator/freeenergyperturbationdata.h
src/gromacs/modularsimulator/simulatoralgorithm.cpp

index 6e80825afef2ddc81dff8a3f7fda9f7d1f473587..a7487ff8cf79ef7f7dd245e223e74edfc8625f6e 100644 (file)
@@ -257,7 +257,8 @@ void gmx::LegacySimulator::do_md()
 
     int*                fep_state = MASTER(cr) ? &state_global->fep_state : nullptr;
     gmx::ArrayRef<real> lambda    = MASTER(cr) ? state_global->lambda : gmx::ArrayRef<real>();
-    initialize_lambdas(fplog, *ir, MASTER(cr), fep_state, lambda);
+    initialize_lambdas(
+            fplog, *ir, gmx::arrayRefFromArray(ir->opts.ref_t, ir->opts.ngtc), MASTER(cr), fep_state, lambda);
     Update upd(*ir, deform);
     bool   doSimulatedAnnealing = false;
     {
index ce23edc526ff30eb75d9c175c38209dd5ca70f1a..a434de0126b78830281683731e4fc81f90615ab6 100644 (file)
@@ -228,7 +228,12 @@ void gmx::LegacySimulator::do_mimic()
         nonConstGlobalTopology->intermolecularExclusionGroup = genQmmmIndices(top_global);
     }
 
-    initialize_lambdas(fplog, *ir, MASTER(cr), &state_global->fep_state, state_global->lambda);
+    initialize_lambdas(fplog,
+                       *ir,
+                       gmx::arrayRefFromArray(ir->opts.ref_t, ir->opts.ngtc),
+                       MASTER(cr),
+                       &state_global->fep_state,
+                       state_global->lambda);
 
     const bool        simulationsShareState = false;
     gmx_mdoutf*       outf                  = init_mdoutf(fplog,
index 2083308bf01eaa385f6f9abe57e3e258706b90b3..c7da246c7cf15394c9c086e0386c08c74d119f51 100644 (file)
@@ -394,7 +394,8 @@ static void init_em(FILE*                fplog,
     }
     int*                fep_state = MASTER(cr) ? &state_global->fep_state : nullptr;
     gmx::ArrayRef<real> lambda    = MASTER(cr) ? state_global->lambda : gmx::ArrayRef<real>();
-    initialize_lambdas(fplog, *ir, MASTER(cr), fep_state, lambda);
+    initialize_lambdas(
+            fplog, *ir, gmx::arrayRefFromArray(ir->opts.ref_t, ir->opts.ngtc), MASTER(cr), fep_state, lambda);
 
     if (ir->eI == IntegrationAlgorithm::NM)
     {
index d3202d64fdd8f791a0527ad93bc33aeb601c14b4..ed4ede713415fbcef6297b86b9f560ee0e67f053 100644 (file)
@@ -275,7 +275,8 @@ void gmx::LegacySimulator::do_rerun()
     }
     int*                fep_state = MASTER(cr) ? &state_global->fep_state : nullptr;
     gmx::ArrayRef<real> lambda    = MASTER(cr) ? state_global->lambda : gmx::ArrayRef<real>();
-    initialize_lambdas(fplog, *ir, MASTER(cr), fep_state, lambda);
+    initialize_lambdas(
+            fplog, *ir, gmx::arrayRefFromArray(ir->opts.ref_t, ir->opts.ngtc), MASTER(cr), fep_state, lambda);
     const bool        simulationsShareState = false;
     gmx_mdoutf*       outf                  = init_mdoutf(fplog,
                                    nfile,
index 79a98f08455eb70017c5147bf25d4ca22ed555b6..4722245ff8570a13f4ec7b342b11258929953689 100644 (file)
@@ -352,7 +352,12 @@ void printLambdaStateToLog(FILE* fplog, gmx::ArrayRef<const real> lambda, const
     }
 }
 
-void initialize_lambdas(FILE* fplog, const t_inputrec& ir, bool isMaster, int* fep_state, gmx::ArrayRef<real> lambda)
+void initialize_lambdas(FILE*               fplog,
+                        const t_inputrec&   ir,
+                        gmx::ArrayRef<real> ref_t,
+                        bool                isMaster,
+                        int*                fep_state,
+                        gmx::ArrayRef<real> lambda)
 {
     /* TODO: Clean up initialization of fep_state and lambda in
        t_state.  This function works, but could probably use a logic
@@ -391,11 +396,11 @@ void initialize_lambdas(FILE* fplog, const t_inputrec& ir, bool isMaster, int* f
     if (ir.bSimTemp)
     {
         /* need to rescale control temperatures to match current state */
-        for (int i = 0; i < ir.opts.ngtc; i++)
+        for (int i = 0; i < ref_t.ssize(); i++)
         {
-            if (ir.opts.ref_t[i] > 0)
+            if (ref_t[i] > 0)
             {
-                ir.opts.ref_t[i] = ir.simtempvals->temperatures[fep->init_fep_state];
+                ref_t[i] = ir.simtempvals->temperatures[fep->init_fep_state];
             }
         }
     }
index 8998b75abcc86683c37af1d99ea4a0b774788405..faa51d21a45082bfeb4553cfaeaf46ea12010352 100644 (file)
@@ -293,17 +293,6 @@ struct t_extmass
     tensor              Winvm; /* inverse pressure mass tensor, computed       */
 };
 
-
-typedef struct
-{
-    real    veta;
-    double  rscale;
-    double  vscale;
-    double  rvscale;
-    double  alpha;
-    double* vscale_nhc;
-} t_vetavars;
-
 #endif // DOXYGEN
 
 //! Resizes the T- and P-coupling state variables
@@ -376,6 +365,11 @@ void printLambdaStateToLog(FILE* fplog, gmx::ArrayRef<const real> lambda, bool i
  * and lambda on master rank.
  *
  * Reports the initial lambda state to the log file. */
-void initialize_lambdas(FILE* fplog, const t_inputrec& ir, bool isMaster, int* fep_state, gmx::ArrayRef<real> lambda);
+void initialize_lambdas(FILE*               fplog,
+                        const t_inputrec&   ir,
+                        gmx::ArrayRef<real> ref_t,
+                        bool                isMaster,
+                        int*                fep_state,
+                        gmx::ArrayRef<real> lambda);
 
 #endif
index d780e69216eb38303ddadc503d75b77a04bfa11c..6671c20d85783871dfdec2700849df9f27f2a353 100644 (file)
@@ -60,8 +60,8 @@
 namespace gmx
 {
 
-FreeEnergyPerturbationData::FreeEnergyPerturbationData(FILE* fplog, const t_inputrec* inputrec, MDAtoms* mdAtoms) :
-    element_(std::make_unique<Element>(this, inputrec->fepvals->delta_lambda)),
+FreeEnergyPerturbationData::FreeEnergyPerturbationData(FILE* fplog, const t_inputrec& inputrec, MDAtoms* mdAtoms) :
+    element_(std::make_unique<Element>(this, inputrec.fepvals->delta_lambda)),
     lambda_(),
     currentFEPState_(0),
     fplog_(fplog),
@@ -72,7 +72,12 @@ FreeEnergyPerturbationData::FreeEnergyPerturbationData(FILE* fplog, const t_inpu
     // The legacy implementation only filled the lambda vector in state_global, which is only
     // available on master. We have the lambda vector available everywhere, so we pass a `true`
     // for isMaster on all ranks. See #3647.
-    initialize_lambdas(fplog_, *inputrec_, true, &currentFEPState_, lambda_);
+    initialize_lambdas(fplog_,
+                       inputrec_,
+                       gmx::arrayRefFromArray(inputrec_.opts.ref_t, inputrec_.opts.ngtc),
+                       true,
+                       &currentFEPState_,
+                       lambda_);
 }
 
 void FreeEnergyPerturbationData::Element::scheduleTask(Step step,
@@ -88,7 +93,7 @@ void FreeEnergyPerturbationData::Element::scheduleTask(Step step,
 void FreeEnergyPerturbationData::updateLambdas(Step step)
 {
     // at beginning of step (if lambdas change...)
-    lambda_ = currentLambdas(step, *(inputrec_->fepvals), currentFEPState_);
+    lambda_ = currentLambdas(step, *(inputrec_.fepvals), currentFEPState_);
     updateMDAtoms();
 }
 
@@ -202,7 +207,7 @@ void FreeEnergyPerturbationData::readCheckpointToTrxFrame(t_trxframe* trxFrame,
 {
     if (readCheckpointData)
     {
-        FreeEnergyPerturbationData freeEnergyPerturbationData;
+        FreeEnergyPerturbationData freeEnergyPerturbationData(nullptr, t_inputrec(), nullptr);
         freeEnergyPerturbationData.doCheckpointData(&readCheckpointData.value());
         trxFrame->lambda = freeEnergyPerturbationData.lambda_[FreeEnergyPerturbationCouplingType::Fep];
         trxFrame->fep_state = freeEnergyPerturbationData.currentFEPState_;
index 7ef19d86f1e26fc21a8b4203dbb0081903356d97..98d9cfc80e5155a170783d0bf01288dfb9e583b7 100644 (file)
@@ -77,7 +77,7 @@ class FreeEnergyPerturbationData final
 {
 public:
     //! Constructor
-    FreeEnergyPerturbationData(FILE* fplog, const t_inputrec* inputrec, MDAtoms* mdAtoms);
+    FreeEnergyPerturbationData(FILE* fplog, const t_inputrec& inputrec, MDAtoms* mdAtoms);
 
     //! Get a view of the current lambda vector
     ArrayRef<real> lambdaView();
@@ -100,8 +100,6 @@ public:
     static const std::string& checkpointID();
 
 private:
-    //! Default constructor - only used internally
-    FreeEnergyPerturbationData() = default;
     //! Update the lambda values
     void updateLambdas(Step step);
     //! Helper function to read from / write to CheckpointData
@@ -119,7 +117,7 @@ private:
     //! Handles logging.
     FILE* fplog_;
     //! Contains user input mdp options.
-    const t_inputrec* inputrec_;
+    const t_inputrec& inputrec_;
     //! Atom parameters for this domain.
     MDAtoms* mdAtoms_;
 };
index 423c3e463c2ce2f4a29c3254c58e3c0d48fdf71f..5c1106f1fdebc87732e7b052aa359b629db7d8af 100644 (file)
@@ -407,7 +407,7 @@ ModularSimulatorAlgorithmBuilder::ModularSimulatorAlgorithmBuilder(
     if (legacySimulatorData->inputrec->efep != FreeEnergyPerturbationType::No)
     {
         freeEnergyPerturbationData_ = std::make_unique<FreeEnergyPerturbationData>(
-                legacySimulatorData->fplog, legacySimulatorData->inputrec, legacySimulatorData->mdAtoms);
+                legacySimulatorData->fplog, *legacySimulatorData->inputrec, legacySimulatorData->mdAtoms);
     }
 
     statePropagatorData_ = std::make_unique<StatePropagatorData>(