From: Berk Hess Date: Wed, 31 Mar 2021 13:54:14 +0000 (+0000) Subject: Use EnumerationArray for pull isAngleType X-Git-Url: http://biod.pnpi.spb.ru/gitweb/?a=commitdiff_plain;h=0d6f1071431ea5388cb3dca21cc8bd157506fd82;p=alexxy%2Fgromacs.git Use EnumerationArray for pull isAngleType --- diff --git a/src/gromacs/applied_forces/awh/awh.cpp b/src/gromacs/applied_forces/awh/awh.cpp index 4a5e5aeeb4..51948de02b 100644 --- a/src/gromacs/applied_forces/awh/awh.cpp +++ b/src/gromacs/applied_forces/awh/awh.cpp @@ -238,7 +238,7 @@ Awh::Awh(FILE* fplog, GMX_THROW(InvalidInputError( "Pull geometry 'direction-periodic' is not supported by AWH")); } - double conversionFactor = pull_coordinate_is_angletype(&pullCoord) ? gmx::c_deg2Rad : 1; + double conversionFactor = pull_conversion_factor_userinput2internal(pullCoord); pullCoordIndex.push_back(awhDimParam.coordinateIndex()); dimParams.emplace_back(DimParams::pullDimParams( conversionFactor, awhDimParam.forceConstant(), beta)); diff --git a/src/gromacs/applied_forces/awh/read_params.cpp b/src/gromacs/applied_forces/awh/read_params.cpp index c2b8e33d7d..3ddb5a2cb5 100644 --- a/src/gromacs/applied_forces/awh/read_params.cpp +++ b/src/gromacs/applied_forces/awh/read_params.cpp @@ -1233,7 +1233,7 @@ static void setStateDependentAwhPullDimParams(AwhDimParams* dimParams, /* The initial coordinate value, converted to external user units. */ double initialCoordinate = get_pull_coord_value(pull_work, dimParams->coordinateIndex(), &pbc); - initialCoordinate *= pull_conversion_factor_internal2userinput(&pullCoordParams); + initialCoordinate *= pull_conversion_factor_internal2userinput(pullCoordParams); dimParams->setInitialCoordinate(initialCoordinate); } diff --git a/src/gromacs/gmxana/gmx_wham.cpp b/src/gromacs/gmxana/gmx_wham.cpp index 9f263dc3a0..55be0a50c6 100644 --- a/src/gromacs/gmxana/gmx_wham.cpp +++ b/src/gromacs/gmxana/gmx_wham.cpp @@ -1580,16 +1580,15 @@ static void read_tpr_header(const char* fn, t_UmbrellaHeader* header, t_Umbrella * so we need to multiply with the internal units (radians for angle) * to user units (degrees for an angle) with the same power. */ - header->pcrd[i].k = - ir->pull->coord[i].k - / gmx::square(pull_conversion_factor_internal2userinput(&ir->pull->coord[i])); + header->pcrd[i].k = ir->pull->coord[i].k + / gmx::square(pull_conversion_factor_internal2userinput(ir->pull->coord[i])); header->pcrd[i].init_dist = ir->pull->coord[i].init; copy_ivec(ir->pull->coord[i].dim, header->pcrd[i].dim); header->pcrd[i].ndim = header->pcrd[i].dim[XX] + header->pcrd[i].dim[YY] + header->pcrd[i].dim[ZZ]; - std::strcpy(header->pcrd[i].coord_unit, pull_coordinate_units(&ir->pull->coord[i])); + std::strcpy(header->pcrd[i].coord_unit, pull_coordinate_units(ir->pull->coord[i])); if (ir->efep != FreeEnergyPerturbationType::No && ir->pull->coord[i].k != ir->pull->coord[i].kB) { diff --git a/src/gromacs/gmxpreprocess/readpull.cpp b/src/gromacs/gmxpreprocess/readpull.cpp index 2c1cae68f2..dc724a3599 100644 --- a/src/gromacs/gmxpreprocess/readpull.cpp +++ b/src/gromacs/gmxpreprocess/readpull.cpp @@ -662,8 +662,8 @@ pull_t* set_pull_init(t_inputrec* ir, value = get_pull_coord_value(pull_work, c, &pbc); - value *= pull_conversion_factor_internal2userinput(pcrd); - fprintf(stderr, " %10.3f %s", value, pull_coordinate_units(pcrd)); + value *= pull_conversion_factor_internal2userinput(*pcrd); + fprintf(stderr, " %10.3f %s", value, pull_coordinate_units(*pcrd)); if (pcrd->bStart) { @@ -712,7 +712,7 @@ pull_t* set_pull_init(t_inputrec* ir, } - fprintf(stderr, " %10.3f %s\n", pcrd->init, pull_coordinate_units(pcrd)); + fprintf(stderr, " %10.3f %s\n", pcrd->init, pull_coordinate_units(*pcrd)); } return pull_work; diff --git a/src/gromacs/pulling/output.cpp b/src/gromacs/pulling/output.cpp index 220b2d6063..af388dd9a6 100644 --- a/src/gromacs/pulling/output.cpp +++ b/src/gromacs/pulling/output.cpp @@ -174,7 +174,7 @@ static void pull_print_coord_dr(FILE* out, double referenceValue, const int numValuesInSum) { - const double unit_factor = pull_conversion_factor_internal2userinput(&coordParams); + const double unit_factor = pull_conversion_factor_internal2userinput(coordParams); fprintf(out, "\t%g", pcrdData.value * unit_factor / numValuesInSum); diff --git a/src/gromacs/pulling/pull.cpp b/src/gromacs/pulling/pull.cpp index ba88e3aed5..e56fb9761a 100644 --- a/src/gromacs/pulling/pull.cpp +++ b/src/gromacs/pulling/pull.cpp @@ -72,6 +72,7 @@ #include "gromacs/topology/topology.h" #include "gromacs/utility/arrayref.h" #include "gromacs/utility/cstringutil.h" +#include "gromacs/utility/enumerationhelpers.h" #include "gromacs/utility/exceptions.h" #include "gromacs/utility/fatalerror.h" #include "gromacs/utility/futil.h" @@ -92,6 +93,11 @@ extern template LocalAtomSet LocalAtomSetManager::add(ArrayRef sc_isAngleType = { + { false, false, false, false, false, true, true, true } +}; + static int groupPbcFromParams(const t_pull_group& params, bool setPbcRefToPrevStepCOM) { if (params.ind.size() <= 1) @@ -138,12 +144,6 @@ pull_group_work_t::pull_group_work_t(const t_pull_group& params, clear_dvec(xp); }; -bool pull_coordinate_is_angletype(const t_pull_coord* pcrd) -{ - return (pcrd->eGeom == PullGroupGeometry::Angle || pcrd->eGeom == PullGroupGeometry::Dihedral - || pcrd->eGeom == PullGroupGeometry::AngleAxis); -} - static bool pull_coordinate_is_directional(const t_pull_coord* pcrd) { return (pcrd->eGeom == PullGroupGeometry::Direction || pcrd->eGeom == PullGroupGeometry::DirectionPBC @@ -151,14 +151,14 @@ static bool pull_coordinate_is_directional(const t_pull_coord* pcrd) || pcrd->eGeom == PullGroupGeometry::Cylinder); } -const char* pull_coordinate_units(const t_pull_coord* pcrd) +const char* pull_coordinate_units(const t_pull_coord& pcrd) { - return pull_coordinate_is_angletype(pcrd) ? "deg" : "nm"; + return sc_isAngleType[pcrd.eGeom] ? "deg" : "nm"; } -double pull_conversion_factor_userinput2internal(const t_pull_coord* pcrd) +double pull_conversion_factor_userinput2internal(const t_pull_coord& pcrd) { - if (pull_coordinate_is_angletype(pcrd)) + if (sc_isAngleType[pcrd.eGeom]) { return gmx::c_deg2Rad; } @@ -168,9 +168,9 @@ double pull_conversion_factor_userinput2internal(const t_pull_coord* pcrd) } } -double pull_conversion_factor_internal2userinput(const t_pull_coord* pcrd) +double pull_conversion_factor_internal2userinput(const t_pull_coord& pcrd) { - if (pull_coordinate_is_angletype(pcrd)) + if (sc_isAngleType[pcrd.eGeom]) { return gmx::c_rad2Deg; } @@ -665,7 +665,7 @@ static double sanitizePullCoordReferenceValue(const t_pull_coord& pcrdParams, do "Pull reference angle for coordinate %d (%f) needs to be in the allowed " "interval [0,180] deg", pcrdParams.coordIndex + 1, - value_ref * pull_conversion_factor_internal2userinput(&pcrdParams)); + value_ref * pull_conversion_factor_internal2userinput(pcrdParams)); } } else if (pcrdParams.eGeom == PullGroupGeometry::Dihedral) @@ -685,7 +685,7 @@ static void updatePullCoordReferenceValue(double* referenceValue, const t_pull_c if (pcrdParams.rate != 0) { const double inputValue = (pcrdParams.init + pcrdParams.rate * t) - * pull_conversion_factor_userinput2internal(&pcrdParams); + * pull_conversion_factor_userinput2internal(pcrdParams); *referenceValue = sanitizePullCoordReferenceValue(pcrdParams, inputValue); } } @@ -2128,7 +2128,7 @@ struct pull_t* init_pull(FILE* fplog, { pcrd->value_ref = sanitizePullCoordReferenceValue( pcrd->params, - pcrd->params.init * pull_conversion_factor_userinput2internal(&pcrd->params)); + pcrd->params.init * pull_conversion_factor_userinput2internal(pcrd->params)); } else { diff --git a/src/gromacs/pulling/pull.h b/src/gromacs/pulling/pull.h index 4016bed9fa..a4e0cf0551 100644 --- a/src/gromacs/pulling/pull.h +++ b/src/gromacs/pulling/pull.h @@ -78,33 +78,26 @@ class ForceWithVirial; class LocalAtomSetManager; } // namespace gmx -/*! \brief Returns if the pull coordinate is an angle - * - * \param[in] pcrd The pull coordinate to query the type for. - * \returns a boolean telling if the coordinate is of angle type. - */ -bool pull_coordinate_is_angletype(const t_pull_coord* pcrd); - /*! \brief Returns the units of the pull coordinate. * * \param[in] pcrd The pull coordinate to query the units for. * \returns a string with the units of the coordinate. */ -const char* pull_coordinate_units(const t_pull_coord* pcrd); +const char* pull_coordinate_units(const t_pull_coord& pcrd); /*! \brief Returns the conversion factor from the pull coord init/rate unit to internal value unit. * * \param[in] pcrd The pull coordinate to get the conversion factor for. * \returns the conversion factor. */ -double pull_conversion_factor_userinput2internal(const t_pull_coord* pcrd); +double pull_conversion_factor_userinput2internal(const t_pull_coord& pcrd); /*! \brief Returns the conversion factor from the pull coord internal value unit to the init/rate unit. * * \param[in] pcrd The pull coordinate to get the conversion factor for. * \returns the conversion factor. */ -double pull_conversion_factor_internal2userinput(const t_pull_coord* pcrd); +double pull_conversion_factor_internal2userinput(const t_pull_coord& pcrd); /*! \brief Get the value for pull coord coord_ind. *