From: Artem Zhmurov Date: Fri, 22 May 2020 07:40:50 +0000 (+0000) Subject: Remove dependence of constraints on t_mdatoms X-Git-Url: http://biod.pnpi.spb.ru/gitweb/?a=commitdiff_plain;h=10c31b140c8791266c7c0940265686ca93864ec4;p=alexxy%2Fgromacs.git Remove dependence of constraints on t_mdatoms This changes the signatures of the constrants set/compute functions by removing the t_mdatoms from the set of arguments and replacing it with the parts that are actually used. Refs. #3535 --- diff --git a/src/gromacs/mdlib/constr.cpp b/src/gromacs/mdlib/constr.cpp index 16b5c6eba0..a158fd6c05 100644 --- a/src/gromacs/mdlib/constr.cpp +++ b/src/gromacs/mdlib/constr.cpp @@ -457,9 +457,10 @@ bool Constraints::Impl::apply(bool bLog, if (lincsd != nullptr) { - bOK = constrain_lincs(bLog || bEner, ir, step, lincsd, md, cr, ms, x, xprime, min_proj, box, - pbc_null, lambda, dvdlambda, invdt, v.unpaddedArrayRef(), - vir != nullptr, vir_r_m_dr, econq, nrnb, maxwarn, &warncount_lincs); + bOK = constrain_lincs(bLog || bEner, ir, step, lincsd, md.invmass, cr, ms, x, xprime, + min_proj, box, pbc_null, md.nMassPerturbed != 0, lambda, dvdlambda, + invdt, v.unpaddedArrayRef(), vir != nullptr, vir_r_m_dr, econq, nrnb, + maxwarn, &warncount_lincs); if (!bOK && maxwarn < INT_MAX) { if (log != nullptr) @@ -699,7 +700,7 @@ bool Constraints::Impl::apply(bool bLog, } return bOK; -} +} // namespace gmx ArrayRef Constraints::rmsdData() const { @@ -884,7 +885,7 @@ void Constraints::Impl::setConstraints(gmx_localtop_t* top, const t_mdatoms& md) */ if (ir.eConstrAlg == econtLINCS) { - set_lincs(*idef, md, EI_DYNAMICS(ir.eI), cr, lincsd); + set_lincs(*idef, md.nr, md.invmass, md.lambda, EI_DYNAMICS(ir.eI), cr, lincsd); } if (ir.eConstrAlg == econtSHAKE) { @@ -896,7 +897,7 @@ void Constraints::Impl::setConstraints(gmx_localtop_t* top, const t_mdatoms& md) } else { - make_shake_sblock_serial(shaked.get(), &top->idef, md); + make_shake_sblock_serial(shaked.get(), &top->idef, md.nr); } } } diff --git a/src/gromacs/mdlib/lincs.cpp b/src/gromacs/mdlib/lincs.cpp index 8b13772cdf..b9fff7fe05 100644 --- a/src/gromacs/mdlib/lincs.cpp +++ b/src/gromacs/mdlib/lincs.cpp @@ -68,7 +68,6 @@ #include "gromacs/mdtypes/commrec.h" #include "gromacs/mdtypes/inputrec.h" #include "gromacs/mdtypes/md_enums.h" -#include "gromacs/mdtypes/mdatom.h" #include "gromacs/pbcutil/pbc.h" #include "gromacs/pbcutil/pbc_simd.h" #include "gromacs/simd/simd.h" @@ -575,7 +574,7 @@ static void do_lincsp(ArrayRefWithPadding xPadded, t_pbc* pbc, Lincs* lincsd, int th, - real* invmass, + const real* invmass, ConstraintVariable econq, bool bCalcDHDL, bool bCalcVir, @@ -1302,7 +1301,7 @@ static void set_lincs_matrix_task(Lincs* li, Task* li_task, const real* invmass, } /*! \brief Sets the elements in the LINCS matrix. */ -static void set_lincs_matrix(Lincs* li, real* invmass, real lambda) +static void set_lincs_matrix(Lincs* li, const real* invmass, real lambda) { const real invsqrt2 = 0.7071067811865475244; @@ -1850,7 +1849,13 @@ static void set_matrix_indices(Lincs* li, const Task& li_task, const ListOfLists } } -void set_lincs(const InteractionDefinitions& idef, const t_mdatoms& md, bool bDynamics, const t_commrec* cr, Lincs* li) +void set_lincs(const InteractionDefinitions& idef, + const int numAtoms, + const real* invmass, + const real lambda, + bool bDynamics, + const t_commrec* cr, + Lincs* li) { li->nc_real = 0; li->nc = 0; @@ -1899,7 +1904,7 @@ void set_lincs(const InteractionDefinitions& idef, const t_mdatoms& md, bool bDy } else { - natoms = md.homenr; + natoms = numAtoms; } const ListOfLists at2con = @@ -2123,10 +2128,10 @@ void set_lincs(const InteractionDefinitions& idef, const t_mdatoms& md, bool bDy if (li->ntask > 1) { - lincs_thread_setup(li, md.nr); + lincs_thread_setup(li, numAtoms); } - set_lincs_matrix(li, md.invmass, md.lambda); + set_lincs_matrix(li, invmass, lambda); } //! Issues a warning when LINCS constraints cannot be satisfied. @@ -2252,7 +2257,7 @@ bool constrain_lincs(bool computeRmsd, const t_inputrec& ir, int64_t step, Lincs* lincsd, - const t_mdatoms& md, + const real* invmass, const t_commrec* cr, const gmx_multisim_t* ms, ArrayRefWithPadding xPadded, @@ -2260,6 +2265,7 @@ bool constrain_lincs(bool computeRmsd, ArrayRef min_proj, const matrix box, t_pbc* pbc, + const bool hasMassPerturbed, real lambda, real* dvdlambda, real invdt, @@ -2299,9 +2305,9 @@ bool constrain_lincs(bool computeRmsd, */ if (ir.efep != efepNO) { - if (md.nMassPerturbed && lincsd->matlam != md.lambda) + if (hasMassPerturbed && lincsd->matlam != lambda) { - set_lincs_matrix(lincsd, md.invmass, md.lambda); + set_lincs_matrix(lincsd, invmass, lambda); } for (int i = 0; i < lincsd->nc; i++) @@ -2365,7 +2371,7 @@ bool constrain_lincs(bool computeRmsd, clear_mat(lincsd->task[th].vir_r_m_dr); - do_lincs(xPadded, xprimePadded, box, pbc, lincsd, th, md.invmass, cr, bCalcDHDL, + do_lincs(xPadded, xprimePadded, box, pbc, lincsd, th, invmass, cr, bCalcDHDL, ir.LincsWarnAngle, &bWarn, invdt, v, bCalcVir, th == 0 ? vir_r_m_dr : lincsd->task[th].vir_r_m_dr); } @@ -2444,7 +2450,7 @@ bool constrain_lincs(bool computeRmsd, { int th = gmx_omp_get_thread_num(); - do_lincsp(xPadded, xprimePadded, min_proj, pbc, lincsd, th, md.invmass, econq, + do_lincsp(xPadded, xprimePadded, min_proj, pbc, lincsd, th, invmass, econq, bCalcDHDL, bCalcVir, th == 0 ? vir_r_m_dr : lincsd->task[th].vir_r_m_dr); } GMX_CATCH_ALL_AND_EXIT_WITH_FATAL_ERROR diff --git a/src/gromacs/mdlib/lincs.h b/src/gromacs/mdlib/lincs.h index 2ab8ee8062..0b00370c01 100644 --- a/src/gromacs/mdlib/lincs.h +++ b/src/gromacs/mdlib/lincs.h @@ -56,7 +56,6 @@ struct gmx_multisim_t; class InteractionDefinitions; struct t_commrec; struct t_inputrec; -struct t_mdatoms; struct t_nrnb; struct t_pbc; @@ -91,7 +90,9 @@ void done_lincs(Lincs* li); /*! \brief Initialize lincs stuff */ void set_lincs(const InteractionDefinitions& idef, - const t_mdatoms& md, + int numAtoms, + const real* invmass, + real lambda, bool bDynamics, const t_commrec* cr, Lincs* li); @@ -103,7 +104,7 @@ bool constrain_lincs(bool computeRmsd, const t_inputrec& ir, int64_t step, Lincs* lincsd, - const t_mdatoms& md, + const real* invmass, const t_commrec* cr, const gmx_multisim_t* ms, ArrayRefWithPadding x, @@ -111,6 +112,7 @@ bool constrain_lincs(bool computeRmsd, ArrayRef min_proj, const matrix box, t_pbc* pbc, + bool hasMassPerturbed, real lambda, real* dvdlambda, real invdt, diff --git a/src/gromacs/mdlib/lincs_gpu.cu b/src/gromacs/mdlib/lincs_gpu.cu index 1d6bf1ca71..7b3c067d46 100644 --- a/src/gromacs/mdlib/lincs_gpu.cu +++ b/src/gromacs/mdlib/lincs_gpu.cu @@ -725,9 +725,8 @@ bool LincsGpu::isNumCoupledConstraintsSupported(const gmx_mtop_t& mtop) return true; } -void LincsGpu::set(const InteractionDefinitions& idef, const t_mdatoms& md) +void LincsGpu::set(const InteractionDefinitions& idef, const int numAtoms, const real* invmass) { - int numAtoms = md.nr; // List of constrained atoms (CPU memory) std::vector constraintsHost; // Equilibrium distances for the constraints (CPU) @@ -861,10 +860,10 @@ void LincsGpu::set(const InteractionDefinitions& idef, const t_mdatoms& md) int center = c1a1; - float sqrtmu1 = 1.0 / sqrt(md.invmass[c1a1] + md.invmass[c1a2]); - float sqrtmu2 = 1.0 / sqrt(md.invmass[c2a1] + md.invmass[c2a2]); + float sqrtmu1 = 1.0 / sqrt(invmass[c1a1] + invmass[c1a2]); + float sqrtmu2 = 1.0 / sqrt(invmass[c2a1] + invmass[c2a2]); - massFactorsHost.at(index) = -sign * md.invmass[center] * sqrtmu1 * sqrtmu2; + massFactorsHost.at(index) = -sign * invmass[center] * sqrtmu1 * sqrtmu2; coupledConstraintsCountsHost.at(splitMap.at(c1))++; } @@ -888,10 +887,10 @@ void LincsGpu::set(const InteractionDefinitions& idef, const t_mdatoms& md) int center = c1a2; - float sqrtmu1 = 1.0 / sqrt(md.invmass[c1a1] + md.invmass[c1a2]); - float sqrtmu2 = 1.0 / sqrt(md.invmass[c2a1] + md.invmass[c2a2]); + float sqrtmu1 = 1.0 / sqrt(invmass[c1a1] + invmass[c1a2]); + float sqrtmu2 = 1.0 / sqrt(invmass[c2a1] + invmass[c2a2]); - massFactorsHost.at(index) = sign * md.invmass[center] * sqrtmu1 * sqrtmu2; + massFactorsHost.at(index) = sign * invmass[center] * sqrtmu1 * sqrtmu2; coupledConstraintsCountsHost.at(splitMap.at(c1))++; } @@ -958,8 +957,8 @@ void LincsGpu::set(const InteractionDefinitions& idef, const t_mdatoms& md) maxCoupledConstraints * kernelParams_.numConstraintsThreads, deviceStream_, GpuApiCallBehavior::Sync, nullptr); - GMX_RELEASE_ASSERT(md.invmass != nullptr, "Masses of atoms should be specified.\n"); - copyToDeviceBuffer(&kernelParams_.d_inverseMasses, md.invmass, 0, numAtoms, deviceStream_, + GMX_RELEASE_ASSERT(invmass != nullptr, "Masses of atoms should be specified.\n"); + copyToDeviceBuffer(&kernelParams_.d_inverseMasses, invmass, 0, numAtoms, deviceStream_, GpuApiCallBehavior::Sync, nullptr); } diff --git a/src/gromacs/mdlib/lincs_gpu.cuh b/src/gromacs/mdlib/lincs_gpu.cuh index ef03516431..15a962bac1 100644 --- a/src/gromacs/mdlib/lincs_gpu.cuh +++ b/src/gromacs/mdlib/lincs_gpu.cuh @@ -47,7 +47,6 @@ #include "gromacs/gpu_utils/device_context.h" #include "gromacs/gpu_utils/gputraits.cuh" #include "gromacs/mdlib/constr.h" -#include "gromacs/mdtypes/mdatom.h" #include "gromacs/pbcutil/pbc_aiuc.h" #include "gromacs/utility/classhelpers.h" @@ -155,13 +154,12 @@ public: * Information about constraints is taken from: * idef.il[F_CONSTR].iatoms --- type (T) of constraint and two atom indexes (i1, i2) * idef.iparams[T].constr.dA --- target length for constraint of type T - * From t_mdatom, the code takes: - * md.invmass --- array of inverse square root of masses for each atom in the system. * - * \param[in] idef Local topology data to get information on constraints from. - * \param[in] md Atoms data to get atom masses from. + * \param[in] idef Local topology data to get information on constraints from. + * \param[in] numAtoms Number of atoms. + * \param[in] invmass Inverse masses of atoms. */ - void set(const InteractionDefinitions& idef, const t_mdatoms& md); + void set(const InteractionDefinitions& idef, int numAtoms, const real* invmass); /*! \brief * Returns whether the maximum number of coupled constraints is supported diff --git a/src/gromacs/mdlib/shake.cpp b/src/gromacs/mdlib/shake.cpp index 62e351dfe9..fd72911b77 100644 --- a/src/gromacs/mdlib/shake.cpp +++ b/src/gromacs/mdlib/shake.cpp @@ -58,7 +58,6 @@ #include "gromacs/mdlib/splitter.h" #include "gromacs/mdtypes/inputrec.h" #include "gromacs/mdtypes/md_enums.h" -#include "gromacs/mdtypes/mdatom.h" #include "gromacs/pbcutil/pbc.h" #include "gromacs/topology/invblock.h" #include "gromacs/utility/fatalerror.h" @@ -122,7 +121,7 @@ static void resizeLagrangianData(shakedata* shaked, int ncons) shaked->scaled_lagrange_multiplier.resize(ncons); } -void make_shake_sblock_serial(shakedata* shaked, InteractionDefinitions* idef, const t_mdatoms& md) +void make_shake_sblock_serial(shakedata* shaked, InteractionDefinitions* idef, const int numAtoms) { int i, m, ncons; int bstart, bnr; @@ -137,7 +136,7 @@ void make_shake_sblock_serial(shakedata* shaked, InteractionDefinitions* idef, c init_blocka(&sblocks); sfree(sblocks.index); // To solve memory leak - gen_sblocks(nullptr, md.homenr, *idef, &sblocks, FALSE); + gen_sblocks(nullptr, numAtoms, *idef, &sblocks, FALSE); /* bstart=(idef->nodeid > 0) ? blocks->multinr[idef->nodeid-1] : 0; @@ -150,7 +149,7 @@ void make_shake_sblock_serial(shakedata* shaked, InteractionDefinitions* idef, c } /* Calculate block number for each atom */ - inv_sblock = make_invblocka(&sblocks, md.nr); + inv_sblock = make_invblocka(&sblocks, numAtoms); done_blocka(&sblocks); diff --git a/src/gromacs/mdlib/shake.h b/src/gromacs/mdlib/shake.h index 34db77879c..8559c57553 100644 --- a/src/gromacs/mdlib/shake.h +++ b/src/gromacs/mdlib/shake.h @@ -52,7 +52,6 @@ struct InteractionList; class InteractionDefinitions; struct t_inputrec; -struct t_mdatoms; struct t_nrnb; struct t_pbc; @@ -96,7 +95,7 @@ struct shakedata }; //! Make SHAKE blocks when not using DD. -void make_shake_sblock_serial(shakedata* shaked, InteractionDefinitions* idef, const t_mdatoms& md); +void make_shake_sblock_serial(shakedata* shaked, InteractionDefinitions* idef, int numAtoms); //! Make SHAKE blocks when using DD. void make_shake_sblock_dd(shakedata* shaked, const InteractionList& ilcon); diff --git a/src/gromacs/mdlib/tests/constrtestdata.cpp b/src/gromacs/mdlib/tests/constrtestdata.cpp index b410504f8d..3de0c11b89 100644 --- a/src/gromacs/mdlib/tests/constrtestdata.cpp +++ b/src/gromacs/mdlib/tests/constrtestdata.cpp @@ -101,13 +101,6 @@ ConstraintsTestData::ConstraintsTestData(const std::string& title, ir_.delta_t = timestep; ir_.eI = 0; - // MD atoms data - md_.nMassPerturbed = 0; - md_.lambda = 0.0; - md_.invmass = invmass_.data(); - md_.nr = numAtoms; - md_.homenr = numAtoms; - // Virial evaluation computeVirial_ = computeVirial; if (computeVirial) @@ -122,8 +115,9 @@ ConstraintsTestData::ConstraintsTestData(const std::string& title, } } - // Free energy evaluation + hasMassPerturbed_ = false; + lambda_ = 0.0; compute_dHdLambda_ = compute_dHdLambda; dHdLambda_ = 0; if (compute_dHdLambda_) diff --git a/src/gromacs/mdlib/tests/constrtestdata.h b/src/gromacs/mdlib/tests/constrtestdata.h index 7374fa1924..65d6e03f56 100644 --- a/src/gromacs/mdlib/tests/constrtestdata.h +++ b/src/gromacs/mdlib/tests/constrtestdata.h @@ -58,7 +58,6 @@ #include "gromacs/mdlib/lincs.h" #include "gromacs/mdlib/shake.h" #include "gromacs/mdtypes/inputrec.h" -#include "gromacs/mdtypes/mdatom.h" #include "gromacs/pbcutil/pbc.h" #include "gromacs/topology/idef.h" #include "gromacs/topology/ifunc.h" @@ -92,8 +91,6 @@ public: t_inputrec ir_; //! Local topology std::unique_ptr idef_; - //! MD atoms - t_mdatoms md_; //! Computational time array (normally used to benchmark performance) t_nrnb nrnb_; @@ -109,6 +106,10 @@ public: tensor virialScaledRef_; //! If the free energy is computed bool compute_dHdLambda_; + //! If there are atoms with perturbed mass + bool hasMassPerturbed_ = false; + //! Lambda value + real lambda_ = 0.0; //! For free energy computation real dHdLambda_; //! For free energy computation (reference value) diff --git a/src/gromacs/mdlib/tests/constrtestrunners.cpp b/src/gromacs/mdlib/tests/constrtestrunners.cpp index 9ef51bdee7..1700e6e90b 100644 --- a/src/gromacs/mdlib/tests/constrtestrunners.cpp +++ b/src/gromacs/mdlib/tests/constrtestrunners.cpp @@ -66,7 +66,6 @@ #include "gromacs/mdrunutility/multisim.h" #include "gromacs/mdtypes/commrec.h" #include "gromacs/mdtypes/inputrec.h" -#include "gromacs/mdtypes/mdatom.h" #include "gromacs/pbcutil/pbc.h" #include "gromacs/topology/idef.h" #include "gromacs/topology/ifunc.h" @@ -90,10 +89,10 @@ namespace test void applyShake(ConstraintsTestData* testData, t_pbc gmx_unused pbc) { shakedata shaked; - make_shake_sblock_serial(&shaked, testData->idef_.get(), testData->md_); + make_shake_sblock_serial(&shaked, testData->idef_.get(), testData->numAtoms_); bool success = constrain_shake( nullptr, &shaked, testData->invmass_.data(), *testData->idef_, testData->ir_, testData->x_, - testData->xPrime_, testData->xPrime2_, nullptr, &testData->nrnb_, testData->md_.lambda, + testData->xPrime_, testData->xPrime2_, nullptr, &testData->nrnb_, testData->lambda_, &testData->dHdLambda_, testData->invdt_, testData->v_, testData->computeVirial_, testData->virialScaled_, false, gmx::ConstraintVariable::Positions); EXPECT_TRUE(success) << "Test failed with a false return value in SHAKE."; @@ -133,14 +132,15 @@ void applyLincs(ConstraintsTestData* testData, t_pbc pbc) // Initialize LINCS lincsd = init_lincs(nullptr, testData->mtop_, testData->nflexcon_, at2con_mt, false, testData->ir_.nLincsIter, testData->ir_.nProjOrder); - set_lincs(*testData->idef_, testData->md_, EI_DYNAMICS(testData->ir_.eI), &cr, lincsd); + set_lincs(*testData->idef_, testData->numAtoms_, testData->invmass_.data(), testData->lambda_, + EI_DYNAMICS(testData->ir_.eI), &cr, lincsd); // Evaluate constraints bool success = constrain_lincs( - false, testData->ir_, 0, lincsd, testData->md_, &cr, &ms, + false, testData->ir_, 0, lincsd, testData->invmass_.data(), &cr, &ms, testData->x_.arrayRefWithPadding(), testData->xPrime_.arrayRefWithPadding(), testData->xPrime2_.arrayRefWithPadding().unpaddedArrayRef(), pbc.box, &pbc, - testData->md_.lambda, &testData->dHdLambda_, testData->invdt_, + testData->hasMassPerturbed_, testData->lambda_, &testData->dHdLambda_, testData->invdt_, testData->v_.arrayRefWithPadding().unpaddedArrayRef(), testData->computeVirial_, testData->virialScaled_, gmx::ConstraintVariable::Positions, &testData->nrnb_, maxwarn, &warncount_lincs); diff --git a/src/gromacs/mdlib/tests/constrtestrunners.cu b/src/gromacs/mdlib/tests/constrtestrunners.cu index 8c2385daab..161c63b5d5 100644 --- a/src/gromacs/mdlib/tests/constrtestrunners.cu +++ b/src/gromacs/mdlib/tests/constrtestrunners.cu @@ -81,7 +81,7 @@ void applyLincsGpu(ConstraintsTestData* testData, t_pbc pbc) int numAtoms = testData->numAtoms_; float3 *d_x, *d_xp, *d_v; - lincsGpu->set(*testData->idef_, testData->md_); + lincsGpu->set(*testData->idef_, testData->numAtoms_, testData->invmass_.data()); PbcAiuc pbcAiuc; setPbcAiuc(pbc.ndim_ePBC, pbc.box, &pbcAiuc); diff --git a/src/gromacs/mdlib/update_constrain_gpu_impl.cu b/src/gromacs/mdlib/update_constrain_gpu_impl.cu index 76899bbd82..b40978bf90 100644 --- a/src/gromacs/mdlib/update_constrain_gpu_impl.cu +++ b/src/gromacs/mdlib/update_constrain_gpu_impl.cu @@ -214,7 +214,7 @@ void UpdateConstrainGpu::Impl::set(DeviceBuffer d_x, // Integrator should also update something, but it does not even have a method yet integrator_->set(md, numTempScaleValues, md.cTC); - lincsGpu_->set(idef, md); + lincsGpu_->set(idef, md.nr, md.invmass); settleGpu_->set(idef, md); coordinateScalingKernelLaunchConfig_.gridSize[0] =