Make EnergyOutput::addDataAtEnergyStep independent of t_state

EnergyOutput::addDataAtEnergyStep took a pointer to a full t_state,
although it only uses a small part of the object. This makes the
dependencies explicit. This is a necessary step to make EnergyData
independent of t_state.

Refs #3419

diff --git a/src/gromacs/mdlib/energyoutput.cpp b/src/gromacs/mdlib/energyoutput.cpp
index 120591f0a74af11cb7e8bbd97ea0be92f50729ec..f658ea377f5cfe768f30c019ea30b484b2d34274 100644 (file)
--- a/src/gromacs/mdlib/energyoutput.cpp
+++ b/src/gromacs/mdlib/energyoutput.cpp
@@ -855,10 +855,11 @@ void EnergyOutput::addDataAtEnergyStep(bool                    bDoDHDL,
                                        double                  time,
                                        real                    tmass,
                                        const gmx_enerdata_t*   enerd,
-                                       const t_state*          state,
                                        const t_lambda*         fep,
                                        const t_expanded*       expand,
                                        const matrix            box,
+                                       PTCouplingArrays        ptCouplingArrays,
+                                       int                     fep_state,
                                        const tensor            svir,
                                        const tensor            fvir,
                                        const tensor            vir,
@@ -936,12 +937,12 @@ void EnergyOutput::addDataAtEnergyStep(bool                    bDoDHDL,
     }
     if (epc_ == epcPARRINELLORAHMAN || epc_ == epcMTTK)
     {
-        tmp6[0] = state->boxv[XX][XX];
-        tmp6[1] = state->boxv[YY][YY];
-        tmp6[2] = state->boxv[ZZ][ZZ];
-        tmp6[3] = state->boxv[YY][XX];
-        tmp6[4] = state->boxv[ZZ][XX];
-        tmp6[5] = state->boxv[ZZ][YY];
+        tmp6[0] = ptCouplingArrays.boxv[XX][XX];
+        tmp6[1] = ptCouplingArrays.boxv[YY][YY];
+        tmp6[2] = ptCouplingArrays.boxv[ZZ][ZZ];
+        tmp6[3] = ptCouplingArrays.boxv[YY][XX];
+        tmp6[4] = ptCouplingArrays.boxv[ZZ][XX];
+        tmp6[5] = ptCouplingArrays.boxv[ZZ][YY];
         add_ebin(ebin_, ipc_, bTricl_ ? 6 : 3, tmp6, bSum);
     }
     if (bMu_)
@@ -1000,8 +1001,8 @@ void EnergyOutput::addDataAtEnergyStep(bool                    bDoDHDL,
                         for (j = 0; j < nNHC_; j++)
                         {
                             k                 = i * nNHC_ + j;
-                            tmp_r_[2 * k]     = state->nosehoover_xi[k];
-                            tmp_r_[2 * k + 1] = state->nosehoover_vxi[k];
+                            tmp_r_[2 * k]     = ptCouplingArrays.nosehoover_xi[k];
+                            tmp_r_[2 * k + 1] = ptCouplingArrays.nosehoover_vxi[k];
                         }
                     }
                     add_ebin(ebin_, itc_, mde_n_, tmp_r_, bSum);
@@ -1013,8 +1014,8 @@ void EnergyOutput::addDataAtEnergyStep(bool                    bDoDHDL,
                             for (j = 0; j < nNHC_; j++)
                             {
                                 k                 = i * nNHC_ + j;
-                                tmp_r_[2 * k]     = state->nhpres_xi[k];
-                                tmp_r_[2 * k + 1] = state->nhpres_vxi[k];
+                                tmp_r_[2 * k]     = ptCouplingArrays.nhpres_xi[k];
+                                tmp_r_[2 * k + 1] = ptCouplingArrays.nhpres_vxi[k];
                             }
                         }
                         add_ebin(ebin_, itcb_, mdeb_n_, tmp_r_, bSum);
@@ -1024,8 +1025,8 @@ void EnergyOutput::addDataAtEnergyStep(bool                    bDoDHDL,
                 {
                     for (int i = 0; (i < nTC_); i++)
                     {
-                        tmp_r_[2 * i]     = state->nosehoover_xi[i];
-                        tmp_r_[2 * i + 1] = state->nosehoover_vxi[i];
+                        tmp_r_[2 * i]     = ptCouplingArrays.nosehoover_xi[i];
+                        tmp_r_[2 * i + 1] = ptCouplingArrays.nosehoover_vxi[i];
                     }
                     add_ebin(ebin_, itc_, mde_n_, tmp_r_, bSum);
                 }
@@ -1062,14 +1063,14 @@ void EnergyOutput::addDataAtEnergyStep(bool                    bDoDHDL,
             dE_[i] = foreignTerms.deltaH(i);
             if (numTemperatures_ > 0)
             {
-                GMX_RELEASE_ASSERT(numTemperatures_ > state->fep_state,
+                GMX_RELEASE_ASSERT(numTemperatures_ > fep_state,
                                    "Number of lambdas in state is bigger then in input record");
                 GMX_RELEASE_ASSERT(
                         numTemperatures_ >= foreignTerms.numLambdas(),
                         "Number of lambdas in energy data is bigger then in input record");
                 /* MRS: is this right, given the way we have defined the exchange probabilities? */
                 /* is this even useful to have at all? */
-                dE_[i] += (temperatures_[i] / temperatures_[state->fep_state] - 1.0) * enerd->term[F_EKIN];
+                dE_[i] += (temperatures_[i] / temperatures_[fep_state] - 1.0) * enerd->term[F_EKIN];
             }
         }
 
@@ -1081,7 +1082,7 @@ void EnergyOutput::addDataAtEnergyStep(bool                    bDoDHDL,
             /* print the current state if we are doing expanded ensemble */
             if (expand->elmcmove > elmcmoveNO)
             {
-                fprintf(fp_dhdl_, " %4d", state->fep_state);
+                fprintf(fp_dhdl_, " %4d", fep_state);
             }
             /* total energy (for if the temperature changes */
 
@@ -1137,7 +1138,7 @@ void EnergyOutput::addDataAtEnergyStep(bool                    bDoDHDL,
             }
             store_energy = enerd->term[F_ETOT];
             /* store_dh is dE */
-            mde_delta_h_coll_add_dh(dhc_, static_cast<double>(state->fep_state), store_energy, pv,
+            mde_delta_h_coll_add_dh(dhc_, static_cast<double>(fep_state), store_energy, pv,
                                     store_dhdl, dE_ + fep->lambda_start_n, time);
         }
     }

diff --git a/src/gromacs/mdlib/energyoutput.h b/src/gromacs/mdlib/energyoutput.h
index 6bfff750de79b9613eac34e42a72f545bba1e4bc..18bfa98de9e12900150a7bb191014403ba5459cc 100644 (file)
--- a/src/gromacs/mdlib/energyoutput.h
+++ b/src/gromacs/mdlib/energyoutput.h
@@ -98,6 +98,23 @@ enum
 namespace gmx
 {
 
+/*! \internal
+ * \brief Arrays connected to Pressure and Temperature coupling
+ */
+struct PTCouplingArrays
+{
+    //! Box velocities for Parrinello-Rahman P-coupling.
+    const rvec* boxv;
+    //! Nose-Hoover coordinates.
+    ArrayRef<const double> nosehoover_xi;
+    //! Nose-Hoover velocities.
+    ArrayRef<const double> nosehoover_vxi;
+    //! Pressure Nose-Hoover coordinates.
+    ArrayRef<const double> nhpres_xi;
+    //! Pressure Nose-Hoover velocities.
+    ArrayRef<const double> nhpres_vxi;
+};
+
 /* Energy output class
  *
  * This is the collection of energy averages collected during mdrun, and to
@@ -136,32 +153,34 @@ public:
      *
      * Called every step on which the thermodynamic values are evaluated.
      *
-     * \param[in] bDoDHDL  Whether the FEP is enabled.
-     * \param[in] bSum     If this stepshould be recorded to compute sums and averaes.
-     * \param[in] time     Current simulation time.
-     * \param[in] tmass    Total mass
-     * \param[in] enerd    Energy data object.
-     * \param[in] state    System state.
-     * \param[in] fep      FEP data.
-     * \param[in] expand   Expanded ensemble (for FEP).
-     * \param[in] lastbox  PBC data.
-     * \param[in] svir     Constraint virial.
-     * \param[in] fvir     Force virial.
-     * \param[in] vir      Total virial.
-     * \param[in] pres     Pressure.
-     * \param[in] ekind    Kinetic energy data.
-     * \param[in] mu_tot   Total dipole.
-     * \param[in] constr   Constraints object to get RMSD from (for LINCS only).
+     * \param[in] bDoDHDL           Whether the FEP is enabled.
+     * \param[in] bSum              If this stepshould be recorded to compute sums and averages.
+     * \param[in] time              Current simulation time.
+     * \param[in] tmass             Total mass
+     * \param[in] enerd             Energy data object.
+     * \param[in] fep               FEP data.
+     * \param[in] expand            Expanded ensemble (for FEP).
+     * \param[in] lastbox           PBC data.
+     * \param[in] ptCouplingArrays  Arrays connected to pressure and temperature coupling.
+     * \param[in] fep_state         The current alchemical state we are in.
+     * \param[in] svir              Constraint virial.
+     * \param[in] fvir              Force virial.
+     * \param[in] vir               Total virial.
+     * \param[in] pres              Pressure.
+     * \param[in] ekind             Kinetic energy data.
+     * \param[in] mu_tot            Total dipole.
+     * \param[in] constr            Constraints object to get RMSD from (for LINCS only).
      */
     void addDataAtEnergyStep(bool                    bDoDHDL,
                              bool                    bSum,
                              double                  time,
                              real                    tmass,
                              const gmx_enerdata_t*   enerd,
-                             const t_state*          state,
                              const t_lambda*         fep,
                              const t_expanded*       expand,
                              const matrix            lastbox,
+                             PTCouplingArrays        ptCouplingArrays,
+                             int                     fep_state,
                              const tensor            svir,
                              const tensor            fvir,
                              const tensor            vir,

diff --git a/src/gromacs/mdlib/tests/energyoutput.cpp b/src/gromacs/mdlib/tests/energyoutput.cpp
index 422213bb854ea75e7e1e840a90672e915d20f725..bb4ce4ce00b914871d5dfdc9eabddb5141f793e0 100644 (file)
--- a/src/gromacs/mdlib/tests/energyoutput.cpp
+++ b/src/gromacs/mdlib/tests/energyoutput.cpp
@@ -634,9 +634,12 @@ TEST_P(EnergyOutputTest, CheckOutput)
     for (int frame = 0; frame < parameters.numFrames; frame++)
     {
         setStepData(&testValue);
-        energyOutput->addDataAtEnergyStep(false, true, time_, tmass_, enerdata_.get(), &state_, nullptr,
-                                          nullptr, box_, constraintsVirial_, forceVirial_, totalVirial_,
-                                          pressure_, &ekindata_, muTotal_, constraints_.get());
+        energyOutput->addDataAtEnergyStep(
+                false, true, time_, tmass_, enerdata_.get(), nullptr, nullptr, box_,
+                PTCouplingArrays({ state_.boxv, state_.nosehoover_xi, state_.nosehoover_vxi,
+                                   state_.nhpres_xi, state_.nhpres_vxi }),
+                state_.fep_state, constraintsVirial_, forceVirial_, totalVirial_, pressure_,
+                &ekindata_, muTotal_, constraints_.get());
 
         energyOutput->printAnnealingTemperatures(log_, &mtop_.groups, &inputrec_.opts);
         energyOutput->printStepToEnergyFile(energyFile_, true, false, false, log_, 100 * frame,

diff --git a/src/gromacs/mdrun/md.cpp b/src/gromacs/mdrun/md.cpp
index 3d1c95a411a8b2e951db9efa242174ce97623f8c..a97e0eeb75703f174c08ff5c2f848494ff2ad7a5 100644 (file)
--- a/src/gromacs/mdrun/md.cpp
+++ b/src/gromacs/mdrun/md.cpp
@@ -1509,9 +1509,12 @@ void gmx::LegacySimulator::do_md()
             }
             if (bCalcEner)
             {
-                energyOutput.addDataAtEnergyStep(bDoDHDL, bCalcEnerStep, t, mdatoms->tmass, enerd, state,
-                                                 ir->fepvals, ir->expandedvals, lastbox, shake_vir,
-                                                 force_vir, total_vir, pres, ekind, mu_tot, constr);
+                energyOutput.addDataAtEnergyStep(
+                        bDoDHDL, bCalcEnerStep, t, mdatoms->tmass, enerd, ir->fepvals,
+                        ir->expandedvals, lastbox,
+                        PTCouplingArrays{ state->boxv, state->nosehoover_xi, state->nosehoover_vxi,
+                                          state->nhpres_xi, state->nhpres_vxi },
+                        state->fep_state, shake_vir, force_vir, total_vir, pres, ekind, mu_tot, constr);
             }
             else
             {

diff --git a/src/gromacs/mdrun/mimic.cpp b/src/gromacs/mdrun/mimic.cpp
index e6c924c2eb3b8f70a1e79400f556a750da40a533..edd2814e5897fa0d2e91a01ee0be2eb2408c6eea 100644 (file)
--- a/src/gromacs/mdrun/mimic.cpp
+++ b/src/gromacs/mdrun/mimic.cpp
@@ -515,10 +515,12 @@ void gmx::LegacySimulator::do_mimic()
         if (MASTER(cr))
         {
             const bool bCalcEnerStep = true;
-            energyOutput.addDataAtEnergyStep(doFreeEnergyPerturbation, bCalcEnerStep, t,
-                                             mdatoms->tmass, enerd, state, ir->fepvals,
-                                             ir->expandedvals, state->box, shake_vir, force_vir,
-                                             total_vir, pres, ekind, mu_tot, constr);
+            energyOutput.addDataAtEnergyStep(
+                    doFreeEnergyPerturbation, bCalcEnerStep, t, mdatoms->tmass, enerd, ir->fepvals,
+                    ir->expandedvals, state->box,
+                    PTCouplingArrays({ state->boxv, state->nosehoover_xi, state->nosehoover_vxi,
+                                       state->nhpres_xi, state->nhpres_vxi }),
+                    state->fep_state, shake_vir, force_vir, total_vir, pres, ekind, mu_tot, constr);
 
             const bool do_ene = true;
             const bool do_log = true;

diff --git a/src/gromacs/mdrun/minimize.cpp b/src/gromacs/mdrun/minimize.cpp
index 95815d24b7af85284d152793e38ec01272a6286e..9c2b75760f55f502de822d36e24ff4eca47fa5e1 100644 (file)
--- a/src/gromacs/mdrun/minimize.cpp
+++ b/src/gromacs/mdrun/minimize.cpp
@@ -1126,8 +1126,8 @@ void LegacySimulator::do_cg()
         /* Copy stuff to the energy bin for easy printing etc. */
         matrix nullBox = {};
         energyOutput.addDataAtEnergyStep(false, false, static_cast<double>(step), mdatoms->tmass,
-                                         enerd, nullptr, nullptr, nullptr, nullBox, nullptr,
-                                         nullptr, vir, pres, nullptr, mu_tot, constr);
+                                         enerd, nullptr, nullptr, nullBox, PTCouplingArrays(), 0,
+                                         nullptr, nullptr, vir, pres, nullptr, mu_tot, constr);
 
         EnergyOutput::printHeader(fplog, step, step);
         energyOutput.printStepToEnergyFile(mdoutf_get_fp_ene(outf), TRUE, FALSE, FALSE, fplog, step,
@@ -1535,8 +1535,8 @@ void LegacySimulator::do_cg()
             /* Store the new (lower) energies */
             matrix nullBox = {};
             energyOutput.addDataAtEnergyStep(false, false, static_cast<double>(step), mdatoms->tmass,
-                                             enerd, nullptr, nullptr, nullptr, nullBox, nullptr,
-                                             nullptr, vir, pres, nullptr, mu_tot, constr);
+                                             enerd, nullptr, nullptr, nullBox, PTCouplingArrays(), 0,
+                                             nullptr, nullptr, vir, pres, nullptr, mu_tot, constr);
 
             do_log = do_per_step(step, inputrec->nstlog);
             do_ene = do_per_step(step, inputrec->nstenergy);
@@ -1776,8 +1776,8 @@ void LegacySimulator::do_lbfgs()
         /* Copy stuff to the energy bin for easy printing etc. */
         matrix nullBox = {};
         energyOutput.addDataAtEnergyStep(false, false, static_cast<double>(step), mdatoms->tmass,
-                                         enerd, nullptr, nullptr, nullptr, nullBox, nullptr,
-                                         nullptr, vir, pres, nullptr, mu_tot, constr);
+                                         enerd, nullptr, nullptr, nullBox, PTCouplingArrays(), 0,
+                                         nullptr, nullptr, vir, pres, nullptr, mu_tot, constr);
 
         EnergyOutput::printHeader(fplog, step, step);
         energyOutput.printStepToEnergyFile(mdoutf_get_fp_ene(outf), TRUE, FALSE, FALSE, fplog, step,
@@ -2261,8 +2261,8 @@ void LegacySimulator::do_lbfgs()
             /* Store the new (lower) energies */
             matrix nullBox = {};
             energyOutput.addDataAtEnergyStep(false, false, static_cast<double>(step), mdatoms->tmass,
-                                             enerd, nullptr, nullptr, nullptr, nullBox, nullptr,
-                                             nullptr, vir, pres, nullptr, mu_tot, constr);
+                                             enerd, nullptr, nullptr, nullBox, PTCouplingArrays(), 0,
+                                             nullptr, nullptr, vir, pres, nullptr, mu_tot, constr);
 
             do_log = do_per_step(step, inputrec->nstlog);
             do_ene = do_per_step(step, inputrec->nstenergy);
@@ -2473,9 +2473,10 @@ void LegacySimulator::do_steep()
             {
                 /* Store the new (lower) energies  */
                 matrix nullBox = {};
-                energyOutput.addDataAtEnergyStep(false, false, static_cast<double>(count), mdatoms->tmass,
-                                                 enerd, nullptr, nullptr, nullptr, nullBox, nullptr,
-                                                 nullptr, vir, pres, nullptr, mu_tot, constr);
+                energyOutput.addDataAtEnergyStep(false, false, static_cast<double>(count),
+                                                 mdatoms->tmass, enerd, nullptr, nullptr, nullBox,
+                                                 PTCouplingArrays(), 0, nullptr, nullptr, vir, pres,
+                                                 nullptr, mu_tot, constr);
 
                 imdSession->fillEnergyRecord(count, TRUE);
 

diff --git a/src/gromacs/mdrun/rerun.cpp b/src/gromacs/mdrun/rerun.cpp
index 12dcb98b7b49845d389be1122d11e95f6a7df52f..8230157ded92d002fd10ea995d2f99484102e7b0 100644 (file)
--- a/src/gromacs/mdrun/rerun.cpp
+++ b/src/gromacs/mdrun/rerun.cpp
@@ -594,10 +594,12 @@ void gmx::LegacySimulator::do_rerun()
         if (MASTER(cr))
         {
             const bool bCalcEnerStep = true;
-            energyOutput.addDataAtEnergyStep(doFreeEnergyPerturbation, bCalcEnerStep, t,
-                                             mdatoms->tmass, enerd, state, ir->fepvals,
-                                             ir->expandedvals, state->box, shake_vir, force_vir,
-                                             total_vir, pres, ekind, mu_tot, constr);
+            energyOutput.addDataAtEnergyStep(
+                    doFreeEnergyPerturbation, bCalcEnerStep, t, mdatoms->tmass, enerd, ir->fepvals,
+                    ir->expandedvals, state->box,
+                    PTCouplingArrays({ state->boxv, state->nosehoover_xi, state->nosehoover_vxi,
+                                       state->nhpres_xi, state->nhpres_vxi }),
+                    state->fep_state, shake_vir, force_vir, total_vir, pres, ekind, mu_tot, constr);
 
             const bool do_ene = true;
             const bool do_log = true;

diff --git a/src/gromacs/modularsimulator/energydata.cpp b/src/gromacs/modularsimulator/energydata.cpp
index 88383cc0f329621b94afa8fe5427318d1e44d158..d8ec87cb3d536185648faebaa5a1a79a1e370164 100644 (file)
--- a/src/gromacs/modularsimulator/energydata.cpp
+++ b/src/gromacs/modularsimulator/energydata.cpp
@@ -241,7 +241,6 @@ void EnergyData::doStep(Time time, bool isEnergyCalculationStep, bool isFreeEner
     {
         accumulateKineticLambdaComponents(enerd_, freeEnergyPerturbationData_->constLambdaView(),
                                           *inputrec_->fepvals);
-        dummyLegacyState_.fep_state = freeEnergyPerturbationData_->currentFEPState();
     }
     if (parrinelloRahmanBarostat_)
     {
@@ -253,11 +252,17 @@ void EnergyData::doStep(Time time, bool isEnergyCalculationStep, bool isFreeEner
         enerd_->term[F_ECONSERVED] =
                 enerd_->term[F_ETOT] + NPT_energy(inputrec_, &dummyLegacyState_, nullptr);
     }
-    energyOutput_->addDataAtEnergyStep(isFreeEnergyCalculationStep, isEnergyCalculationStep, time,
-                                       mdAtoms_->mdatoms()->tmass, enerd_, &dummyLegacyState_,
-                                       inputrec_->fepvals, inputrec_->expandedvals,
-                                       statePropagatorData_->constPreviousBox(), shakeVirial_,
-                                       forceVirial_, totalVirial_, pressure_, ekind_, muTot_, constr_);
+    matrix nullMatrix = {};
+    energyOutput_->addDataAtEnergyStep(
+            isFreeEnergyCalculationStep, isEnergyCalculationStep, time, mdAtoms_->mdatoms()->tmass, enerd_,
+            inputrec_->fepvals, inputrec_->expandedvals, statePropagatorData_->constPreviousBox(),
+            PTCouplingArrays({ parrinelloRahmanBarostat_ ? parrinelloRahmanBarostat_->boxVelocities() : nullMatrix,
+                               {},
+                               {},
+                               {},
+                               {} }),
+            freeEnergyPerturbationData_ ? freeEnergyPerturbationData_->currentFEPState() : 0,
+            shakeVirial_, forceVirial_, totalVirial_, pressure_, ekind_, muTot_, constr_);
 }
 
 void EnergyData::write(gmx_mdoutf* outf, Step step, Time time, bool writeTrajectory, bool writeLog)