From cd9941cd2185866809dfa6c3f953d56bf373b4b7 Mon Sep 17 00:00:00 2001 From: Berk Hess Date: Thu, 28 Oct 2021 13:50:39 +0000 Subject: [PATCH] Document and clean up do_numbering() --- src/gromacs/gmxpreprocess/readir.cpp | 81 +++++++++++++++------------- src/gromacs/gmxpreprocess/readir.h | 12 +++-- 2 files changed, 51 insertions(+), 42 deletions(-) diff --git a/src/gromacs/gmxpreprocess/readir.cpp b/src/gromacs/gmxpreprocess/readir.cpp index 4be0d276fa..ff47ab5a35 100644 --- a/src/gromacs/gmxpreprocess/readir.cpp +++ b/src/gromacs/gmxpreprocess/readir.cpp @@ -139,15 +139,13 @@ void done_inputrec_strings() } -enum +//! How to treat coverage of the whole system for a set of atom groupsx +enum class GroupCoverage { - egrptpALL, /* All particles have to be a member of a group. */ - egrptpALL_GENREST, /* A rest group with name is generated for particles * - * that are not part of any group. */ - egrptpPART, /* As egrptpALL_GENREST, but no name is generated * - * for the rest group. */ - egrptpONE /* Merge all selected groups into one group, * - * make a rest group for the remaining particles. */ + All, //!< All particles have to be a member of a group + AllGenerateRest, // groupsFromMdpFile, - t_blocka* block, - char* gnames[], - SimulationAtomGroupType gtype, - int restnm, - int grptp, - bool bVerbose, - warninp_t wi) +/*! Creates the groups of atom indices for group type \p gtype + * + * \param[in] natoms The total number of atoms in the system + * \param[in,out] groups Index \p gtype in this list of list of groups will be set + * \param[in] groupsFromMdpFile The list of group names set for \p gtype in the mdp file + * \param[in] block The list of atom indices for all available index groups + * \param[in] gnames The list of names for all available index groups + * \param[in] gtype The group type to creates groups for + * \param[in] restnm The index of rest group name in \p gnames + * \param[in] coverage How to treat coverage of all atoms in the system + * \param[in] bVerbose Whether to print when we make a rest group + * \param[in,out] wi List of warnings + */ +static void do_numbering(const int natoms, + SimulationGroups* groups, + gmx::ArrayRef groupsFromMdpFile, + const t_blocka* block, + char* const gnames[], + const SimulationAtomGroupType gtype, + const int restnm, + const GroupCoverage coverage, + const bool bVerbose, + warninp_t wi) { unsigned short* cbuf; AtomGroupIndices* grps = &(groups->groups[gtype]); @@ -2890,7 +2897,7 @@ static void do_numbering(int natoms, { /* Lookup the group name in the block structure */ const int gid = search_string(groupsFromMdpFile[i].c_str(), block->nr, gnames); - if ((grptp != egrptpONE) || (i == 0)) + if ((coverage != GroupCoverage::OneGroup) || (i == 0)) { grps->emplace_back(gid); } @@ -2909,7 +2916,7 @@ static void do_numbering(int natoms, else { /* Store the group number in buffer */ - if (grptp == egrptpONE) + if (coverage == GroupCoverage::OneGroup) { cbuf[aj] = 0; } @@ -2925,11 +2932,11 @@ static void do_numbering(int natoms, /* Now check whether we have done all atoms */ if (ntot != natoms) { - if (grptp == egrptpALL) + if (coverage == GroupCoverage::All) { gmx_fatal(FARGS, "%d atoms are not part of any of the %s groups", natoms - ntot, title); } - else if (grptp == egrptpPART) + else if (coverage == GroupCoverage::Partial) { sprintf(warn_buf, "%d atoms are not part of any of the %s groups", natoms - ntot, title); warning_note(wi, warn_buf); @@ -2942,7 +2949,7 @@ static void do_numbering(int natoms, cbuf[j] = grps->size(); } } - if (grptp != egrptpPART) + if (coverage != GroupCoverage::Partial) { if (bVerbose) { @@ -3555,7 +3562,7 @@ void do_index(const char* mdparin, gnames, SimulationAtomGroupType::TemperatureCoupling, restnm, - useReferenceTemperature ? egrptpALL : egrptpALL_GENREST, + useReferenceTemperature ? GroupCoverage::All : GroupCoverage::AllGenerateRest, bVerbose, wi); nr = groups->groups[SimulationAtomGroupType::TemperatureCoupling].size(); @@ -3916,7 +3923,7 @@ void do_index(const char* mdparin, gnames, SimulationAtomGroupType::Freeze, restnm, - egrptpALL_GENREST, + GroupCoverage::AllGenerateRest, bVerbose, wi); nr = groups->groups[SimulationAtomGroupType::Freeze].size(); @@ -3956,7 +3963,7 @@ void do_index(const char* mdparin, gnames, SimulationAtomGroupType::EnergyOutput, restnm, - egrptpALL_GENREST, + GroupCoverage::AllGenerateRest, bVerbose, wi); add_wall_energrps(groups, ir->nwall, symtab); @@ -3969,7 +3976,7 @@ void do_index(const char* mdparin, gnames, SimulationAtomGroupType::MassCenterVelocityRemoval, restnm, - vcmGroupNames.empty() ? egrptpALL_GENREST : egrptpPART, + vcmGroupNames.empty() ? GroupCoverage::AllGenerateRest : GroupCoverage::Partial, bVerbose, wi); @@ -3989,7 +3996,7 @@ void do_index(const char* mdparin, gnames, SimulationAtomGroupType::User1, restnm, - egrptpALL_GENREST, + GroupCoverage::AllGenerateRest, bVerbose, wi); auto user2GroupNames = gmx::splitString(inputrecStrings->user2); @@ -4000,7 +4007,7 @@ void do_index(const char* mdparin, gnames, SimulationAtomGroupType::User2, restnm, - egrptpALL_GENREST, + GroupCoverage::AllGenerateRest, bVerbose, wi); auto compressedXGroupNames = gmx::splitString(inputrecStrings->x_compressed_groups); @@ -4011,7 +4018,7 @@ void do_index(const char* mdparin, gnames, SimulationAtomGroupType::CompressedPositionOutput, restnm, - egrptpONE, + GroupCoverage::OneGroup, bVerbose, wi); auto orirefFitGroupNames = gmx::splitString(inputrecStrings->orirefitgrp); @@ -4022,7 +4029,7 @@ void do_index(const char* mdparin, gnames, SimulationAtomGroupType::OrientationRestraintsFit, restnm, - egrptpALL_GENREST, + GroupCoverage::AllGenerateRest, bVerbose, wi); @@ -4040,7 +4047,7 @@ void do_index(const char* mdparin, gnames, SimulationAtomGroupType::QuantumMechanics, restnm, - egrptpALL_GENREST, + GroupCoverage::AllGenerateRest, bVerbose, wi); ir->opts.ngQM = qmGroupNames.size(); diff --git a/src/gromacs/gmxpreprocess/readir.h b/src/gromacs/gmxpreprocess/readir.h index 9b7aa8194d..fef3384e11 100644 --- a/src/gromacs/gmxpreprocess/readir.h +++ b/src/gromacs/gmxpreprocess/readir.h @@ -112,16 +112,18 @@ void init_inputrec_strings(); /*! \brief Clean up object that holds strings parsed from an .mdp file */ void done_inputrec_strings(); +/*! \brief Performs all validation on \p ir that can be done without index groups and topology + * + * Any errors, warnings or notes are added to \p wi + */ void check_ir(const char* mdparin, const gmx::MDModulesNotifiers& mdModulesNotifiers, t_inputrec* ir, t_gromppopts* opts, warninp_t wi); -/* Validate inputrec data. - * Fatal errors will be added to nerror. - */ -int search_string(const char* s, int ng, char* gn[]); -/* Returns the index of string s in the index groups */ + +//! Returns the index of string \p s in \p gn or exit with a verbose fatal error when not found +int search_string(const char* s, int ng, char* const gn[]); void double_check(t_inputrec* ir, matrix box, bool bHasNormalConstraints, bool bHasAnyConstraints, warninp_t wi); /* Do more checks */ -- 2.22.0