#include "gromacs/mdlib/mdatoms.h"
#include "gromacs/mdlib/vsite.h"
#include "gromacs/mdtypes/commrec.h"
+#include "gromacs/mdtypes/forcebuffers.h"
#include "gromacs/mdtypes/forcerec.h"
#include "gromacs/mdtypes/inputrec.h"
#include "gromacs/mdtypes/interaction_const.h"
* The final solution should be an MD algorithm base class with methods
* for initialization and atom-data setup.
*/
-void mdAlgorithmsSetupAtomData(const t_commrec* cr,
- const t_inputrec* ir,
- const gmx_mtop_t& top_global,
- gmx_localtop_t* top,
- t_forcerec* fr,
- PaddedHostVector<RVec>* force,
- MDAtoms* mdAtoms,
- Constraints* constr,
- VirtualSitesHandler* vsite,
- gmx_shellfc_t* shellfc)
+void mdAlgorithmsSetupAtomData(const t_commrec* cr,
+ const t_inputrec* ir,
+ const gmx_mtop_t& top_global,
+ gmx_localtop_t* top,
+ t_forcerec* fr,
+ ForceBuffers* force,
+ MDAtoms* mdAtoms,
+ Constraints* constr,
+ VirtualSitesHandler* vsite,
+ gmx_shellfc_t* shellfc)
{
bool usingDomDec = DOMAINDECOMP(cr);
if (force != nullptr)
{
- /* We need to allocate one element extra, since we might use
- * (unaligned) 4-wide SIMD loads to access rvec entries.
- */
- force->resizeWithPadding(numTotalAtoms);
+ force->resize(numTotalAtoms);
}
atoms2md(&top_global, ir, numAtomIndex,
#ifndef GMX_DOMDEC_MDSETUP_H
#define GMX_DOMDEC_MDSETUP_H
-#include "gromacs/gpu_utils/hostallocator.h"
#include "gromacs/math/vectypes.h"
struct bonded_threading_t;
namespace gmx
{
class Constraints;
+class ForceBuffers;
class MDAtoms;
class VirtualSitesHandler;
* \param[in,out] vsite The virtual site data, can be NULL
* \param[in,out] shellfc The shell/flexible-constraint data, can be NULL
*/
-void mdAlgorithmsSetupAtomData(const t_commrec* cr,
- const t_inputrec* ir,
- const gmx_mtop_t& top_global,
- gmx_localtop_t* top,
- t_forcerec* fr,
- PaddedHostVector<RVec>* force,
- MDAtoms* mdAtoms,
- Constraints* constr,
- VirtualSitesHandler* vsite,
- gmx_shellfc_t* shellfc);
+void mdAlgorithmsSetupAtomData(const t_commrec* cr,
+ const t_inputrec* ir,
+ const gmx_mtop_t& top_global,
+ gmx_localtop_t* top,
+ t_forcerec* fr,
+ ForceBuffers* force,
+ MDAtoms* mdAtoms,
+ Constraints* constr,
+ VirtualSitesHandler* vsite,
+ gmx_shellfc_t* shellfc);
} // namespace gmx
}
//!\brief TODO Remove fplog when group scheme and charge groups are gone
-void dd_partition_system(FILE* fplog,
- const gmx::MDLogger& mdlog,
- int64_t step,
- const t_commrec* cr,
- gmx_bool bMasterState,
- int nstglobalcomm,
- t_state* state_global,
- const gmx_mtop_t& top_global,
- const t_inputrec* ir,
- gmx::ImdSession* imdSession,
- pull_t* pull_work,
- t_state* state_local,
- PaddedHostVector<gmx::RVec>* f,
- gmx::MDAtoms* mdAtoms,
- gmx_localtop_t* top_local,
- t_forcerec* fr,
- gmx::VirtualSitesHandler* vsite,
- gmx::Constraints* constr,
- t_nrnb* nrnb,
- gmx_wallcycle* wcycle,
- gmx_bool bVerbose)
+void dd_partition_system(FILE* fplog,
+ const gmx::MDLogger& mdlog,
+ int64_t step,
+ const t_commrec* cr,
+ gmx_bool bMasterState,
+ int nstglobalcomm,
+ t_state* state_global,
+ const gmx_mtop_t& top_global,
+ const t_inputrec* ir,
+ gmx::ImdSession* imdSession,
+ pull_t* pull_work,
+ t_state* state_local,
+ gmx::ForceBuffers* f,
+ gmx::MDAtoms* mdAtoms,
+ gmx_localtop_t* top_local,
+ t_forcerec* fr,
+ gmx::VirtualSitesHandler* vsite,
+ gmx::Constraints* constr,
+ t_nrnb* nrnb,
+ gmx_wallcycle* wcycle,
+ gmx_bool bVerbose)
{
gmx_domdec_t* dd;
gmx_domdec_comm_t* comm;
#ifndef GMX_DOMDEC_PARTITION_H
#define GMX_DOMDEC_PARTITION_H
-#include "gromacs/gpu_utils/hostallocator.h"
+#include <cstdio>
+
#include "gromacs/math/vectypes.h"
#include "gromacs/utility/basedefinitions.h"
#include "gromacs/utility/real.h"
namespace gmx
{
+template<typename>
+class ArrayRef;
class Constraints;
+class ForceBuffers;
class ImdSession;
class MDAtoms;
class MDLogger;
* \param[in] pull_work Pulling data
* \param[in] state_local Local state
* \param[in] f Force buffer
- * \param[in] mdatoms MD atoms
+ * \param[in] mdAtoms MD atoms
* \param[in] top_local Local topology
* \param[in] fr Force record
* \param[in] vsite Virtual sites handler
* \param[in] wcycle Timers
* \param[in] bVerbose Be verbose
*/
-void dd_partition_system(FILE* fplog,
- const gmx::MDLogger& mdlog,
- int64_t step,
- const t_commrec* cr,
- gmx_bool bMasterState,
- int nstglobalcomm,
- t_state* state_global,
- const gmx_mtop_t& top_global,
- const t_inputrec* ir,
- gmx::ImdSession* imdSession,
- pull_t* pull_work,
- t_state* state_local,
- gmx::PaddedHostVector<gmx::RVec>* f,
- gmx::MDAtoms* mdatoms,
- gmx_localtop_t* top_local,
- t_forcerec* fr,
- gmx::VirtualSitesHandler* vsite,
- gmx::Constraints* constr,
- t_nrnb* nrnb,
- gmx_wallcycle* wcycle,
- gmx_bool bVerbose);
+void dd_partition_system(FILE* fplog,
+ const gmx::MDLogger& mdlog,
+ int64_t step,
+ const t_commrec* cr,
+ gmx_bool bMasterState,
+ int nstglobalcomm,
+ t_state* state_global,
+ const gmx_mtop_t& top_global,
+ const t_inputrec* ir,
+ gmx::ImdSession* imdSession,
+ pull_t* pull_work,
+ t_state* state_local,
+ gmx::ForceBuffers* f,
+ gmx::MDAtoms* mdAtoms,
+ gmx_localtop_t* top_local,
+ t_forcerec* fr,
+ gmx::VirtualSitesHandler* vsite,
+ gmx::Constraints* constr,
+ t_nrnb* nrnb,
+ gmx_wallcycle* wcycle,
+ gmx_bool bVerbose);
/*! \brief Check whether bonded interactions are missing, if appropriate
*
namespace gmx
{
class Awh;
+class ForceBuffersView;
class ForceOutputs;
class ForceWithVirial;
class ImdSession;
const matrix box,
gmx::ArrayRefWithPadding<gmx::RVec> coordinates,
history_t* hist,
- gmx::ArrayRefWithPadding<gmx::RVec> force,
+ gmx::ForceBuffersView* force,
tensor vir_force,
const t_mdatoms* mdatoms,
gmx_enerdata_t* enerd,
#endif /* end GMX_FAHCORE block */
}
-void mdoutf_write_to_trajectory_files(FILE* fplog,
- const t_commrec* cr,
- gmx_mdoutf_t of,
- int mdof_flags,
- int natoms,
- int64_t step,
- double t,
- t_state* state_local,
- t_state* state_global,
- ObservablesHistory* observablesHistory,
- gmx::ArrayRef<gmx::RVec> f_local)
+void mdoutf_write_to_trajectory_files(FILE* fplog,
+ const t_commrec* cr,
+ gmx_mdoutf_t of,
+ int mdof_flags,
+ int natoms,
+ int64_t step,
+ double t,
+ t_state* state_local,
+ t_state* state_global,
+ ObservablesHistory* observablesHistory,
+ gmx::ArrayRef<const gmx::RVec> f_local)
{
- rvec* f_global;
+ const rvec* f_global;
if (DOMAINDECOMP(cr))
{
f_global = of->f_global;
if (mdof_flags & MDOF_F)
{
- dd_collect_vec(cr->dd, state_local, f_local,
- gmx::arrayRefFromArray(reinterpret_cast<gmx::RVec*>(f_global), f_local.size()));
+ dd_collect_vec(
+ cr->dd, state_local, f_local,
+ gmx::arrayRefFromArray(reinterpret_cast<gmx::RVec*>(of->f_global), f_local.size()));
}
}
else
* \param[in] observablesHistory Pointer to the ObservableHistory object.
* \param[in] f_local The local forces.
*/
-void mdoutf_write_to_trajectory_files(FILE* fplog,
- const t_commrec* cr,
- gmx_mdoutf_t of,
- int mdof_flags,
- int natoms,
- int64_t step,
- double t,
- t_state* state_local,
- t_state* state_global,
- ObservablesHistory* observablesHistory,
- gmx::ArrayRef<gmx::RVec> f_local);
+void mdoutf_write_to_trajectory_files(FILE* fplog,
+ const t_commrec* cr,
+ gmx_mdoutf_t of,
+ int mdof_flags,
+ int natoms,
+ int64_t step,
+ double t,
+ t_state* state_local,
+ t_state* state_global,
+ ObservablesHistory* observablesHistory,
+ gmx::ArrayRef<const gmx::RVec> f_local);
/*! \brief Get the output interval of box size of uncompressed TNG output.
* Returns 0 if no uncompressed TNG file is open.
#include "gromacs/mdlib/wholemoleculetransform.h"
#include "gromacs/mdtypes/commrec.h"
#include "gromacs/mdtypes/enerdata.h"
+#include "gromacs/mdtypes/forcebuffers.h"
#include "gromacs/mdtypes/forceoutput.h"
#include "gromacs/mdtypes/forcerec.h"
#include "gromacs/mdtypes/iforceprovider.h"
const matrix box,
gmx::ArrayRefWithPadding<gmx::RVec> x,
history_t* hist,
- gmx::ArrayRefWithPadding<gmx::RVec> force,
+ gmx::ForceBuffersView* forceView,
tensor vir_force,
const t_mdatoms* mdatoms,
gmx_enerdata_t* enerd,
int legacyFlags,
const DDBalanceRegionHandler& ddBalanceRegionHandler)
{
+ auto force = forceView->forceWithPadding();
GMX_ASSERT(force.unpaddedArrayRef().ssize() >= fr->natoms_force_constr,
"The size of the force buffer should be at least the number of atoms to compute "
"forces for");
#include "gromacs/topology/topology.h"
#include "gromacs/utility/smalloc.h"
-void do_md_trajectory_writing(FILE* fplog,
- t_commrec* cr,
- int nfile,
- const t_filenm fnm[],
- int64_t step,
- int64_t step_rel,
- double t,
- t_inputrec* ir,
- t_state* state,
- t_state* state_global,
- ObservablesHistory* observablesHistory,
- const gmx_mtop_t* top_global,
- t_forcerec* fr,
- gmx_mdoutf_t outf,
- const gmx::EnergyOutput& energyOutput,
- gmx_ekindata_t* ekind,
- gmx::ArrayRef<gmx::RVec> f,
- gmx_bool bCPT,
- gmx_bool bRerunMD,
- gmx_bool bLastStep,
- gmx_bool bDoConfOut,
- gmx_bool bSumEkinhOld)
+void do_md_trajectory_writing(FILE* fplog,
+ t_commrec* cr,
+ int nfile,
+ const t_filenm fnm[],
+ int64_t step,
+ int64_t step_rel,
+ double t,
+ t_inputrec* ir,
+ t_state* state,
+ t_state* state_global,
+ ObservablesHistory* observablesHistory,
+ const gmx_mtop_t* top_global,
+ t_forcerec* fr,
+ gmx_mdoutf_t outf,
+ const gmx::EnergyOutput& energyOutput,
+ gmx_ekindata_t* ekind,
+ gmx::ArrayRef<const gmx::RVec> f,
+ gmx_bool bCPT,
+ gmx_bool bRerunMD,
+ gmx_bool bLastStep,
+ gmx_bool bDoConfOut,
+ gmx_bool bSumEkinhOld)
{
int mdof_flags;
rvec* x_for_confout = nullptr;
*
* This routine does communication (e.g. collecting distributed coordinates)
*/
-void do_md_trajectory_writing(FILE* fplog,
- struct t_commrec* cr,
- int nfile,
- const t_filenm fnm[],
- int64_t step,
- int64_t step_rel,
- double t,
- t_inputrec* ir,
- t_state* state,
- t_state* state_global,
- ObservablesHistory* observablesHistory,
- const gmx_mtop_t* top_global,
- t_forcerec* fr,
- gmx_mdoutf_t outf,
- const gmx::EnergyOutput& energyOutput,
- gmx_ekindata_t* ekind,
- gmx::ArrayRef<gmx::RVec> f,
- gmx_bool bCPT,
- gmx_bool bRerunMD,
- gmx_bool bLastStep,
- gmx_bool bDoConfOut,
- gmx_bool bSumEkinhOld);
+void do_md_trajectory_writing(FILE* fplog,
+ struct t_commrec* cr,
+ int nfile,
+ const t_filenm fnm[],
+ int64_t step,
+ int64_t step_rel,
+ double t,
+ t_inputrec* ir,
+ t_state* state,
+ t_state* state_global,
+ ObservablesHistory* observablesHistory,
+ const gmx_mtop_t* top_global,
+ t_forcerec* fr,
+ gmx_mdoutf_t outf,
+ const gmx::EnergyOutput& energyOutput,
+ gmx_ekindata_t* ekind,
+ gmx::ArrayRef<const gmx::RVec> f,
+ gmx_bool bCPT,
+ gmx_bool bRerunMD,
+ gmx_bool bLastStep,
+ gmx_bool bDoConfOut,
+ gmx_bool bSumEkinhOld);
#endif
#include "gromacs/mdtypes/df_history.h"
#include "gromacs/mdtypes/energyhistory.h"
#include "gromacs/mdtypes/fcdata.h"
+#include "gromacs/mdtypes/forcebuffers.h"
#include "gromacs/mdtypes/forcerec.h"
#include "gromacs/mdtypes/group.h"
#include "gromacs/mdtypes/inputrec.h"
int i, m;
rvec mu_tot;
matrix pressureCouplingMu, M;
- gmx_repl_ex_t repl_ex = nullptr;
- PaddedHostVector<gmx::RVec> f{};
- gmx_global_stat_t gstat;
- gmx_shellfc_t* shellfc;
- gmx_bool bSumEkinhOld, bDoReplEx, bExchanged, bNeedRepartition;
- gmx_bool bTemp, bPres, bTrotter;
- real dvdl_constr;
- std::vector<RVec> cbuf;
- matrix lastbox;
- int lamnew = 0;
+ gmx_repl_ex_t repl_ex = nullptr;
+ gmx_global_stat_t gstat;
+ gmx_shellfc_t* shellfc;
+ gmx_bool bSumEkinhOld, bDoReplEx, bExchanged, bNeedRepartition;
+ gmx_bool bTemp, bPres, bTrotter;
+ real dvdl_constr;
+ std::vector<RVec> cbuf;
+ matrix lastbox;
+ int lamnew = 0;
/* for FEP */
int nstfep = 0;
double cycles;
auto mdatoms = mdAtoms->mdatoms();
- std::unique_ptr<UpdateConstrainGpu> integrator;
+ const auto& simulationWork = runScheduleWork->simulationWork;
+ const bool useGpuForPme = simulationWork.useGpuPme;
+ const bool useGpuForNonbonded = simulationWork.useGpuNonbonded;
+ const bool useGpuForBufferOps = simulationWork.useGpuBufferOps;
+ const bool useGpuForUpdate = simulationWork.useGpuUpdate;
+ ForceBuffers f(((useGpuForNonbonded && useGpuForBufferOps) || useGpuForUpdate)
+ ? PinningPolicy::PinnedIfSupported
+ : PinningPolicy::CannotBePinned);
if (DOMAINDECOMP(cr))
{
stateInstance = std::make_unique<t_state>();
upd.setNumAtoms(state->natoms);
}
- const auto& simulationWork = runScheduleWork->simulationWork;
- const bool useGpuForPme = simulationWork.useGpuPme;
- const bool useGpuForNonbonded = simulationWork.useGpuNonbonded;
- const bool useGpuForBufferOps = simulationWork.useGpuBufferOps;
- const bool useGpuForUpdate = simulationWork.useGpuUpdate;
+ std::unique_ptr<UpdateConstrainGpu> integrator;
StatePropagatorDataGpu* stateGpu = fr->stateGpu;
{
changePinningPolicy(&state->x, PinningPolicy::PinnedIfSupported);
}
- if ((useGpuForNonbonded && useGpuForBufferOps) || useGpuForUpdate)
- {
- changePinningPolicy(&f, PinningPolicy::PinnedIfSupported);
- }
if (useGpuForUpdate)
{
changePinningPolicy(&state->v, PinningPolicy::PinnedIfSupported);
relax_shell_flexcon(fplog, cr, ms, mdrunOptions.verbose, enforcedRotation, step, ir,
imdSession, pull_work, bNS, force_flags, &top, constr, enerd,
state->natoms, state->x.arrayRefWithPadding(),
- state->v.arrayRefWithPadding(), state->box, state->lambda, &state->hist,
- f.arrayRefWithPadding(), force_vir, mdatoms, nrnb, wcycle, shellfc,
+ state->v.arrayRefWithPadding(), state->box, state->lambda,
+ &state->hist, &f.view(), force_vir, mdatoms, nrnb, wcycle, shellfc,
fr, runScheduleWork, t, mu_tot, vsite, ddBalanceRegionHandler);
}
else
*/
do_force(fplog, cr, ms, ir, awh.get(), enforcedRotation, imdSession, pull_work, step,
nrnb, wcycle, &top, state->box, state->x.arrayRefWithPadding(), &state->hist,
- f.arrayRefWithPadding(), force_vir, mdatoms, enerd, state->lambda, fr,
- runScheduleWork, vsite, mu_tot, t, ed ? ed->getLegacyED() : nullptr,
+ &f.view(), force_vir, mdatoms, enerd, state->lambda, fr, runScheduleWork,
+ vsite, mu_tot, t, ed ? ed->getLegacyED() : nullptr,
(bNS ? GMX_FORCE_NS : 0) | force_flags, ddBalanceRegionHandler);
}
trotter_seq, ettTSEQ1);
}
- upd.update_coords(*ir, step, mdatoms, state, f.arrayRefWithPadding(), fcdata, ekind, M,
- etrtVELOCITY1, cr, constr != nullptr);
+ upd.update_coords(*ir, step, mdatoms, state, f.view().forceWithPadding(), fcdata, ekind,
+ M, etrtVELOCITY1, cr, constr != nullptr);
wallcycle_stop(wcycle, ewcUPDATE);
constrain_velocities(constr, do_log, do_ene, step, state, nullptr, bCalcVir, shake_vir);
if (runScheduleWork->stepWork.useGpuFBufferOps && (simulationWork.useGpuUpdate && !vsite)
&& do_per_step(step, ir->nstfout))
{
- stateGpu->copyForcesFromGpu(ArrayRef<RVec>(f), AtomLocality::Local);
+ stateGpu->copyForcesFromGpu(f.view().force(), AtomLocality::Local);
stateGpu->waitForcesReadyOnHost(AtomLocality::Local);
}
/* Now we have the energies and forces corresponding to the
* the update.
*/
do_md_trajectory_writing(fplog, cr, nfile, fnm, step, step_rel, t, ir, state, state_global,
- observablesHistory, top_global, fr, outf, energyOutput, ekind, f,
- checkpointHandler->isCheckpointingStep(), bRerunMD, bLastStep,
- mdrunOptions.writeConfout, bSumEkinhOld);
+ observablesHistory, top_global, fr, outf, energyOutput, ekind,
+ f.view().force(), checkpointHandler->isCheckpointingStep(),
+ bRerunMD, bLastStep, mdrunOptions.writeConfout, bSumEkinhOld);
/* Check if IMD step and do IMD communication, if bIMD is TRUE. */
bInteractiveMDstep = imdSession->run(step, bNS, state->box, state->x.rvec_array(), t);
if (EI_VV(ir->eI))
{
/* velocity half-step update */
- upd.update_coords(*ir, step, mdatoms, state, f.arrayRefWithPadding(), fcdata, ekind, M,
- etrtVELOCITY2, cr, constr != nullptr);
+ upd.update_coords(*ir, step, mdatoms, state, f.view().forceWithPadding(), fcdata, ekind,
+ M, etrtVELOCITY2, cr, constr != nullptr);
}
/* Above, initialize just copies ekinh into ekin,
// If the buffer ops were not offloaded this step, the forces are on the host and have to be copied
if (!runScheduleWork->stepWork.useGpuFBufferOps)
{
- stateGpu->copyForcesToGpu(ArrayRef<RVec>(f), AtomLocality::Local);
+ stateGpu->copyForcesToGpu(f.view().force(), AtomLocality::Local);
}
const bool doTemperatureScaling =
}
else
{
- upd.update_coords(*ir, step, mdatoms, state, f.arrayRefWithPadding(), fcdata, ekind, M,
- etrtPOSITION, cr, constr != nullptr);
+ upd.update_coords(*ir, step, mdatoms, state, f.view().forceWithPadding(), fcdata, ekind,
+ M, etrtPOSITION, cr, constr != nullptr);
wallcycle_stop(wcycle, ewcUPDATE);
/* now we know the scaling, we can compute the positions again */
std::copy(cbuf.begin(), cbuf.end(), state->x.begin());
- upd.update_coords(*ir, step, mdatoms, state, f.arrayRefWithPadding(), fcdata, ekind, M,
- etrtPOSITION, cr, constr != nullptr);
+ upd.update_coords(*ir, step, mdatoms, state, f.view().forceWithPadding(), fcdata, ekind,
+ M, etrtPOSITION, cr, constr != nullptr);
wallcycle_stop(wcycle, ewcUPDATE);
/* do we need an extra constraint here? just need to copy out of as_rvec_array(state->v.data()) to upd->xp? */
#include "gromacs/mdtypes/df_history.h"
#include "gromacs/mdtypes/enerdata.h"
#include "gromacs/mdtypes/energyhistory.h"
+#include "gromacs/mdtypes/forcebuffers.h"
#include "gromacs/mdtypes/forcerec.h"
#include "gromacs/mdtypes/group.h"
#include "gromacs/mdtypes/inputrec.h"
void gmx::LegacySimulator::do_mimic()
{
- t_inputrec* ir = inputrec;
- int64_t step, step_rel;
- double t, lam0[efptNR];
- bool isLastStep = false;
- bool doFreeEnergyPerturbation = false;
- unsigned int force_flags;
- tensor force_vir, shake_vir, total_vir, pres;
- rvec mu_tot;
- PaddedHostVector<gmx::RVec> f{};
- gmx_global_stat_t gstat;
- gmx_shellfc_t* shellfc;
+ t_inputrec* ir = inputrec;
+ int64_t step, step_rel;
+ double t, lam0[efptNR];
+ bool isLastStep = false;
+ bool doFreeEnergyPerturbation = false;
+ unsigned int force_flags;
+ tensor force_vir, shake_vir, total_vir, pres;
+ rvec mu_tot;
+ ForceBuffers f;
+ gmx_global_stat_t gstat;
+ gmx_shellfc_t* shellfc;
double cycles;
relax_shell_flexcon(fplog, cr, ms, mdrunOptions.verbose, enforcedRotation, step, ir,
imdSession, pull_work, bNS, force_flags, &top, constr, enerd,
state->natoms, state->x.arrayRefWithPadding(),
- state->v.arrayRefWithPadding(), state->box, state->lambda, &state->hist,
- f.arrayRefWithPadding(), force_vir, mdatoms, nrnb, wcycle, shellfc,
+ state->v.arrayRefWithPadding(), state->box, state->lambda,
+ &state->hist, &f.view(), force_vir, mdatoms, nrnb, wcycle, shellfc,
fr, runScheduleWork, t, mu_tot, vsite, ddBalanceRegionHandler);
}
else
gmx_edsam* ed = nullptr;
do_force(fplog, cr, ms, ir, awh, enforcedRotation, imdSession, pull_work, step, nrnb,
wcycle, &top, state->box, state->x.arrayRefWithPadding(), &state->hist,
- f.arrayRefWithPadding(), force_vir, mdatoms, enerd, state->lambda, fr, runScheduleWork,
+ &f.view(), force_vir, mdatoms, enerd, state->lambda, fr, runScheduleWork,
vsite, mu_tot, t, ed, GMX_FORCE_NS | force_flags, ddBalanceRegionHandler);
}
const bool bSumEkinhOld = false;
do_md_trajectory_writing(fplog, cr, nfile, fnm, step, step_rel, t, ir, state,
state_global, observablesHistory, top_global, fr, outf,
- energyOutput, ekind, f, isCheckpointingStep, doRerun,
- isLastStep, mdrunOptions.writeConfout, bSumEkinhOld);
+ energyOutput, ekind, f.view().force(), isCheckpointingStep,
+ doRerun, isLastStep, mdrunOptions.writeConfout, bSumEkinhOld);
}
stopHandler->setSignal();
{
gmx::HostVector<gmx::RVec> fglobal(top_global->natoms);
gmx::ArrayRef<gmx::RVec> ftemp;
- gmx::ArrayRef<const gmx::RVec> flocal = gmx::makeArrayRef(f);
+ gmx::ArrayRef<const gmx::RVec> flocal = f.view().force();
if (DOMAINDECOMP(cr))
{
ftemp = gmx::makeArrayRef(fglobal);
}
else
{
- ftemp = gmx::makeArrayRef(f);
+ ftemp = f.view().force();
}
if (MASTER(cr))
#include "gromacs/mdrunutility/handlerestart.h"
#include "gromacs/mdrunutility/printtime.h"
#include "gromacs/mdtypes/commrec.h"
+#include "gromacs/mdtypes/forcebuffers.h"
#include "gromacs/mdtypes/forcerec.h"
#include "gromacs/mdtypes/inputrec.h"
#include "gromacs/mdtypes/interaction_const.h"
//! Copy of the global state
t_state s;
//! Force array
- PaddedHostVector<gmx::RVec> f;
+ gmx::ForceBuffers f;
//! Potential energy
real epot;
//! Norm of the force
}
//! Compute the norm and max of the force array in parallel
-static void get_f_norm_max(const t_commrec* cr,
- t_grpopts* opts,
- t_mdatoms* mdatoms,
- const rvec* f,
- real* fnorm,
- real* fmax,
- int* a_fmax)
+static void get_f_norm_max(const t_commrec* cr,
+ t_grpopts* opts,
+ t_mdatoms* mdatoms,
+ gmx::ArrayRef<const gmx::RVec> f,
+ real* fnorm,
+ real* fmax,
+ int* a_fmax)
{
double fnorm2, *sum;
real fmax2, fam;
//! Compute the norm of the force
static void get_state_f_norm_max(const t_commrec* cr, t_grpopts* opts, t_mdatoms* mdatoms, em_state_t* ems)
{
- get_f_norm_max(cr, opts, mdatoms, ems->f.rvec_array(), &ems->fnorm, &ems->fmax, &ems->a_fmax);
+ get_f_norm_max(cr, opts, mdatoms, ems->f.view().force(), &ems->fnorm, &ems->fmax, &ems->a_fmax);
}
//! Initialize the energy minimization
mdoutf_write_to_trajectory_files(fplog, cr, outf, mdof_flags, top_global->natoms, step,
static_cast<double>(step), &state->s, state_global,
- observablesHistory, state->f);
+ observablesHistory, state->f.view().force());
if (confout != nullptr)
{
//! \brief Do one minimization step
//
// \returns true when the step succeeded, false when a constraint error occurred
-static bool do_em_step(const t_commrec* cr,
- t_inputrec* ir,
- t_mdatoms* md,
- em_state_t* ems1,
- real a,
- const PaddedHostVector<gmx::RVec>* force,
- em_state_t* ems2,
- gmx::Constraints* constr,
- int64_t count)
+static bool do_em_step(const t_commrec* cr,
+ t_inputrec* ir,
+ t_mdatoms* md,
+ em_state_t* ems1,
+ real a,
+ gmx::ArrayRefWithPadding<const gmx::RVec> force,
+ em_state_t* ems2,
+ gmx::Constraints* constr,
+ int64_t count)
{
t_state *s1, *s2;
if (s2->natoms != s1->natoms)
{
state_change_natoms(s2, s1->natoms);
- ems2->f.resizeWithPadding(s2->natoms);
+ ems2->f.resize(s2->natoms);
}
if (DOMAINDECOMP(cr) && s2->cg_gl.size() != s1->cg_gl.size())
{
{
const rvec* x1 = s1->x.rvec_array();
rvec* x2 = s2->x.rvec_array();
- const rvec* f = force->rvec_array();
+ const rvec* f = as_rvec_array(force.unpaddedArrayRef().data());
int gf = 0;
#pragma omp for schedule(static) nowait
* We do not unshift, so molecules are always whole in congrad.c
*/
do_force(fplog, cr, ms, inputrec, nullptr, nullptr, imdSession, pull_work, count, nrnb, wcycle,
- top, ems->s.box, ems->s.x.arrayRefWithPadding(), &ems->s.hist,
- ems->f.arrayRefWithPadding(), force_vir, mdAtoms->mdatoms(), enerd, ems->s.lambda, fr,
- runScheduleWork, vsite, mu_tot, t, nullptr,
+ top, ems->s.box, ems->s.x.arrayRefWithPadding(), &ems->s.hist, &ems->f.view(), force_vir,
+ mdAtoms->mdatoms(), enerd, ems->s.lambda, fr, runScheduleWork, vsite, mu_tot, t, nullptr,
GMX_FORCE_STATECHANGED | GMX_FORCE_ALLFORCES | GMX_FORCE_VIRIAL | GMX_FORCE_ENERGY
| (bNS ? GMX_FORCE_NS : 0),
DDBalanceRegionHandler(cr));
bool computeEnergy = false;
bool computeVirial = true;
dvdl_constr = 0;
- auto f = ems->f.arrayRefWithPadding();
+ auto f = ems->f.view().forceWithPadding();
constr->apply(needsLogging, computeEnergy, count, 0, 1.0, ems->s.x.arrayRefWithPadding(), f,
f.unpaddedArrayRef(), ems->s.box, ems->s.lambda[efptBONDED], &dvdl_constr,
gmx::ArrayRefWithPadding<RVec>(), computeVirial, shake_vir,
static double reorder_partsum(const t_commrec* cr,
t_grpopts* opts,
const gmx_mtop_t* top_global,
- em_state_t* s_min,
- em_state_t* s_b)
+ const em_state_t* s_min,
+ const em_state_t* s_b)
{
if (debug)
{
fprintf(debug, "Doing reorder_partsum\n");
}
- const rvec* fm = s_min->f.rvec_array();
- const rvec* fb = s_b->f.rvec_array();
+ auto fm = s_min->f.view().force();
+ auto fb = s_b->f.view().force();
/* Collect fm in a global vector fmg.
* This conflicts with the spirit of domain decomposition,
t_grpopts* opts,
t_mdatoms* mdatoms,
const gmx_mtop_t* top_global,
- em_state_t* s_min,
- em_state_t* s_b)
+ const em_state_t* s_min,
+ const em_state_t* s_b)
{
double sum;
if (!DOMAINDECOMP(cr)
|| (s_min->s.ddp_count == cr->dd->ddp_count && s_b->s.ddp_count == cr->dd->ddp_count))
{
- const rvec* fm = s_min->f.rvec_array();
- const rvec* fb = s_b->f.rvec_array();
- sum = 0;
- int gf = 0;
+ auto fm = s_min->f.view().force();
+ auto fb = s_b->f.view().force();
+ sum = 0;
+ int gf = 0;
/* This part of code can be incorrect with DD,
* since the atom ordering in s_b and s_min might differ.
*/
*/
/* Calculate the new direction in p, and the gradient in this direction, gpa */
- rvec* pm = s_min->s.cg_p.rvec_array();
- const rvec* sfm = s_min->f.rvec_array();
- double gpa = 0;
- int gf = 0;
+ gmx::ArrayRef<gmx::RVec> pm = s_min->s.cg_p;
+ gmx::ArrayRef<const gmx::RVec> sfm = s_min->f.view().force();
+ double gpa = 0;
+ int gf = 0;
for (int i = 0; i < mdatoms->homenr; i++)
{
if (mdatoms->cFREEZE)
}
/* Take a trial step (new coords in s_c) */
- do_em_step(cr, inputrec, mdatoms, s_min, c, &s_min->s.cg_p, s_c, constr, -1);
+ do_em_step(cr, inputrec, mdatoms, s_min, c, s_min->s.cg_p.constArrayRefWithPadding(), s_c,
+ constr, -1);
neval++;
/* Calculate energy for the trial step */
energyEvaluator.run(s_c, mu_tot, vir, pres, -1, FALSE);
/* Calc derivative along line */
- const rvec* pc = s_c->s.cg_p.rvec_array();
- const rvec* sfc = s_c->f.rvec_array();
- double gpc = 0;
+ const rvec* pc = s_c->s.cg_p.rvec_array();
+ gmx::ArrayRef<const gmx::RVec> sfc = s_c->f.view().force();
+ double gpc = 0;
for (int i = 0; i < mdatoms->homenr; i++)
{
for (m = 0; m < DIM; m++)
}
/* Take a trial step to this new point - new coords in s_b */
- do_em_step(cr, inputrec, mdatoms, s_min, b, &s_min->s.cg_p, s_b, constr, -1);
+ do_em_step(cr, inputrec, mdatoms, s_min, b,
+ s_min->s.cg_p.constArrayRefWithPadding(), s_b, constr, -1);
neval++;
/* Calculate energy for the trial step */
/* p does not change within a step, but since the domain decomposition
* might change, we have to use cg_p of s_b here.
*/
- const rvec* pb = s_b->s.cg_p.rvec_array();
- const rvec* sfb = s_b->f.rvec_array();
- gpb = 0;
+ const rvec* pb = s_b->s.cg_p.rvec_array();
+ gmx::ArrayRef<const gmx::RVec> sfb = s_b->f.view().force();
+ gpb = 0;
for (int i = 0; i < mdatoms->homenr; i++)
{
for (m = 0; m < DIM; m++)
point = 0;
// Set initial search direction to the force (-gradient), or 0 for frozen particles.
- real* fInit = static_cast<real*>(ems.f.rvec_array()[0]);
+ real* fInit = static_cast<real*>(ems.f.view().force().data()[0]);
for (i = 0; i < n; i++)
{
if (!frozen[i])
mdoutf_write_to_trajectory_files(fplog, cr, outf, mdof_flags, top_global->natoms, step,
static_cast<real>(step), &ems.s, state_global,
- observablesHistory, ems.f);
+ observablesHistory, ems.f.view().force());
/* Do the linesearching in the direction dx[point][0..(n-1)] */
s = dx[point];
real* xx = static_cast<real*>(ems.s.x.rvec_array()[0]);
- real* ff = static_cast<real*>(ems.f.rvec_array()[0]);
+ real* ff = static_cast<real*>(ems.f.view().force().data()[0]);
// calculate line gradient in position A
for (gpa = 0, i = 0; i < n; i++)
// Before taking any steps along the line, store the old position
*last = ems;
real* lastx = static_cast<real*>(last->s.x.data()[0]);
- real* lastf = static_cast<real*>(last->f.data()[0]);
+ real* lastf = static_cast<real*>(last->f.view().force().data()[0]);
Epot0 = ems.epot;
*sa = ems;
energyEvaluator.run(sc, mu_tot, vir, pres, step, FALSE);
// Calc line gradient in position C
- real* fc = static_cast<real*>(sc->f.rvec_array()[0]);
+ real* fc = static_cast<real*>(sc->f.view().force()[0]);
for (gpc = 0, i = 0; i < n; i++)
{
gpc -= s[i] * fc[i]; /* f is negative gradient, thus the sign */
fnorm = sb->fnorm;
// Calculate gradient in point B
- real* fb = static_cast<real*>(sb->f.rvec_array()[0]);
+ real* fb = static_cast<real*>(sb->f.view().force()[0]);
for (gpb = 0, i = 0; i < n; i++)
{
gpb -= s[i] * fb[i]; /* f is negative gradient, thus the sign */
bool validStep = true;
if (count > 0)
{
- validStep = do_em_step(cr, inputrec, mdatoms, s_min, stepsize, &s_min->f, s_try, constr, count);
+ validStep = do_em_step(cr, inputrec, mdatoms, s_min, stepsize,
+ s_min->f.view().forceWithPadding(), s_try, constr, count);
}
if (validStep)
/* Steps are divided one by one over the nodes */
bool bNS = true;
auto state_work_x = makeArrayRef(state_work.s.x);
- auto state_work_f = makeArrayRef(state_work.f);
+ auto state_work_f = state_work.f.view().force();
for (index aid = cr->nodeid; aid < ssize(atom_index); aid += nnodes)
{
size_t atom = atom_index[aid];
if (shellfc)
{
/* Now is the time to relax the shells */
- relax_shell_flexcon(
- fplog, cr, ms, mdrunOptions.verbose, nullptr, step, inputrec, imdSession,
- pull_work, bNS, force_flags, &top, constr, enerd, state_work.s.natoms,
- state_work.s.x.arrayRefWithPadding(), state_work.s.v.arrayRefWithPadding(),
- state_work.s.box, state_work.s.lambda, &state_work.s.hist,
- state_work.f.arrayRefWithPadding(), vir, mdatoms, nrnb, wcycle, shellfc,
- fr, runScheduleWork, t, mu_tot, vsite, DDBalanceRegionHandler(nullptr));
+ relax_shell_flexcon(fplog, cr, ms, mdrunOptions.verbose, nullptr, step, inputrec,
+ imdSession, pull_work, bNS, force_flags, &top, constr, enerd,
+ state_work.s.natoms, state_work.s.x.arrayRefWithPadding(),
+ state_work.s.v.arrayRefWithPadding(), state_work.s.box,
+ state_work.s.lambda, &state_work.s.hist, &state_work.f.view(),
+ vir, mdatoms, nrnb, wcycle, shellfc, fr, runScheduleWork, t,
+ mu_tot, vsite, DDBalanceRegionHandler(nullptr));
bNS = false;
step++;
}
#include "gromacs/mdtypes/commrec.h"
#include "gromacs/mdtypes/df_history.h"
#include "gromacs/mdtypes/energyhistory.h"
+#include "gromacs/mdtypes/forcebuffers.h"
#include "gromacs/mdtypes/forcerec.h"
#include "gromacs/mdtypes/group.h"
#include "gromacs/mdtypes/inputrec.h"
// alias to avoid a large ripple of nearly useless changes.
// t_inputrec is being replaced by IMdpOptionsProvider, so this
// will go away eventually.
- t_inputrec* ir = inputrec;
- int64_t step, step_rel;
- double t, lam0[efptNR];
- bool isLastStep = false;
- bool doFreeEnergyPerturbation = false;
- unsigned int force_flags;
- tensor force_vir, shake_vir, total_vir, pres;
- t_trxstatus* status = nullptr;
- rvec mu_tot;
- t_trxframe rerun_fr;
- gmx_localtop_t top(top_global->ffparams);
- PaddedHostVector<gmx::RVec> f{};
- gmx_global_stat_t gstat;
- gmx_shellfc_t* shellfc;
+ t_inputrec* ir = inputrec;
+ int64_t step, step_rel;
+ double t, lam0[efptNR];
+ bool isLastStep = false;
+ bool doFreeEnergyPerturbation = false;
+ unsigned int force_flags;
+ tensor force_vir, shake_vir, total_vir, pres;
+ t_trxstatus* status = nullptr;
+ rvec mu_tot;
+ t_trxframe rerun_fr;
+ gmx_localtop_t top(top_global->ffparams);
+ ForceBuffers f;
+ gmx_global_stat_t gstat;
+ gmx_shellfc_t* shellfc;
double cycles;
relax_shell_flexcon(fplog, cr, ms, mdrunOptions.verbose, enforcedRotation, step, ir,
imdSession, pull_work, bNS, force_flags, &top, constr, enerd,
state->natoms, state->x.arrayRefWithPadding(),
- state->v.arrayRefWithPadding(), state->box, state->lambda, &state->hist,
- f.arrayRefWithPadding(), force_vir, mdatoms, nrnb, wcycle, shellfc,
+ state->v.arrayRefWithPadding(), state->box, state->lambda,
+ &state->hist, &f.view(), force_vir, mdatoms, nrnb, wcycle, shellfc,
fr, runScheduleWork, t, mu_tot, vsite, ddBalanceRegionHandler);
}
else
gmx_edsam* ed = nullptr;
do_force(fplog, cr, ms, ir, awh, enforcedRotation, imdSession, pull_work, step, nrnb,
wcycle, &top, state->box, state->x.arrayRefWithPadding(), &state->hist,
- f.arrayRefWithPadding(), force_vir, mdatoms, enerd, state->lambda, fr, runScheduleWork,
+ &f.view(), force_vir, mdatoms, enerd, state->lambda, fr, runScheduleWork,
vsite, mu_tot, t, ed, GMX_FORCE_NS | force_flags, ddBalanceRegionHandler);
}
const bool bSumEkinhOld = false;
do_md_trajectory_writing(fplog, cr, nfile, fnm, step, step_rel, t, ir, state,
state_global, observablesHistory, top_global, fr, outf,
- energyOutput, ekind, f, isCheckpointingStep, doRerun,
- isLastStep, mdrunOptions.writeConfout, bSumEkinhOld);
+ energyOutput, ekind, f.view().force(), isCheckpointingStep,
+ doRerun, isLastStep, mdrunOptions.writeConfout, bSumEkinhOld);
}
stopHandler->setSignal();
#include "gromacs/mdlib/vsite.h"
#include "gromacs/mdtypes/commrec.h"
#include "gromacs/mdtypes/enerdata.h"
+#include "gromacs/mdtypes/forcebuffers.h"
#include "gromacs/mdtypes/forcerec.h"
#include "gromacs/mdtypes/inputrec.h"
#include "gromacs/mdtypes/md_enums.h"
const matrix box,
ArrayRef<real> lambda,
history_t* hist,
- ArrayRefWithPadding<RVec> f,
+ gmx::ForceBuffersView* f,
tensor force_vir,
const t_mdatoms* md,
t_nrnb* nrnb,
{
pr_rvecs(debug, 0, "x b4 do_force", as_rvec_array(x.data()), homenr);
}
- int shellfc_flags = force_flags | (bVerbose ? GMX_FORCE_ENERGY : 0);
+ int shellfc_flags = force_flags | (bVerbose ? GMX_FORCE_ENERGY : 0);
+ gmx::ForceBuffersView forceViewInit = gmx::ForceBuffersView(forceWithPadding[Min]);
do_force(fplog, cr, ms, inputrec, nullptr, enforcedRotation, imdSession, pull_work, mdstep,
- nrnb, wcycle, top, box, xPadded, hist, forceWithPadding[Min], force_vir, md, enerd,
- lambda, fr, runScheduleWork, vsite, mu_tot, t, nullptr,
+ nrnb, wcycle, top, box, xPadded, hist, &forceViewInit, force_vir, md, enerd, lambda,
+ fr, runScheduleWork, vsite, mu_tot, t, nullptr,
(bDoNS ? GMX_FORCE_NS : 0) | shellfc_flags, ddBalanceRegionHandler);
sf_dir = 0;
pr_rvecs(debug, 0, "RELAX: pos[Try] ", as_rvec_array(pos[Try].data()), homenr);
}
/* Try the new positions */
- do_force(fplog, cr, ms, inputrec, nullptr, enforcedRotation, imdSession, pull_work, 1, nrnb,
- wcycle, top, box, posWithPadding[Try], hist, forceWithPadding[Try], force_vir, md,
- enerd, lambda, fr, runScheduleWork, vsite, mu_tot, t, nullptr, shellfc_flags,
- ddBalanceRegionHandler);
+ gmx::ForceBuffersView forceViewTry = gmx::ForceBuffersView(forceWithPadding[Try]);
+ do_force(fplog, cr, ms, inputrec, nullptr, enforcedRotation, imdSession, pull_work, 1, nrnb, wcycle,
+ top, box, posWithPadding[Try], hist, &forceViewTry, force_vir, md, enerd, lambda, fr,
+ runScheduleWork, vsite, mu_tot, t, nullptr, shellfc_flags, ddBalanceRegionHandler);
accumulatePotentialEnergies(enerd, lambda, inputrec->fepvals);
if (gmx_debug_at)
{
/* Copy back the coordinates and the forces */
std::copy(pos[Min].begin(), pos[Min].end(), x.data());
- std::copy(force[Min].begin(), force[Min].end(), f.unpaddedArrayRef().begin());
+ std::copy(force[Min].begin(), force[Min].end(), f->force().begin());
}
void done_shellfc(FILE* fplog, gmx_shellfc_t* shfc, int64_t numSteps)
template<typename>
class ArrayRefWithPadding;
class Constraints;
+class ForceBuffersView;
class ImdSession;
class MdrunScheduleWorkload;
class VirtualSitesHandler;
const matrix box,
gmx::ArrayRef<real> lambda,
history_t* hist,
- gmx::ArrayRefWithPadding<gmx::RVec> f,
+ gmx::ForceBuffersView* f,
tensor force_vir,
const t_mdatoms* md,
t_nrnb* nrnb,
#include "gromacs/mdlib/vsite.h"
#include "gromacs/mdrunutility/printtime.h"
#include "gromacs/mdtypes/commrec.h"
+#include "gromacs/mdtypes/forcebuffers.h"
#include "gromacs/mdtypes/forcerec.h"
#include "gromacs/mdtypes/group.h"
#include "gromacs/mdtypes/inputrec.h"
{
GMX_RELEASE_ASSERT(gmx_omp_nthreads_get(emntDefault) == 1, "TPI does not support OpenMP");
- gmx_localtop_t top(top_global->ffparams);
- PaddedHostVector<gmx::RVec> f{};
- real lambda, t, temp, beta, drmax, epot;
- double embU, sum_embU, *sum_UgembU, V, V_all, VembU_all;
- t_trxstatus* status;
- t_trxframe rerun_fr;
- gmx_bool bDispCorr, bCharge, bRFExcl, bNotLastFrame, bStateChanged, bNS;
- tensor force_vir, shake_vir, vir, pres;
- int a_tp0, a_tp1, ngid, gid_tp, nener, e;
- rvec* x_mol;
- rvec mu_tot, x_init, dx;
- int nnodes, frame;
- int64_t frame_step_prev, frame_step;
- int64_t nsteps, stepblocksize = 0, step;
- int64_t seed;
- int i;
- FILE* fp_tpi = nullptr;
- char * ptr, *dump_pdb, **leg, str[STRLEN], str2[STRLEN];
- double dbl, dump_ener;
- gmx_bool bCavity;
- int nat_cavity = 0, d;
- real * mass_cavity = nullptr, mass_tot;
- int nbin;
- double invbinw, *bin, refvolshift, logV, bUlogV;
- gmx_bool bEnergyOutOfBounds;
- const char* tpid_leg[2] = { "direct", "reweighted" };
- auto mdatoms = mdAtoms->mdatoms();
+ gmx_localtop_t top(top_global->ffparams);
+ gmx::ForceBuffers f;
+ real lambda, t, temp, beta, drmax, epot;
+ double embU, sum_embU, *sum_UgembU, V, V_all, VembU_all;
+ t_trxstatus* status;
+ t_trxframe rerun_fr;
+ gmx_bool bDispCorr, bCharge, bRFExcl, bNotLastFrame, bStateChanged, bNS;
+ tensor force_vir, shake_vir, vir, pres;
+ int a_tp0, a_tp1, ngid, gid_tp, nener, e;
+ rvec* x_mol;
+ rvec mu_tot, x_init, dx;
+ int nnodes, frame;
+ int64_t frame_step_prev, frame_step;
+ int64_t nsteps, stepblocksize = 0, step;
+ int64_t seed;
+ int i;
+ FILE* fp_tpi = nullptr;
+ char * ptr, *dump_pdb, **leg, str[STRLEN], str2[STRLEN];
+ double dbl, dump_ener;
+ gmx_bool bCavity;
+ int nat_cavity = 0, d;
+ real * mass_cavity = nullptr, mass_tot;
+ int nbin;
+ double invbinw, *bin, refvolshift, logV, bUlogV;
+ gmx_bool bEnergyOutOfBounds;
+ const char* tpid_leg[2] = { "direct", "reweighted" };
+ auto mdatoms = mdAtoms->mdatoms();
GMX_UNUSED_VALUE(outputProvider);
atoms2md(top_global, inputrec, -1, {}, top_global->natoms, mdAtoms);
update_mdatoms(mdatoms, inputrec->fepvals->init_lambda);
- f.resizeWithPadding(top_global->natoms);
+ f.resize(top_global->natoms);
/* Print to log file */
walltime_accounting_start_time(walltime_accounting);
std::feholdexcept(&floatingPointEnvironment);
do_force(fplog, cr, ms, inputrec, nullptr, nullptr, imdSession, pull_work, step, nrnb,
wcycle, &top, state_global->box, state_global->x.arrayRefWithPadding(),
- &state_global->hist, f.arrayRefWithPadding(), force_vir, mdatoms, enerd,
+ &state_global->hist, &f.view(), force_vir, mdatoms, enerd,
state_global->lambda, fr, runScheduleWork, nullptr, mu_tot, t, nullptr,
GMX_FORCE_NONBONDED | GMX_FORCE_ENERGY | (bStateChanged ? GMX_FORCE_STATECHANGED : 0),
DDBalanceRegionHandler(nullptr));
file(GLOB MDTYPES_SOURCES
df_history.cpp
+ forcebuffers.cpp
group.cpp
iforceprovider.cpp
inputrec.cpp
)
endif()
+if (BUILD_TESTING)
+ add_subdirectory(tests)
+endif()
set(LIBGROMACS_SOURCES ${LIBGROMACS_SOURCES} ${MDTYPES_SOURCES} PARENT_SCOPE)
--- /dev/null
+/*
+ * This file is part of the GROMACS molecular simulation package.
+ *
+ * Copyright (c) 2020, 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.
+ *
+ * GROMACS is free software; you can redistribute it and/or
+ * modify it under the terms of the GNU Lesser General Public License
+ * as published by the Free Software Foundation; either version 2.1
+ * of the License, or (at your option) any later version.
+ *
+ * GROMACS is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ * Lesser General Public License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public
+ * License along with GROMACS; if not, see
+ * http://www.gnu.org/licenses, or write to the Free Software Foundation,
+ * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
+ *
+ * If you want to redistribute modifications to GROMACS, please
+ * consider that scientific software is very special. Version
+ * control is crucial - bugs must be traceable. We will be happy to
+ * consider code for inclusion in the official distribution, but
+ * derived work must not be called official GROMACS. Details are found
+ * in the README & COPYING files - if they are missing, get the
+ * official version at http://www.gromacs.org.
+ *
+ * To help us fund GROMACS development, we humbly ask that you cite
+ * the research papers on the package. Check out http://www.gromacs.org.
+ */
+/*! \internal \file
+ * \brief
+ * Implements the ForceBuffers class
+ *
+ * \author Berk Hess <hess@kth.se>
+ * \ingroup module_mdtypes
+ */
+
+#include "gmxpre.h"
+
+#include "forcebuffers.h"
+
+#include "gromacs/gpu_utils/hostallocator.h"
+
+namespace gmx
+{
+
+ForceBuffers::ForceBuffers(PinningPolicy pinningPolicy) : force_({}, { pinningPolicy }), view_({})
+{
+}
+
+ForceBuffers::~ForceBuffers() = default;
+
+ForceBuffers& ForceBuffers::operator=(ForceBuffers const& o)
+{
+ auto oForce = o.view().force();
+ resize(oForce.size());
+ std::copy(oForce.begin(), oForce.end(), view().force().begin());
+
+ return *this;
+}
+
+PinningPolicy ForceBuffers::pinningPolicy() const
+{
+ return force_.get_allocator().pinningPolicy();
+}
+
+void ForceBuffers::resize(int numAtoms)
+{
+ force_.resizeWithPadding(numAtoms);
+ view_ = ForceBuffersView(force_.arrayRefWithPadding());
+}
+
+} // namespace gmx
--- /dev/null
+/*
+ * This file is part of the GROMACS molecular simulation package.
+ *
+ * Copyright (c) 2020, 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.
+ *
+ * GROMACS is free software; you can redistribute it and/or
+ * modify it under the terms of the GNU Lesser General Public License
+ * as published by the Free Software Foundation; either version 2.1
+ * of the License, or (at your option) any later version.
+ *
+ * GROMACS is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ * Lesser General Public License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public
+ * License along with GROMACS; if not, see
+ * http://www.gnu.org/licenses, or write to the Free Software Foundation,
+ * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
+ *
+ * If you want to redistribute modifications to GROMACS, please
+ * consider that scientific software is very special. Version
+ * control is crucial - bugs must be traceable. We will be happy to
+ * consider code for inclusion in the official distribution, but
+ * derived work must not be called official GROMACS. Details are found
+ * in the README & COPYING files - if they are missing, get the
+ * official version at http://www.gromacs.org.
+ *
+ * To help us fund GROMACS development, we humbly ask that you cite
+ * the research papers on the package. Check out http://www.gromacs.org.
+ */
+
+/*! \libinternal \file
+ *
+ * \brief
+ * This file contains the definition of a container for force buffers.
+ *
+ * \author Berk Hess
+ *
+ * \inlibraryapi
+ * \ingroup module_mdtypes
+ */
+
+#ifndef GMX_MDTYPES_FORCEBUFFERS_H
+#define GMX_MDTYPES_FORCEBUFFERS_H
+
+#include "gromacs/gpu_utils/hostallocator.h"
+#include "gromacs/math/arrayrefwithpadding.h"
+#include "gromacs/math/vectypes.h"
+#include "gromacs/utility/arrayref.h"
+#include "gromacs/utility/classhelpers.h"
+
+namespace gmx
+{
+
+enum class PinningPolicy : int;
+
+/*! \libinternal \brief A view of the force buffer
+ *
+ * This is used for read and/or write access to the force buffer.
+ */
+class ForceBuffersView
+{
+public:
+ //! Constructor, creates a view to \p force
+ ForceBuffersView(const ArrayRefWithPadding<RVec>& force) : force_(force) {}
+
+ //! Copy constructor deleted to avoid creating non-const from const
+ ForceBuffersView(const ForceBuffersView& o) = delete;
+
+ //! Move constructor deleted, but could be implemented
+ ForceBuffersView(ForceBuffersView&& o) = delete;
+
+ ~ForceBuffersView() = default;
+
+ //! Copy assignment deleted to avoid creating non-const from const
+ ForceBuffersView& operator=(const ForceBuffersView& v) = delete;
+
+ //! Move assignment, moves the view
+ ForceBuffersView& operator=(ForceBuffersView&& v) = default;
+
+ //! Returns a const arrayref to the force buffer without padding
+ ArrayRef<const RVec> force() const { return force_.unpaddedConstArrayRef(); }
+
+ //! Returns an arrayref to the force buffer without padding
+ ArrayRef<RVec> force() { return force_.unpaddedArrayRef(); }
+
+ //! Returns an ArrayRefWithPadding to the force buffer
+ ArrayRefWithPadding<RVec> forceWithPadding() { return force_; }
+
+private:
+ //! The force buffer
+ ArrayRefWithPadding<RVec> force_;
+};
+
+/*! \libinternal \brief Object that holds the force buffers
+ *
+ * More buffers can be added when needed. Those should also be added
+ * to ForceBuffersView.
+ * Can be pinned for efficient transfer to/from GPUs.
+ * All access happens through the ForceBuffersView object.
+ */
+class ForceBuffers
+{
+public:
+ //! Constructor, creates an empty force buffer with pinning not active
+ ForceBuffers(PinningPolicy pinningPolicy = PinningPolicy::CannotBePinned);
+
+ //! Copy constructor deleted, but could be implemented
+ ForceBuffers(const ForceBuffers& o) = delete;
+
+ //! Move constructor deleted, but could be implemented
+ ForceBuffers(ForceBuffers&& o) = delete;
+
+ ~ForceBuffers();
+
+ //! Copy assignment operator, sets the pinning policy to CannotBePinned
+ ForceBuffers& operator=(ForceBuffers const& o);
+
+ //! Move assignment operator, deleted but could be implemented
+ ForceBuffers& operator=(ForceBuffers&& o) = delete;
+
+ //! Returns a const view to the force buffers
+ const ForceBuffersView& view() const { return view_; }
+
+ //! Returns a view to the force buffer
+ ForceBuffersView& view() { return view_; }
+
+ //! Resize the force buffer
+ void resize(int numAtoms);
+
+ /*! \brief Returns the active pinning policy.
+ *
+ * Does not throw.
+ */
+ PinningPolicy pinningPolicy() const;
+
+private:
+ //! The force buffer
+ PaddedHostVector<RVec> force_;
+ //! The view to the force buffer
+ ForceBuffersView view_;
+};
+
+} // namespace gmx
+
+#endif
--- /dev/null
+#
+# This file is part of the GROMACS molecular simulation package.
+#
+# Copyright (c) 2020, 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.
+#
+# GROMACS is free software; you can redistribute it and/or
+# modify it under the terms of the GNU Lesser General Public License
+# as published by the Free Software Foundation; either version 2.1
+# of the License, or (at your option) any later version.
+#
+# GROMACS is distributed in the hope that it will be useful,
+# but WITHOUT ANY WARRANTY; without even the implied warranty of
+# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+# Lesser General Public License for more details.
+#
+# You should have received a copy of the GNU Lesser General Public
+# License along with GROMACS; if not, see
+# http://www.gnu.org/licenses, or write to the Free Software Foundation,
+# Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
+#
+# If you want to redistribute modifications to GROMACS, please
+# consider that scientific software is very special. Version
+# control is crucial - bugs must be traceable. We will be happy to
+# consider code for inclusion in the official distribution, but
+# derived work must not be called official GROMACS. Details are found
+# in the README & COPYING files - if they are missing, get the
+# official version at http://www.gromacs.org.
+#
+# To help us fund GROMACS development, we humbly ask that you cite
+# the research papers on the package. Check out http://www.gromacs.org.
+
+gmx_add_unit_test(MdtypesUnitTest mdtypes-test
+ CPP_SOURCE_FILES
+ forcebuffers.cpp
+ )
--- /dev/null
+/*
+ * This file is part of the GROMACS molecular simulation package.
+ *
+ * Copyright (c) 2020, 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.
+ *
+ * GROMACS is free software; you can redistribute it and/or
+ * modify it under the terms of the GNU Lesser General Public License
+ * as published by the Free Software Foundation; either version 2.1
+ * of the License, or (at your option) any later version.
+ *
+ * GROMACS is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ * Lesser General Public License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public
+ * License along with GROMACS; if not, see
+ * http://www.gnu.org/licenses, or write to the Free Software Foundation,
+ * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
+ *
+ * If you want to redistribute modifications to GROMACS, please
+ * consider that scientific software is very special. Version
+ * control is crucial - bugs must be traceable. We will be happy to
+ * consider code for inclusion in the official distribution, but
+ * derived work must not be called official GROMACS. Details are found
+ * in the README & COPYING files - if they are missing, get the
+ * official version at http://www.gromacs.org.
+ *
+ * To help us fund GROMACS development, we humbly ask that you cite
+ * the research papers on the package. Check out http://www.gromacs.org.
+ */
+/*! \internal \file
+ * \brief
+ * Tests for the ForceBuffers and ForceBuffersView classes.
+ *
+ * \author berk Hess <hess@kth.se>
+ * \ingroup module_mdtypes
+ */
+#include "gmxpre.h"
+
+#include "gromacs/mdtypes/forcebuffers.h"
+
+#include <array>
+
+#include <gmock/gmock.h>
+#include <gtest/gtest.h>
+
+#include "testutils/testasserts.h"
+
+namespace gmx
+{
+
+const std::array<RVec, 2> c_forces = { { { 0.5, 0.1, 1.2 }, { -2.1, 0.2, 0.3 } } };
+
+TEST(ForceBuffers, ConstructsUnpinned)
+{
+ ForceBuffers forceBuffers;
+
+ EXPECT_EQ(forceBuffers.pinningPolicy(), PinningPolicy::CannotBePinned);
+}
+
+TEST(ForceBuffers, ConstructsPinned)
+{
+ ForceBuffers forceBuffers(PinningPolicy::PinnedIfSupported);
+
+ EXPECT_EQ(forceBuffers.pinningPolicy(), PinningPolicy::PinnedIfSupported);
+}
+
+TEST(ForceBuffers, ConstructsEmpty)
+{
+ ForceBuffers forceBuffers;
+
+ EXPECT_EQ(forceBuffers.view().force().size(), 0);
+}
+
+TEST(ForceBuffers, ResizeWorks)
+{
+ ForceBuffers forceBuffers;
+
+ forceBuffers.resize(2);
+ EXPECT_EQ(forceBuffers.view().force().size(), 2);
+}
+
+TEST(ForceBuffers, PaddingWorks)
+{
+ ForceBuffers forceBuffers;
+
+ forceBuffers.resize(2);
+ auto paddedRef = forceBuffers.view().forceWithPadding();
+ EXPECT_EQ(paddedRef.unpaddedArrayRef().size(), 2);
+ EXPECT_GT(paddedRef.size(), 2);
+}
+
+
+TEST(ForceBuffers, CopyWorks)
+{
+ ForceBuffers forceBuffers;
+
+ forceBuffers.resize(2);
+ auto force = forceBuffers.view().force();
+ index i = 0;
+ for (RVec& v : force)
+ {
+ v = c_forces[i];
+ i++;
+ }
+
+ ForceBuffers forceBuffersCopy;
+ forceBuffersCopy = forceBuffers;
+ auto forceCopy = forceBuffersCopy.view().force();
+ EXPECT_EQ(forceBuffersCopy.view().force().size(), 2);
+ for (index i = 0; i < ssize(forceCopy); i++)
+ {
+ for (int d = 0; d < DIM; d++)
+ {
+ EXPECT_EQ(forceCopy[i][d], force[i][d]);
+ }
+ }
+}
+
+TEST(ForceBuffers, CopyDoesNotPin)
+{
+ ForceBuffers forceBuffers(PinningPolicy::PinnedIfSupported);
+
+ ForceBuffers forceBuffersCopy;
+ forceBuffersCopy = forceBuffers;
+ EXPECT_EQ(forceBuffersCopy.pinningPolicy(), PinningPolicy::CannotBePinned);
+}
+
+} // namespace gmx
{
std::unique_ptr<t_state> localState = statePropagatorData_->localState();
t_state* globalState = statePropagatorData_->globalState();
- PaddedHostVector<RVec>* forcePointer = statePropagatorData_->forcePointer();
+ ForceBuffers* forcePointer = statePropagatorData_->forcePointer();
// constant choices for this call to dd_partition_system
const bool verbose = false;
}
std::unique_ptr<t_state> localState = statePropagatorData_->localState();
t_state* globalState = statePropagatorData_->globalState();
- PaddedHostVector<RVec>* forcePointer = statePropagatorData_->forcePointer();
+ ForceBuffers* forcePointer = statePropagatorData_->forcePointer();
// constant choices for this call to dd_partition_system
const bool verbose = isVerbose_ && (step % verbosePrintInterval_ == 0 || step == inputrec_->init_step);
#include "gromacs/mdlib/force_flags.h"
#include "gromacs/mdlib/mdatoms.h"
#include "gromacs/mdrun/shellfc.h"
+#include "gromacs/mdtypes/forcebuffers.h"
#include "gromacs/mdtypes/inputrec.h"
#include "gromacs/mdtypes/mdatom.h"
#include "gromacs/mdtypes/mdrunoptions.h"
* Check comments in sim_util.c
*/
auto x = statePropagatorData_->positionsView();
- auto forces = statePropagatorData_->forcesView();
+ auto& forces = statePropagatorData_->forcesView();
auto box = statePropagatorData_->constBox();
history_t* hist = nullptr; // disabled
fplog_, cr_, ms, isVerbose_, enforcedRotation_, step, inputrec_, imdSession_,
pull_work_, step == nextNSStep_, static_cast<int>(flags), localTopology_, constr_,
energyData_->enerdata(), statePropagatorData_->localNumAtoms(), x, v, box, lambda,
- hist, forces, force_vir, mdAtoms_->mdatoms(), nrnb_, wcycle_, shellfc_, fr_,
+ hist, &forces, force_vir, mdAtoms_->mdatoms(), nrnb_, wcycle_, shellfc_, fr_,
runScheduleWork_, time, energyData_->muTot(), vsite_, ddBalanceRegionHandler_);
nShellRelaxationSteps_++;
}
gmx_edsam* ed = nullptr;
do_force(fplog_, cr_, ms, inputrec_, awh, enforcedRotation_, imdSession_, pull_work_, step,
- nrnb_, wcycle_, localTopology_, box, x, hist, forces, force_vir,
+ nrnb_, wcycle_, localTopology_, box, x, hist, &forces, force_vir,
mdAtoms_->mdatoms(), energyData_->enerdata(), lambda, fr_, runScheduleWork_, vsite_,
energyData_->muTot(), time, ed, static_cast<int>(flags), ddBalanceRegionHandler_);
}
wallcycle_start(wcycle_, ewcUPDATE);
auto v = as_rvec_array(statePropagatorData_->velocitiesView().paddedArrayRef().data());
- auto f = as_rvec_array(statePropagatorData_->constForcesView().paddedArrayRef().data());
+ auto f = as_rvec_array(statePropagatorData_->constForcesView().force().data());
auto invMassPerDim = mdAtoms_->mdatoms()->invMassPerDim;
const real lambda =
auto xp = as_rvec_array(statePropagatorData_->positionsView().paddedArrayRef().data());
auto x = as_rvec_array(statePropagatorData_->constPreviousPositionsView().paddedArrayRef().data());
auto v = as_rvec_array(statePropagatorData_->velocitiesView().paddedArrayRef().data());
- auto f = as_rvec_array(statePropagatorData_->constForcesView().paddedArrayRef().data());
+ auto f = as_rvec_array(statePropagatorData_->constForcesView().force().data());
auto invMassPerDim = mdAtoms_->mdatoms()->invMassPerDim;
const real lambda =
auto xp = as_rvec_array(statePropagatorData_->positionsView().paddedArrayRef().data());
auto x = as_rvec_array(statePropagatorData_->constPreviousPositionsView().paddedArrayRef().data());
auto v = as_rvec_array(statePropagatorData_->velocitiesView().paddedArrayRef().data());
- auto f = as_rvec_array(statePropagatorData_->constForcesView().paddedArrayRef().data());
+ auto f = as_rvec_array(statePropagatorData_->constForcesView().force().data());
auto invMassPerDim = mdAtoms_->mdatoms()->invMassPerDim;
int nth = gmx_omp_nthreads_get(emntUpdate);
#include "gromacs/mdlib/stat.h"
#include "gromacs/mdlib/update.h"
#include "gromacs/mdtypes/commrec.h"
+#include "gromacs/mdtypes/forcebuffers.h"
#include "gromacs/mdtypes/forcerec.h"
#include "gromacs/mdtypes/inputrec.h"
#include "gromacs/mdtypes/mdatom.h"
else
{
state_change_natoms(globalState, globalState->natoms);
- f_.resizeWithPadding(globalState->natoms);
+ f_.resize(globalState->natoms);
localNAtoms_ = globalState->natoms;
x_ = globalState->x;
v_ = globalState->v;
return v_.constArrayRefWithPadding();
}
-ArrayRefWithPadding<RVec> StatePropagatorData::forcesView()
+ForceBuffersView& StatePropagatorData::forcesView()
{
- return f_.arrayRefWithPadding();
+ return f_.view();
}
-ArrayRefWithPadding<const RVec> StatePropagatorData::constForcesView() const
+const ForceBuffersView& StatePropagatorData::constForcesView() const
{
- return f_.constArrayRefWithPadding();
+ return f_.view();
}
rvec* StatePropagatorData::box()
return globalState_;
}
-PaddedHostVector<RVec>* StatePropagatorData::forcePointer()
+ForceBuffers* StatePropagatorData::forcePointer()
{
return &f_;
}
mdoutf_write_to_trajectory_files(fplog_, cr_, outf, static_cast<int>(mdof_flags),
statePropagatorData_->totalNumAtoms_, currentStep, currentTime,
localStateBackup_.get(), statePropagatorData_->globalState_,
- observablesHistory, statePropagatorData_->f_);
+ observablesHistory, statePropagatorData_->f_.view().force());
if (currentStep != lastStep_ || !isRegularSimulationEnd_)
{
#include "gromacs/gpu_utils/hostallocator.h"
#include "gromacs/math/paddedvector.h"
#include "gromacs/math/vectypes.h"
+#include "gromacs/mdtypes/forcebuffers.h"
#include "modularsimulatorinterfaces.h"
#include "topologyholder.h"
//! Get read access to velocity vector
ArrayRefWithPadding<const RVec> constVelocitiesView() const;
//! Get write access to force vector
- ArrayRefWithPadding<RVec> forcesView();
+ ForceBuffersView& forcesView();
//! Get read access to force vector
- ArrayRefWithPadding<const RVec> constForcesView() const;
+ const ForceBuffersView& constForcesView() const;
//! Get pointer to box
rvec* box();
//! Get const pointer to box
//! The velocity vector
PaddedHostVector<RVec> v_;
//! The force vector
- PaddedHostVector<RVec> f_;
+ ForceBuffers f_;
//! The box matrix
matrix box_;
//! The box matrix of the previous step
//! Get a pointer to the global state
t_state* globalState();
//! Get a force pointer
- PaddedHostVector<gmx::RVec>* forcePointer();
+ ForceBuffers* forcePointer();
//! Whether we're doing VV and need to reset velocities after the first half step
bool vvResetVelocities_;