#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_;
}
}