From: Artem Zhmurov Date: Thu, 4 Feb 2021 09:56:25 +0000 (+0000) Subject: Make LJ combination rule enum into an enum class X-Git-Url: http://biod.pnpi.spb.ru/gitweb/?a=commitdiff_plain;h=f3dbd7334bdec032738d7bc7342fb62e06c30a01;p=alexxy%2Fgromacs.git Make LJ combination rule enum into an enum class To introduce type-safety, the combination rule for LJ is madde into enum class. This also fixes: 1. Out-of-order string assignment to enum values 2. Unsafe utilization of the enum for LJ PME combination rule (the values of this enum coincided so no bug by luck) --- diff --git a/src/gromacs/nbnxm/atomdata.cpp b/src/gromacs/nbnxm/atomdata.cpp index 35e675b847..315b05f3f2 100644 --- a/src/gromacs/nbnxm/atomdata.cpp +++ b/src/gromacs/nbnxm/atomdata.cpp @@ -72,7 +72,14 @@ using namespace gmx; // TODO: Remove when this file is moved into gmx namespace -const char* const c_ljcrNames[ljcrNR + 1] = { "none", "geometric", "Lorentz-Berthelot", nullptr }; + +const char* enumValueToString(LJCombinationRule enumValue) +{ + static constexpr gmx::EnumerationArray s_ljCombinationRuleNames = { + "Geometric", "Lorentz-Berthelot", "None" + }; + return s_ljCombinationRuleNames[enumValue]; +} void nbnxn_atomdata_t::resizeCoordinateBuffer(int numAtoms) { @@ -292,11 +299,9 @@ static void set_lj_parameter_data(nbnxn_atomdata_t::Params* params, gmx_bool bSI * not per pair of atom types. */ params->nbfp_comb.resize(nt * 2); - switch (params->comb_rule) + switch (params->ljCombinationRule) { - case ljcrGEOM: - params->comb_rule = ljcrGEOM; - + case LJCombinationRule::Geometric: for (int i = 0; i < nt; i++) { /* Store the sqrt of the diagonal from the nbfp matrix */ @@ -304,7 +309,7 @@ static void set_lj_parameter_data(nbnxn_atomdata_t::Params* params, gmx_bool bSI params->nbfp_comb[i * 2 + 1] = std::sqrt(params->nbfp[(i * nt + i) * 2 + 1]); } break; - case ljcrLB: + case LJCombinationRule::LorentzBerthelot: for (int i = 0; i < nt; i++) { /* Get 6*C6 and 12*C12 from the diagonal of the nbfp matrix */ @@ -325,7 +330,7 @@ static void set_lj_parameter_data(nbnxn_atomdata_t::Params* params, gmx_bool bSI } } break; - case ljcrNONE: + case LJCombinationRule::None: /* We always store the full matrix (see code above) */ break; default: gmx_incons("Unknown combination rule"); @@ -563,37 +568,41 @@ static void nbnxn_atomdata_params_init(const gmx::MDLogger& mdlog, */ if (bCombGeom) { - params->comb_rule = ljcrGEOM; + params->ljCombinationRule = LJCombinationRule::Geometric; } else if (bCombLB) { - params->comb_rule = ljcrLB; + params->ljCombinationRule = LJCombinationRule::LorentzBerthelot; } else { - params->comb_rule = ljcrNONE; + params->ljCombinationRule = LJCombinationRule::None; params->nbfp_comb.clear(); } { std::string mesg; - if (params->comb_rule == ljcrNONE) + if (params->ljCombinationRule == LJCombinationRule::None) { mesg = "Using full Lennard-Jones parameter combination matrix"; } else { mesg = gmx::formatString("Using %s Lennard-Jones combination rule", - enum_name(params->comb_rule, ljcrNR, c_ljcrNames)); + enumValueToString(params->ljCombinationRule)); } GMX_LOG(mdlog.info).asParagraph().appendText(mesg); } break; - case enbnxninitcombruleGEOM: params->comb_rule = ljcrGEOM; break; - case enbnxninitcombruleLB: params->comb_rule = ljcrLB; break; + case enbnxninitcombruleGEOM: + params->ljCombinationRule = LJCombinationRule::Geometric; + break; + case enbnxninitcombruleLB: + params->ljCombinationRule = LJCombinationRule::LorentzBerthelot; + break; case enbnxninitcombruleNONE: - params->comb_rule = ljcrNONE; + params->ljCombinationRule = LJCombinationRule::None; params->nbfp_comb.clear(); break; @@ -748,7 +757,7 @@ static void nbnxn_atomdata_set_ljcombparams(nbnxn_atomdata_t::Params* params, { params->lj_comb.resize(gridSet.numGridAtomsTotal() * 2); - if (params->comb_rule != ljcrNONE) + if (params->ljCombinationRule != LJCombinationRule::None) { for (const Nbnxm::Grid& grid : gridSet.grids()) { diff --git a/src/gromacs/nbnxm/atomdata.h b/src/gromacs/nbnxm/atomdata.h index 6f988eb14f..0a269c280b 100644 --- a/src/gromacs/nbnxm/atomdata.h +++ b/src/gromacs/nbnxm/atomdata.h @@ -2,7 +2,7 @@ * This file is part of the GROMACS molecular simulation package. * * Copyright (c) 2012,2013,2014,2015,2016 by the GROMACS development team. - * Copyright (c) 2017,2018,2019,2020, by the GROMACS development team, led by + * Copyright (c) 2017,2018,2019,2020,2021, 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. @@ -155,17 +155,21 @@ struct nbnxn_atomdata_output_t #define NBNXN_BUFFERFLAG_MAX_THREADS (BITMASK_SIZE) -/*! \brief LJ combination rules: geometric, Lorentz-Berthelot, none */ -enum +//! LJ combination rules +enum class LJCombinationRule : int { - ljcrGEOM, - ljcrLB, - ljcrNONE, - ljcrNR + //! Geometric + Geometric, + //! Lorentz-Berthelot + LorentzBerthelot, + //! No rule + None, + //! Size of the enum + Count }; //! String corresponding to LJ combination rule -extern const char* const c_ljcrNames[ljcrNR + 1]; +const char* enumValueToString(LJCombinationRule enumValue); /*! \internal * \brief Struct that stores atom related data for the nbnxn module @@ -190,7 +194,7 @@ struct nbnxn_atomdata_t //! Lennard-Jone 6*C6 and 12*C12 parameters, size numTypes*2*2 gmx::HostVector nbfp; //! Combination rule, see enum defined above - int comb_rule; + LJCombinationRule ljCombinationRule; //! LJ parameters per atom type, size numTypes*2 gmx::HostVector nbfp_comb; //! As nbfp, but with a stride for the present SIMD architecture diff --git a/src/gromacs/nbnxm/cuda/nbnxm_cuda_data_mgmt.cu b/src/gromacs/nbnxm/cuda/nbnxm_cuda_data_mgmt.cu index 5fa2e00eec..f9c590e9db 100644 --- a/src/gromacs/nbnxm/cuda/nbnxm_cuda_data_mgmt.cu +++ b/src/gromacs/nbnxm/cuda/nbnxm_cuda_data_mgmt.cu @@ -137,7 +137,7 @@ static void init_nbparam(NBParamGpu* nbp, * combination is rarely used. LJ force-switch with LB rule is more common, * but gives only 1% speed-up. */ - nbp->vdwType = nbnxmGpuPickVdwKernelType(ic, nbatParams.comb_rule); + nbp->vdwType = nbnxmGpuPickVdwKernelType(ic, nbatParams.ljCombinationRule); nbp->elecType = nbnxmGpuPickElectrostaticsKernelType(ic, deviceContext.deviceInfo()); /* generate table for PME */ diff --git a/src/gromacs/nbnxm/gpu_data_mgmt.h b/src/gromacs/nbnxm/gpu_data_mgmt.h index 1bcc676879..5b0c22085c 100644 --- a/src/gromacs/nbnxm/gpu_data_mgmt.h +++ b/src/gromacs/nbnxm/gpu_data_mgmt.h @@ -135,7 +135,8 @@ enum ElecType nbnxmGpuPickElectrostaticsKernelType(const interaction_const_t gmx /** Return the enum value of VdW kernel type for given \p ic and \p combRule. */ GPU_FUNC_QUALIFIER -enum VdwType nbnxmGpuPickVdwKernelType(const interaction_const_t gmx_unused* ic, int gmx_unused combRule) +enum VdwType nbnxmGpuPickVdwKernelType(const interaction_const_t gmx_unused* ic, + LJCombinationRule gmx_unused ljCombinationRule) GPU_FUNC_TERM_WITH_RETURN(VdwType::Count); /** Returns an opaque pointer to the GPU command stream diff --git a/src/gromacs/nbnxm/kerneldispatch.cpp b/src/gromacs/nbnxm/kerneldispatch.cpp index 5b4c5d7906..c840addd4c 100644 --- a/src/gromacs/nbnxm/kerneldispatch.cpp +++ b/src/gromacs/nbnxm/kerneldispatch.cpp @@ -2,7 +2,7 @@ * This file is part of the GROMACS molecular simulation package. * * Copyright (c) 2012,2013,2014,2015,2016 by the GROMACS development team. - * Copyright (c) 2017,2018,2019,2020, by the GROMACS development team, led by + * Copyright (c) 2017,2018,2019,2020,2021, 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. @@ -203,11 +203,11 @@ static void nbnxn_kernel_cpu(const PairlistSet& pairlistSet, { case eintmodNONE: case eintmodPOTSHIFT: - switch (nbatParams.comb_rule) + switch (nbatParams.ljCombinationRule) { - case ljcrGEOM: vdwkt = vdwktLJCUT_COMBGEOM; break; - case ljcrLB: vdwkt = vdwktLJCUT_COMBLB; break; - case ljcrNONE: vdwkt = vdwktLJCUT_COMBNONE; break; + case LJCombinationRule::Geometric: vdwkt = vdwktLJCUT_COMBGEOM; break; + case LJCombinationRule::LorentzBerthelot: vdwkt = vdwktLJCUT_COMBLB; break; + case LJCombinationRule::None: vdwkt = vdwktLJCUT_COMBNONE; break; default: GMX_RELEASE_ASSERT(false, "Unknown combination rule"); } break; diff --git a/src/gromacs/nbnxm/nbnxm.h b/src/gromacs/nbnxm/nbnxm.h index bf94b022e7..9a212ce929 100644 --- a/src/gromacs/nbnxm/nbnxm.h +++ b/src/gromacs/nbnxm/nbnxm.h @@ -2,7 +2,7 @@ * This file is part of the GROMACS molecular simulation package. * * Copyright (c) 2012,2013,2014,2015,2017 by the GROMACS development team. - * Copyright (c) 2018,2019,2020, by the GROMACS development team, led by + * Copyright (c) 2018,2019,2020,2021, 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. @@ -128,6 +128,7 @@ struct gmx_mtop_t; struct NbnxmGpu; struct gmx_wallcycle; struct interaction_const_t; +enum class LJCombinationRule; struct nbnxn_atomdata_t; struct nonbonded_verlet_t; class PairSearch; diff --git a/src/gromacs/nbnxm/nbnxm_gpu_data_mgmt.cpp b/src/gromacs/nbnxm/nbnxm_gpu_data_mgmt.cpp index 114d758353..87ec7b5b21 100644 --- a/src/gromacs/nbnxm/nbnxm_gpu_data_mgmt.cpp +++ b/src/gromacs/nbnxm/nbnxm_gpu_data_mgmt.cpp @@ -357,7 +357,7 @@ enum ElecType nbnxmGpuPickElectrostaticsKernelType(const interaction_const_t* ic } -enum VdwType nbnxmGpuPickVdwKernelType(const interaction_const_t* ic, int combRule) +enum VdwType nbnxmGpuPickVdwKernelType(const interaction_const_t* ic, LJCombinationRule ljCombinationRule) { if (ic->vdwtype == evdwCUT) { @@ -365,17 +365,16 @@ enum VdwType nbnxmGpuPickVdwKernelType(const interaction_const_t* ic, int combRu { case eintmodNONE: case eintmodPOTSHIFT: - switch (combRule) + switch (ljCombinationRule) { - case ljcrNONE: return VdwType::Cut; - case ljcrGEOM: return VdwType::CutCombGeom; - case ljcrLB: return VdwType::CutCombLB; + case LJCombinationRule::None: return VdwType::Cut; + case LJCombinationRule::Geometric: return VdwType::CutCombGeom; + case LJCombinationRule::LorentzBerthelot: return VdwType::CutCombLB; default: GMX_THROW(gmx::InconsistentInputError(gmx::formatString( - "The requested LJ combination rule %s (%d) is not implemented in " + "The requested LJ combination rule %s is not implemented in " "the GPU accelerated kernels!", - enum_name(combRule, ljcrNR, c_ljcrNames), - combRule))); + enumValueToString(ljCombinationRule)))); } case eintmodFORCESWITCH: return VdwType::FSwitch; case eintmodPOTSWITCH: return VdwType::PSwitch; @@ -389,14 +388,14 @@ enum VdwType nbnxmGpuPickVdwKernelType(const interaction_const_t* ic, int combRu } else if (ic->vdwtype == evdwPME) { - if (ic->ljpme_comb_rule == ljcrGEOM) + if (ic->ljpme_comb_rule == eljpmeGEOM) { - assert(combRule == ljcrGEOM); + assert(ljCombinationRule == LJCombinationRule::Geometric); return VdwType::EwaldGeom; } else { - assert(combRule == ljcrLB); + assert(ljCombinationRule == LJCombinationRule::LorentzBerthelot); return VdwType::EwaldLB; } } diff --git a/src/gromacs/nbnxm/opencl/nbnxm_ocl_data_mgmt.cpp b/src/gromacs/nbnxm/opencl/nbnxm_ocl_data_mgmt.cpp index 290d5f669a..1badb5c04a 100644 --- a/src/gromacs/nbnxm/opencl/nbnxm_ocl_data_mgmt.cpp +++ b/src/gromacs/nbnxm/opencl/nbnxm_ocl_data_mgmt.cpp @@ -135,18 +135,20 @@ static void init_nbparam(NBParamGpu* nbp, { set_cutoff_parameters(nbp, ic, listParams); - nbp->vdwType = nbnxmGpuPickVdwKernelType(ic, nbatParams.comb_rule); + nbp->vdwType = nbnxmGpuPickVdwKernelType(ic, nbatParams.ljCombinationRule); nbp->elecType = nbnxmGpuPickElectrostaticsKernelType(ic, deviceContext.deviceInfo()); if (ic->vdwtype == evdwPME) { - if (ic->ljpme_comb_rule == ljcrGEOM) + if (ic->ljpme_comb_rule == eljpmeGEOM) { - GMX_ASSERT(nbatParams.comb_rule == ljcrGEOM, "Combination rule mismatch!"); + GMX_ASSERT(nbatParams.ljCombinationRule == LJCombinationRule::Geometric, + "Combination rule mismatch!"); } else { - GMX_ASSERT(nbatParams.comb_rule == ljcrLB, "Combination rule mismatch!"); + GMX_ASSERT(nbatParams.ljCombinationRule == LJCombinationRule::LorentzBerthelot, + "Combination rule mismatch!"); } } /* generate table for PME */