{
do_mtop(serializer, mtop, tpx->fileVersion);
set_disres_npair(mtop);
- gmx_mtop_finalize(mtop);
+ mtop->finalize();
}
else
{
/* Copy structures from msys to sys */
molinfo2mtop(*mi, sys);
- gmx_mtop_finalize(sys);
+ sys->finalize();
/* Discourage using the, previously common, setup of applying constraints
* to all bonds with force fields that have been parametrized with H-bond
// Re-establish the invariants of the derived values within
// mtop.
- gmx_mtop_finalize(mtop);
+ mtop->finalize();
if (ftp2bSet(efTOP, nfile, fnm))
{
/*
* This file is part of the GROMACS molecular simulation package.
*
- * Copyright (c) 2018,2019, by the GROMACS development team, led by
+ * Copyright (c) 2018,2019,2020, by the GROMACS development team, led by
* Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
* and including many others, as listed in the AUTHORS file in the
* top-level source directory and at http://www.gromacs.org.
mtop.molblock[0].type = 0;
mtop.molblock[0].nmol = numMolecules;
mtop.natoms = numAtoms;
- gmx_mtop_finalize(&mtop);
+ mtop.finalize();
// Set up the SETTLE parameters.
const real dOH = 0.1;
molblock.nmol = 2;
molblock.type = 0;
testTopology_.natoms = moltype.atoms.nr * molblock.nmol;
- gmx_mtop_finalize(&testTopology_);
+ testTopology_.finalize();
FloatingPointTolerance tolerance(
FloatingPointTolerance(1.0e-6, 1.0e-6, 1.0e-8, 1.0e-12, 10000, 100, false));
checker_.setDefaultTolerance(tolerance);
mtop_->molblock[0].type = 0;
mtop_->molblock[0].nmol = 1;
mtop_->natoms = count;
- mtop_->maxres_renum = 0;
- gmx_mtop_finalize(mtop_.get());
- GMX_RELEASE_ASSERT(mtop_->maxres_renum == 0,
+ mtop_->finalize();
+ GMX_RELEASE_ASSERT(mtop_->maxResiduesPerMoleculeToTriggerRenumber() == 0,
"maxres_renum in mtop can be modified by an env.var., that is not supported "
"in this test");
t_atoms& atoms = this->atoms();
{
GMX_RELEASE_ASSERT(mtop_ != nullptr, "Topology not initialized");
- mtop_->maxres_renum = 0;
mtop_->haveMoleculeIndices = true;
- gmx_mtop_finalize(mtop_.get());
+ mtop_->finalize();
}
void TopologyManager::initUniformResidues(int residueSize)
"The residues should break at molecule boundaries");
atoms.nres = nres;
mtop_->haveMoleculeIndices = true;
- gmx_mtop_finalize(mtop_.get());
+ mtop_->finalize();
}
void TopologyManager::initFrameIndices(const ArrayRef<const int>& index)
/*
* This file is part of the GROMACS molecular simulation package.
*
- * Copyright (c) 2016,2017,2018,2019, by the GROMACS development team, led by
+ * Copyright (c) 2016,2017,2018,2019,2020, by the GROMACS development team, led by
* Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
* and including many others, as listed in the AUTHORS file in the
* top-level source directory and at http://www.gromacs.org.
}
if (residueNumber != nullptr)
{
- if (atoms.nres > mtop->maxres_renum)
+ if (atoms.nres > mtop->maxResiduesPerMoleculeToTriggerRenumber())
{
*residueNumber = atoms.resinfo[atoms.atom[atomIndexInMolecule].resind].nr;
}
#include "gromacs/utility/real.h"
#include "gromacs/utility/smalloc.h"
-static int gmx_mtop_maxresnr(const gmx_mtop_t* mtop, int maxres_renum)
-{
- int maxresnr = 0;
-
- for (const gmx_moltype_t& moltype : mtop->moltype)
- {
- const t_atoms& atoms = moltype.atoms;
- if (atoms.nres > maxres_renum)
- {
- for (int r = 0; r < atoms.nres; r++)
- {
- if (atoms.resinfo[r].nr > maxresnr)
- {
- maxresnr = atoms.resinfo[r].nr;
- }
- }
- }
- }
-
- return maxresnr;
-}
-
-static void buildMolblockIndices(gmx_mtop_t* mtop)
-{
- mtop->moleculeBlockIndices.resize(mtop->molblock.size());
-
- int atomIndex = 0;
- int residueIndex = 0;
- int residueNumberStart = mtop->maxresnr + 1;
- int moleculeIndexStart = 0;
- for (size_t mb = 0; mb < mtop->molblock.size(); mb++)
- {
- const gmx_molblock_t& molb = mtop->molblock[mb];
- MoleculeBlockIndices& indices = mtop->moleculeBlockIndices[mb];
- const int numResPerMol = mtop->moltype[molb.type].atoms.nres;
-
- indices.numAtomsPerMolecule = mtop->moltype[molb.type].atoms.nr;
- indices.globalAtomStart = atomIndex;
- indices.globalResidueStart = residueIndex;
- atomIndex += molb.nmol * indices.numAtomsPerMolecule;
- residueIndex += molb.nmol * numResPerMol;
- indices.globalAtomEnd = atomIndex;
- indices.residueNumberStart = residueNumberStart;
- if (numResPerMol <= mtop->maxres_renum)
- {
- residueNumberStart += molb.nmol * numResPerMol;
- }
- indices.moleculeIndexStart = moleculeIndexStart;
- moleculeIndexStart += molb.nmol;
- }
-}
-
-void gmx_mtop_finalize(gmx_mtop_t* mtop)
-{
- char* env;
-
- if (mtop->molblock.size() == 1 && mtop->molblock[0].nmol == 1)
- {
- /* We have a single molecule only, no renumbering needed.
- * This case also covers an mtop converted from pdb/gro/... input,
- * so we retain the original residue numbering.
- */
- mtop->maxres_renum = 0;
- }
- else
- {
- /* We only renumber single residue molecules. Their intra-molecular
- * residue numbering is anyhow irrelevant.
- */
- mtop->maxres_renum = 1;
- }
-
- env = getenv("GMX_MAXRESRENUM");
- if (env != nullptr)
- {
- sscanf(env, "%d", &mtop->maxres_renum);
- }
- if (mtop->maxres_renum == -1)
- {
- /* -1 signals renumber residues in all molecules */
- mtop->maxres_renum = INT_MAX;
- }
-
- mtop->maxresnr = gmx_mtop_maxresnr(mtop, mtop->maxres_renum);
-
- buildMolblockIndices(mtop);
-}
-
void gmx_mtop_count_atomtypes(const gmx_mtop_t* mtop, int state, int typecount[])
{
for (int i = 0; i < mtop->ffparams.atnr; ++i)
mblock_(0),
atoms_(&mtop.moltype[mtop.molblock[0].type].atoms),
currentMolecule_(0),
- highestResidueNumber_(mtop.maxresnr),
+ highestResidueNumber_(mtop.maxResNumberNotRenumbered()),
localAtomNumber_(0),
globalAtomNumber_(globalAtomNumber)
{
if (localAtomNumber_ >= atoms_->nr)
{
- if (atoms_->nres <= mtop_->maxresnr)
+ if (atoms_->nres <= mtop_->maxResiduesPerMoleculeToTriggerRenumber())
{
/* Single residue molecule, increase the count with one */
highestResidueNumber_ += atoms_->nres;
int AtomProxy::residueNumber() const
{
int residueIndexInMolecule = it_->atoms_->atom[it_->localAtomNumber_].resind;
- if (it_->atoms_->nres <= it_->mtop_->maxres_renum)
+ if (it_->atoms_->nres <= it_->mtop_->maxResiduesPerMoleculeToTriggerRenumber())
{
return it_->highestResidueNumber_ + 1 + residueIndexInMolecule;
}
init_t_atoms(&atoms, 0, FALSE);
- int maxresnr = mtop->maxresnr;
+ int maxresnr = mtop->maxResNumberNotRenumbered();
for (const gmx_molblock_t& molb : mtop->molblock)
{
- atomcat(&atoms, &mtop->moltype[molb.type].atoms, molb.nmol, mtop->maxres_renum, &maxresnr);
+ atomcat(&atoms, &mtop->moltype[molb.type].atoms, molb.nmol,
+ mtop->maxResiduesPerMoleculeToTriggerRenumber(), &maxresnr);
}
return atoms;
mtop->haveMoleculeIndices = false;
- gmx_mtop_finalize(mtop);
+ mtop->finalize();
}
bool haveFepPerturbedNBInteractions(const gmx_mtop_t* mtop)
// and should be replaced by versions taking const gmx_mtop & when
// their callers are refactored similarly.
-/* Should be called after generating or reading mtop,
- * to set some compute intesive variables to avoid
- * N^2 operations later on.
- */
-void gmx_mtop_finalize(gmx_mtop_t* mtop);
-
/* Counts the number of atoms of each type. State should be 0 for
* state A and 1 for state B types. typecount should have at
* least mtop->ffparams.atnr elements.
mtop->molblock[0].type = 0;
mtop->molblock[0].nmol = 3;
mtop->natoms = moltype.atoms.nr * mtop->molblock[0].nmol;
- gmx_mtop_finalize(mtop);
+ mtop->finalize();
}
/*! \brief
mtop->molblock[0].type = 0;
mtop->molblock[0].nmol = 1;
mtop->natoms = moltype.atoms.nr * mtop->molblock[0].nmol;
- gmx_mtop_finalize(mtop);
+ mtop->finalize();
std::vector<gmx::Range<int>> residueRange;
residueRange.emplace_back(0, residueOneSize);
residueRange.emplace_back(residueOneSize, residueOneSize + residueTwoSize);
done_atom(&atoms);
}
+static int gmx_mtop_maxresnr(const gmx::ArrayRef<const gmx_moltype_t> moltypes, int maxres_renum)
+{
+ int maxresnr = 0;
+
+ for (const gmx_moltype_t& moltype : moltypes)
+ {
+ const t_atoms& atoms = moltype.atoms;
+ if (atoms.nres > maxres_renum)
+ {
+ for (int r = 0; r < atoms.nres; r++)
+ {
+ if (atoms.resinfo[r].nr > maxresnr)
+ {
+ maxresnr = atoms.resinfo[r].nr;
+ }
+ }
+ }
+ }
+
+ return maxresnr;
+}
+
gmx_mtop_t::gmx_mtop_t()
{
init_atomtypes(&atomtypes);
done_atomtypes(&atomtypes);
}
+void gmx_mtop_t::finalize()
+{
+ if (molblock.size() == 1 && molblock[0].nmol == 1)
+ {
+ /* We have a single molecule only, no renumbering needed.
+ * This case also covers an mtop converted from pdb/gro/... input,
+ * so we retain the original residue numbering.
+ */
+ maxResiduesPerMoleculeToTriggerRenumber_ = 0;
+ }
+ else
+ {
+ /* We only renumber single residue molecules. Their intra-molecular
+ * residue numbering is anyhow irrelevant.
+ */
+ maxResiduesPerMoleculeToTriggerRenumber_ = 1;
+ }
+
+ const char* env = getenv("GMX_MAXRESRENUM");
+ if (env != nullptr)
+ {
+ sscanf(env, "%d", &maxResiduesPerMoleculeToTriggerRenumber_);
+ }
+ if (maxResiduesPerMoleculeToTriggerRenumber_ == -1)
+ {
+ /* -1 signals renumber residues in all molecules */
+ maxResiduesPerMoleculeToTriggerRenumber_ = std::numeric_limits<int>::max();
+ }
+
+ maxResNumberNotRenumbered_ = gmx_mtop_maxresnr(moltype, maxResiduesPerMoleculeToTriggerRenumber_);
+
+ buildMolblockIndices();
+}
+
+void gmx_mtop_t::buildMolblockIndices()
+{
+ moleculeBlockIndices.resize(molblock.size());
+
+ int atomIndex = 0;
+ int residueIndex = 0;
+ int residueNumberStart = maxResNumberNotRenumbered_ + 1;
+ int moleculeIndexStart = 0;
+ for (size_t mb = 0; mb < molblock.size(); mb++)
+ {
+ const gmx_molblock_t& molb = molblock[mb];
+ MoleculeBlockIndices& indices = moleculeBlockIndices[mb];
+ const int numResPerMol = moltype[molb.type].atoms.nres;
+
+ indices.numAtomsPerMolecule = moltype[molb.type].atoms.nr;
+ indices.globalAtomStart = atomIndex;
+ indices.globalResidueStart = residueIndex;
+ atomIndex += molb.nmol * indices.numAtomsPerMolecule;
+ residueIndex += molb.nmol * numResPerMol;
+ indices.globalAtomEnd = atomIndex;
+ indices.residueNumberStart = residueNumberStart;
+ if (numResPerMol <= maxResiduesPerMoleculeToTriggerRenumber_)
+ {
+ residueNumberStart += molb.nmol * numResPerMol;
+ }
+ indices.moleculeIndexStart = moleculeIndexStart;
+ moleculeIndexStart += molb.nmol;
+ }
+}
+
void done_top(t_topology* top)
{
done_idef(&top->idef);
fprintf(fp, "comparing mtop topology\n");
cmp_str(fp, "Name", -1, *mtop1.name, *mtop2.name);
cmp_int(fp, "natoms", -1, mtop1.natoms, mtop2.natoms);
- cmp_int(fp, "maxres_renum", -1, mtop1.maxres_renum, mtop2.maxres_renum);
- cmp_int(fp, "maxresnr", -1, mtop1.maxresnr, mtop2.maxresnr);
+ cmp_int(fp, "maxres_renum", -1, mtop1.maxResiduesPerMoleculeToTriggerRenumber(),
+ mtop2.maxResiduesPerMoleculeToTriggerRenumber());
+ cmp_int(fp, "maxresnr", -1, mtop1.maxResNumberNotRenumbered(), mtop2.maxResNumberNotRenumbered());
cmp_bool(fp, "bIntermolecularInteractions", -1, mtop1.bIntermolecularInteractions,
mtop2.bIntermolecularInteractions);
cmp_bool(fp, "haveMoleculeIndices", -1, mtop1.haveMoleculeIndices, mtop2.haveMoleculeIndices);
std::unique_ptr<InteractionLists> intermolecular_ilist = nullptr;
//! Number of global atoms.
int natoms = 0;
- //! Parameter for residue numbering.
- int maxres_renum = 0;
- //! The maximum residue number in moltype
- int maxresnr = -1;
//! Atomtype properties
t_atomtypes atomtypes;
//! Groups of atoms for different purposes
*/
std::vector<int> intermolecularExclusionGroup;
- /* Derived data below */
+ //! Maximum number of residues in molecule to trigger renumbering of residues
+ int maxResiduesPerMoleculeToTriggerRenumber() const
+ {
+ return maxResiduesPerMoleculeToTriggerRenumber_;
+ }
+ //! Maximum residue number that is not renumbered.
+ int maxResNumberNotRenumbered() const { return maxResNumberNotRenumbered_; }
+ /*! \brief Finalize this data structure.
+ *
+ * Should be called after generating or reading mtop, to set some compute
+ * intesive variables to avoid N^2 operations later on.
+ *
+ * \todo Move into a builder class, once available.
+ */
+ void finalize();
+
+ /* Derived data below */
//! Indices for each molblock entry for fast lookup of atom properties
std::vector<MoleculeBlockIndices> moleculeBlockIndices;
+
+private:
+ //! Build the molblock indices
+ void buildMolblockIndices();
+ //! Maximum number of residues in molecule to trigger renumbering of residues
+ int maxResiduesPerMoleculeToTriggerRenumber_ = 0;
+ //! The maximum residue number in moltype that is not renumbered
+ int maxResNumberNotRenumbered_ = -1;
};
/*! \brief