From: ejjordan Date: Tue, 6 Apr 2021 09:21:22 +0000 (+0200) Subject: Minor forcerec cleanup X-Git-Url: http://biod.pnpi.spb.ru/gitweb/?a=commitdiff_plain;h=93a1ce64588b4db3b7b0520c38871dd18e29cb68;p=alexxy%2Fgromacs.git Minor forcerec cleanup * gmx_bool -> bool * rename a variable * rvec -> RVec * more array and vector * now has default destructor --- diff --git a/src/gromacs/mdlib/force.cpp b/src/gromacs/mdlib/force.cpp index 8c964b9945..081d603dc1 100644 --- a/src/gromacs/mdlib/force.cpp +++ b/src/gromacs/mdlib/force.cpp @@ -83,7 +83,7 @@ static void clearEwaldThreadOutput(ewald_corr_thread_t* ewc_t) clear_mat(ewc_t->vir_lj); } -static void reduceEwaldThreadOuput(int nthreads, ewald_corr_thread_t* ewc_t) +static void reduceEwaldThreadOuput(int nthreads, gmx::ArrayRef ewc_t) { ewald_corr_thread_t& dest = ewc_t[0]; diff --git a/src/gromacs/mdlib/forcerec.cpp b/src/gromacs/mdlib/forcerec.cpp index c2004aa45e..d2a3bd8692 100644 --- a/src/gromacs/mdlib/forcerec.cpp +++ b/src/gromacs/mdlib/forcerec.cpp @@ -213,7 +213,7 @@ static std::vector init_cginfo_mb(const gmx_mtop_t& mtop, const t_f type_VDW[ai] = FALSE; for (int j = 0; j < fr->ntype; j++) { - type_VDW[ai] = type_VDW[ai] || fr->bBHAM || C6(fr->nbfp, fr->ntype, ai, j) != 0 + type_VDW[ai] = type_VDW[ai] || fr->haveBuckingham || C6(fr->nbfp, fr->ntype, ai, j) != 0 || C12(fr->nbfp, fr->ntype, ai, j) != 0; } } @@ -1040,7 +1040,7 @@ void init_forcerec(FILE* fplog, } } - forcerec->bBHAM = (mtop.ffparams.functype[0] == F_BHAM); + forcerec->haveBuckingham = (mtop.ffparams.functype[0] == F_BHAM); /* Neighbour searching stuff */ forcerec->pbcType = inputrec.pbcType; @@ -1151,7 +1151,7 @@ void init_forcerec(FILE* fplog, switch (interactionConst->vdwtype) { case VanDerWaalsType::Cut: - if (forcerec->bBHAM) + if (forcerec->haveBuckingham) { forcerec->nbkernel_vdw_interaction = NbkernelVdwType::Buckingham; } @@ -1189,9 +1189,6 @@ void init_forcerec(FILE* fplog, enumValueToString(inputrec.coulombtype)); } - forcerec->bvdwtab = FALSE; - forcerec->bcoultab = FALSE; - /* 1-4 interaction electrostatics */ forcerec->fudgeQQ = mtop.ffparams.fudgeQQ; @@ -1226,7 +1223,7 @@ void init_forcerec(FILE* fplog, if (forcerec->nbfp.empty()) { forcerec->ntype = mtop.ffparams.atnr; - forcerec->nbfp = makeNonBondedParameterLists(mtop.ffparams, forcerec->bBHAM); + forcerec->nbfp = makeNonBondedParameterLists(mtop.ffparams, forcerec->haveBuckingham); if (EVDW_PME(interactionConst->vdwtype)) { forcerec->ljpme_c6grid = makeLJPmeC6GridCorrectionParameters(mtop.ffparams, *forcerec); @@ -1238,7 +1235,7 @@ void init_forcerec(FILE* fplog, /* Van der Waals stuff */ if ((interactionConst->vdwtype != VanDerWaalsType::Cut) - && (interactionConst->vdwtype != VanDerWaalsType::User) && !forcerec->bBHAM) + && (interactionConst->vdwtype != VanDerWaalsType::User) && !forcerec->haveBuckingham) { if (interactionConst->rvdw_switch >= interactionConst->rvdw) { @@ -1257,19 +1254,19 @@ void init_forcerec(FILE* fplog, } } - if (forcerec->bBHAM && EVDW_PME(interactionConst->vdwtype)) + if (forcerec->haveBuckingham && EVDW_PME(interactionConst->vdwtype)) { gmx_fatal(FARGS, "LJ PME not supported with Buckingham"); } - if (forcerec->bBHAM + if (forcerec->haveBuckingham && (interactionConst->vdwtype == VanDerWaalsType::Shift || interactionConst->vdwtype == VanDerWaalsType::Switch)) { gmx_fatal(FARGS, "Switch/shift interaction not supported with Buckingham"); } - if (forcerec->bBHAM) + if (forcerec->haveBuckingham) { gmx_fatal(FARGS, "The Verlet cutoff-scheme does not (yet) support Buckingham"); } @@ -1394,12 +1391,12 @@ void init_forcerec(FILE* fplog, forcerec->print_force = print_force; forcerec->nthread_ewc = gmx_omp_nthreads_get(emntBonded); - snew(forcerec->ewc_t, forcerec->nthread_ewc); + forcerec->ewc_t.resize(forcerec->nthread_ewc); if (inputrec.eDispCorr != DispersionCorrectionType::No) { forcerec->dispersionCorrection = std::make_unique( - mtop, inputrec, forcerec->bBHAM, forcerec->ntype, forcerec->nbfp, *forcerec->ic, tabfn); + mtop, inputrec, forcerec->haveBuckingham, forcerec->ntype, forcerec->nbfp, *forcerec->ic, tabfn); forcerec->dispersionCorrection->print(mdlog); } @@ -1415,8 +1412,4 @@ void init_forcerec(FILE* fplog, t_forcerec::t_forcerec() = default; -t_forcerec::~t_forcerec() -{ - /* Note: This code will disappear when types are converted to C++ */ - sfree(ewc_t); -} +t_forcerec::~t_forcerec() = default; diff --git a/src/gromacs/mdrun/tpi.cpp b/src/gromacs/mdrun/tpi.cpp index 819532068e..056cbb29d4 100644 --- a/src/gromacs/mdrun/tpi.cpp +++ b/src/gromacs/mdrun/tpi.cpp @@ -862,7 +862,7 @@ void LegacySimulator::do_tpi() /* Determine the weighted energy contributions of each energy group */ e = 0; sum_UgembU[e++] += epot * embU; - if (fr->bBHAM) + if (fr->haveBuckingham) { for (i = 0; i < ngid; i++) { diff --git a/src/gromacs/mdtypes/forcerec.h b/src/gromacs/mdtypes/forcerec.h index ee8425551f..de1c67659f 100644 --- a/src/gromacs/mdtypes/forcerec.h +++ b/src/gromacs/mdtypes/forcerec.h @@ -107,16 +107,6 @@ class WholeMoleculeTransform; */ #define GMX_CUTOFF_INF 1E+18 -/* enums for the neighborlist type */ -enum -{ - enbvdwNONE, - enbvdwLJ, - enbvdwBHAM, - enbvdwTAB, - enbvdwNR -}; - struct cginfo_mb_t { int cg_start = 0; @@ -190,12 +180,12 @@ struct t_forcerec //! Tells whether atoms inside a molecule can be in different periodic images, // i.e. whether we need to take into account PBC when computing distances inside molecules. // This determines whether PBC must be considered for e.g. bonded interactions. - gmx_bool bMolPBC = FALSE; + bool bMolPBC = false; RefCoordScaling rc_scaling = RefCoordScaling::No; - rvec posres_com = { 0 }; - rvec posres_comB = { 0 }; + gmx::RVec posres_com = { 0, 0, 0 }; + gmx::RVec posres_comB = { 0, 0, 0 }; - gmx_bool use_simd_kernels = FALSE; + bool use_simd_kernels = false; /* Interaction for calculated in kernels. In many cases this is similar to * the electrostatics settings in the inputrecord, but the difference is that @@ -217,9 +207,9 @@ struct t_forcerec real rlist = 0; /* Charge sum for topology A/B ([0]/[1]) for Ewald corrections */ - double qsum[2] = { 0 }; - double q2sum[2] = { 0 }; - double c6sum[2] = { 0 }; + std::array qsum = { 0 }; + std::array q2sum = { 0 }; + std::array c6sum = { 0 }; /* Dispersion correction stuff */ std::unique_ptr dispersionCorrection; @@ -227,10 +217,6 @@ struct t_forcerec /* Fudge factors */ real fudgeQQ = 0; - /* Table stuff */ - gmx_bool bcoultab = FALSE; - gmx_bool bvdwtab = FALSE; - std::unique_ptr pairsTable; /* for 1-4 interactions, [pairs] and [pairs_nb] */ /* Free energy */ @@ -261,15 +247,15 @@ struct t_forcerec std::vector forceHelperBuffers; /* Data for PPPM/PME/Ewald */ - struct gmx_pme_t* pmedata = nullptr; - LongRangeVdW ljpme_combination_rule = LongRangeVdW::Geom; + gmx_pme_t* pmedata = nullptr; + LongRangeVdW ljpme_combination_rule = LongRangeVdW::Geom; /* PME/Ewald stuff */ std::unique_ptr ewald_table; /* Non bonded Parameter lists */ - int ntype = 0; /* Number of atom types */ - gmx_bool bBHAM = FALSE; + int ntype = 0; /* Number of atom types */ + bool haveBuckingham = false; std::vector nbfp; std::vector ljpme_c6grid; /* C6-values used on grid in LJPME */ @@ -312,8 +298,8 @@ struct t_forcerec gmx::GpuBonded* gpuBonded = nullptr; /* Ewald correction thread local virial and energy data */ - int nthread_ewc = 0; - struct ewald_corr_thread_t* ewc_t = nullptr; + int nthread_ewc = 0; + std::vector ewc_t; gmx::ForceProviders* forceProviders = nullptr; diff --git a/src/gromacs/nbnxm/kerneldispatch.cpp b/src/gromacs/nbnxm/kerneldispatch.cpp index 24608cf456..9b026403c9 100644 --- a/src/gromacs/nbnxm/kerneldispatch.cpp +++ b/src/gromacs/nbnxm/kerneldispatch.cpp @@ -466,8 +466,9 @@ void nonbonded_verlet_t::dispatchNonbondedKernel(gmx::InteractionLocality iLoc stepWork, clearF, enerd->grpp.energyGroupPairTerms[NonBondedEnergyTerms::CoulombSR].data(), - fr.bBHAM ? enerd->grpp.energyGroupPairTerms[NonBondedEnergyTerms::BuckinghamSR].data() - : enerd->grpp.energyGroupPairTerms[NonBondedEnergyTerms::LJSR].data(), + fr.haveBuckingham + ? enerd->grpp.energyGroupPairTerms[NonBondedEnergyTerms::BuckinghamSR].data() + : enerd->grpp.energyGroupPairTerms[NonBondedEnergyTerms::LJSR].data(), wcycle_); break; @@ -486,8 +487,9 @@ void nonbonded_verlet_t::dispatchNonbondedKernel(gmx::InteractionLocality iLoc nbat->out[0].f, nbat->out[0].fshift.data(), enerd->grpp.energyGroupPairTerms[NonBondedEnergyTerms::CoulombSR].data(), - fr.bBHAM ? enerd->grpp.energyGroupPairTerms[NonBondedEnergyTerms::BuckinghamSR].data() - : enerd->grpp.energyGroupPairTerms[NonBondedEnergyTerms::LJSR].data()); + fr.haveBuckingham + ? enerd->grpp.energyGroupPairTerms[NonBondedEnergyTerms::BuckinghamSR].data() + : enerd->grpp.energyGroupPairTerms[NonBondedEnergyTerms::LJSR].data()); break; default: GMX_RELEASE_ASSERT(false, "Invalid nonbonded kernel type passed!");