From f10bc6258033cefe70f29f32d40d749ffb21fdb9 Mon Sep 17 00:00:00 2001 From: Berk Hess Date: Wed, 20 May 2020 16:56:32 +0200 Subject: [PATCH] Remove t_mdatoms from Nbnxm setAtomProperties() --- src/gromacs/mdlib/sim_util.cpp | 3 +- src/gromacs/mdrun/tpi.cpp | 4 ++- src/gromacs/nbnxm/atomdata.cpp | 31 +++++++++++---------- src/gromacs/nbnxm/atomdata.h | 10 +++---- src/gromacs/nbnxm/benchmark/bench_setup.cpp | 6 +--- src/gromacs/nbnxm/nbnxm.cpp | 6 ++-- src/gromacs/nbnxm/nbnxm.h | 4 ++- 7 files changed, 35 insertions(+), 29 deletions(-) diff --git a/src/gromacs/mdlib/sim_util.cpp b/src/gromacs/mdlib/sim_util.cpp index 31b0955b5b..ca42600eda 100644 --- a/src/gromacs/mdlib/sim_util.cpp +++ b/src/gromacs/mdlib/sim_util.cpp @@ -1213,7 +1213,8 @@ void do_force(FILE* fplog, wallcycle_sub_stop(wcycle, ewcsNBS_GRID_NONLOCAL); } - nbv->setAtomProperties(*mdatoms, fr->cginfo); + nbv->setAtomProperties(gmx::constArrayRefFromArray(mdatoms->typeA, mdatoms->nr), + gmx::constArrayRefFromArray(mdatoms->chargeA, mdatoms->nr), fr->cginfo); wallcycle_stop(wcycle, ewcNS); diff --git a/src/gromacs/mdrun/tpi.cpp b/src/gromacs/mdrun/tpi.cpp index ac2340ce17..280d1f22cb 100644 --- a/src/gromacs/mdrun/tpi.cpp +++ b/src/gromacs/mdrun/tpi.cpp @@ -677,7 +677,9 @@ void LegacySimulator::do_tpi() -1, fr->cginfo, x, 0, nullptr); /* TODO: Avoid updating all atoms at every bNS step */ - fr->nbv->setAtomProperties(*mdatoms, fr->cginfo); + fr->nbv->setAtomProperties(gmx::constArrayRefFromArray(mdatoms->typeA, mdatoms->nr), + gmx::constArrayRefFromArray(mdatoms->chargeA, mdatoms->nr), + fr->cginfo); fr->nbv->constructPairlist(InteractionLocality::Local, top.excls, step, nrnb); diff --git a/src/gromacs/nbnxm/atomdata.cpp b/src/gromacs/nbnxm/atomdata.cpp index a64fe1a381..35901f40f1 100644 --- a/src/gromacs/nbnxm/atomdata.cpp +++ b/src/gromacs/nbnxm/atomdata.cpp @@ -52,7 +52,6 @@ #include "gromacs/math/vec.h" #include "gromacs/mdlib/gmx_omp_nthreads.h" #include "gromacs/mdtypes/forcerec.h" // only for GET_CGINFO_* -#include "gromacs/mdtypes/mdatom.h" #include "gromacs/nbnxm/nbnxm.h" #include "gromacs/pbcutil/ishift.h" #include "gromacs/simd/simd.h" @@ -724,7 +723,7 @@ static void copy_lj_to_nbat_lj_comb(gmx::ArrayRef ljparam_type, cons /* Sets the atom type in nbnxn_atomdata_t */ static void nbnxn_atomdata_set_atomtypes(nbnxn_atomdata_t::Params* params, const Nbnxm::GridSet& gridSet, - const int* type) + ArrayRef atomTypes) { params->type.resize(gridSet.numGridAtomsTotal()); @@ -736,8 +735,9 @@ static void nbnxn_atomdata_set_atomtypes(nbnxn_atomdata_t::Params* params, const int numAtoms = grid.paddedNumAtomsInColumn(i); const int atomOffset = grid.firstAtomInColumn(i); - copy_int_to_nbat_int(gridSet.atomIndices().data() + atomOffset, grid.numAtomsInColumn(i), - numAtoms, type, params->numTypes - 1, params->type.data() + atomOffset); + copy_int_to_nbat_int(gridSet.atomIndices().data() + atomOffset, + grid.numAtomsInColumn(i), numAtoms, atomTypes.data(), + params->numTypes - 1, params->type.data() + atomOffset); } } } @@ -782,7 +782,9 @@ static void nbnxn_atomdata_set_ljcombparams(nbnxn_atomdata_t::Params* params, } /* Sets the charges in nbnxn_atomdata_t *nbat */ -static void nbnxn_atomdata_set_charges(nbnxn_atomdata_t* nbat, const Nbnxm::GridSet& gridSet, const real* charge) +static void nbnxn_atomdata_set_charges(nbnxn_atomdata_t* nbat, + const Nbnxm::GridSet& gridSet, + ArrayRef charges) { if (nbat->XFormat != nbatXYZQ) { @@ -804,7 +806,7 @@ static void nbnxn_atomdata_set_charges(nbnxn_atomdata_t* nbat, const Nbnxm::Grid int i; for (i = 0; i < numAtoms; i++) { - *q = charge[gridSet.atomIndices()[atomOffset + i]]; + *q = charges[gridSet.atomIndices()[atomOffset + i]]; q += STRIDE_XYZQ; } /* Complete the partially filled last cell with zeros */ @@ -820,7 +822,7 @@ static void nbnxn_atomdata_set_charges(nbnxn_atomdata_t* nbat, const Nbnxm::Grid int i; for (i = 0; i < numAtoms; i++) { - *q = charge[gridSet.atomIndices()[atomOffset + i]]; + *q = charges[gridSet.atomIndices()[atomOffset + i]]; q++; } /* Complete the partially filled last cell with zeros */ @@ -926,7 +928,7 @@ copy_egp_to_nbat_egps(const int* a, int na, int na_round, int na_c, int bit_shif /* Set the energy group indices for atoms in nbnxn_atomdata_t */ static void nbnxn_atomdata_set_energygroups(nbnxn_atomdata_t::Params* params, const Nbnxm::GridSet& gridSet, - const int* atinfo) + ArrayRef atomInfo) { if (params->nenergrp == 1) { @@ -944,7 +946,7 @@ static void nbnxn_atomdata_set_energygroups(nbnxn_atomdata_t::Params* params, const int atomOffset = grid.firstAtomInColumn(i); copy_egp_to_nbat_egps(gridSet.atomIndices().data() + atomOffset, grid.numAtomsInColumn(i), - numAtoms, c_nbnxnCpuIClusterSize, params->neg_2log, atinfo, + numAtoms, c_nbnxnCpuIClusterSize, params->neg_2log, atomInfo.data(), params->energrp.data() + grid.atomToCluster(atomOffset)); } } @@ -953,14 +955,15 @@ static void nbnxn_atomdata_set_energygroups(nbnxn_atomdata_t::Params* params, /* Sets all required atom parameter data in nbnxn_atomdata_t */ void nbnxn_atomdata_set(nbnxn_atomdata_t* nbat, const Nbnxm::GridSet& gridSet, - const t_mdatoms* mdatoms, - const int* atinfo) + ArrayRef atomTypes, + ArrayRef atomCharges, + ArrayRef atomInfo) { nbnxn_atomdata_t::Params& params = nbat->paramsDeprecated(); - nbnxn_atomdata_set_atomtypes(¶ms, gridSet, mdatoms->typeA); + nbnxn_atomdata_set_atomtypes(¶ms, gridSet, atomTypes); - nbnxn_atomdata_set_charges(nbat, gridSet, mdatoms->chargeA); + nbnxn_atomdata_set_charges(nbat, gridSet, atomCharges); if (gridSet.haveFep()) { @@ -970,7 +973,7 @@ void nbnxn_atomdata_set(nbnxn_atomdata_t* nbat, /* This must be done after masking types for FEP */ nbnxn_atomdata_set_ljcombparams(¶ms, nbat->XFormat, gridSet); - nbnxn_atomdata_set_energygroups(¶ms, gridSet, atinfo); + nbnxn_atomdata_set_energygroups(¶ms, gridSet, atomInfo); } /* Copies the shift vector array to nbnxn_atomdata_t */ diff --git a/src/gromacs/nbnxm/atomdata.h b/src/gromacs/nbnxm/atomdata.h index 22d115706e..4b3a398562 100644 --- a/src/gromacs/nbnxm/atomdata.h +++ b/src/gromacs/nbnxm/atomdata.h @@ -64,7 +64,6 @@ class MDLogger; struct NbnxmGpu; struct nbnxn_atomdata_t; struct nonbonded_verlet_t; -struct t_mdatoms; struct tMPI_Atomic; class GpuEventSynchronizer; @@ -340,10 +339,11 @@ void nbnxn_atomdata_init(const gmx::MDLogger& mdlog, int nout); //! Sets the atomdata after pair search -void nbnxn_atomdata_set(nbnxn_atomdata_t* nbat, - const Nbnxm::GridSet& gridSet, - const t_mdatoms* mdatoms, - const int* atinfo); +void nbnxn_atomdata_set(nbnxn_atomdata_t* nbat, + const Nbnxm::GridSet& gridSet, + gmx::ArrayRef atomTypes, + gmx::ArrayRef atomCharges, + gmx::ArrayRef atomInfo); //! Copy the shift vectors to nbat void nbnxn_atomdata_copy_shiftvec(gmx_bool dynamic_box, rvec* shift_vec, nbnxn_atomdata_t* nbat); diff --git a/src/gromacs/nbnxm/benchmark/bench_setup.cpp b/src/gromacs/nbnxm/benchmark/bench_setup.cpp index b933cee1b2..d5d4ab9d22 100644 --- a/src/gromacs/nbnxm/benchmark/bench_setup.cpp +++ b/src/gromacs/nbnxm/benchmark/bench_setup.cpp @@ -227,11 +227,7 @@ static std::unique_ptr setupNbnxmForBenchInstance(const Kern nbv->constructPairlist(gmx::InteractionLocality::Local, system.excls, 0, &nrnb); - t_mdatoms mdatoms; - // We only use (read) the atom type and charge from mdatoms - mdatoms.typeA = const_cast(system.atomTypes.data()); - mdatoms.chargeA = const_cast(system.charges.data()); - nbv->setAtomProperties(mdatoms, atomInfo); + nbv->setAtomProperties(system.atomTypes, system.charges, atomInfo); return nbv; } diff --git a/src/gromacs/nbnxm/nbnxm.cpp b/src/gromacs/nbnxm/nbnxm.cpp index ab79a8cf5b..20495969e0 100644 --- a/src/gromacs/nbnxm/nbnxm.cpp +++ b/src/gromacs/nbnxm/nbnxm.cpp @@ -119,9 +119,11 @@ void nonbonded_verlet_t::setLocalAtomOrder() pairSearch_->setLocalAtomOrder(); } -void nonbonded_verlet_t::setAtomProperties(const t_mdatoms& mdatoms, gmx::ArrayRef atomInfo) +void nonbonded_verlet_t::setAtomProperties(gmx::ArrayRef atomTypes, + gmx::ArrayRef atomCharges, + gmx::ArrayRef atomInfo) { - nbnxn_atomdata_set(nbat.get(), pairSearch_->gridSet(), &mdatoms, atomInfo.data()); + nbnxn_atomdata_set(nbat.get(), pairSearch_->gridSet(), atomTypes, atomCharges, atomInfo); } void nonbonded_verlet_t::convertCoordinates(const gmx::AtomLocality locality, diff --git a/src/gromacs/nbnxm/nbnxm.h b/src/gromacs/nbnxm/nbnxm.h index d5ced753f7..f1d7bed568 100644 --- a/src/gromacs/nbnxm/nbnxm.h +++ b/src/gromacs/nbnxm/nbnxm.h @@ -266,7 +266,9 @@ public: t_nrnb* nrnb); //! Updates all the atom properties in Nbnxm - void setAtomProperties(const t_mdatoms& mdatoms, gmx::ArrayRef atomInfo); + void setAtomProperties(gmx::ArrayRef atomTypes, + gmx::ArrayRef atomCharges, + gmx::ArrayRef atomInfo); /*!\brief Convert the coordinates to NBNXM format for the given locality. * -- 2.22.0