assign_param(t_functype ftype, t_iparams *newparam,
real old[MAXFORCEPARAM], int comb, double reppow)
{
- int i, j;
- bool all_param_zero = TRUE;
+ bool all_param_zero = true;
/* Set to zero */
- for (j = 0; (j < MAXFORCEPARAM); j++)
+ for (int j = 0; (j < MAXFORCEPARAM); j++)
{
newparam->generic.buf[j] = 0.0;
/* If all parameters are zero we might not add some interaction types (selected below).
break;
case F_QUARTIC_ANGLES:
newparam->qangle.theta = old[0];
- for (i = 0; i < 5; i++)
+ for (int i = 0; i < 5; i++)
{
newparam->qangle.c[i] = old[i+1];
}
newparam->dihres.kfacB = old[5];
break;
case F_RBDIHS:
- for (i = 0; (i < NR_RBDIHS); i++)
+ for (int i = 0; (i < NR_RBDIHS); i++)
{
newparam->rbdihs.rbcA[i] = old[i];
newparam->rbdihs.rbcB[i] = old[NR_RBDIHS+i];
}
break;
case F_CBTDIHS:
- for (i = 0; (i < NR_CBTDIHS); i++)
+ for (int i = 0; (i < NR_CBTDIHS); i++)
{
newparam->cbtdihs.cbtcA[i] = old[i];
}
return type;
}
-static int vsite_nm2type(const char *name, gpp_atomtype *atype)
+static int vsite_nm2type(const char *name, PreprocessingAtomTypes *atype)
{
int tp;
- tp = get_atomtype_type(name, atype);
+ tp = atype->atomTypeFromName(name);
if (tp == NOTSET)
{
gmx_fatal(FARGS, "Dummy mass type (%s) not found in atom type database",
}
-static int gen_vsites_trp(gpp_atomtype *atype,
+static int gen_vsites_trp(PreprocessingAtomTypes *atype,
std::vector<gmx::RVec> *newx,
t_atom *newatom[], char ***newatomname[],
int *o2n[], int *newvsite_type[], int *newcgnr[],
}
-static int gen_vsites_tyr(gpp_atomtype *atype,
+static int gen_vsites_tyr(PreprocessingAtomTypes *atype,
std::vector<gmx::RVec> *newx,
t_atom *newatom[], char ***newatomname[],
int *o2n[], int *newvsite_type[], int *newcgnr[],
static char atomnamesuffix[] = "1234";
-void do_vsites(gmx::ArrayRef<const PreprocessResidue> rtpFFDB, gpp_atomtype *atype,
+void do_vsites(gmx::ArrayRef<const PreprocessResidue> rtpFFDB, PreprocessingAtomTypes *atype,
t_atoms *at, t_symtab *symtab,
std::vector<gmx::RVec> *x,
gmx::ArrayRef<InteractionTypeParameters> plist, int *vsite_type[], int *cgnr[],
&nrbonds, &nrHatoms, Hatoms, &Heavy, &nrheavies, heavies);
/* get Heavy atom type */
tpHeavy = get_atype(Heavy, at, rtpFFDB, &rt);
- strcpy(tpname, get_atomtype_name(tpHeavy, atype));
+ strcpy(tpname, atype->atomNameFromAtomType(tpHeavy));
bWARNING = FALSE;
bAddVsiteParam = TRUE;
}
/* get dummy mass type from first char of heavy atom type (N or C) */
- strcpy(nexttpname, get_atomtype_name(get_atype(heavies[0], at, rtpFFDB, &rt), atype));
+ strcpy(nexttpname, atype->atomNameFromAtomType(get_atype(heavies[0], at, rtpFFDB, &rt)));
std::string ch = get_dummymass_name(vsiteconflist, tpname, nexttpname);
std::string name;
if (ch.empty())
#include "gromacs/utility/arrayref.h"
#include "gromacs/utility/real.h"
-struct gpp_atomtype;
+class PreprocessingAtomTypes;
struct t_atoms;
struct InteractionTypeParameters;
struct PreprocessResidue;
/* stuff for pdb2gmx */
-void do_vsites(gmx::ArrayRef<const PreprocessResidue> rtpFFDB, gpp_atomtype *atype,
+void do_vsites(gmx::ArrayRef<const PreprocessResidue> rtpFFDB, PreprocessingAtomTypes *atype,
t_atoms *at, t_symtab *symtab, std::vector<gmx::RVec> *x,
gmx::ArrayRef<InteractionTypeParameters> plist, int *dummy_type[], int *cgnr[],
real mHmult, bool bVSiteAromatics,
#include <cmath>
#include <cstring>
+#include <algorithm>
+
#include "gromacs/gmxpreprocess/grompp_impl.h"
#include "gromacs/gmxpreprocess/notset.h"
#include "gromacs/gmxpreprocess/topdirs.h"
#include "gromacs/utility/fatalerror.h"
#include "gromacs/utility/smalloc.h"
-struct gpp_atomtype
+struct AtomTypeData
{
- int nr; /* The number of atomtypes */
- t_atom *atom; /* Array of atoms */
- char ***atomname; /* Names of the atomtypes */
- t_param *nb; /* Nonbonded force default params */
- int *bondatomtype; /* The bond_atomtype for each atomtype */
- int *atomnumber; /* Atomic number, used for QM/MM */
+ //! Constructor initializes all fields.
+ AtomTypeData(const t_atom &a,
+ char **name,
+ const t_param *nb,
+ const int bondAtomType,
+ const int atomNumber) :
+ atom_(a), name_(name),
+ bondAtomType_(bondAtomType),
+ atomNumber_(atomNumber)
+ {
+ cp_param(&nb_, nb);
+ }
+ //! Actual atom data.
+ t_atom atom_;
+ //! Atom name. The pointer is the result of a symtab operation and can be shallow copied.
+ char **name_;
+ //! Nonbonded data.
+ t_param nb_;
+ //! Bonded atomtype for the type.
+ int bondAtomType_;
+ //! Atom number for the atom type.
+ int atomNumber_;
};
-int get_atomtype_type(const char *str, gpp_atomtype *ga)
+class PreprocessingAtomTypes::Impl
{
- int i;
+ public:
+ //! The number for currently loaded entries.
+ size_t size() const { return types.size(); }
+ //! The actual atom type data.
+ std::vector<AtomTypeData> types;
+};
+bool PreprocessingAtomTypes::isSet(int nt) const
+{
+ return ((nt >= 0) || (nt < gmx::ssize(*this)));
+}
+
+int PreprocessingAtomTypes::atomTypeFromName(const std::string &str) const
+{
/* Atom types are always case sensitive */
- for (i = 0; (i < ga->nr); i++)
+ auto found = std::find_if(impl_->types.begin(), impl_->types.end(),
+ [&str](const auto &type)
+ { return str == *type.name_; });
+ if (found == impl_->types.end())
{
- if (strcmp(str, *(ga->atomname[i])) == 0)
- {
- return i;
- }
+ return NOTSET;
+ }
+ else
+ {
+ return std::distance(impl_->types.begin(), found);
}
-
- return NOTSET;
}
-int get_atomtype_ntypes(gpp_atomtype *ga)
+size_t PreprocessingAtomTypes::size() const
{
- return ga->nr;
+ return impl_->size();
}
-char *get_atomtype_name(int nt, gpp_atomtype *ga)
+const char *PreprocessingAtomTypes::atomNameFromAtomType(int nt) const
{
- if ((nt < 0) || (nt >= ga->nr))
- {
- return nullptr;
- }
-
- return *(ga->atomname[nt]);
+ return isSet(nt) ? *(impl_->types[nt].name_) : nullptr;
}
-real get_atomtype_massA(int nt, gpp_atomtype *ga)
+real PreprocessingAtomTypes::atomMassAFromAtomType(int nt) const
{
- if ((nt < 0) || (nt >= ga->nr))
- {
- return NOTSET;
- }
-
- return ga->atom[nt].m;
+ return isSet(nt) ? impl_->types[nt].atom_.m : NOTSET;
}
-real get_atomtype_massB(int nt, gpp_atomtype *ga)
+real PreprocessingAtomTypes::atomMassBFromAtomType(int nt) const
{
- if ((nt < 0) || (nt >= ga->nr))
- {
- return NOTSET;
- }
-
- return ga->atom[nt].mB;
+ return isSet(nt) ? impl_->types[nt].atom_.mB : NOTSET;
}
-real get_atomtype_qA(int nt, gpp_atomtype *ga)
+real PreprocessingAtomTypes::atomChargeAFromAtomType(int nt) const
{
- if ((nt < 0) || (nt >= ga->nr))
- {
- return NOTSET;
- }
-
- return ga->atom[nt].q;
+ return isSet(nt) ? impl_->types[nt].atom_.q : NOTSET;
}
-real get_atomtype_qB(int nt, gpp_atomtype *ga)
+real PreprocessingAtomTypes::atomChargeBFromAtomType(int nt) const
{
- if ((nt < 0) || (nt >= ga->nr))
- {
- return NOTSET;
- }
-
- return ga->atom[nt].qB;
+ return isSet(nt) ? impl_->types[nt].atom_.qB : NOTSET;
}
-int get_atomtype_ptype(int nt, gpp_atomtype *ga)
+int PreprocessingAtomTypes::atomParticleTypeFromAtomType(int nt) const
{
- if ((nt < 0) || (nt >= ga->nr))
- {
- return NOTSET;
- }
-
- return ga->atom[nt].ptype;
+ return isSet(nt) ? impl_->types[nt].atom_.ptype : NOTSET;
}
-int get_atomtype_batype(int nt, const gpp_atomtype* ga)
+int PreprocessingAtomTypes::bondAtomTypeFromAtomType(int nt) const
{
- if ((nt < 0) || (nt >= ga->nr))
- {
- return NOTSET;
- }
-
- return ga->bondatomtype[nt];
+ return isSet(nt) ? impl_->types[nt].bondAtomType_ : NOTSET;
}
-int get_atomtype_atomnumber(int nt, gpp_atomtype *ga)
+int PreprocessingAtomTypes::atomNumberFromAtomType(int nt) const
{
- if ((nt < 0) || (nt >= ga->nr))
- {
- return NOTSET;
- }
-
- return ga->atomnumber[nt];
+ return isSet(nt) ? impl_->types[nt].atomNumber_ : NOTSET;
}
-real get_atomtype_nbparam(int nt, int param, gpp_atomtype *ga)
+real PreprocessingAtomTypes::atomNonBondedParamFromAtomType(int nt, int param) const
{
- if ((nt < 0) || (nt >= ga->nr))
+ if (!isSet(nt))
{
return NOTSET;
}
{
return NOTSET;
}
- return ga->nb[nt].c[param];
+ return impl_->types[nt].nb_.c[param];
}
-gpp_atomtype *init_atomtype()
-{
- gpp_atomtype *ga;
-
- snew(ga, 1);
+PreprocessingAtomTypes::PreprocessingAtomTypes()
+ : impl_(new Impl)
+{}
- ga->nr = 0;
- ga->atom = nullptr;
- ga->atomname = nullptr;
- ga->nb = nullptr;
- ga->bondatomtype = nullptr;
- ga->atomnumber = nullptr;
-
- return ga;
-}
+PreprocessingAtomTypes::PreprocessingAtomTypes(PreprocessingAtomTypes &&old) noexcept
+ : impl_(std::move(old.impl_))
+{}
-int set_atomtype(int nt, gpp_atomtype *ga, t_symtab *tab,
- t_atom *a, const char *name, t_param *nb,
- int bondatomtype, int atomnumber)
+PreprocessingAtomTypes &PreprocessingAtomTypes::operator=(PreprocessingAtomTypes &&old) noexcept
{
- if ((nt < 0) || (nt >= ga->nr))
- {
- return NOTSET;
- }
-
- ga->atom[nt] = *a;
- ga->atomname[nt] = put_symtab(tab, name);
- ga->nb[nt] = *nb;
- ga->bondatomtype[nt] = bondatomtype;
- ga->atomnumber[nt] = atomnumber;
-
- return nt;
+ impl_ = std::move(old.impl_);
+ return *this;
}
-int add_atomtype(gpp_atomtype *ga, t_symtab *tab,
- t_atom *a, const char *name, t_param *nb,
- int bondatomtype, int atomnumber)
+PreprocessingAtomTypes::~PreprocessingAtomTypes()
+{}
+
+int PreprocessingAtomTypes::addType(t_symtab *tab,
+ const t_atom &a,
+ const char *name,
+ const t_param *nb,
+ int bondAtomType,
+ int atomNumber)
{
- int i;
+ auto found = std::find_if(impl_->types.begin(), impl_->types.end(),
+ [&name](const AtomTypeData &data)
+ { return strcmp(name, *data.name_) == 0; });
- for (i = 0; (i < ga->nr); i++)
- {
- if (strcmp(*ga->atomname[i], name) == 0)
- {
- break;
- }
- }
- if (i == ga->nr)
+ if (found == impl_->types.end())
{
- ga->nr++;
- srenew(ga->atom, ga->nr);
- srenew(ga->atomname, ga->nr);
- srenew(ga->nb, ga->nr);
- srenew(ga->bondatomtype, ga->nr);
- srenew(ga->atomnumber, ga->nr);
-
- return set_atomtype(ga->nr-1, ga, tab, a, name, nb, bondatomtype, atomnumber);
+ impl_->types.emplace_back(a,
+ put_symtab(tab, name),
+ nb,
+ bondAtomType,
+ atomNumber);
+ return size() - 1;
}
else
{
- return i;
+ return std::distance(impl_->types.begin(), found);
}
}
-void print_at (FILE * out, gpp_atomtype *ga)
+int PreprocessingAtomTypes::setType(int nt,
+ t_symtab *tab,
+ const t_atom &a,
+ const char *name,
+ const t_param *nb,
+ int bondAtomType,
+ int atomNumber)
{
- int i;
- t_atom *atom = ga->atom;
- t_param *nb = ga->nb;
+ if (!isSet(nt))
+ {
+ return NOTSET;
+ }
+
+ impl_->types[nt].atom_ = a;
+ impl_->types[nt].name_ = put_symtab(tab, name);
+ cp_param(&impl_->types[nt].nb_, nb);
+ impl_->types[nt].bondAtomType_ = bondAtomType;
+ impl_->types[nt].atomNumber_ = atomNumber;
+
+ return nt;
+}
+void PreprocessingAtomTypes::printTypes(FILE * out)
+{
fprintf (out, "[ %s ]\n", dir2str(Directive::d_atomtypes));
fprintf (out, "; %6s %8s %8s %8s %12s %12s\n",
"type", "mass", "charge", "particle", "c6", "c12");
- for (i = 0; (i < ga->nr); i++)
+ for (auto &entry : impl_->types)
{
fprintf(out, "%8s %8.3f %8.3f %8s %12e %12e\n",
- *(ga->atomname[i]), atom[i].m, atom[i].q, "A",
- nb[i].c0(), nb[i].c1());
+ *(entry.name_), entry.atom_.m, entry.atom_.q, "A",
+ entry.nb_.c0(), entry.nb_.c1());
}
fprintf (out, "\n");
}
-void done_atomtype(gpp_atomtype *ga)
+static int search_atomtypes(const PreprocessingAtomTypes *ga,
+ int *n,
+ gmx::ArrayRef<int> typelist,
+ int thistype,
+ t_param param[],
+ int ftype)
{
- sfree(ga->atom);
- sfree(ga->atomname);
- sfree(ga->nb);
- sfree(ga->bondatomtype);
- sfree(ga->atomnumber);
- ga->nr = 0;
- sfree(ga);
-}
-
-static int search_atomtypes(gpp_atomtype *ga, int *n, int typelist[],
- int thistype,
- t_param param[], int ftype)
-{
- int i, nn, nrfp, j, k, ntype, tli;
- bool bFound = FALSE;
+ int i, nn, nrfp, ntype;
nn = *n;
nrfp = NRFP(ftype);
- ntype = get_atomtype_ntypes(ga);
+ ntype = ga->size();
for (i = 0; (i < nn); i++)
{
/* Otherwise, check if the parameters are identical to any previously added type */
- bFound = TRUE;
- for (j = 0; j < ntype && bFound; j++)
+ bool bFound = true;
+ for (int j = 0; j < ntype && bFound; j++)
{
/* Check nonbonded parameters */
- for (k = 0; k < nrfp && bFound; k++)
+ for (int k = 0; k < nrfp && bFound; k++)
{
bFound = (param[ntype*typelist[i]+j].c[k] == param[ntype*thistype+j].c[k]);
}
/* Check atomnumber */
- tli = typelist[i];
+ int tli = typelist[i];
bFound = bFound &&
- (get_atomtype_atomnumber(tli, ga) == get_atomtype_atomnumber(thistype, ga));
+ (ga->atomNumberFromAtomType(tli) == ga->atomNumberFromAtomType(thistype));
}
if (bFound)
{
return i;
}
-void renum_atype(gmx::ArrayRef<InteractionTypeParameters> plist, gmx_mtop_t *mtop,
- int *wall_atomtype,
- gpp_atomtype *ga, bool bVerbose)
+void PreprocessingAtomTypes::renumberTypes(gmx::ArrayRef<InteractionTypeParameters> plist,
+ gmx_mtop_t *mtop,
+ int *wall_atomtype,
+ bool bVerbose)
{
- int i, j, k, l, mi, mj, nat, nrfp, ftype, ntype;
- t_atoms *atoms;
- t_param *nbsnew;
- int *typelist;
- int *new_atomnumber;
- char ***new_atomname;
+ int nat, ftype, ntype;
- ntype = get_atomtype_ntypes(ga);
- snew(typelist, ntype);
+ ntype = size();
+ std::vector<int> typelist(ntype);
if (bVerbose)
{
* can determine if two types should be merged.
*/
nat = 0;
- for (gmx_moltype_t &moltype : mtop->moltype)
+ for (const gmx_moltype_t &moltype : mtop->moltype)
{
- atoms = &moltype.atoms;
- for (i = 0; (i < atoms->nr); i++)
+ const t_atoms *atoms = &moltype.atoms;
+ for (int i = 0; (i < atoms->nr); i++)
{
atoms->atom[i].type =
- search_atomtypes(ga, &nat, typelist, atoms->atom[i].type,
+ search_atomtypes(this, &nat, typelist, atoms->atom[i].type,
plist[ftype].param, ftype);
atoms->atom[i].typeB =
- search_atomtypes(ga, &nat, typelist, atoms->atom[i].typeB,
+ search_atomtypes(this, &nat, typelist, atoms->atom[i].typeB,
plist[ftype].param, ftype);
}
}
- for (i = 0; i < 2; i++)
+ for (int i = 0; i < 2; i++)
{
if (wall_atomtype[i] >= 0)
{
- wall_atomtype[i] = search_atomtypes(ga, &nat, typelist, wall_atomtype[i],
+ wall_atomtype[i] = search_atomtypes(this, &nat, typelist, wall_atomtype[i],
plist[ftype].param, ftype);
}
}
- snew(new_atomnumber, nat);
- snew(new_atomname, nat);
+ std::vector<AtomTypeData> new_types;
/* We now have a list of unique atomtypes in typelist */
/* Renumber nlist */
- nbsnew = nullptr;
+ /* Renumber nlist */
+ t_param *nbsnew = nullptr;
snew(nbsnew, plist[ftype].nr);
- nrfp = NRFP(ftype);
+ int nrfp = NRFP(ftype);
- for (i = k = 0; (i < nat); i++)
+ int k = 0;
+ for (int i = 0; (i < nat); i++)
{
- mi = typelist[i];
- for (j = 0; (j < nat); j++, k++)
+ int mi = typelist[i];
+ for (int j = 0; (j < nat); j++, k++)
{
- mj = typelist[j];
- for (l = 0; (l < nrfp); l++)
+ int mj = typelist[j];
+ for (int l = 0; (l < nrfp); l++)
{
nbsnew[k].c[l] = plist[ftype].param[ntype*mi+mj].c[l];
}
}
- new_atomnumber[i] = get_atomtype_atomnumber(mi, ga);
- new_atomname[i] = ga->atomname[mi];
+ new_types.push_back(impl_->types[mi]);
}
+ int i;
for (i = 0; (i < nat*nat); i++)
{
- for (l = 0; (l < nrfp); l++)
+ for (int l = 0; (l < nrfp); l++)
{
plist[ftype].param[i].c[l] = nbsnew[i].c[l];
}
plist[ftype].nr = i;
mtop->ffparams.atnr = nat;
- sfree(ga->atomnumber);
- /* Dangling atomname pointers ? */
- sfree(ga->atomname);
-
- ga->atomnumber = new_atomnumber;
- ga->atomname = new_atomname;
-
- ga->nr = nat;
+ impl_->types = new_types;
sfree(nbsnew);
- sfree(typelist);
}
-void copy_atomtype_atomtypes(gpp_atomtype *ga, t_atomtypes *atomtypes)
+void PreprocessingAtomTypes::copyTot_atomtypes(t_atomtypes *atomtypes) const
{
- int i, ntype;
-
/* Copy the atomtype data to the topology atomtype list */
- ntype = get_atomtype_ntypes(ga);
+ int ntype = size();
atomtypes->nr = ntype;
snew(atomtypes->atomnumber, ntype);
- for (i = 0; i < ntype; i++)
+ for (int i = 0; i < ntype; i++)
{
- atomtypes->atomnumber[i] = ga->atomnumber[i];
+ atomtypes->atomnumber[i] = impl_->types[i].atomNumber_;
}
}
#include <cstdio>
#include "gromacs/utility/arrayref.h"
+#include "gromacs/utility/classhelpers.h"
#include "gromacs/utility/real.h"
struct gmx_mtop_t;
-struct gpp_atomtype;
struct t_atom;
struct t_atomtypes;
struct t_param;
struct InteractionTypeParameters;
struct t_symtab;
-int get_atomtype_type(const char *str, gpp_atomtype *at);
-/* Return atomtype corresponding to case-insensitive str
- or NOTSET if not found */
-
-int get_atomtype_ntypes(gpp_atomtype *at);
-/* Return number of atomtypes */
-
-char *get_atomtype_name(int nt, gpp_atomtype *at);
-/* Return name corresponding to atomtype nt, or NULL if not found */
-
-real get_atomtype_massA(int nt, gpp_atomtype *at);
-real get_atomtype_massB(int nt, gpp_atomtype *at);
-real get_atomtype_qA(int nt, gpp_atomtype *at);
-real get_atomtype_qB(int nt, gpp_atomtype *at);
-int get_atomtype_ptype(int nt, gpp_atomtype *at);
-int get_atomtype_batype(int nt, const gpp_atomtype* at);
-int get_atomtype_atomnumber(int nt, gpp_atomtype *at);
-
-/* Return the above variable for atomtype nt, or NOTSET if not found */
-
-real get_atomtype_nbparam(int nt, int param, gpp_atomtype *at);
-/* Similar to the previous but returns the paramth parameter or NOTSET */
-
-gpp_atomtype *init_atomtype();
-/* Return a new atomtype structure */
-
-void done_atomtype(gpp_atomtype *at);
-/* Free the memory in the structure */
-
-int set_atomtype(int nt, gpp_atomtype *at, t_symtab *tab,
- t_atom *a, const char *name, t_param *nb,
- int bondatomtype, int atomnumber);
-/* Set the values of an existing atom type nt. Returns nt on success or
- NOTSET on error. */
-
-int add_atomtype(gpp_atomtype *at, t_symtab *tab,
- t_atom *a, const char *name, t_param *nb,
- int bondatomtype, int atomnumber);
-/* Add a complete new atom type to an existing atomtype structure. Returns
- the number of the atom type. */
-
-void print_at (FILE * out, gpp_atomtype *at);
-/* Print an atomtype record to a text file */
-
-void renum_atype(gmx::ArrayRef<InteractionTypeParameters> plist, gmx_mtop_t *mtop,
- int *wall_atomtype,
- gpp_atomtype *at, bool bVerbose);
-
-void copy_atomtype_atomtypes(gpp_atomtype *atype, t_atomtypes *atypes);
-/* Copy from one structure to another */
+/*! \libinternal \brief
+ * Storage of all atom types used during preprocessing of a simulation
+ * input.
+ */
+class PreprocessingAtomTypes
+{
+ public:
+ PreprocessingAtomTypes();
+ //! Move constructor.
+ PreprocessingAtomTypes(PreprocessingAtomTypes &&old) noexcept;
+ //! Move assignment constructor.
+ PreprocessingAtomTypes &operator=(PreprocessingAtomTypes &&old) noexcept;
+
+ ~PreprocessingAtomTypes();
+
+ /*! \brief
+ * Get atom type index for atom type name if present in the database, or NOTSET.
+ *
+ * \todo The code should be changed to instead use a gmx::compat version
+ * of std::optional to return an iterator to the element being searched,
+ * or an empty optional construct if the entry has not been found.
+ *
+ * \param[in] str Input string to search type for.
+ * \returns Atomtype as integer.
+ */
+ int atomTypeFromName(const std::string &str) const;
+
+ //! Get number of defined atom types.
+ size_t size() const;
+
+ /*! \brief
+ * Get name of atom from internal atom type number.
+ *
+ * \param[in] nt Internal number of atom type.
+ * \returns The type name.
+ */
+ const char *atomNameFromAtomType(int nt) const;
+
+ /*! \brief
+ * Get normal mass of atom from internal atom type number.
+ *
+ * \param[in] nt Internal number of atom type.
+ * \returns The mass for the atom in state A or NOTSET.
+ */
+ real atomMassAFromAtomType(int nt) const;
+
+ /*! \brief
+ * Get mass for B state of atom from internal atom type number.
+ *
+ * \param[in] nt Internal number of atom type.
+ * \returns The mass for the atom in state B or NOTSET.
+ */
+ real atomMassBFromAtomType(int nt) const;
+
+ /*! \brief
+ * Get normal charge of atom from internal atom type number.
+ *
+ * \param[in] nt Internal number of atom type.
+ * \returns The charge for the atom in state A or NOTSET.
+ */
+ real atomChargeAFromAtomType(int nt) const;
+
+ /*! \brief
+ * Get charge for B state of atom from internal atom type number.
+ *
+ * \param[in] nt Internal number of atom type.
+ * \returns The charge for the atom in state B or NOTSET.
+ */
+ real atomChargeBFromAtomType(int nt) const;
+
+ /*! \brief
+ * Get particle type for atom type \p nt
+ *
+ * \param[in] nt Internal number of atom type.
+ * \returns The particle type or NOTSET.
+ */
+ int atomParticleTypeFromAtomType(int nt) const;
+
+ /*! \brief
+ * Get bond atom parameter of atom from internal atom type number.
+ *
+ * \param[in] nt Internal number of atom type.
+ * \returns The bond atom parameter or NOTSET.
+ */
+ int bondAtomTypeFromAtomType(int nt) const;
+
+ /*! \brief
+ * Get atomic number of atom from internal atom type number.
+ *
+ * \param[in] nt Internal number of atom type.
+ * \returns The atomic number type or NOTSET.
+ */
+ int atomNumberFromAtomType(int nt) const;
+
+ /*! \brief
+ * Get the value of \p param of type \p nt.
+ *
+ * \param[in] param The parameter value to find.
+ * \param[in] nt The number of the type.
+ * \returns The value of the parameter or NOTSET.
+ */
+ real atomNonBondedParamFromAtomType(int nt, int param) const;
+
+ /*! \brief
+ * If a value is within the range of the current types or not.
+ *
+ * \param[in] nt Value to check.
+ * \returns True if value is in range.
+ */
+ bool isSet(int nt) const;
+
+ /*! \brief
+ * Print data to file.
+ *
+ * \param[in] out File pointer.
+ */
+ void printTypes(FILE *out);
+
+ /*! \brief
+ * Set the values of an existing atom type \p nt.
+ *
+ * \param[in] nt Type that should be set.
+ * \param[in] tab Symbol table.
+ * \param[in] a Atom information.
+ * \param[in] name Atom name.
+ * \param[in] nb Nonbonded parameters.
+ * \param[in] bondAtomType What kind of bonded interactions are there.
+ * \param[in] atomNumber Atomic number of the entry.
+ * \returns Number of the type set or NOTSET
+ */
+ int setType(int nt,
+ t_symtab *tab,
+ const t_atom &a,
+ const char *name,
+ const t_param *nb,
+ int bondAtomType,
+ int atomNumber);
+
+ /*! \brief
+ * Add new atom type to database.
+ *
+ * \param[in] tab Symbol table.
+ * \param[in] a Atom information.
+ * \param[in] name Atom name.
+ * \param[in] nb Nonbonded parameters.
+ * \param[in] bondAtomType What kind of bonded interactions are there.
+ * \param[in] atomNumber Atomic number of the entry.
+ * \returns Number of entries in database.
+ */
+ int addType(t_symtab *tab,
+ const t_atom &a,
+ const char *name,
+ const t_param *nb,
+ int bondAtomType,
+ int atomNumber);
+
+ /*! \brief
+ * Renumber existing atom type entries.
+ *
+ * \param[in] plist List of parameters.
+ * \param[in] mtop Global topology.
+ * \param[inout] wallAtomType Atom types of wall atoms, which may also be renumbered
+ * \param[in] verbose If we want to print additional info.
+ */
+ void renumberTypes(gmx::ArrayRef<InteractionTypeParameters> plist,
+ gmx_mtop_t *mtop,
+ int *wallAtomType,
+ bool verbose);
+
+ /*! \brief
+ * Copy information to other structure.
+ *
+ * \param[inout] atypes Other datastructure to copy to.
+ */
+ void copyTot_atomtypes(t_atomtypes *atypes) const;
+ private:
+ class Impl;
+ gmx::PrivateImplPointer<Impl> impl_;
+};
#endif
new_status(const char *topfile, const char *topppfile, const char *confin,
t_gromppopts *opts, t_inputrec *ir, gmx_bool bZero,
bool bGenVel, bool bVerbose, t_state *state,
- gpp_atomtype *atype, gmx_mtop_t *sys,
+ PreprocessingAtomTypes *atypes, gmx_mtop_t *sys,
std::vector<MoleculeInformation> *mi,
std::unique_ptr<MoleculeInformation> *intermolecular_interactions,
gmx::ArrayRef<InteractionTypeParameters> plist,
warninp *wi)
{
std::vector<gmx_molblock_t> molblock;
- int i, nrmols, nmismatch;
+ int i, nmismatch;
bool ffParametrizedWithHBondConstraints;
char buf[STRLEN];
char warn_buf[STRLEN];
/* TOPOLOGY processing */
sys->name = do_top(bVerbose, topfile, topppfile, opts, bZero, &(sys->symtab),
plist, comb, reppow, fudgeQQ,
- atype, &nrmols, mi, intermolecular_interactions,
+ atypes, mi, intermolecular_interactions,
ir,
&molblock,
&ffParametrizedWithHBondConstraints,
if (bMorse)
{
- convert_harmonics(*mi, atype);
+ convert_harmonics(*mi, atypes);
}
if (ir->eDisre == edrNone)
read_posres(mtop, mi, TRUE, fnB, rc_scaling, ePBC, comB, wi);
}
-static void set_wall_atomtype(gpp_atomtype *at, t_gromppopts *opts,
+static void set_wall_atomtype(PreprocessingAtomTypes *at, t_gromppopts *opts,
t_inputrec *ir, warninp *wi)
{
int i;
}
for (i = 0; i < ir->nwall; i++)
{
- ir->wall_atomtype[i] = get_atomtype_type(opts->wall_atomtype[i], at);
+ ir->wall_atomtype[i] = at->atomTypeFromName(opts->wall_atomtype[i]);
if (ir->wall_atomtype[i] == NOTSET)
{
sprintf(warn_buf, "Specified wall atom type %s is not defined", opts->wall_atomtype[i]);
t_gromppopts *opts;
std::vector<MoleculeInformation> mi;
std::unique_ptr<MoleculeInformation> intermolecular_interactions;
- gpp_atomtype *atype;
int nvsite, comb;
real fudgeQQ;
double reppow;
const char *mdparin;
- int ntype;
bool bNeedVel, bGenVel;
gmx_bool have_atomnumber;
gmx_output_env_t *oenv;
std::array<InteractionTypeParameters, F_NRE> plist;
init_plist(plist);
- gmx_mtop_t sys;
- atype = init_atomtype();
+ gmx_mtop_t sys;
+ PreprocessingAtomTypes atypes;
if (debug)
{
pr_symtab(debug, 0, "Just opened", &sys.symtab);
t_state state;
new_status(fn, opt2fn_null("-pp", NFILE, fnm), opt2fn("-c", NFILE, fnm),
opts, ir, bZero, bGenVel, bVerbose, &state,
- atype, &sys, &mi, &intermolecular_interactions,
+ &atypes, &sys, &mi, &intermolecular_interactions,
plist, &comb, &reppow, &fudgeQQ,
opts->bMorse,
wi);
for (size_t mt = 0; mt < sys.moltype.size(); mt++)
{
nvsite +=
- set_vsites(bVerbose, &sys.moltype[mt].atoms, atype, mi[mt].plist);
+ set_vsites(bVerbose, &sys.moltype[mt].atoms, &atypes, mi[mt].plist);
}
/* now throw away all obsolete bonds, angles and dihedrals: */
/* note: constraints are ALWAYS removed */
/* If we are doing QM/MM, check that we got the atom numbers */
have_atomnumber = TRUE;
- for (i = 0; i < get_atomtype_ntypes(atype); i++)
+ for (int i = 0; i < gmx::ssize(atypes); i++)
{
- have_atomnumber = have_atomnumber && (get_atomtype_atomnumber(i, atype) >= 0);
+ have_atomnumber = have_atomnumber && (atypes.atomNumberFromAtomType(i) >= 0);
}
if (!have_atomnumber && ir->bQMMM)
{
setup_cmap(plist[F_CMAP].cmakeGridSpacing, plist[F_CMAP].cmapAngles, plist[F_CMAP].cmap, &sys.ffparams.cmap_grid);
}
- set_wall_atomtype(atype, opts, ir, wi);
+ set_wall_atomtype(&atypes, opts, ir, wi);
if (bRenum)
{
- renum_atype(plist, &sys, ir->wall_atomtype, atype, bVerbose);
- get_atomtype_ntypes(atype);
+ atypes.renumberTypes(plist, &sys, ir->wall_atomtype, bVerbose);
}
if (ir->implicit_solvent)
}
/* PELA: Copy the atomtype data to the topology atomtype list */
- copy_atomtype_atomtypes(atype, &(sys.atomtypes));
+ atypes.copyTot_atomtypes(&(sys.atomtypes));
if (debug)
{
- pr_symtab(debug, 0, "After renum_atype", &sys.symtab);
+ pr_symtab(debug, 0, "After atype.renumberTypes", &sys.symtab);
}
if (bVerbose)
fprintf(stderr, "converting bonded parameters...\n");
}
- ntype = get_atomtype_ntypes(atype);
+ const int ntype = atypes.size();
convertInteractionTypeParameters(ntype, plist, mi, intermolecular_interactions.get(),
comb, reppow, fudgeQQ, &sys);
// fullCleanUp can't be called.
mol.partialCleanUp();
}
- done_atomtype(atype);
done_inputrec_strings();
output_env_done(oenv);
}
int nm2type(int nnm, t_nm2type nm2t[], t_symtab *tab, t_atoms *atoms,
- gpp_atomtype *atype, int *nbonds, InteractionTypeParameters *bonds)
+ PreprocessingAtomTypes *atype, int *nbonds, InteractionTypeParameters *bonds)
{
int cur = 0;
#define prev (1-cur)
- int i, j, k, m, n, nresolved, nb, maxbond, ai, aj, best, im, nqual[2][ematchNR];
+ int nresolved, nb, maxbond, nqual[2][ematchNR];
int *bbb, *n_mask, *m_mask, **match;
- char *aname_i, *aname_m, *aname_n, *type;
- double qq, mm;
+ char *aname_i, *aname_n;
t_param *param;
snew(param, 1);
maxbond = 0;
- for (i = 0; (i < atoms->nr); i++)
+ for (int i = 0; (i < atoms->nr); i++)
{
maxbond = std::max(maxbond, nbonds[i]);
}
snew(n_mask, maxbond);
snew(m_mask, maxbond);
snew(match, maxbond);
- for (i = 0; (i < maxbond); i++)
+ for (int i = 0; (i < maxbond); i++)
{
snew(match[i], maxbond);
}
nresolved = 0;
- for (i = 0; (i < atoms->nr); i++)
+ for (int i = 0; (i < atoms->nr); i++)
{
aname_i = *atoms->atomname[i];
nb = 0;
- for (j = 0; (j < bonds->nr); j++)
+ for (int j = 0; (j < bonds->nr); j++)
{
- ai = bonds->param[j].ai();
- aj = bonds->param[j].aj();
+ int ai = bonds->param[j].ai();
+ int aj = bonds->param[j].aj();
if (ai == i)
{
bbb[nb++] = aj;
if (debug)
{
fprintf(debug, "%4s has bonds to", aname_i);
- for (j = 0; (j < nb); j++)
+ for (int j = 0; (j < nb); j++)
{
fprintf(debug, " %4s", *atoms->atomname[bbb[j]]);
}
fprintf(debug, "\n");
}
- best = -1;
- for (k = 0; (k < ematchNR); k++)
+ int best = -1;
+ for (int k = 0; (k < ematchNR); k++)
{
nqual[prev][k] = 0;
}
/* First check for names */
- for (k = 0; (k < nnm); k++)
+ for (int k = 0; (k < nnm); k++)
{
if (nm2t[k].nbonds == nb)
{
- im = match_str(*atoms->atomname[i], nm2t[k].elem);
+ int im = match_str(*atoms->atomname[i], nm2t[k].elem);
if (im > ematchWild)
{
- for (j = 0; (j < ematchNR); j++)
+ for (int j = 0; (j < ematchNR); j++)
{
nqual[cur][j] = 0;
}
/* Fill a matrix with matching quality */
- for (m = 0; (m < nb); m++)
+ for (int m = 0; (m < nb); m++)
{
- aname_m = *atoms->atomname[bbb[m]];
- for (n = 0; (n < nb); n++)
+ const char *aname_m = *atoms->atomname[bbb[m]];
+ for (int n = 0; (n < nb); n++)
{
aname_n = nm2t[k].bond[n];
match[m][n] = match_str(aname_m, aname_n);
}
}
/* Now pick the best matches */
- for (m = 0; (m < nb); m++)
+ for (int m = 0; (m < nb); m++)
{
n_mask[m] = 0;
m_mask[m] = 0;
}
- for (j = ematchNR-1; (j > 0); j--)
+ for (int j = ematchNR-1; (j > 0); j--)
{
- for (m = 0; (m < nb); m++)
+ for (int m = 0; (m < nb); m++)
{
- for (n = 0; (n < nb); n++)
+ for (int n = 0; (n < nb); n++)
{
if ((n_mask[n] == 0) &&
(m_mask[m] == 0) &&
}
if (best != -1)
{
- int atomnr = 0;
- real alpha = 0;
+ int atomnr = 0;
+ real alpha = 0;
- qq = nm2t[best].q;
- mm = nm2t[best].m;
- type = nm2t[best].type;
+ double qq = nm2t[best].q;
+ double mm = nm2t[best].m;
+ const char *type = nm2t[best].type;
- if ((k = get_atomtype_type(type, atype)) == NOTSET)
+ int k;
+ if ((k = atype->atomTypeFromName(type)) == NOTSET)
{
atoms->atom[i].qB = alpha;
atoms->atom[i].m = atoms->atom[i].mB = mm;
- k = add_atomtype(atype, tab, &(atoms->atom[i]), type, param,
- atoms->atom[i].type, atomnr);
+ k = atype->addType(tab, atoms->atom[i], type, param,
+ atoms->atom[i].type, atomnr);
}
atoms->atom[i].type = k;
atoms->atom[i].typeB = k;
#include <cstdio>
-struct gpp_atomtype;
+class PreprocessingAtomTypes;
struct t_atoms;
struct InteractionTypeParameters;
struct t_symtab;
/* Dump the database for debugging. Can be reread by the program */
int nm2type(int nnm, t_nm2type nm2t[], t_symtab *tab, t_atoms *atoms,
- gpp_atomtype *atype, int *nbonds, InteractionTypeParameters *bond);
+ PreprocessingAtomTypes *atype, int *nbonds, InteractionTypeParameters *bond);
/* Try to determine the atomtype (force field dependent) for the atoms
* with help of the bond list
*/
check_occupancy(&pdba_all, inputConfFile_.c_str(), bVerbose_);
/* Read atomtypes... */
- gpp_atomtype *atype = read_atype(ffdir_, &symtab);
+ PreprocessingAtomTypes atype = read_atype(ffdir_, &symtab);
/* read residue database */
printf("Reading residue database... (%s)\n", forcefield_);
std::vector<PreprocessResidue> rtpFFDB;
for (const auto &filename : rtpf)
{
- readResidueDatabase(filename, &rtpFFDB, atype, &symtab, false);
+ readResidueDatabase(filename, &rtpFFDB, &atype, &symtab, false);
}
if (bNewRTP_)
{
std::vector<MoleculePatchDatabase> ntdb;
std::vector<MoleculePatchDatabase> ctdb;
std::vector<MoleculePatchDatabase *> tdblist;
- int nNtdb = read_ter_db(ffdir_, 'n', &ntdb, atype);
- int nCtdb = read_ter_db(ffdir_, 'c', &ctdb, atype);
+ int nNtdb = read_ter_db(ffdir_, 'n', &ntdb, &atype);
+ int nCtdb = read_ter_db(ffdir_, 'c', &ctdb, &atype);
FILE *top_file = gmx_fio_fopen(topologyFile_.c_str(), "w");
top_file2 = top_file;
}
- pdb2top(top_file2, posre_fn.c_str(), molname.c_str(), pdba, &x, atype, &symtab,
+ pdb2top(top_file2, posre_fn.c_str(), molname.c_str(), pdba, &x, &atype, &symtab,
rtpFFDB,
restp_chain, hb_chain,
bAllowMissing_,
cc->pdba = pdba;
cc->x = x;
}
- done_atomtype(atype);
if (watermodel_ == nullptr)
{
void write_top(FILE *out, const char *pr, const char *molname,
t_atoms *at, bool bRTPresname,
int bts[], gmx::ArrayRef<const InteractionTypeParameters> plist, t_excls excls[],
- gpp_atomtype *atype, int *cgnr, int nrexcl)
+ PreprocessingAtomTypes *atype, int *cgnr, int nrexcl)
/* NOTE: nrexcl is not the size of *excl! */
{
if (at && atype && cgnr)
void pdb2top(FILE *top_file, const char *posre_fn, const char *molname,
t_atoms *atoms,
- std::vector<gmx::RVec> *x, gpp_atomtype *atype, t_symtab *tab,
+ std::vector<gmx::RVec> *x, PreprocessingAtomTypes *atype, t_symtab *tab,
gmx::ArrayRef<const PreprocessResidue> rtpFFDB,
gmx::ArrayRef<PreprocessResidue> usedPpResidues,
gmx::ArrayRef<MoleculePatchDatabase> globalPatches,
#include "gromacs/math/vectypes.h"
#include "gromacs/utility/arrayref.h"
-struct gpp_atomtype;
+class PreprocessingAtomTypes;
struct t_atoms;
struct t_excls;
struct MoleculePatchDatabase;
void write_top(FILE *out, const char *pr, const char *molname,
t_atoms *at, bool bRTPresname,
int bts[], gmx::ArrayRef<const InteractionTypeParameters> plist, t_excls excls[],
- gpp_atomtype *atype, int *cgnr, int nrexcl);
+ PreprocessingAtomTypes *atype, int *cgnr, int nrexcl);
/* NOTE: nrexcl is not the size of *excl! */
void pdb2top(FILE *top_file, const char *posre_fn, const char *molname,
t_atoms *atoms,
std::vector<gmx::RVec> *x,
- gpp_atomtype *atype, t_symtab *tab,
+ PreprocessingAtomTypes *atype, t_symtab *tab,
gmx::ArrayRef<const PreprocessResidue> rtpFFDB,
gmx::ArrayRef<PreprocessResidue> usedPpResidues,
gmx::ArrayRef<MoleculePatchDatabase> globalPatches,
#include "hackblock.h"
-gpp_atomtype *read_atype(const char *ffdir, t_symtab *tab)
+PreprocessingAtomTypes read_atype(const char *ffdir, t_symtab *tab)
{
FILE *in;
char buf[STRLEN], name[STRLEN];
double m;
int nratt = 0;
- gpp_atomtype *at;
t_atom *a;
t_param *nb;
std::vector<std::string> files = fflib_search_file_end(ffdir, ".atp", TRUE);
- at = init_atomtype();
snew(a, 1);
snew(nb, 1);
+ PreprocessingAtomTypes at;
for (const auto &filename : files)
{
if (sscanf(buf, "%s%lf", name, &m) == 2)
{
a->m = m;
- add_atomtype(at, tab, a, name, nb, 0, 0);
+ at.addType(tab, *a, name, nb, 0, 0);
fprintf(stderr, "\rAtomtype %d", ++nratt);
fflush(stderr);
}
}
fprintf(stderr, "\n");
sfree(a);
- sfree(nb);
-
return at;
}
-static void print_resatoms(FILE *out, gpp_atomtype *atype, const PreprocessResidue &rtpDBEntry)
+static void print_resatoms(FILE *out, const PreprocessingAtomTypes &atype, const PreprocessResidue &rtpDBEntry)
{
fprintf(out, "[ %s ]\n", rtpDBEntry.resname.c_str());
fprintf(out, " [ atoms ]\n");
for (int j = 0; (j < rtpDBEntry.natom()); j++)
{
int tp = rtpDBEntry.atom[j].type;
- const char *tpnm = get_atomtype_name(tp, atype);
+ const char *tpnm = atype.atomNameFromAtomType(tp);
if (tpnm == nullptr)
{
gmx_fatal(FARGS, "Incorrect atomtype (%d)", tp);
}
static bool read_atoms(FILE *in, char *line,
- PreprocessResidue *r0, t_symtab *tab, gpp_atomtype *atype)
+ PreprocessResidue *r0, t_symtab *tab, PreprocessingAtomTypes *atype)
{
int cg;
char buf[256], buf1[256];
r0->atom.emplace_back();
r0->atom.back().q = q;
r0->cgnr.push_back(cg);
- int j = get_atomtype_type(buf1, atype);
+ int j = atype->atomTypeFromName(buf1);
if (j == NOTSET)
{
gmx_fatal(FARGS, "Atom type %s (residue %s) not found in atomtype "
"database", buf1, r0->resname.c_str());
}
r0->atom.back().type = j;
- r0->atom.back().m = get_atomtype_massA(j, atype);
+ r0->atom.back().m = atype->atomMassAFromAtomType(j);
}
return TRUE;
}
void print_resall(FILE *out, gmx::ArrayRef<const PreprocessResidue> rtpDBEntry,
- gpp_atomtype *atype)
+ const PreprocessingAtomTypes &atype)
{
if (rtpDBEntry.empty())
{
}
void readResidueDatabase(const std::string &rrdb, std::vector<PreprocessResidue> *rtpDBEntry,
- gpp_atomtype *atype, t_symtab *tab,
+ PreprocessingAtomTypes *atype, t_symtab *tab,
bool bAllowOverrideRTP)
{
FILE *in;
#include "gromacs/utility/arrayref.h"
-struct gpp_atomtype;
+class PreprocessingAtomTypes;
struct PreprocessResidue;
struct t_symtab;
* \param[in] tab Symbol table for names.
* \returns Atom type database.
*/
-gpp_atomtype *read_atype(const char *ffdir, t_symtab *tab);
+PreprocessingAtomTypes read_atype(const char *ffdir, t_symtab *tab);
/*! \brief
* Read in database, append to exisiting.
*
* \param[in] resdb Name of database file.
* \param[inout] rtpDBEntry Database to populate.
- * \param[in] atype Atomtype information.
+ * \param[inout] atype Atomtype information.
* \param[inout] tab Symbol table for names.
* \param[in] bAllowOverrideRTP If entries can be overwritten in the database.
*/
-void readResidueDatabase(const std::string &resdb,
- std::vector<PreprocessResidue> *rtpDBEntry,
- gpp_atomtype *atype,
- t_symtab *tab,
- bool bAllowOverrideRTP);
+void readResidueDatabase(const std::string &resdb,
+ std::vector<PreprocessResidue> *rtpDBEntry,
+ PreprocessingAtomTypes *atype,
+ t_symtab *tab,
+ bool bAllowOverrideRTP);
/*! \brief
* Print out database.
* \param[in] atype Atom type information.
*/
void print_resall(FILE *out, gmx::ArrayRef<const PreprocessResidue> rtpDBEntry,
- gpp_atomtype *atype);
+ const PreprocessingAtomTypes &atype);
#endif
#define FATAL() gmx_fatal(FARGS, "Reading Termini Database: not enough items on line\n%s", line)
static void read_atom(char *line, bool bAdd,
- std::string *nname, t_atom *a, gpp_atomtype *atype, int *cgnr)
+ std::string *nname, t_atom *a, PreprocessingAtomTypes *atype, int *cgnr)
{
int nr, i;
char buf[5][30];
*nname = "";
}
}
- a->type = get_atomtype_type(buf[i++], atype);
+ a->type = atype->atomTypeFromName(buf[i++]);
sscanf(buf[i++], "%lf", &m);
a->m = m;
sscanf(buf[i++], "%lf", &q);
}
}
-static void print_atom(FILE *out, const t_atom *a, gpp_atomtype *atype)
+static void print_atom(FILE *out, const t_atom &a, PreprocessingAtomTypes *atype)
{
fprintf(out, "\t%s\t%g\t%g\n",
- get_atomtype_name(a->type, atype), a->m, a->q);
+ atype->atomNameFromAtomType(a.type), a.m, a.q);
}
static void print_ter_db(const char *ff, char C, gmx::ArrayRef<const MoleculePatchDatabase> tb,
- gpp_atomtype *atype)
+ PreprocessingAtomTypes *atype)
{
std::string buf = gmx::formatString("%s-%c.tdb", ff, C);
FILE *out = gmx_fio_fopen(buf.c_str(), "w");
if (hack.type() == MoleculePatchType::Replace)
{
fprintf(out, "%s\t", hack.oname.c_str());
- print_atom(out, &hack.atom.back(), atype);
+ print_atom(out, hack.atom.back(), atype);
}
}
}
if (hack.type() == MoleculePatchType::Add)
{
print_ab(out, hack, hack.nname.c_str());
- print_atom(out, &hack.atom.back(), atype);
+ print_atom(out, hack.atom.back(), atype);
}
}
}
gmx_fio_fclose(out);
}
-static void read_ter_db_file(const char *fn,
- std::vector<MoleculePatchDatabase> *tbptr,
- gpp_atomtype *atype)
+static void read_ter_db_file(const char *fn,
+ std::vector<MoleculePatchDatabase> *tbptr,
+ PreprocessingAtomTypes *atype)
{
char filebase[STRLEN];
char header[STRLEN], buf[STRLEN], line[STRLEN];
}
int read_ter_db(const char *ffdir, char ter,
- std::vector<MoleculePatchDatabase> *tbptr, gpp_atomtype *atype)
+ std::vector<MoleculePatchDatabase> *tbptr, PreprocessingAtomTypes *atype)
{
std::string ext = gmx::formatString(".%c.tdb", ter);
#include "gromacs/utility/arrayref.h"
-struct gpp_atomtype;
+class PreprocessingAtomTypes;
struct MoleculePatchDatabase;
/*! \brief
* \returns Number of entries entered into database.
*/
int read_ter_db(const char *ffdir, char ter,
- std::vector<MoleculePatchDatabase> *tbptr, gpp_atomtype *atype);
+ std::vector<MoleculePatchDatabase> *tbptr, PreprocessingAtomTypes *atype);
/*! \brief
* Return entries for modification blocks that match a residue name.
return t2m;
}
-static int nequal(char *a1, char *a2)
+static int nequal(const char *a1, const char *a2)
{
int i;
return i;
}
-static real search_e_diss(int n2m, t_2morse t2m[], char *ai, char *aj)
+static real search_e_diss(int n2m, t_2morse t2m[], const char *ai, const char *aj)
{
int i;
int ibest = -1;
}
}
-void convert_harmonics(gmx::ArrayRef<MoleculeInformation> mols, gpp_atomtype *atype)
+void convert_harmonics(gmx::ArrayRef<MoleculeInformation> mols, PreprocessingAtomTypes *atype)
{
int n2m;
t_2morse *t2m;
int nj = mol.plist[bb].param[j].aj();
real edis =
search_e_diss(n2m, t2m,
- get_atomtype_name(mol.atoms.atom[ni].type, atype),
- get_atomtype_name(mol.atoms.atom[nj].type, atype));
+ atype->atomNameFromAtomType(mol.atoms.atom[ni].type),
+ atype->atomNameFromAtomType(mol.atoms.atom[nj].type));
if (edis != 0)
{
bRemoveHarm[j] = true;
#include "gromacs/utility/arrayref.h"
-struct gpp_atomtype;
+class PreprocessingAtomTypes;
struct MoleculeInformation;
-void convert_harmonics(gmx::ArrayRef<MoleculeInformation> mols, gpp_atomtype *atype);
+void convert_harmonics(gmx::ArrayRef<MoleculeInformation> mols, PreprocessingAtomTypes *atype);
#endif
#define OPENDIR '[' /* starting sign for directive */
#define CLOSEDIR ']' /* ending sign for directive */
-static void gen_pairs(InteractionTypeParameters *nbs, InteractionTypeParameters *pairs, real fudge, int comb)
+static void gen_pairs(const InteractionTypeParameters &nbs, InteractionTypeParameters *pairs, real fudge, int comb)
{
- int i, j, ntp, nrfp, nrfpA, nrfpB, nnn;
real scaling;
- ntp = nbs->nr;
- nnn = static_cast<int>(std::sqrt(static_cast<double>(ntp)));
+ int ntp = nbs.nr;
+ int nnn = static_cast<int>(std::sqrt(static_cast<double>(ntp)));
GMX_ASSERT(nnn * nnn == ntp, "Number of pairs of generated non-bonded parameters should be a perfect square");
- nrfp = NRFP(F_LJ);
- nrfpA = interaction_function[F_LJ14].nrfpA;
- nrfpB = interaction_function[F_LJ14].nrfpB;
+ int nrfp = NRFP(F_LJ);
+ int nrfpA = interaction_function[F_LJ14].nrfpA;
+ int nrfpB = interaction_function[F_LJ14].nrfpB;
pairs->nr = ntp;
if ((nrfp != nrfpA) || (nrfpA != nrfpB))
fprintf(stderr, "Generating 1-4 interactions: fudge = %g\n", fudge);
snew(pairs->param, pairs->nr);
- for (i = 0; (i < ntp); i++)
+ for (int i = 0; (i < ntp); i++)
{
/* Copy param.a */
pairs->param[i].a[0] = i / nnn;
- for (j = 0; (j < nrfp); j++)
+ for (int j = 0; (j < nrfp); j++)
{
/* If we are using sigma/epsilon values, only the epsilon values
* should be scaled, but not sigma.
scaling = fudge;
}
- pairs->param[i].c[j] = scaling*nbs->param[i].c[j];
+ pairs->param[i].c[j] = scaling*nbs.param[i].c[j];
/* NOTE: this should be clear to the compiler, but some gcc 5.2 versions
* issue false positive warnings for the pairs->param.c[] indexing below.
*/
assert(2*nrfp <= MAXFORCEPARAM);
- pairs->param[i].c[nrfp+j] = scaling*nbs->param[i].c[j];
+ pairs->param[i].c[nrfp+j] = scaling*nbs.param[i].c[j];
}
}
}
static char **read_topol(const char *infile, const char *outfile,
const char *define, const char *include,
t_symtab *symtab,
- gpp_atomtype *atype,
- int *nrmols,
+ PreprocessingAtomTypes *atypes,
std::vector<MoleculeInformation> *molinfo,
std::unique_ptr<MoleculeInformation> *intermolecular_interactions,
gmx::ArrayRef<InteractionTypeParameters> plist,
char line[STRLEN], errbuf[256], comb_str[256], nb_str[256];
char genpairs[32];
char *dirstr, *dummy2;
- int nrcopies, nmol, nscan, ncombs, ncopy;
+ int nrcopies, nscan, ncombs, ncopy;
double fLJ, fQQ, fPOW;
MoleculeInformation *mi0 = nullptr;
DirStack *DS;
/* some local variables */
DS_Init(&DS); /* directive stack */
- nmol = 0; /* no molecules yet... */
d = Directive::d_invalid; /* first thing should be a directive */
nbparam = nullptr; /* The temporary non-bonded matrix */
pair = nullptr; /* The temporary pair interaction matrix */
break;
case Directive::d_atomtypes:
- push_at(symtab, atype, batype, pline, nb_funct,
+ push_at(symtab, atypes, batype, pline, nb_funct,
&nbparam, bGenPairs ? &pair : nullptr, wi);
break;
case Directive::d_pairtypes:
if (bGenPairs)
{
- push_nbt(d, pair, atype, pline, F_LJ14, wi);
+ push_nbt(d, pair, atypes, pline, F_LJ14, wi);
}
else
{
- push_bt(d, plist, 2, atype, nullptr, pline, wi);
+ push_bt(d, plist, 2, atypes, nullptr, pline, wi);
}
break;
case Directive::d_angletypes:
break;
case Directive::d_nonbond_params:
- push_nbt(d, nbparam, atype, pline, nb_funct, wi);
+ push_nbt(d, nbparam, atypes, pline, nb_funct, wi);
break;
case Directive::d_implicit_genborn_params:
break;
case Directive::d_cmaptypes:
- push_cmaptype(d, plist, 5, atype, batype, pline, wi);
+ push_cmaptype(d, plist, 5, atypes, batype, pline, wi);
break;
case Directive::d_moleculetype:
opts->couple_lam1 == ecouplamNONE ||
opts->couple_lam1 == ecouplamQ))
{
- dcatt = add_atomtype_decoupled(symtab, atype,
+ dcatt = add_atomtype_decoupled(symtab, atypes,
&nbparam, bGenPairs ? &pair : nullptr);
}
- ntype = get_atomtype_ntypes(atype);
+ ntype = atypes->size();
ncombs = (ntype*(ntype+1))/2;
- generate_nbparams(*combination_rule, nb_funct, &(plist[nb_funct]), atype, wi);
+ generate_nbparams(*combination_rule, nb_funct, &(plist[nb_funct]), atypes, wi);
ncopy = copy_nbparams(nbparam, nb_funct, &(plist[nb_funct]),
ntype);
fprintf(stderr, "Generated %d of the %d non-bonded parameter combinations\n", ncombs-ncopy, ncombs);
free_nbparam(nbparam, ntype);
if (bGenPairs)
{
- gen_pairs(&(plist[nb_funct]), &(plist[F_LJ14]), fudgeLJ, *combination_rule);
+ gen_pairs((plist[nb_funct]), &(plist[F_LJ14]), fudgeLJ, *combination_rule);
ncopy = copy_nbparams(pair, nb_funct, &(plist[F_LJ14]),
ntype);
fprintf(stderr, "Generated %d of the %d 1-4 parameter combinations\n", ncombs-ncopy, ncombs);
}
push_molt(symtab, molinfo, pline, wi);
- nmol = molinfo->size();
exclusionBlocks.emplace_back();
mi0 = &molinfo->back();
mi0->atoms.haveMass = TRUE;
break;
}
case Directive::d_atoms:
- push_atom(symtab, &(mi0->cgs), &(mi0->atoms), atype, pline, &lastcg, wi);
+ push_atom(symtab, &(mi0->cgs), &(mi0->atoms), atypes, pline, &lastcg, wi);
break;
case Directive::d_pairs:
GMX_RELEASE_ASSERT(mi0, "Need to have a valid MoleculeInformation object to work on");
- push_bond(d, plist, mi0->plist, &(mi0->atoms), atype, pline, FALSE,
+ push_bond(d, plist, mi0->plist, &(mi0->atoms), atypes, pline, FALSE,
bGenPairs, *fudgeQQ, bZero, &bWarn_copy_A_B, wi);
break;
case Directive::d_pairs_nb:
GMX_RELEASE_ASSERT(mi0, "Need to have a valid MoleculeInformation object to work on");
- push_bond(d, plist, mi0->plist, &(mi0->atoms), atype, pline, FALSE,
+ push_bond(d, plist, mi0->plist, &(mi0->atoms), atypes, pline, FALSE,
FALSE, 1.0, bZero, &bWarn_copy_A_B, wi);
break;
case Directive::d_water_polarization:
case Directive::d_thole_polarization:
GMX_RELEASE_ASSERT(mi0, "Need to have a valid MoleculeInformation object to work on");
- push_bond(d, plist, mi0->plist, &(mi0->atoms), atype, pline, TRUE,
+ push_bond(d, plist, mi0->plist, &(mi0->atoms), atypes, pline, TRUE,
bGenPairs, *fudgeQQ, bZero, &bWarn_copy_A_B, wi);
break;
case Directive::d_cmap:
GMX_RELEASE_ASSERT(mi0, "Need to have a valid MoleculeInformation object to work on");
- push_cmap(d, plist, mi0->plist, &(mi0->atoms), atype, pline, wi);
+ push_cmap(d, plist, mi0->plist, &(mi0->atoms), atypes, pline, wi);
break;
case Directive::d_vsitesn:
sfree(intermolecular_interactions->get()->atoms.atom);
}
- *nrmols = nmol;
-
return title;
}
-char **do_top(bool bVerbose,
- const char *topfile,
- const char *topppfile,
- t_gromppopts *opts,
- bool bZero,
- t_symtab *symtab,
- gmx::ArrayRef<InteractionTypeParameters> plist,
- int *combination_rule,
- double *repulsion_power,
- real *fudgeQQ,
- gpp_atomtype *atype,
- int *nrmols,
- std::vector<MoleculeInformation> *molinfo,
- std::unique_ptr<MoleculeInformation> *intermolecular_interactions,
- const t_inputrec *ir,
- std::vector<gmx_molblock_t> *molblock,
- bool *ffParametrizedWithHBondConstraints,
- warninp *wi)
+char **do_top(bool bVerbose,
+ const char *topfile,
+ const char *topppfile,
+ t_gromppopts *opts,
+ bool bZero,
+ t_symtab *symtab,
+ gmx::ArrayRef<InteractionTypeParameters> plist,
+ int *combination_rule,
+ double *repulsion_power,
+ real *fudgeQQ,
+ PreprocessingAtomTypes *atypes,
+ std::vector<MoleculeInformation> *molinfo,
+ std::unique_ptr<MoleculeInformation> *intermolecular_interactions,
+ const t_inputrec *ir,
+ std::vector<gmx_molblock_t> *molblock,
+ bool *ffParametrizedWithHBondConstraints,
+ warninp *wi)
{
/* Tmpfile might contain a long path */
const char *tmpfile;
printf("processing topology...\n");
}
title = read_topol(topfile, tmpfile, opts->define, opts->include,
- symtab, atype,
- nrmols, molinfo, intermolecular_interactions,
+ symtab, atypes,
+ molinfo, intermolecular_interactions,
plist, combination_rule, repulsion_power,
opts, fudgeQQ, molblock,
ffParametrizedWithHBondConstraints,
struct gmx_molblock_t;
struct gmx_mtop_t;
-struct gpp_atomtype;
+class PreprocessingAtomTypes;
struct t_gromppopts;
struct t_inputrec;
struct MoleculeInformation;
double check_mol(const gmx_mtop_t *mtop, warninp_t wi);
/* Check mass and charge */
-char **do_top(bool bVerbose,
- const char *topfile,
- const char *topppfile,
- t_gromppopts *opts,
- bool bZero,
- t_symtab *symtab,
- gmx::ArrayRef<InteractionTypeParameters> plist,
- int *combination_rule,
- double *repulsion_power,
- real *fudgeQQ,
- gpp_atomtype *atype,
- int *nrmols,
- std::vector<MoleculeInformation> *molinfo,
- std::unique_ptr<MoleculeInformation> *intermolecular_interactions,
- const t_inputrec *ir,
- std::vector<gmx_molblock_t> *molblock,
- bool *ffParametrizedWithHBondConstraints,
- warninp_t wi);
+char **do_top(bool bVerbose,
+ const char *topfile,
+ const char *topppfile,
+ t_gromppopts *opts,
+ bool bZero,
+ t_symtab *symtab,
+ gmx::ArrayRef<InteractionTypeParameters> plist,
+ int *combination_rule,
+ double *repulsion_power,
+ real *fudgeQQ,
+ PreprocessingAtomTypes *atype,
+ std::vector<MoleculeInformation> *molinfo,
+ std::unique_ptr<MoleculeInformation> *intermolecular_interactions,
+ const t_inputrec *ir,
+ std::vector<gmx_molblock_t> *molblock,
+ bool *ffParametrizedWithHBondConstraints,
+ warninp_t wi);
/* This routine expects sys->molt[m].ilist to be of size F_NRE and ordered. */
void generate_qmexcl(gmx_mtop_t *sys, t_inputrec *ir, warninp_t wi, GmxQmmmMode qmmmMode);
#include "gromacs/utility/smalloc.h"
#include "gromacs/utility/stringutil.h"
-void generate_nbparams(int comb, int ftype, InteractionTypeParameters *plist, gpp_atomtype *atype,
- warninp *wi)
+void generate_nbparams(int comb,
+ int ftype,
+ InteractionTypeParameters *plist,
+ PreprocessingAtomTypes *atypes,
+ warninp *wi)
{
- int i, j, k = -1, nf;
+ int k = -1;
int nr, nrfp;
real c, bi, bj, ci, cj, ci0, ci1, ci2, cj0, cj1, cj2;
/* Lean mean shortcuts */
- nr = get_atomtype_ntypes(atype);
+ nr = atypes->size();
nrfp = NRFP(ftype);
snew(plist->param, nr*nr);
plist->nr = nr*nr;
{
case eCOMB_GEOMETRIC:
/* Gromos rules */
- for (i = k = 0; (i < nr); i++)
+ for (int i = k = 0; (i < nr); i++)
{
- for (j = 0; (j < nr); j++, k++)
+ for (int j = 0; (j < nr); j++, k++)
{
- for (nf = 0; (nf < nrfp); nf++)
+ for (int nf = 0; (nf < nrfp); nf++)
{
- ci = get_atomtype_nbparam(i, nf, atype);
- cj = get_atomtype_nbparam(j, nf, atype);
+ ci = atypes->atomNonBondedParamFromAtomType(i, nf);
+ cj = atypes->atomNonBondedParamFromAtomType(j, nf);
c = std::sqrt(ci * cj);
plist->param[k].c[nf] = c;
}
case eCOMB_ARITHMETIC:
/* c0 and c1 are sigma and epsilon */
- for (i = k = 0; (i < nr); i++)
+ for (int i = k = 0; (i < nr); i++)
{
- for (j = 0; (j < nr); j++, k++)
+ for (int j = 0; (j < nr); j++, k++)
{
- ci0 = get_atomtype_nbparam(i, 0, atype);
- cj0 = get_atomtype_nbparam(j, 0, atype);
- ci1 = get_atomtype_nbparam(i, 1, atype);
- cj1 = get_atomtype_nbparam(j, 1, atype);
+ ci0 = atypes->atomNonBondedParamFromAtomType(i, 0);
+ cj0 = atypes->atomNonBondedParamFromAtomType(j, 0);
+ ci1 = atypes->atomNonBondedParamFromAtomType(i, 1);
+ cj1 = atypes->atomNonBondedParamFromAtomType(j, 1);
plist->param[k].c[0] = (fabs(ci0) + fabs(cj0))*0.5;
/* Negative sigma signals that c6 should be set to zero later,
* so we need to propagate that through the combination rules.
break;
case eCOMB_GEOM_SIG_EPS:
/* c0 and c1 are sigma and epsilon */
- for (i = k = 0; (i < nr); i++)
+ for (int i = k = 0; (i < nr); i++)
{
- for (j = 0; (j < nr); j++, k++)
+ for (int j = 0; (j < nr); j++, k++)
{
- ci0 = get_atomtype_nbparam(i, 0, atype);
- cj0 = get_atomtype_nbparam(j, 0, atype);
- ci1 = get_atomtype_nbparam(i, 1, atype);
- cj1 = get_atomtype_nbparam(j, 1, atype);
+ ci0 = atypes->atomNonBondedParamFromAtomType(i, 0);
+ cj0 = atypes->atomNonBondedParamFromAtomType(j, 0);
+ ci1 = atypes->atomNonBondedParamFromAtomType(i, 1);
+ cj1 = atypes->atomNonBondedParamFromAtomType(j, 1);
plist->param[k].c[0] = std::sqrt(std::fabs(ci0*cj0));
/* Negative sigma signals that c6 should be set to zero later,
* so we need to propagate that through the combination rules.
case F_BHAM:
/* Buckingham rules */
- for (i = k = 0; (i < nr); i++)
+ for (int i = k = 0; (i < nr); i++)
{
- for (j = 0; (j < nr); j++, k++)
+ for (int j = 0; (j < nr); j++, k++)
{
- ci0 = get_atomtype_nbparam(i, 0, atype);
- cj0 = get_atomtype_nbparam(j, 0, atype);
- ci2 = get_atomtype_nbparam(i, 2, atype);
- cj2 = get_atomtype_nbparam(j, 2, atype);
- bi = get_atomtype_nbparam(i, 1, atype);
- bj = get_atomtype_nbparam(j, 1, atype);
+ ci0 = atypes->atomNonBondedParamFromAtomType(i, 0);
+ cj0 = atypes->atomNonBondedParamFromAtomType(j, 0);
+ ci2 = atypes->atomNonBondedParamFromAtomType(i, 2);
+ cj2 = atypes->atomNonBondedParamFromAtomType(j, 2);
+ bi = atypes->atomNonBondedParamFromAtomType(i, 1);
+ bj = atypes->atomNonBondedParamFromAtomType(j, 1);
plist->param[k].c[0] = std::sqrt(ci0 * cj0);
if ((bi == 0) || (bj == 0))
{
real c[4];
};
-static void realloc_nb_params(gpp_atomtype *at,
+static void realloc_nb_params(PreprocessingAtomTypes *atypes,
t_nbparam ***nbparam, t_nbparam ***pair)
{
/* Add space in the non-bonded parameters matrix */
- int atnr = get_atomtype_ntypes(at);
+ int atnr = atypes->size();
srenew(*nbparam, atnr);
snew((*nbparam)[atnr-1], atnr);
if (pair)
}
}
-void push_at (t_symtab *symtab, gpp_atomtype *at, gpp_bond_atomtype *bat,
+void push_at (t_symtab *symtab, PreprocessingAtomTypes *at, gpp_bond_atomtype *bat,
char *line, int nb_funct,
t_nbparam ***nbparam, t_nbparam ***pair,
warninp *wi)
auto message = gmx::formatString("Invalid function type %d in push_at", nb_funct);
warning_error_and_exit(wi, message, FARGS);
}
- for (j = nfp0; (j < MAXFORCEPARAM); j++)
+ for (int j = nfp0; (j < MAXFORCEPARAM); j++)
{
c[j] = 0.0;
}
}
batype_nr = get_bond_atomtype_type(btype, bat);
- if ((nr = get_atomtype_type(type, at)) != NOTSET)
+ if ((nr = at->atomTypeFromName(type)) != NOTSET)
{
auto message = gmx::formatString("Overriding atomtype %s", type);
warning(wi, message);
- if ((nr = set_atomtype(nr, at, symtab, atom, type, param, batype_nr,
- atomnr)) == NOTSET)
+ if ((nr = at->setType(nr, symtab, *atom, type, param, batype_nr,
+ atomnr)) == NOTSET)
{
auto message = gmx::formatString("Replacing atomtype %s failed", type);
warning_error_and_exit(wi, message, FARGS);
}
}
- else if ((add_atomtype(at, symtab, atom, type, param,
- batype_nr, atomnr)) == NOTSET)
+ else if ((at->addType(symtab, *atom, type, param,
+ batype_nr, atomnr)) == NOTSET)
{
auto message = gmx::formatString("Adding atomtype %s failed", type);
warning_error_and_exit(wi, message, FARGS);
}
}
-void push_bt(Directive d, gmx::ArrayRef<InteractionTypeParameters> bt, int nral,
- gpp_atomtype *at,
- gpp_bond_atomtype *bat, char *line,
- warninp *wi)
+void push_bt(Directive d,
+ gmx::ArrayRef<InteractionTypeParameters> bt,
+ int nral,
+ PreprocessingAtomTypes *at,
+ gpp_bond_atomtype *bat,
+ char *line,
+ warninp *wi)
{
const char *formal[MAXATOMLIST+1] = {
"%s",
}
for (i = 0; (i < nral); i++)
{
- if (at && ((p.a[i] = get_atomtype_type(alc[i], at)) == NOTSET))
+ if ((at != nullptr) && ((p.a[i] = at->atomTypeFromName(alc[i])) == NOTSET))
{
auto message = gmx::formatString("Unknown atomtype %s\n", alc[i]);
warning_error_and_exit(wi, message, FARGS);
}
-void push_nbt(Directive d, t_nbparam **nbt, gpp_atomtype *atype,
+void push_nbt(Directive d, t_nbparam **nbt, PreprocessingAtomTypes *atypes,
char *pline, int nb_funct,
warninp *wi)
{
}
/* Put the parameters in the matrix */
- if ((ai = get_atomtype_type (a0, atype)) == NOTSET)
+ if ((ai = atypes->atomTypeFromName(a0)) == NOTSET)
{
auto message = gmx::formatString("Atomtype %s not found", a0);
warning_error_and_exit(wi, message, FARGS);
}
- if ((aj = get_atomtype_type (a1, atype)) == NOTSET)
+ if ((aj = atypes->atomTypeFromName(a1)) == NOTSET)
{
auto message = gmx::formatString("Atomtype %s not found", a1);
warning_error_and_exit(wi, message, FARGS);
}
void
-push_cmaptype(Directive d, gmx::ArrayRef<InteractionTypeParameters> bt, int nral, gpp_atomtype *at,
- gpp_bond_atomtype *bat, char *line,
- warninp *wi)
+push_cmaptype(Directive d,
+ gmx::ArrayRef<InteractionTypeParameters> bt,
+ int nral,
+ PreprocessingAtomTypes *at,
+ gpp_bond_atomtype *bat,
+ char *line,
+ warninp *wi)
{
const char *formal = "%s%s%s%s%s%s%s%s%n";
- int i, ft, ftype, nn, nrfp, nrfpA, nrfpB;
+ int ft, ftype, nn, nrfp, nrfpA, nrfpB;
int start, nchar_consumed;
int nxcmap, nycmap, ncmap, read_cmap, sl, nct;
char s[20], alc[MAXATOMLIST+2][20];
/* Read in CMAP parameters */
sl = 0;
- for (i = 0; i < ncmap; i++)
+ for (int i = 0; i < ncmap; i++)
{
while (isspace(*(line+start+sl)))
{
/* Check do that we got the number of parameters we expected */
if (read_cmap == nrfpA)
{
- for (i = 0; i < ncmap; i++)
+ for (int i = 0; i < ncmap; i++)
{
bt[F_CMAP].cmap.emplace_back(bt[F_CMAP].cmap[i]);
}
bt[F_CMAP].cmapAngles++; /* Since we are incrementing here, we need to subtract later, see (*****) */
nct = (nral+1) * bt[F_CMAP].cmapAngles;
- for (i = 0; (i < nral); i++)
+ for (int i = 0; (i < nral); i++)
{
if (at && ((p.a[i] = get_bond_atomtype_type(alc[i], bat)) == NOTSET))
{
}
/* Is this correct?? */
- for (i = 0; i < MAXFORCEPARAM; i++)
+ for (int i = 0; i < MAXFORCEPARAM; i++)
{
p.c[i] = NOTSET;
}
}
void push_atom(t_symtab *symtab, t_block *cgs,
- t_atoms *at, gpp_atomtype *atype, char *line, int *lastcg,
+ t_atoms *at, PreprocessingAtomTypes *atypes, char *line, int *lastcg,
warninp *wi)
{
int nr, ptype;
return;
}
sscanf(id, "%d", &atomnr);
- if ((type = get_atomtype_type(ctype, atype)) == NOTSET)
+ if ((type = atypes->atomTypeFromName(ctype)) == NOTSET)
{
auto message = gmx::formatString("Atomtype %s not found", ctype);
warning_error_and_exit(wi, message, FARGS);
}
- ptype = get_atomtype_ptype(type, atype);
+ ptype = atypes->atomParticleTypeFromAtomType(type);
/* Set default from type */
- q0 = get_atomtype_qA(type, atype);
- m0 = get_atomtype_massA(type, atype);
+ q0 = atypes->atomChargeAFromAtomType(type);
+ m0 = atypes->atomMassAFromAtomType(type);
typeB = type;
qB = q0;
mB = m0;
m0 = mB = m;
if (nscan > 2)
{
- if ((typeB = get_atomtype_type(ctypeB, atype)) == NOTSET)
+ if ((typeB = atypes->atomTypeFromName(ctypeB)) == NOTSET)
{
auto message = gmx::formatString("Atomtype %s not found", ctypeB);
warning_error_and_exit(wi, message, FARGS);
}
- qB = get_atomtype_qA(typeB, atype);
- mB = get_atomtype_massA(typeB, atype);
+ qB = atypes->atomChargeAFromAtomType(typeB);
+ mB = atypes->atomMassAFromAtomType(typeB);
if (nscan > 3)
{
qB = qb;
push_cg(cgs, lastcg, cgnumber, nr);
- push_atom_now(symtab, at, atomnr, get_atomtype_atomnumber(type, atype),
+ push_atom_now(symtab, at, atomnr, atypes->atomNumberFromAtomType(type),
type, ctype, ptype, resnumberic,
resname, name, m0, q0, typeB,
typeB == type ? ctype : ctypeB, mB, qB, wi);
}
static bool default_cmap_params(gmx::ArrayRef<InteractionTypeParameters> bondtype,
- t_atoms *at, gpp_atomtype *atype,
+ t_atoms *at, PreprocessingAtomTypes *atypes,
t_param *p, bool bB,
int *cmap_type, int *nparam_def,
warninp *wi)
{
- int i, nparam_found;
+ int nparam_found;
int ct;
- bool bFound = FALSE;
+ bool bFound = false;
nparam_found = 0;
ct = 0;
/* Match the current cmap angle against the list of cmap_types */
- for (i = 0; i < bondtype[F_CMAP].nct() && !bFound; i += 6)
+ for (int i = 0; i < bondtype[F_CMAP].nct() && !bFound; i += 6)
{
if (bB)
{
else
{
if (
- (get_atomtype_batype(at->atom[p->a[0]].type, atype) == bondtype[F_CMAP].cmapAtomTypes[i]) &&
- (get_atomtype_batype(at->atom[p->a[1]].type, atype) == bondtype[F_CMAP].cmapAtomTypes[i+1]) &&
- (get_atomtype_batype(at->atom[p->a[2]].type, atype) == bondtype[F_CMAP].cmapAtomTypes[i+2]) &&
- (get_atomtype_batype(at->atom[p->a[3]].type, atype) == bondtype[F_CMAP].cmapAtomTypes[i+3]) &&
- (get_atomtype_batype(at->atom[p->a[4]].type, atype) == bondtype[F_CMAP].cmapAtomTypes[i+4]))
+ (atypes->bondAtomTypeFromAtomType(at->atom[p->a[0]].type) == bondtype[F_CMAP].cmapAtomTypes[i]) &&
+ (atypes->bondAtomTypeFromAtomType(at->atom[p->a[1]].type) == bondtype[F_CMAP].cmapAtomTypes[i+1]) &&
+ (atypes->bondAtomTypeFromAtomType(at->atom[p->a[2]].type) == bondtype[F_CMAP].cmapAtomTypes[i+2]) &&
+ (atypes->bondAtomTypeFromAtomType(at->atom[p->a[3]].type) == bondtype[F_CMAP].cmapAtomTypes[i+3]) &&
+ (atypes->bondAtomTypeFromAtomType(at->atom[p->a[4]].type) == bondtype[F_CMAP].cmapAtomTypes[i+4]))
{
/* Found cmap torsion */
- bFound = TRUE;
+ bFound = true;
ct = bondtype[F_CMAP].cmapAtomTypes[i+5];
nparam_found = 1;
}
*/
static int natom_match(t_param *pi,
int type_i, int type_j, int type_k, int type_l,
- const gpp_atomtype* atype)
+ const PreprocessingAtomTypes* atypes)
{
- if ((pi->ai() == -1 || get_atomtype_batype(type_i, atype) == pi->ai()) &&
- (pi->aj() == -1 || get_atomtype_batype(type_j, atype) == pi->aj()) &&
- (pi->ak() == -1 || get_atomtype_batype(type_k, atype) == pi->ak()) &&
- (pi->al() == -1 || get_atomtype_batype(type_l, atype) == pi->al()))
+ if ((pi->ai() == -1 || atypes->bondAtomTypeFromAtomType(type_i) == pi->ai()) &&
+ (pi->aj() == -1 || atypes->bondAtomTypeFromAtomType(type_j) == pi->aj()) &&
+ (pi->ak() == -1 || atypes->bondAtomTypeFromAtomType(type_k) == pi->ak()) &&
+ (pi->al() == -1 || atypes->bondAtomTypeFromAtomType(type_l) == pi->al()))
{
return
(pi->ai() == -1 ? 0 : 1) +
}
static bool defaultInteractionTypeParameters(int ftype, gmx::ArrayRef<InteractionTypeParameters> bt,
- t_atoms *at, gpp_atomtype *atype,
+ t_atoms *at, PreprocessingAtomTypes *atypes,
t_param *p, bool bB,
t_param **param_def,
int *nparam_def)
pt = &(bt[ftype].param[t]);
if (bB)
{
- nmatch = natom_match(pt, at->atom[p->ai()].typeB, at->atom[p->aj()].typeB, at->atom[p->ak()].typeB, at->atom[p->al()].typeB, atype);
+ nmatch = natom_match(pt, at->atom[p->ai()].typeB, at->atom[p->aj()].typeB, at->atom[p->ak()].typeB, at->atom[p->al()].typeB, atypes);
}
else
{
- nmatch = natom_match(pt, at->atom[p->ai()].type, at->atom[p->aj()].type, at->atom[p->ak()].type, at->atom[p->al()].type, atype);
+ nmatch = natom_match(pt, at->atom[p->ai()].type, at->atom[p->aj()].type, at->atom[p->ak()].type, at->atom[p->al()].type, atypes);
}
if (nmatch > nmatch_max)
{
if (bB)
{
for (j = 0; ((j < nral) &&
- (get_atomtype_batype(at->atom[p->a[j]].typeB, atype) == pi->a[j])); j++)
+ (atypes->bondAtomTypeFromAtomType(at->atom[p->a[j]].typeB) == pi->a[j])); j++)
{
;
}
else
{
for (j = 0; ((j < nral) &&
- (get_atomtype_batype(at->atom[p->a[j]].type, atype) == pi->a[j])); j++)
+ (atypes->bondAtomTypeFromAtomType(at->atom[p->a[j]].type) == pi->a[j])); j++)
{
;
}
void push_bond(Directive d, gmx::ArrayRef<InteractionTypeParameters> bondtype,
gmx::ArrayRef<InteractionTypeParameters> bond,
- t_atoms *at, gpp_atomtype *atype, char *line,
+ t_atoms *at, PreprocessingAtomTypes *atypes, char *line,
bool bBonded, bool bGenPairs, real fudgeQQ,
bool bZero, bool *bWarn_copy_A_B,
warninp *wi)
if (bBonded)
{
- bFoundA = defaultInteractionTypeParameters(ftype, bondtype, at, atype, ¶m, FALSE, ¶m_defA, &nparam_defA);
+ bFoundA = defaultInteractionTypeParameters(ftype, bondtype, at, atypes, ¶m, FALSE, ¶m_defA, &nparam_defA);
if (bFoundA)
{
/* Copy the A-state and B-state default parameters. */
param.c[j] = param_defA->c[j];
}
}
- bFoundB = defaultInteractionTypeParameters(ftype, bondtype, at, atype, ¶m, TRUE, ¶m_defB, &nparam_defB);
+ bFoundB = defaultInteractionTypeParameters(ftype, bondtype, at, atypes, ¶m, TRUE, ¶m_defB, &nparam_defB);
if (bFoundB)
{
/* Copy only the B-state default parameters */
void push_cmap(Directive d, gmx::ArrayRef<InteractionTypeParameters> bondtype,
gmx::ArrayRef<InteractionTypeParameters> bond,
- t_atoms *at, gpp_atomtype *atype, char *line,
+ t_atoms *at, PreprocessingAtomTypes *atypes, char *line,
warninp *wi)
{
const char *aaformat[MAXATOMLIST+1] =
}
/* Get the cmap type for this cmap angle */
- bFound = default_cmap_params(bondtype, at, atype, ¶m, FALSE, &cmap_type, &ncmap_params, wi);
+ bFound = default_cmap_params(bondtype, at, atypes, ¶m, FALSE, &cmap_type, &ncmap_params, wi);
/* We want exactly one parameter (the cmap type in state A (currently no state B) back */
if (bFound && ncmap_params == 1)
while (n == 1);
}
-int add_atomtype_decoupled(t_symtab *symtab, gpp_atomtype *at,
+int add_atomtype_decoupled(t_symtab *symtab, PreprocessingAtomTypes *at,
t_nbparam ***nbparam, t_nbparam ***pair)
{
t_atom atom;
param.c[i] = 0.0;
}
- nr = add_atomtype(at, symtab, &atom, "decoupled", ¶m, -1, 0);
+ nr = at->addType(symtab, atom, "decoupled", ¶m, -1, 0);
/* Add space in the non-bonded parameters matrix */
realloc_nb_params(at, nbparam, pair);
#include "gromacs/utility/real.h"
enum class Directive : int;
-struct gpp_atomtype;
+class PreprocessingAtomTypes;
struct gpp_bond_atomtype;
struct t_atoms;
struct t_block;
} // namespace gmx
void generate_nbparams(int comb, int funct, InteractionTypeParameters *plist,
- gpp_atomtype *atype,
+ PreprocessingAtomTypes *atype,
warninp *wi);
-void push_at (struct t_symtab *symtab, gpp_atomtype *at,
+void push_at (struct t_symtab *symtab, PreprocessingAtomTypes *at,
gpp_bond_atomtype *bat, char *line, int nb_funct,
t_nbparam ***nbparam, t_nbparam ***pair,
warninp *wi);
void push_bt(Directive d, gmx::ArrayRef<InteractionTypeParameters> bt, int nral,
- gpp_atomtype *at, gpp_bond_atomtype *bat, char *line,
+ PreprocessingAtomTypes *at, gpp_bond_atomtype *bat, char *line,
warninp *wi);
void push_dihedraltype(Directive d, gmx::ArrayRef<InteractionTypeParameters> bt,
gpp_bond_atomtype *bat, char *line,
warninp *wi);
-void push_cmaptype(Directive d, gmx::ArrayRef<InteractionTypeParameters> bt, int nral, gpp_atomtype *at,
+void push_cmaptype(Directive d, gmx::ArrayRef<InteractionTypeParameters> bt, int nral, PreprocessingAtomTypes *at,
gpp_bond_atomtype *bat, char *line,
warninp *wi);
-void push_nbt(Directive d, t_nbparam **nbt, gpp_atomtype *atype,
+void push_nbt(Directive d, t_nbparam **nbt, PreprocessingAtomTypes *atype,
char *plines, int nb_funct,
warninp *wi);
-void push_atom(struct t_symtab *symtab,
- t_block *cgs,
- t_atoms *at,
- gpp_atomtype *atype,
- char *line,
- int *lastcg,
- warninp *wi);
+void push_atom(struct t_symtab *symtab,
+ t_block *cgs,
+ t_atoms *at,
+ PreprocessingAtomTypes *atype,
+ char *line,
+ int *lastcg,
+ warninp *wi);
void push_bond(Directive d, gmx::ArrayRef<InteractionTypeParameters> bondtype,
gmx::ArrayRef<InteractionTypeParameters> bond,
- t_atoms *at, gpp_atomtype *atype, char *line,
+ t_atoms *at, PreprocessingAtomTypes *atype, char *line,
bool bBonded, bool bGenPairs, real fudgeQQ,
bool bZero, bool *bWarn_copy_A_B,
warninp *wi);
void push_cmap(Directive d,
gmx::ArrayRef<InteractionTypeParameters> bondtype,
gmx::ArrayRef<InteractionTypeParameters> bond,
- t_atoms *at, gpp_atomtype *atype, char *line,
+ t_atoms *at, PreprocessingAtomTypes *atype, char *line,
warninp *wi);
void push_vsitesn(Directive d, gmx::ArrayRef<InteractionTypeParameters> bond,
void free_nbparam(t_nbparam **param, int nr);
-int add_atomtype_decoupled(struct t_symtab *symtab, gpp_atomtype *at,
+int add_atomtype_decoupled(struct t_symtab *symtab, PreprocessingAtomTypes *at,
t_nbparam ***nbparam, t_nbparam ***pair);
/* Add an atom type with all parameters set to zero (no interactions).
* Returns the atom type number.
static int count_hydrogens (char ***atomname, int nra, const int a[])
{
- int i, nh;
+ int nh;
if (!atomname)
{
}
nh = 0;
- for (i = 0; (i < nra); i++)
+ for (int i = 0; (i < nra); i++)
{
if (toupper(**(atomname[a[i]])) == 'H')
{
void make_shake(gmx::ArrayRef<InteractionTypeParameters> plist, t_atoms *atoms, int nshake)
{
- char ***info = atoms->atomname;
- InteractionTypeParameters *pr;
- InteractionTypeParameters *bonds;
- t_param p, *bond, *ang;
- real b_ij, b_jk;
- int i, j, ftype, ftype_a;
- bool bFound;
+ char ***info = atoms->atomname;
+ real b_ij, b_jk;
+ int i, j;
if (nshake != eshNONE)
{
/* Add all the angles with hydrogens to the shake list
* and remove them from the bond list
*/
- for (ftype = 0; (ftype < F_NRE); ftype++)
+ for (int ftype = 0; (ftype < F_NRE); ftype++)
{
if (interaction_function[ftype].flags & IF_BTYPE)
{
- bonds = &(plist[ftype]);
+ InteractionTypeParameters *bonds = &(plist[ftype]);
- for (ftype_a = 0; (bonds->nr > 0 && ftype_a < F_NRE); ftype_a++)
+ for (int ftype_a = 0; (bonds->nr > 0 && ftype_a < F_NRE); ftype_a++)
{
if (interaction_function[ftype_a].flags & IF_ATYPE)
{
- pr = &(plist[ftype_a]);
+ InteractionTypeParameters *pr = &(plist[ftype_a]);
- for (i = 0; (i < pr->nr); )
+ for (int i = 0; (i < pr->nr); )
{
- int numhydrogens;
-
- ang = &(pr->param[i]);
+ t_param *ang = &(pr->param[i]);
#ifdef DEBUG
printf("Angle: %d-%d-%d\n", ang->ai(), ang->aj(), ang->ak());
#endif
- numhydrogens = count_hydrogens(info, 3, ang->a);
+ int numhydrogens = count_hydrogens(info, 3, ang->a);
if ((nshake == eshALLANGLES) ||
(numhydrogens > 1) ||
(numhydrogens == 1 && toupper(**(info[ang->a[1]])) == 'O'))
* are constrained.
* append this angle to the shake list
*/
+ t_param p;
p.ai() = ang->ai();
p.aj() = ang->ak();
/* Calculate length of constraint */
- bFound = FALSE;
+ bool bFound = false;
b_ij = b_jk = 0.0;
- for (j = 0; !bFound && (j < bonds->nr); j++)
+ for (int j = 0; !bFound && (j < bonds->nr); j++)
{
- bond = &(bonds->param[j]);
+ t_param *bond = &(bonds->param[j]);
if (((bond->ai() == ang->ai()) &&
(bond->aj() == ang->aj())) ||
((bond->ai() == ang->aj()) &&
/* Add all the bonds with hydrogens to the shake list
* and remove them from the bond list
*/
- for (ftype = 0; (ftype < F_NRE); ftype++)
+ for (int ftype = 0; (ftype < F_NRE); ftype++)
{
if (interaction_function[ftype].flags & IF_BTYPE)
{
- pr = &(plist[ftype]);
+ InteractionTypeParameters *pr = &(plist[ftype]);
j = 0;
for (i = 0; i < pr->nr; i++)
{
(count_hydrogens (info, 2, pr->param[i].a) > 0) )
{
/* append this bond to the shake list */
+ t_param p;
p.ai() = pr->param[i].ai();
p.aj() = pr->param[i].aj();
p.c0() = pr->param[i].c0();
}
}
-void cp_param(t_param *dest, t_param *src)
+void cp_param(t_param *dest, const t_param *src)
{
int j;
void add_param_to_list(InteractionTypeParameters *list, t_param *b)
{
- int j;
-
/* allocate one position extra */
pr_alloc (1, list);
/* fill the arrays */
- for (j = 0; (j < MAXFORCEPARAM); j++)
+ for (int j = 0; (j < MAXFORCEPARAM); j++)
{
list->param[list->nr].c[j] = b->c[j];
}
- for (j = 0; (j < MAXATOMLIST); j++)
+ for (int j = 0; (j < MAXATOMLIST); j++)
{
list->param[list->nr].a[j] = b->a[j];
}
/* PRINTING STRUCTURES */
-static void print_bt(FILE *out, Directive d, gpp_atomtype *at,
+static void print_bt(FILE *out, Directive d, PreprocessingAtomTypes *at,
int ftype, int fsubtype, gmx::ArrayRef<const InteractionTypeParameters> plist,
bool bFullDih)
{
* all four atoms to determine the type.
*/
const int dihp[2][2] = { { 1, 2 }, { 0, 3 } };
- int i, j, f, nral, nrfp;
- bool bDih = FALSE, bSwapParity;
+ int nral, nrfp;
+ bool bDih = false, bSwapParity;
const InteractionTypeParameters *bt = &(plist[ftype]);
- if (!bt->nr)
+ if (bt->nr == 0)
{
return;
}
- f = 0;
+ int f = 0;
switch (ftype)
{
case F_G96ANGLES:
if (!bDih)
{
fprintf (out, "%3s %4s", "ai", "aj");
- for (j = 2; (j < nral); j++)
+ for (int j = 2; (j < nral); j++)
{
fprintf (out, " %3c%c", 'a', 'i'+j);
}
}
else
{
- for (j = 0; (j < 2); j++)
+ for (int j = 0; (j < 2); j++)
{
fprintf (out, "%3c%c", 'a', 'i'+dihp[f][j]);
}
}
fprintf (out, " funct");
- for (j = 0; (j < nrfp); j++)
+ for (int j = 0; (j < nrfp); j++)
{
fprintf (out, " %12c%1d", 'c', j);
}
fprintf (out, "\n");
/* print bondtypes */
- for (i = 0; (i < bt->nr); i++)
+ for (int i = 0; (i < bt->nr); i++)
{
bSwapParity = (bt->param[i].c0() == NOTSET) && (bt->param[i].c1() == -1);
if (!bDih)
{
- for (j = 0; (j < nral); j++)
+ for (int j = 0; (j < nral); j++)
{
- fprintf (out, "%5s ", get_atomtype_name(bt->param[i].a[j], at));
+ fprintf (out, "%5s ", at->atomNameFromAtomType(bt->param[i].a[j]));
}
}
else
{
- for (j = 0; (j < 2); j++)
+ for (int j = 0; (j < 2); j++)
{
- fprintf (out, "%5s ", get_atomtype_name(bt->param[i].a[dihp[f][j]], at));
+ fprintf (out, "%5s ", at->atomNameFromAtomType(bt->param[i].a[dihp[f][j]]));
}
}
fprintf (out, "%5d ", bSwapParity ? -f-1 : f+1);
}
else
{
- for (j = 0; (j < nrfp && (bt->param[i].c[j] != NOTSET)); j++)
+ for (int j = 0; (j < nrfp && (bt->param[i].c[j] != NOTSET)); j++)
{
fprintf (out, "%13.6e ", bt->param[i].c[j]);
}
return q;
}
-void print_atoms(FILE *out, gpp_atomtype *atype, t_atoms *at, int *cgnr,
+void print_atoms(FILE *out, PreprocessingAtomTypes *atype, t_atoms *at, int *cgnr,
bool bRTPresname)
{
- int i, ri;
- int tpA, tpB;
- const char *as;
- char *tpnmA, *tpnmB;
- double qres, qtot;
+ int i, ri;
+ int tpA, tpB;
+ const char *as;
+ const char *tpnmA, *tpnmB;
+ double qres, qtot;
as = dir2str(Directive::d_atoms);
fprintf(out, "[ %s ]\n", as);
fprintf(out, "\n");
}
tpA = at->atom[i].type;
- if ((tpnmA = get_atomtype_name(tpA, atype)) == nullptr)
+ if ((tpnmA = atype->atomNameFromAtomType(tpA)) == nullptr)
{
gmx_fatal(FARGS, "tpA = %d, i= %d in print_atoms", tpA, i);
}
if (PERTURBED(at->atom[i]))
{
tpB = at->atom[i].typeB;
- if ((tpnmB = get_atomtype_name(tpB, atype)) == nullptr)
+ if ((tpnmB = atype->atomNameFromAtomType(tpB)) == nullptr)
{
gmx_fatal(FARGS, "tpB = %d, i= %d in print_atoms", tpB, i);
}
void print_bondeds(FILE *out, int natoms, Directive d,
int ftype, int fsubtype, gmx::ArrayRef<const InteractionTypeParameters> plist)
{
- t_symtab stab;
- gpp_atomtype *atype;
- t_param *param;
- t_atom *a;
- int i;
+ t_symtab stab;
+ t_param *param;
+ t_atom *a;
- atype = init_atomtype();
+ PreprocessingAtomTypes atype;
snew(a, 1);
snew(param, 1);
open_symtab(&stab);
- for (i = 0; (i < natoms); i++)
+ for (int i = 0; (i < natoms); i++)
{
char buf[12];
sprintf(buf, "%4d", (i+1));
- add_atomtype(atype, &stab, a, buf, param, 0, 0);
+ atype.addType(&stab, *a, buf, param, 0, 0);
}
- print_bt(out, d, atype, ftype, fsubtype, plist, TRUE);
+ print_bt(out, d, &atype, ftype, fsubtype, plist, TRUE);
done_symtab(&stab);
sfree(a);
sfree(param);
- done_atomtype(atype);
}
#include "gromacs/utility/arrayref.h"
+#include "gromacs/utility/arrayref.h"
+
enum class Directive : int;
-struct gpp_atomtype;
+class PreprocessingAtomTypes;
struct t_atoms;
struct t_blocka;
struct t_excls;
void set_p_string(t_param *p, const char *s);
-void cp_param(t_param *dest, t_param *src);
+void cp_param(t_param *dest, const t_param *src);
void add_param_to_list(InteractionTypeParameters *list, t_param *b);
void print_blocka(FILE *out, const char *szName, const char *szIndex,
const char *szA, t_blocka *block);
-void print_atoms(FILE *out, gpp_atomtype *atype, t_atoms *at, int *cgnr,
+void print_atoms(FILE *out, PreprocessingAtomTypes *atype, t_atoms *at, int *cgnr,
bool bRTPresname);
void print_bondeds(FILE *out, int natoms, Directive d,
return angle;
}
-static char *get_atomtype_name_AB(t_atom *atom, gpp_atomtype *atype)
+static const char *get_atomtype_name_AB(t_atom *atom, PreprocessingAtomTypes *atypes)
{
- char *name;
-
- name = get_atomtype_name(atom->type, atype);
+ const char* name = atypes->atomNameFromAtomType(atom->type);
/* When using the decoupling option, atom types are changed
* to decoupled for the non-bonded interactions, but the virtual
*/
if (strcmp(name, "decoupled") == 0)
{
- name = get_atomtype_name(atom->typeB, atype);
+ name = atypes->atomNameFromAtomType(atom->typeB);
}
return name;
}
-static bool calc_vsite3_param(gpp_atomtype *atype,
+static bool calc_vsite3_param(PreprocessingAtomTypes *atypes,
t_param *param, t_atoms *at,
int nrbond, t_mybonded *bonds,
int nrang, t_mybonded *angles )
/* check if this is part of a NH3 , NH2-umbrella or CH3 group,
* i.e. if atom k and l are dummy masses (MNH* or MCH3*) */
bXH3 =
- ( (gmx_strncasecmp(get_atomtype_name_AB(&at->atom[param->ak()], atype), "MNH", 3) == 0) &&
- (gmx_strncasecmp(get_atomtype_name_AB(&at->atom[param->al()], atype), "MNH", 3) == 0) ) ||
- ( (gmx_strncasecmp(get_atomtype_name_AB(&at->atom[param->ak()], atype), "MCH3", 4) == 0) &&
- (gmx_strncasecmp(get_atomtype_name_AB(&at->atom[param->al()], atype), "MCH3", 4) == 0) );
+ ( (gmx_strncasecmp(get_atomtype_name_AB(&at->atom[param->ak()], atypes), "MNH", 3) == 0) &&
+ (gmx_strncasecmp(get_atomtype_name_AB(&at->atom[param->al()], atypes), "MNH", 3) == 0) ) ||
+ ( (gmx_strncasecmp(get_atomtype_name_AB(&at->atom[param->ak()], atypes), "MCH3", 4) == 0) &&
+ (gmx_strncasecmp(get_atomtype_name_AB(&at->atom[param->al()], atypes), "MCH3", 4) == 0) );
bjk = get_bond_length(nrbond, bonds, param->aj(), param->ak());
bjl = get_bond_length(nrbond, bonds, param->aj(), param->al());
return bError;
}
-static bool calc_vsite3out_param(gpp_atomtype *atype,
+static bool calc_vsite3out_param(PreprocessingAtomTypes *atypes,
t_param *param, t_atoms *at,
int nrbond, t_mybonded *bonds,
int nrang, t_mybonded *angles)
/* check if this is part of a NH2-umbrella, NH3 or CH3 group,
* i.e. if atom k and l are dummy masses (MNH* or MCH3*) */
bXH3 =
- ( (gmx_strncasecmp(get_atomtype_name_AB(&at->atom[param->ak()], atype), "MNH", 3) == 0) &&
- (gmx_strncasecmp(get_atomtype_name_AB(&at->atom[param->al()], atype), "MNH", 3) == 0) ) ||
- ( (gmx_strncasecmp(get_atomtype_name_AB(&at->atom[param->ak()], atype), "MCH3", 4) == 0) &&
- (gmx_strncasecmp(get_atomtype_name_AB(&at->atom[param->al()], atype), "MCH3", 4) == 0) );
+ ( (gmx_strncasecmp(get_atomtype_name_AB(&at->atom[param->ak()], atypes), "MNH", 3) == 0) &&
+ (gmx_strncasecmp(get_atomtype_name_AB(&at->atom[param->al()], atypes), "MNH", 3) == 0) ) ||
+ ( (gmx_strncasecmp(get_atomtype_name_AB(&at->atom[param->ak()], atypes), "MCH3", 4) == 0) &&
+ (gmx_strncasecmp(get_atomtype_name_AB(&at->atom[param->al()], atypes), "MCH3", 4) == 0) );
/* check if construction parity must be swapped */
bSwapParity = ( param->c1() == -1 );
-int set_vsites(bool bVerbose, t_atoms *atoms, gpp_atomtype *atype,
+int set_vsites(bool bVerbose, t_atoms *atoms, PreprocessingAtomTypes *atypes,
gmx::ArrayRef<InteractionTypeParameters> plist)
{
int i, j, ftype;
{
case F_VSITE3:
bERROR =
- calc_vsite3_param(atype, &(plist[ftype].param[i]), atoms,
+ calc_vsite3_param(atypes, &(plist[ftype].param[i]), atoms,
nrbond, bonds, nrang, angles);
break;
case F_VSITE3FD:
break;
case F_VSITE3OUT:
bERROR =
- calc_vsite3out_param(atype, &(plist[ftype].param[i]), atoms,
+ calc_vsite3out_param(atypes, &(plist[ftype].param[i]), atoms,
nrbond, bonds, nrang, angles);
break;
case F_VSITE4FD:
#include "gromacs/utility/arrayref.h"
-struct gpp_atomtype;
+class PreprocessingAtomTypes;
struct gmx_moltype_t;
struct t_atoms;
struct InteractionTypeParameters;
-int set_vsites(bool bVerbose, t_atoms *atoms, gpp_atomtype *atype,
+int set_vsites(bool bVerbose, t_atoms *atoms, PreprocessingAtomTypes *atype,
gmx::ArrayRef<InteractionTypeParameters> plist);
/* set parameters for virtual sites, return number of virtual sites */
return cgnr;
}
-static gpp_atomtype *set_atom_type(t_symtab *tab, t_atoms *atoms, InteractionTypeParameters *bonds,
- int *nbonds, int nnm, t_nm2type nm2t[])
+static void set_atom_type(PreprocessingAtomTypes *atypes,
+ t_symtab *tab,
+ t_atoms *atoms,
+ InteractionTypeParameters *bonds,
+ int *nbonds,
+ int nnm,
+ t_nm2type nm2t[])
{
- gpp_atomtype *atype;
int nresolved;
- atype = init_atomtype();
snew(atoms->atomtype, atoms->nr);
- nresolved = nm2type(nnm, nm2t, tab, atoms, atype, nbonds, bonds);
+ nresolved = nm2type(nnm, nm2t, tab, atoms, atypes, nbonds, bonds);
if (nresolved != atoms->nr)
{
gmx_fatal(FARGS, "Could only find a forcefield type for %d out of %d atoms",
nresolved, atoms->nr);
}
- fprintf(stderr, "There are %d different atom types in your sample\n",
- get_atomtype_ntypes(atype));
-
- return atype;
+ fprintf(stderr, "There are %zu different atom types in your sample\n",
+ atypes->size());
}
static void lo_set_force_const(InteractionTypeParameters *plist, real c[], int nrfp, bool bRound,
static void dump_hybridization(FILE *fp, t_atoms *atoms, int nbonds[])
{
- int i;
-
- for (i = 0; (i < atoms->nr); i++)
+ for (int i = 0; (i < atoms->nr); i++)
{
fprintf(fp, "Atom %5s has %1d bonds\n", *atoms->atomname[i], nbonds[i]);
}
}
}
-static void print_rtp(const char *filenm, const char *title, t_atoms *atoms,
- gmx::ArrayRef<const InteractionTypeParameters> plist, gpp_atomtype *atype, int cgnr[])
+static void print_rtp(const char *filenm,
+ const char *title,
+ t_atoms *atoms,
+ gmx::ArrayRef<const InteractionTypeParameters> plist,
+ PreprocessingAtomTypes *atypes,
+ int cgnr[])
{
- FILE *fp;
- int i, tp;
- char *tpnm;
+ FILE *fp;
+ int i, tp;
+ const char *tpnm;
fp = gmx_fio_fopen(filenm, "w");
fprintf(fp, "; %s\n", title);
for (i = 0; (i < atoms->nr); i++)
{
tp = atoms->atom[i].type;
- if ((tpnm = get_atomtype_name(tp, atype)) == nullptr)
+ if ((tpnm = atypes->atomNameFromAtomType(tp)) == nullptr)
{
gmx_fatal(FARGS, "tp = %d, i = %d in print_rtp", tp, i);
}
FILE *fp;
std::array<InteractionTypeParameters, F_NRE> plist;
t_excls *excls;
- gpp_atomtype *atype;
t_nextnb nnb;
t_nm2type *nm2t;
t_mols mymol;
mk_bonds(nnm, nm2t, atoms, x, &(plist[F_BONDS]), nbonds, bPBC, box);
open_symtab(&symtab);
- atype = set_atom_type(&symtab, atoms, &(plist[F_BONDS]), nbonds, nnm, nm2t);
+ PreprocessingAtomTypes atypes;
+ set_atom_type(&atypes, &symtab, atoms, &(plist[F_BONDS]), nbonds, nnm, nm2t);
/* Make Angles and Dihedrals */
snew(excls, atoms->nr);
fp = ftp2FILE(efTOP, NFILE, fnm, "w");
print_top_header(fp, ftp2fn(efTOP, NFILE, fnm), TRUE, ffdir, 1.0);
- write_top(fp, nullptr, mymol.name.c_str(), atoms, FALSE, bts, plist, excls, atype,
+ write_top(fp, nullptr, mymol.name.c_str(), atoms, FALSE, bts, plist, excls, &atypes,
cgnr, rtp_header_settings.nrexcl);
print_top_mols(fp, mymol.name.c_str(), ffdir, nullptr, {}, gmx::arrayRefFromArray(&mymol, 1));
if (bRTP)
{
print_rtp(ftp2fn(efRTP, NFILE, fnm), "Generated by x2top",
- atoms, plist, atype, cgnr);
+ atoms, plist, &atypes, cgnr);
}
if (debug)
biod.pnpi.spb.ru / alexxy / gromacs.git / commitdiff