}
}
-void dd_init_bondeds(FILE* fplog,
- gmx_domdec_t* dd,
- const gmx_mtop_t* mtop,
- const gmx_vsite_t* vsite,
- const t_inputrec* ir,
- gmx_bool bBCheck,
- cginfo_mb_t* cginfo_mb)
+void dd_init_bondeds(FILE* fplog,
+ gmx_domdec_t* dd,
+ const gmx_mtop_t& mtop,
+ const gmx_vsite_t* vsite,
+ const t_inputrec* ir,
+ gmx_bool bBCheck,
+ gmx::ArrayRef<cginfo_mb_t> cginfo_mb)
{
gmx_domdec_comm_t* comm;
- dd_make_reverse_top(fplog, dd, mtop, vsite, ir, bBCheck);
+ dd_make_reverse_top(fplog, dd, &mtop, vsite, ir, bBCheck);
comm = dd->comm;
bool is1DAnd1PulseDD(const gmx_domdec_t& dd);
/*! \brief Initialize data structures for bonded interactions */
-void dd_init_bondeds(FILE* fplog,
- gmx_domdec_t* dd,
- const gmx_mtop_t* mtop,
- const gmx_vsite_t* vsite,
- const t_inputrec* ir,
- gmx_bool bBCheck,
- cginfo_mb_t* cginfo_mb);
+void dd_init_bondeds(FILE* fplog,
+ gmx_domdec_t* dd,
+ const gmx_mtop_t& mtop,
+ const gmx_vsite_t* vsite,
+ const t_inputrec* ir,
+ gmx_bool bBCheck,
+ gmx::ArrayRef<cginfo_mb_t> cginfo_mb);
/*! \brief Returns whether molecules are always whole, i.e. not broken by PBC */
bool dd_moleculesAreAlwaysWhole(const gmx_domdec_t& dd);
*
* Also stores whether atoms are linked in \p cginfo_mb.
*/
-t_blocka* makeBondedLinks(const gmx_mtop_t* mtop, cginfo_mb_t* cginfo_mb);
+t_blocka* makeBondedLinks(const gmx_mtop_t& mtop, gmx::ArrayRef<cginfo_mb_t> cginfo_mb);
/*! \brief Calculate the maximum distance involved in 2-body and multi-body bonded interactions */
void dd_bonded_cg_distance(const gmx::MDLogger& mdlog,
}
}
-t_blocka* makeBondedLinks(const gmx_mtop_t* mtop, cginfo_mb_t* cginfo_mb)
+t_blocka* makeBondedLinks(const gmx_mtop_t& mtop, gmx::ArrayRef<cginfo_mb_t> cginfo_mb)
{
t_blocka* link;
cginfo_mb_t* cgi_mb;
*/
reverse_ilist_t ril_intermol;
- if (mtop->bIntermolecularInteractions)
+ if (mtop.bIntermolecularInteractions)
{
t_atoms atoms;
- atoms.nr = mtop->natoms;
+ atoms.nr = mtop.natoms;
atoms.atom = nullptr;
- GMX_RELEASE_ASSERT(mtop->intermolecular_ilist,
+ GMX_RELEASE_ASSERT(mtop.intermolecular_ilist,
"We should have an ilist when intermolecular interactions are on");
- make_reverse_ilist(*mtop->intermolecular_ilist, &atoms, FALSE, FALSE, FALSE, TRUE, &ril_intermol);
+ make_reverse_ilist(*mtop.intermolecular_ilist, &atoms, FALSE, FALSE, FALSE, TRUE, &ril_intermol);
}
snew(link, 1);
- snew(link->index, mtop->natoms + 1);
+ snew(link->index, mtop.natoms + 1);
link->nalloc_a = 0;
link->a = nullptr;
link->index[0] = 0;
int cg_offset = 0;
int ncgi = 0;
- for (size_t mb = 0; mb < mtop->molblock.size(); mb++)
+ for (size_t mb = 0; mb < mtop.molblock.size(); mb++)
{
- const gmx_molblock_t& molb = mtop->molblock[mb];
+ const gmx_molblock_t& molb = mtop.molblock[mb];
if (molb.nmol == 0)
{
continue;
}
- const gmx_moltype_t& molt = mtop->moltype[molb.type];
+ const gmx_moltype_t& molt = mtop.moltype[molb.type];
/* Make a reverse ilist in which the interactions are linked
* to all atoms, not only the first atom as in gmx_reverse_top.
* The constraints are discarded here.
cgi_mb = &cginfo_mb[mb];
int mol;
- for (mol = 0; mol < (mtop->bIntermolecularInteractions ? molb.nmol : 1); mol++)
+ for (mol = 0; mol < (mtop.bIntermolecularInteractions ? molb.nmol : 1); mol++)
{
for (int a = 0; a < molt.atoms.nr; a++)
{
i += nral_rt(ftype);
}
- if (mtop->bIntermolecularInteractions)
+ if (mtop.bIntermolecularInteractions)
{
int i = ril_intermol.index[a];
while (i < ril_intermol.index[a + 1])
if (debug)
{
- fprintf(debug, "Of the %d atoms %d are linked via bonded interactions\n", mtop->natoms, ncgi);
+ fprintf(debug, "Of the %d atoms %d are linked via bonded interactions\n", mtop.natoms, ncgi);
}
return link;
{
if (fr != nullptr)
{
- const cginfo_mb_t* cginfo_mb = fr->cginfo_mb;
- gmx::ArrayRef<int> cginfo = fr->cginfo;
+ gmx::ArrayRef<cginfo_mb_t> cginfo_mb = fr->cginfo_mb;
+ gmx::ArrayRef<int> cginfo = fr->cginfo;
for (int cg = cg0; cg < cg1; cg++)
{
const int* recv_i,
gmx::ArrayRef<gmx::RVec> x,
gmx::ArrayRef<const gmx::RVec> recv_vr,
- cginfo_mb_t* cginfo_mb,
+ gmx::ArrayRef<cginfo_mb_t> cginfo_mb,
gmx::ArrayRef<int> cginfo)
{
gmx_domdec_ind_t *ind, *ind_p;
gmx_domdec_comm_t* comm;
gmx_domdec_zones_t* zones;
gmx_domdec_comm_dim_t* cd;
- cginfo_mb_t* cginfo_mb;
gmx_bool bBondComm, bDist2B, bDistMB, bDistBonded;
dd_corners_t corners;
rvec * normal, *v_d, *v_0 = nullptr, *v_1 = nullptr;
v_1 = ddbox->v[dim1];
}
- zone_cg_range = zones->cg_range;
- cginfo_mb = fr->cginfo_mb;
+ zone_cg_range = zones->cg_range;
+ gmx::ArrayRef<cginfo_mb_t> cginfo_mb = fr->cginfo_mb;
zone_cg_range[0] = 0;
zone_cg_range[1] = dd->ncg_home;
/* We reuse the intBuffer without reacquiring since we are in the same scope */
DDBufferAccess<int>& flagBuffer = moveBuffer;
- const cginfo_mb_t* cginfo_mb = fr->cginfo_mb;
+ gmx::ArrayRef<const cginfo_mb_t> cginfo_mb = fr->cginfo_mb;
/* Temporarily store atoms passed to our rank at the end of the range */
int home_pos_cg = dd->ncg_home;
#include "gromacs/gpu_utils/hostallocator.h"
#include "gromacs/mdtypes/forcerec.h"
+#include "gromacs/utility/arrayref.h"
#include "domdec_internal.h"
void check_screw_box(const matrix box);
/*! \brief Return the charge group information flags for charge group cg */
-static inline int ddcginfo(const cginfo_mb_t* cginfo_mb, int cg)
+static inline int ddcginfo(gmx::ArrayRef<const cginfo_mb_t> cginfo_mb, int cg)
{
- while (cg >= cginfo_mb->cg_end)
+ size_t index = 0;
+ while (cg >= cginfo_mb[index].cg_end)
{
- cginfo_mb++;
+ index++;
}
+ const cginfo_mb_t& cgimb = cginfo_mb[index];
- return cginfo_mb->cginfo[(cg - cginfo_mb->cg_start) % cginfo_mb->cg_mod];
+ return cgimb.cginfo[(cg - cgimb.cg_start) % cgimb.cg_mod];
};
/*! \brief Returns the number of MD steps for which load has been recorded */
const int* typeA = mdatoms->typeA;
const int* typeB = mdatoms->typeB;
const int ntype = fr->ntype;
- const real* nbfp = fr->nbfp;
+ const real* nbfp = fr->nbfp.data();
const real* nbfp_grid = fr->ljpme_c6grid;
real* Vv = kernel_data->energygrp_vdw;
const real lambda_coul = kernel_data->lambda[efptCOUL];
static const bool c_enableGpuPmePpComms =
(getenv("GMX_GPU_PME_PP_COMMS") != nullptr) && GMX_THREAD_MPI && (GMX_GPU == GMX_GPU_CUDA);
-static real* mk_nbfp(const gmx_ffparams_t* idef, gmx_bool bBHAM)
+static std::vector<real> mk_nbfp(const gmx_ffparams_t* idef, gmx_bool bBHAM)
{
- real* nbfp;
- int i, j, k, atnr;
+ std::vector<real> nbfp;
+ int atnr;
atnr = idef->atnr;
if (bBHAM)
{
- snew(nbfp, 3 * atnr * atnr);
- for (i = k = 0; (i < atnr); i++)
+ nbfp.resize(3 * atnr * atnr);
+ int k = 0;
+ for (int i = 0; (i < atnr); i++)
{
- for (j = 0; (j < atnr); j++, k++)
+ for (int j = 0; (j < atnr); j++, k++)
{
BHAMA(nbfp, atnr, i, j) = idef->iparams[k].bham.a;
BHAMB(nbfp, atnr, i, j) = idef->iparams[k].bham.b;
}
else
{
- snew(nbfp, 2 * atnr * atnr);
- for (i = k = 0; (i < atnr); i++)
+ nbfp.resize(2 * atnr * atnr);
+ int k = 0;
+ for (int i = 0; (i < atnr); i++)
{
- for (j = 0; (j < atnr); j++, k++)
+ for (int j = 0; (j < atnr); j++, k++)
{
/* nbfp now includes the 6.0/12.0 derivative prefactors */
C6(nbfp, atnr, i, j) = idef->iparams[k].lj.c6 * 6.0;
acSETTLE
};
-static cginfo_mb_t* init_cginfo_mb(const gmx_mtop_t* mtop, const t_forcerec* fr, gmx_bool* bFEP_NonBonded)
+static std::vector<cginfo_mb_t> init_cginfo_mb(const gmx_mtop_t* mtop, const t_forcerec* fr, gmx_bool* bFEP_NonBonded)
{
- cginfo_mb_t* cginfo_mb;
- gmx_bool* type_VDW;
- int* cginfo;
- int* a_con;
-
- snew(cginfo_mb, mtop->molblock.size());
+ gmx_bool* type_VDW;
+ int* a_con;
snew(type_VDW, fr->ntype);
for (int ai = 0; ai < fr->ntype; ai++)
*bFEP_NonBonded = FALSE;
- int a_offset = 0;
+ std::vector<cginfo_mb_t> cginfoPerMolblock;
+ int a_offset = 0;
for (size_t mb = 0; mb < mtop->molblock.size(); mb++)
{
const gmx_molblock_t& molb = mtop->molblock[mb];
}
}
- cginfo_mb[mb].cg_start = a_offset;
- cginfo_mb[mb].cg_end = a_offset + molb.nmol * molt.atoms.nr;
- cginfo_mb[mb].cg_mod = (bId ? 1 : molb.nmol) * molt.atoms.nr;
- snew(cginfo_mb[mb].cginfo, cginfo_mb[mb].cg_mod);
- cginfo = cginfo_mb[mb].cginfo;
+ cginfo_mb_t cginfo_mb;
+ cginfo_mb.cg_start = a_offset;
+ cginfo_mb.cg_end = a_offset + molb.nmol * molt.atoms.nr;
+ cginfo_mb.cg_mod = (bId ? 1 : molb.nmol) * molt.atoms.nr;
+ cginfo_mb.cginfo.resize(cginfo_mb.cg_mod);
+ gmx::ArrayRef<int> cginfo = cginfo_mb.cginfo;
/* Set constraints flags for constrained atoms */
snew(a_con, molt.atoms.nr);
sfree(a_con);
+ cginfoPerMolblock.push_back(cginfo_mb);
+
a_offset += molb.nmol * molt.atoms.nr;
}
sfree(type_VDW);
- return cginfo_mb;
+ return cginfoPerMolblock;
}
-static std::vector<int> cginfo_expand(const int nmb, const cginfo_mb_t* cgi_mb)
+static std::vector<int> cginfo_expand(const int nmb, gmx::ArrayRef<const cginfo_mb_t> cgi_mb)
{
const int ncg = cgi_mb[nmb - 1].cg_end;
return cginfo;
}
-static void done_cginfo_mb(cginfo_mb_t* cginfo_mb, int numMolBlocks)
-{
- if (cginfo_mb == nullptr)
- {
- return;
- }
- for (int mb = 0; mb < numMolBlocks; ++mb)
- {
- sfree(cginfo_mb[mb].cginfo);
- }
- sfree(cginfo_mb);
-}
-
/* Sets the sum of charges (squared) and C6 in the system in fr.
* Returns whether the system has a net charge.
*/
fr->shiftForces.resize(SHIFTS);
- if (fr->nbfp == nullptr)
+ if (fr->nbfp.empty())
{
fr->ntype = mtop->ffparams.atnr;
fr->nbfp = mk_nbfp(&mtop->ffparams, fr->bBHAM);
if (ir->eDispCorr != edispcNO)
{
fr->dispersionCorrection = std::make_unique<DispersionCorrection>(
- *mtop, *ir, fr->bBHAM, fr->ntype,
- gmx::arrayRefFromArray(fr->nbfp, fr->ntype * fr->ntype * 2), *fr->ic, tabfn);
+ *mtop, *ir, fr->bBHAM, fr->ntype, fr->nbfp, *fr->ic, tabfn);
fr->dispersionCorrection->print(mdlog);
}
t_forcerec::t_forcerec() = default;
-t_forcerec::~t_forcerec() = default;
-
-/* Frees GPU memory and sets a tMPI node barrier.
- *
- * Note that this function needs to be called even if GPUs are not used
- * in this run because the PME ranks have no knowledge of whether GPUs
- * are used or not, but all ranks need to enter the barrier below.
- * \todo Remove physical node barrier from this function after making sure
- * that it's not needed anymore (with a shared GPU run).
- */
-void free_gpu_resources(t_forcerec* fr,
- const gmx::PhysicalNodeCommunicator& physicalNodeCommunicator,
- const gmx_gpu_info_t& gpu_info)
+t_forcerec::~t_forcerec()
{
- bool isPPrankUsingGPU = (fr != nullptr) && (fr->nbv != nullptr) && fr->nbv->useGpu();
-
- /* stop the GPU profiler (only CUDA) */
- if (gpu_info.n_dev > 0)
- {
- stopGpuProfiler();
- }
-
- if (isPPrankUsingGPU)
- {
- /* Free data in GPU memory and pinned memory before destroying the GPU context */
- fr->nbv.reset();
-
- delete fr->gpuBonded;
- fr->gpuBonded = nullptr;
- }
-
- /* With tMPI we need to wait for all ranks to finish deallocation before
- * destroying the CUDA context in free_gpu() as some tMPI ranks may be sharing
- * GPU and context.
- *
- * This is not a concern in OpenCL where we use one context per rank which
- * is freed in nbnxn_gpu_free().
- *
- * Note: it is safe to not call the barrier on the ranks which do not use GPU,
- * but it is easier and more futureproof to call it on the whole node.
- */
- if (GMX_THREAD_MPI)
- {
- physicalNodeCommunicator.barrier();
- }
-}
-
-void done_forcerec(t_forcerec* fr, int numMolBlocks)
-{
- if (fr == nullptr)
- {
- // PME-only ranks don't have a forcerec
- return;
- }
- done_cginfo_mb(fr->cginfo_mb, numMolBlocks);
- sfree(fr->nbfp);
- delete fr->ic;
- sfree(fr->shift_vec);
- sfree(fr->ewc_t);
- tear_down_bonded_threading(fr->bondedThreading);
- GMX_RELEASE_ASSERT(fr->gpuBonded == nullptr, "Should have been deleted earlier, when used");
- fr->bondedThreading = nullptr;
- delete fr;
+ /* Note: This code will disappear when types are converted to C++ */
+ sfree(shift_vec);
+ sfree(ewc_t);
+ tear_down_bonded_threading(bondedThreading);
}
class PhysicalNodeCommunicator;
} // namespace gmx
-//! Destroy a forcerec.
-void done_forcerec(t_forcerec* fr, int numMolBlocks);
-
/*! \brief Print the contents of the forcerec to a file
*
* \param[in] fplog The log file to print to
*/
void forcerec_set_excl_load(t_forcerec* fr, const gmx_localtop_t* top);
-void free_gpu_resources(t_forcerec* fr,
- const gmx::PhysicalNodeCommunicator& physicalNodeCommunicator,
- const gmx_gpu_info_t& gpu_info);
-
#endif
const int nwall = ir.nwall;
const int ngid = ir.opts.ngener;
const int ntype = fr.ntype;
- const real* nbfp = fr.nbfp;
+ const real* nbfp = fr.nbfp.data();
const int* egp_flags = fr.egp_flags;
for (int w = 0; w < nwall; w++)
/* This call is not included in init_domain_decomposition mainly
* because fr->cginfo_mb is set later.
*/
- dd_init_bondeds(fplog, cr->dd, &mtop, vsite.get(), inputrec,
+ dd_init_bondeds(fplog, cr->dd, mtop, vsite.get(), inputrec,
domdecOptions.checkBondedInteractions, fr->cginfo_mb);
}
}
// FIXME: this is only here to manually unpin mdAtoms->chargeA_ and state->x,
- // before we destroy the GPU context(s) in free_gpu_resources().
+ // before we destroy the GPU context(s) in free_gpu().
// Pinned buffers are associated with contexts in CUDA.
// As soon as we destroy GPU contexts after mdrunner() exits, these lines should go.
mdAtoms.reset(nullptr);
globalState.reset(nullptr);
mdModules_.reset(nullptr); // destruct force providers here as they might also use the GPU
+ /* Free pinned buffers in *fr */
+ delete fr;
+ fr = nullptr;
+
+ if (hwinfo->gpu_info.n_dev > 0)
+ {
+ /* stop the GPU profiler (only CUDA) */
+ stopGpuProfiler();
+ }
+
+ /* With tMPI we need to wait for all ranks to finish deallocation before
+ * destroying the CUDA context in free_gpu() as some tMPI ranks may be sharing
+ * GPU and context.
+ *
+ * This is not a concern in OpenCL where we use one context per rank which
+ * is freed in nbnxn_gpu_free().
+ *
+ * Note: it is safe to not call the barrier on the ranks which do not use GPU,
+ * but it is easier and more futureproof to call it on the whole node.
+ *
+ * Note that this function needs to be called even if GPUs are not used
+ * in this run because the PME ranks have no knowledge of whether GPUs
+ * are used or not, but all ranks need to enter the barrier below.
+ * \todo Remove this physical node barrier after making sure
+ * that it's not needed anymore (with a shared GPU run).
+ */
+ if (GMX_THREAD_MPI)
+ {
+ physicalNodeComm.barrier();
+ }
- /* Free GPU memory and set a physical node tMPI barrier (which should eventually go away) */
- free_gpu_resources(fr, physicalNodeComm, hwinfo->gpu_info);
free_gpu(nonbondedDeviceInfo);
free_gpu(pmeDeviceInfo);
- done_forcerec(fr, mtop.molblock.size());
sfree(fcd);
if (doMembed)
struct cginfo_mb_t
{
- int cg_start;
- int cg_end;
- int cg_mod;
- int* cginfo;
+ int cg_start = 0;
+ int cg_end = 0;
+ int cg_mod = 0;
+ std::vector<int> cginfo;
};
real sc_sigma6_min = 0;
/* Information about atom properties for the molecule blocks in the system */
- struct cginfo_mb_t* cginfo_mb = nullptr;
+ std::vector<cginfo_mb_t> cginfo_mb;
/* Information about atom properties for local and non-local atoms */
std::vector<int> cginfo;
std::vector<gmx::RVec> shiftForces;
/* Non bonded Parameter lists */
- int ntype = 0; /* Number of atom types */
- gmx_bool bBHAM = FALSE;
- real* nbfp = nullptr;
- real* ljpme_c6grid = nullptr; /* C6-values used on grid in LJPME */
+ int ntype = 0; /* Number of atom types */
+ gmx_bool bBHAM = FALSE;
+ std::vector<real> nbfp;
+ real* ljpme_c6grid = nullptr; /* C6-values used on grid in LJPME */
/* Energy group pair flags */
int* egp_flags = nullptr;
const Nbnxm::KernelType kernelType,
int enbnxninitcombrule,
int ntype,
- const real* nbfp,
+ ArrayRef<const real> nbfp,
int n_energygroups)
{
real c6, c12, tol;
const Nbnxm::KernelType kernelType,
int enbnxninitcombrule,
int ntype,
- const real* nbfp,
+ ArrayRef<const real> nbfp,
int n_energygroups,
int nout)
{
* to the atom data structure.
* enbnxninitcombrule sets what combination rule data gets stored in nbat.
*/
-void nbnxn_atomdata_init(const gmx::MDLogger& mdlog,
- nbnxn_atomdata_t* nbat,
- Nbnxm::KernelType kernelType,
- int enbnxninitcombrule,
- int ntype,
- const real* nbfp,
- int n_energygroups,
- int nout);
+void nbnxn_atomdata_init(const gmx::MDLogger& mdlog,
+ nbnxn_atomdata_t* nbat,
+ Nbnxm::KernelType kernelType,
+ int enbnxninitcombrule,
+ int ntype,
+ gmx::ArrayRef<const real> nbfp,
+ int n_energygroups,
+ int nout);
void nbnxn_atomdata_set(nbnxn_atomdata_t* nbat,
const Nbnxm::GridSet& gridSet,
std::move(atomData), kernelSetup, nullptr, nullptr);
nbnxn_atomdata_init(gmx::MDLogger(), nbv->nbat.get(), kernelSetup.kernelType, combinationRule,
- system.numAtomTypes, system.nonbondedParameters.data(), 1, numThreads);
+ system.numAtomTypes, system.nonbondedParameters, 1, numThreads);
t_nrnb nrnb;
}
forceRec.ntype = numAtomTypes;
- forceRec.nbfp = nonbondedParameters.data();
+ forceRec.nbfp = nonbondedParameters;
snew(forceRec.shift_vec, SHIFTS);
calc_shifts(box, forceRec.shift_vec);
}