From 102da8b05a3e85821c41558b9881742bedd76089 Mon Sep 17 00:00:00 2001 From: ejjordan Date: Fri, 9 Apr 2021 10:57:01 +0200 Subject: [PATCH] Constructor for gmx_ekindata_t Moved init_ekin_data from tgroup header/source to group and made it the constructor for gmx_ekindata_t. Also made gmx_ekindata_t a class since the number of threads is only needed by constructor and destructor and can thus be private. Further refactoring of members of gmx_ekindata_t would eliminate the need for the nthreads data member entirely. Other members we not made private to avoid a ripple of changes adding "()" to numerous lines in several files. --- src/gromacs/mdlib/coupling.h | 2 +- src/gromacs/mdlib/energyoutput.h | 2 +- src/gromacs/mdlib/md_support.h | 2 +- src/gromacs/mdlib/stat.h | 2 +- src/gromacs/mdlib/tests/energyoutput.cpp | 17 +++--- src/gromacs/mdlib/tests/leapfrogtestdata.cpp | 16 ++---- src/gromacs/mdlib/tgroup.cpp | 59 -------------------- src/gromacs/mdlib/tgroup.h | 15 +---- src/gromacs/mdlib/trajectory_writing.h | 2 +- src/gromacs/mdlib/update.h | 2 +- src/gromacs/mdlib/update_vv.h | 2 +- src/gromacs/mdrun/isimulator.h | 2 +- src/gromacs/mdrun/runner.cpp | 3 +- src/gromacs/mdrun/simulatorbuilder.h | 2 +- src/gromacs/mdtypes/group.cpp | 51 ++++++++++++++++- src/gromacs/mdtypes/group.h | 12 ++-- src/gromacs/modularsimulator/energydata.h | 2 +- 17 files changed, 82 insertions(+), 111 deletions(-) diff --git a/src/gromacs/mdlib/coupling.h b/src/gromacs/mdlib/coupling.h index 7f8d3c5a26..affd45bd44 100644 --- a/src/gromacs/mdlib/coupling.h +++ b/src/gromacs/mdlib/coupling.h @@ -45,7 +45,7 @@ #include "gromacs/utility/basedefinitions.h" #include "gromacs/utility/real.h" -struct gmx_ekindata_t; +class gmx_ekindata_t; struct gmx_enerdata_t; struct t_commrec; struct t_extmass; diff --git a/src/gromacs/mdlib/energyoutput.h b/src/gromacs/mdlib/energyoutput.h index 902ad5b811..ea917b6e02 100644 --- a/src/gromacs/mdlib/energyoutput.h +++ b/src/gromacs/mdlib/energyoutput.h @@ -53,7 +53,7 @@ class energyhistory_t; struct ener_file; -struct gmx_ekindata_t; +class gmx_ekindata_t; struct gmx_enerdata_t; struct SimulationGroups; struct gmx_mtop_t; diff --git a/src/gromacs/mdlib/md_support.h b/src/gromacs/mdlib/md_support.h index 1975239425..40e1437c13 100644 --- a/src/gromacs/mdlib/md_support.h +++ b/src/gromacs/mdlib/md_support.h @@ -42,7 +42,7 @@ #include "gromacs/mdtypes/md_enums.h" #include "gromacs/timing/wallcycle.h" -struct gmx_ekindata_t; +class gmx_ekindata_t; struct gmx_enerdata_t; struct gmx_global_stat; struct gmx_signalling_t; diff --git a/src/gromacs/mdlib/stat.h b/src/gromacs/mdlib/stat.h index 03687279e4..8001bf0d39 100644 --- a/src/gromacs/mdlib/stat.h +++ b/src/gromacs/mdlib/stat.h @@ -43,7 +43,7 @@ #include "gromacs/math/vectypes.h" -struct gmx_ekindata_t; +class gmx_ekindata_t; struct gmx_enerdata_t; struct t_vcm; struct t_inputrec; diff --git a/src/gromacs/mdlib/tests/energyoutput.cpp b/src/gromacs/mdlib/tests/energyoutput.cpp index fac276b9f6..78a4d92367 100644 --- a/src/gromacs/mdlib/tests/energyoutput.cpp +++ b/src/gromacs/mdlib/tests/energyoutput.cpp @@ -138,6 +138,9 @@ const EnergyOutputTestParameters parametersSets[] = { */ class EnergyOutputTest : public ::testing::TestWithParam { + int numTempCouplingGroups_ = 3; + real cosAccel_ = 1.0; + public: //! File manager TestFileManager fileManager_; @@ -193,6 +196,7 @@ public: TestReferenceChecker checker_; EnergyOutputTest() : + ekindata_(numTempCouplingGroups_, cosAccel_, 1), logFilename_(fileManager_.getTemporaryFilePath(".log")), edrFilename_(fileManager_.getTemporaryFilePath(".edr")), log_(std::fopen(logFilename_.c_str(), "w")), @@ -312,12 +316,12 @@ public: mtop_.groups.groupNames.emplace_back(&handle); } - mtop_.groups.groups[SimulationAtomGroupType::EnergyOutput].resize(3); + mtop_.groups.groups[SimulationAtomGroupType::EnergyOutput].resize(numTempCouplingGroups_); mtop_.groups.groups[SimulationAtomGroupType::EnergyOutput][0] = 0; mtop_.groups.groups[SimulationAtomGroupType::EnergyOutput][1] = 1; mtop_.groups.groups[SimulationAtomGroupType::EnergyOutput][2] = 2; - mtop_.groups.groups[SimulationAtomGroupType::TemperatureCoupling].resize(3); + mtop_.groups.groups[SimulationAtomGroupType::TemperatureCoupling].resize(numTempCouplingGroups_); mtop_.groups.groups[SimulationAtomGroupType::TemperatureCoupling][0] = 0; mtop_.groups.groups[SimulationAtomGroupType::TemperatureCoupling][1] = 1; mtop_.groups.groups[SimulationAtomGroupType::TemperatureCoupling][2] = 2; @@ -344,15 +348,10 @@ public: ekindata_.tcstat.resize(mtop_.groups.groups[SimulationAtomGroupType::TemperatureCoupling].size()); // This is needed so that the ebin space will be allocated - inputrec_.cos_accel = 1.0; - // This is to keep the destructor happy (otherwise sfree() segfaults) - ekindata_.nthreads = 0; - snew(ekindata_.ekin_work_alloc, 1); - snew(ekindata_.ekin_work, 1); - snew(ekindata_.dekindl_work, 1); + inputrec_.cos_accel = cosAccel_; // Group options for annealing output - inputrec_.opts.ngtc = 3; + inputrec_.opts.ngtc = numTempCouplingGroups_; snew(inputrec_.opts.ref_t, inputrec_.opts.ngtc); snew(inputrec_.opts.annealing, inputrec_.opts.ngtc); inputrec_.opts.annealing[0] = SimulatedAnnealing::No; diff --git a/src/gromacs/mdlib/tests/leapfrogtestdata.cpp b/src/gromacs/mdlib/tests/leapfrogtestdata.cpp index 8ae0a74051..c00fbbd4c0 100644 --- a/src/gromacs/mdlib/tests/leapfrogtestdata.cpp +++ b/src/gromacs/mdlib/tests/leapfrogtestdata.cpp @@ -83,6 +83,7 @@ LeapFrogTestData::LeapFrogTestData(int numAtoms, f_(numAtoms), inverseMasses_(numAtoms), inverseMassesPerDim_(numAtoms), + kineticEnergyData_(numTCoupleGroups == 0 ? 1 : numTCoupleGroups, 0.0, 1), numTCoupleGroups_(numTCoupleGroups) { mdAtoms_.nr = numAtoms_; @@ -128,10 +129,9 @@ LeapFrogTestData::LeapFrogTestData(int numAtoms, { mdAtoms_.cTC[i] = 0; } - kineticEnergyData_.ngtc = 1; t_grp_tcstat temperatureCouplingGroupData; temperatureCouplingGroupData.lambda = 1.0; - kineticEnergyData_.tcstat.emplace_back(temperatureCouplingGroupData); + kineticEnergyData_.tcstat[0] = temperatureCouplingGroupData; } else { @@ -140,13 +140,12 @@ LeapFrogTestData::LeapFrogTestData(int numAtoms, { mdAtoms_.cTC[i] = i % numTCoupleGroups_; } - kineticEnergyData_.ngtc = numTCoupleGroups_; - for (int i = 0; i < numTCoupleGroups; i++) + for (int i = 0; i < numTCoupleGroups_; i++) { real tCoupleLambda = 1.0 - (i + 1.0) / 10.0; t_grp_tcstat temperatureCouplingGroupData; temperatureCouplingGroupData.lambda = tCoupleLambda; - kineticEnergyData_.tcstat.emplace_back(temperatureCouplingGroupData); + kineticEnergyData_.tcstat[i] = temperatureCouplingGroupData; } } @@ -167,13 +166,6 @@ LeapFrogTestData::LeapFrogTestData(int numAtoms, state_.box[ZZ][YY] = 0.0; state_.box[ZZ][ZZ] = 10.0; - kineticEnergyData_.cosacc.cos_accel = 0.0; - - kineticEnergyData_.nthreads = 1; - snew(kineticEnergyData_.ekin_work_alloc, kineticEnergyData_.nthreads); - snew(kineticEnergyData_.ekin_work, kineticEnergyData_.nthreads); - snew(kineticEnergyData_.dekindl_work, kineticEnergyData_.nthreads); - mdAtoms_.homenr = numAtoms_; mdAtoms_.haveVsites = false; mdAtoms_.havePartiallyFrozenAtoms = false; diff --git a/src/gromacs/mdlib/tgroup.cpp b/src/gromacs/mdlib/tgroup.cpp index 1fccc924fe..f86fc068b1 100644 --- a/src/gromacs/mdlib/tgroup.cpp +++ b/src/gromacs/mdlib/tgroup.cpp @@ -40,70 +40,11 @@ #include "tgroup.h" -#include - -#include "gromacs/gmxlib/network.h" #include "gromacs/math/vec.h" #include "gromacs/mdlib/coupling.h" -#include "gromacs/mdlib/gmx_omp_nthreads.h" -#include "gromacs/mdlib/rbin.h" #include "gromacs/mdtypes/group.h" #include "gromacs/mdtypes/inputrec.h" -#include "gromacs/mdtypes/mdatom.h" -#include "gromacs/topology/mtop_util.h" -#include "gromacs/topology/topology.h" -#include "gromacs/utility/exceptions.h" -#include "gromacs/utility/fatalerror.h" -#include "gromacs/utility/futil.h" -#include "gromacs/utility/smalloc.h" - -void init_ekindata(FILE gmx_unused* log, const t_grpopts* opts, gmx_ekindata_t* ekind, real cos_accel) -{ - int i; - - ekind->ngtc = opts->ngtc; - ekind->tcstat.resize(opts->ngtc); - /* Set Berendsen tcoupl lambda's to 1, - * so runs without Berendsen coupling are not affected. - */ - for (i = 0; i < opts->ngtc; i++) - { - ekind->tcstat[i].lambda = 1.0; - ekind->tcstat[i].vscale_nhc = 1.0; - ekind->tcstat[i].ekinscaleh_nhc = 1.0; - ekind->tcstat[i].ekinscalef_nhc = 1.0; - } - int nthread = gmx_omp_nthreads_get(emntUpdate); - ekind->nthreads = nthread; - snew(ekind->ekin_work_alloc, nthread); - snew(ekind->ekin_work, nthread); - snew(ekind->dekindl_work, nthread); -#pragma omp parallel for num_threads(nthread) schedule(static) - for (int thread = 0; thread < nthread; thread++) - { - try - { -#define EKIN_WORK_BUFFER_SIZE 2 - /* Allocate 2 extra elements on both sides, so in single - * precision we have - * EKIN_WORK_BUFFER_SIZE*DIM*DIM*sizeof(real) = 72/144 bytes - * buffer on both sides to avoid cache pollution. - */ - snew(ekind->ekin_work_alloc[thread], ekind->ngtc + 2 * EKIN_WORK_BUFFER_SIZE); - ekind->ekin_work[thread] = ekind->ekin_work_alloc[thread] + EKIN_WORK_BUFFER_SIZE; - /* Nasty hack so we can have the per-thread accumulation - * variable for dekindl in the same thread-local cache lines - * as the per-thread accumulation tensors for ekin[fh], - * because they are accumulated in the same loop. */ - ekind->dekindl_work[thread] = &(ekind->ekin_work[thread][ekind->ngtc][0][0]); -#undef EKIN_WORK_BUFFER_SIZE - } - GMX_CATCH_ALL_AND_EXIT_WITH_FATAL_ERROR - } - - ekind->cosacc.cos_accel = cos_accel; -} real sum_ekin(const t_grpopts* opts, gmx_ekindata_t* ekind, real* dekindlambda, gmx_bool bEkinAveVel, gmx_bool bScaleEkin) { diff --git a/src/gromacs/mdlib/tgroup.h b/src/gromacs/mdlib/tgroup.h index 1f9189414a..9a62f8fe0e 100644 --- a/src/gromacs/mdlib/tgroup.h +++ b/src/gromacs/mdlib/tgroup.h @@ -37,23 +37,12 @@ #ifndef GMX_MDLIB_TGROUP_H #define GMX_MDLIB_TGROUP_H -#include - -#include "gromacs/utility/basedefinitions.h" #include "gromacs/utility/real.h" -struct gmx_ekindata_t; -struct t_commrec; +class gmx_ekindata_t; struct t_grpopts; -struct t_mdatoms; - -void init_ekindata(FILE* log, const t_grpopts* opts, gmx_ekindata_t* ekind, real cos_accel); -/* Allocate memory and set the grpnr array. */ - -void done_ekindata(gmx_ekindata_t* ekind); -/* Free the memory */ -real sum_ekin(const t_grpopts* opts, gmx_ekindata_t* ekind, real* dekindlambda, gmx_bool bEkinFullStep, gmx_bool bScaleEkin); +real sum_ekin(const t_grpopts* opts, gmx_ekindata_t* ekind, real* dekindlambda, bool bEkinFullStep, bool bScaleEkin); /* Sum the group ekins into total ekin and calc temp per group, * return total temperature. */ diff --git a/src/gromacs/mdlib/trajectory_writing.h b/src/gromacs/mdlib/trajectory_writing.h index 5d60396834..220b6056eb 100644 --- a/src/gromacs/mdlib/trajectory_writing.h +++ b/src/gromacs/mdlib/trajectory_writing.h @@ -44,7 +44,7 @@ #include "gromacs/mdlib/mdoutf.h" #include "gromacs/timing/wallcycle.h" -struct gmx_ekindata_t; +class gmx_ekindata_t; struct gmx_mtop_t; struct ObservablesHistory; struct t_commrec; diff --git a/src/gromacs/mdlib/update.h b/src/gromacs/mdlib/update.h index 207685f19e..1e7a5d530f 100644 --- a/src/gromacs/mdlib/update.h +++ b/src/gromacs/mdlib/update.h @@ -49,7 +49,7 @@ #include "gromacs/utility/real.h" class ekinstate_t; -struct gmx_ekindata_t; +class gmx_ekindata_t; struct gmx_enerdata_t; enum class PbcType; struct t_fcdata; diff --git a/src/gromacs/mdlib/update_vv.h b/src/gromacs/mdlib/update_vv.h index b1d7cf00ae..ddc4f8b747 100644 --- a/src/gromacs/mdlib/update_vv.h +++ b/src/gromacs/mdlib/update_vv.h @@ -44,7 +44,7 @@ #include "gromacs/mdrunutility/handlerestart.h" #include "gromacs/mdtypes/md_enums.h" -struct gmx_ekindata_t; +class gmx_ekindata_t; struct gmx_enerdata_t; struct gmx_global_stat; struct gmx_localtop_t; diff --git a/src/gromacs/mdrun/isimulator.h b/src/gromacs/mdrun/isimulator.h index 6e9c6c6149..7ef5a99846 100644 --- a/src/gromacs/mdrun/isimulator.h +++ b/src/gromacs/mdrun/isimulator.h @@ -44,7 +44,7 @@ #include "gromacs/mdlib/stophandler.h" class energyhistory_t; -struct gmx_ekindata_t; +class gmx_ekindata_t; struct gmx_enerdata_t; struct gmx_enfrot; struct gmx_mtop_t; diff --git a/src/gromacs/mdrun/runner.cpp b/src/gromacs/mdrun/runner.cpp index f732c60ede..19fe216747 100644 --- a/src/gromacs/mdrun/runner.cpp +++ b/src/gromacs/mdrun/runner.cpp @@ -1868,8 +1868,7 @@ int Mdrunner::mdrunner() "cos_acceleration is only supported by integrator=md"); /* Kinetic energy data */ - gmx_ekindata_t ekind; - init_ekindata(fplog, &(inputrec->opts), &ekind, inputrec->cos_accel); + gmx_ekindata_t ekind(inputrec->opts.ngtc, inputrec->cos_accel, gmx_omp_nthreads_get(emntUpdate)); /* Set up interactive MD (IMD) */ auto imdSession = makeImdSession(inputrec.get(), diff --git a/src/gromacs/mdrun/simulatorbuilder.h b/src/gromacs/mdrun/simulatorbuilder.h index 3bd22bb769..624f8d8e71 100644 --- a/src/gromacs/mdrun/simulatorbuilder.h +++ b/src/gromacs/mdrun/simulatorbuilder.h @@ -45,7 +45,7 @@ class energyhistory_t; -struct gmx_ekindata_t; +class gmx_ekindata_t; struct gmx_enerdata_t; struct gmx_enfrot; struct gmx_mtop_t; diff --git a/src/gromacs/mdtypes/group.cpp b/src/gromacs/mdtypes/group.cpp index b0c15ca852..8e708c737d 100644 --- a/src/gromacs/mdtypes/group.cpp +++ b/src/gromacs/mdtypes/group.cpp @@ -1,7 +1,7 @@ /* * This file is part of the GROMACS molecular simulation package. * - * Copyright (c) 2019, by the GROMACS development team, led by + * Copyright (c) 2019,2021, by the GROMACS development team, led by * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl, * and including many others, as listed in the AUTHORS file in the * top-level source directory and at http://www.gromacs.org. @@ -43,9 +43,56 @@ #include "group.h" +#include "gromacs/utility/exceptions.h" + +gmx_ekindata_t::gmx_ekindata_t(int numTempCoupleGroups, real cos_accel, int numThreads) : + ngtc(numTempCoupleGroups), + nthreads_(numThreads) +{ + tcstat.resize(ngtc); + /* Set Berendsen tcoupl lambda's to 1, + * so runs without Berendsen coupling are not affected. + */ + for (int i = 0; i < ngtc; i++) + { + tcstat[i].lambda = 1.0; + tcstat[i].vscale_nhc = 1.0; + tcstat[i].ekinscaleh_nhc = 1.0; + tcstat[i].ekinscalef_nhc = 1.0; + } + + snew(ekin_work_alloc, nthreads_); + snew(ekin_work, nthreads_); + snew(dekindl_work, nthreads_); + +#pragma omp parallel for num_threads(nthreads_) schedule(static) + for (int thread = 0; thread < nthreads_; thread++) + { + try + { + constexpr int EKIN_WORK_BUFFER_SIZE = 2; + /* Allocate 2 extra elements on both sides, so in single + * precision we have + * EKIN_WORK_BUFFER_SIZE*DIM*DIM*sizeof(real) = 72/144 bytes + * buffer on both sides to avoid cache pollution. + */ + snew(ekin_work_alloc[thread], ngtc + 2 * EKIN_WORK_BUFFER_SIZE); + ekin_work[thread] = ekin_work_alloc[thread] + EKIN_WORK_BUFFER_SIZE; + /* Nasty hack so we can have the per-thread accumulation + * variable for dekindl in the same thread-local cache lines + * as the per-thread accumulation tensors for ekin[fh], + * because they are accumulated in the same loop. */ + dekindl_work[thread] = &(ekin_work[thread][ngtc][0][0]); + } + GMX_CATCH_ALL_AND_EXIT_WITH_FATAL_ERROR + } + + cosacc.cos_accel = cos_accel; +} + gmx_ekindata_t::~gmx_ekindata_t() { - for (int i = 0; i < nthreads; i++) + for (int i = 0; i < nthreads_; i++) { sfree(ekin_work_alloc[i]); } diff --git a/src/gromacs/mdtypes/group.h b/src/gromacs/mdtypes/group.h index 49ece65476..27ebd8e1be 100644 --- a/src/gromacs/mdtypes/group.h +++ b/src/gromacs/mdtypes/group.h @@ -75,12 +75,12 @@ struct t_cos_acc real vcos = 0; }; -struct gmx_ekindata_t +class gmx_ekindata_t { +public: + gmx_ekindata_t(int numTempCoupleGroups, real cos_accel, int numThreads); //! The number of T-coupling groups - int ngtc = 0; - //! For size of ekin_work - int nthreads = 0; + int ngtc; //! T-coupling data std::vector tcstat; //! Allocated locations for *_work members @@ -101,6 +101,10 @@ struct gmx_ekindata_t t_cos_acc cosacc; ~gmx_ekindata_t(); + +private: + //! For size of ekin_work + int nthreads_ = 0; }; #define GID(igid, jgid, gnr) \ diff --git a/src/gromacs/modularsimulator/energydata.h b/src/gromacs/modularsimulator/energydata.h index bc0c55dc30..0020dd04da 100644 --- a/src/gromacs/modularsimulator/energydata.h +++ b/src/gromacs/modularsimulator/energydata.h @@ -49,7 +49,7 @@ #include "modularsimulatorinterfaces.h" -struct gmx_ekindata_t; +class gmx_ekindata_t; struct gmx_enerdata_t; struct gmx_mtop_t; struct ObservablesHistory; -- 2.22.0