Enables range based for loop. Also renamed to AtomIterator.
Change-Id: I87e0ef50c4e0899a8954d584d5200bb6f5c8dd29
const char *format = get_hconf_format(v != nullptr);
- SystemAtomIterator aloop(*mtop);
- while (aloop.nextAtom())
+ for (const AtomProxy &atomP : AtomRange(*mtop))
{
- int i = aloop.globalAtomNumber();
- int residueNumber = aloop.residueNumber();
- const char *atomName = aloop.atomName();
- const char *residueName = aloop.residueName();
+ int i = atomP.globalAtomNumber();
+ int residueNumber = atomP.residueNumber();
+ const char *atomName = atomP.atomName();
+ const char *residueName = atomP.residueName();
fprintf(out, "%5d%-5.5s%5.5s%5d",
residueNumber%100000, residueName, atomName, (i+1)%100000);
snew(*x, mtop->natoms);
nq = 0;
- SystemAtomIterator aloop(*mtop);
- while (aloop.nextAtom())
+ for (const AtomProxy &atomP : AtomRange(*mtop))
{
- const t_atom &local = aloop.atom();
- int i = aloop.globalAtomNumber();
+ const t_atom &local = atomP.atom();
+ int i = atomP.globalAtomNumber();
if (is_charge(local.q))
{
(*q)[nq] = local.q;
boltz = BOLTZ*tempi;
ekin = 0.0;
nrdf = 0;
- SystemAtomIterator aloop(*mtop);
- while (aloop.nextAtom())
+ for (const AtomProxy &atomP : AtomRange(*mtop))
{
- const t_atom &local = aloop.atom();
- int i = aloop.globalAtomNumber();
+ const t_atom &local = atomP.atom();
+ int i = atomP.globalAtomNumber();
real mass = local.m;
if (mass > 0)
{
static void check_vel(gmx_mtop_t *mtop, rvec v[])
{
- SystemAtomIterator aloop(*mtop);
- while (aloop.nextAtom())
+ for (const AtomProxy &atomP : AtomRange(*mtop))
{
- const t_atom &local = aloop.atom();
- int i = aloop.globalAtomNumber();
+ const t_atom &local = atomP.atom();
+ int i = atomP.globalAtomNumber();
if (local.ptype == eptShell ||
local.ptype == eptBond ||
local.ptype == eptVSite)
int nshells = 0;
char warn_buf[STRLEN];
- SystemAtomIterator aloop(*mtop);
- while (aloop.nextAtom())
+ for (const AtomProxy &atomP : AtomRange(*mtop))
{
- const t_atom &local = aloop.atom();
+ const t_atom &local = atomP.atom();
if (local.ptype == eptShell ||
local.ptype == eptBond)
{
real *mass;
snew(mass, state->natoms);
- SystemAtomIterator aloop(*sys);
- while (aloop.nextAtom())
+ for (const AtomProxy &atomP : AtomRange(*sys))
{
- const t_atom &local = aloop.atom();
- int i = aloop.globalAtomNumber();
+ const t_atom &local = atomP.atom();
+ int i = atomP.globalAtomNumber();
mass[i] = local.m;
}
rvec *v)
{
double sum_mv2 = 0;
- SystemAtomIterator aloop(*mtop);
- while (aloop.nextAtom())
+ for (const AtomProxy &atomP : AtomRange(*mtop))
{
- const t_atom &local = aloop.atom();
- int i = aloop.globalAtomNumber();
+ const t_atom &local = atomP.atom();
+ int i = atomP.globalAtomNumber();
sum_mv2 += local.m*norm2(v[i]);
}
}
snew(nrdf2, natoms);
- SystemAtomIterator aloop(*mtop);
- while (aloop.nextAtom())
+ for (const AtomProxy &atomP : AtomRange(*mtop))
{
- const t_atom &local = aloop.atom();
- int i = aloop.globalAtomNumber();
+ const t_atom &local = atomP.atom();
+ int i = atomP.globalAtomNumber();
nrdf2[i] = 0;
if (local.ptype == eptAtom || local.ptype == eptNucleus)
{
{
clear_rvec(acc);
snew(mgrp, sys->groups.grps[egcACC].nr);
- SystemAtomIterator aloop(*sys);
- while (aloop.nextAtom())
+ for (const AtomProxy &atomP : AtomRange(*sys))
{
- const t_atom &local = aloop.atom();
- int i = aloop.globalAtomNumber();
+ const t_atom &local = atomP.atom();
+ int i = atomP.globalAtomNumber();
mgrp[getGroupType(sys->groups, egcACC, i)] += local.m;
}
mt = 0.0;
rvec com = { 0, 0, 0 };
double mtot = 0.0;
int j = 0;
- SystemAtomIterator aloop(*mtop);
- while (aloop.nextAtom())
+ for (const AtomProxy &atomP : AtomRange(*mtop))
{
- const t_atom &local = aloop.atom();
- int i = aloop.globalAtomNumber();
+ const t_atom &local = atomP.atom();
+ int i = atomP.globalAtomNumber();
if (mtop->groups.grpnr[egcORFIT] == nullptr ||
mtop->groups.grpnr[egcORFIT][i] == 0)
{
std::vector<int> qmmmAtoms;
for (int i = 0; i < numQmmmGroups; i++)
{
- SystemAtomIterator aloop(mtop);
- while (aloop.nextAtom())
+ for (const AtomProxy &atomP : AtomRange(mtop))
{
- int index = aloop.globalAtomNumber();
+ int index = atomP.globalAtomNumber();
if (getGroupType(groups, egcQMMM, index) == i)
{
qmmmAtoms.push_back(index);
/* Global system sized array, this should be avoided */
snew(shell_index, mtop->natoms);
- SystemAtomIterator aloop(*mtop);
nshell = 0;
- while (aloop.nextAtom())
+ for (const AtomProxy &atomP : AtomRange(*mtop))
{
- const t_atom &local = aloop.atom();
- int i = aloop.globalAtomNumber();
+ const t_atom &local = atomP.atom();
+ int i = atomP.globalAtomNumber();
if (local.ptype == eptShell)
{
shell_index[i] = nshell++;
if (ngacc > 0)
{
const gmx_groups_t &groups = mtop->groups;
- SystemAtomIterator aloop(*mtop);
- while (aloop.nextAtom())
+ for (const AtomProxy &atomP : AtomRange(*mtop))
{
- const t_atom &local = aloop.atom();
- int i = aloop.globalAtomNumber();
+ const t_atom &local = atomP.atom();
+ int i = atomP.globalAtomNumber();
int grp = getGroupType(groups, egcACC, i);
if ((grp < 0) && (grp >= ngacc))
{
)
if (BUILD_TESTING)
-# add_subdirectory(tests)
+ add_subdirectory(tests)
endif()
}
}
-SystemAtomIterator::SystemAtomIterator(const gmx_mtop_t &mtop)
+AtomIterator::AtomIterator(const gmx_mtop_t &mtop, int globalAtomNumber)
: mtop_(&mtop), mblock_(0),
atoms_(&mtop.moltype[mtop.molblock[0].type].atoms),
currentMolecule_(0), highestResidueNumber_(mtop.maxresnr),
- localAtomNumber_(-1), globalAtomNumber_(-1)
-{}
+ localAtomNumber_(0), globalAtomNumber_(globalAtomNumber)
+{
+ GMX_ASSERT(globalAtomNumber == 0 || globalAtomNumber == mtop.natoms,
+ "Starting at other atoms not implemented yet");
+}
-bool SystemAtomIterator::nextAtom()
+AtomIterator &AtomIterator::operator++()
{
localAtomNumber_++;
globalAtomNumber_++;
mblock_++;
if (mblock_ >= mtop_->molblock.size())
{
- return false;
+ return *this;
}
atoms_ = &mtop_->moltype[mtop_->molblock[mblock_].type].atoms;
currentMolecule_ = 0;
}
}
- return true;
+ return *this;
+}
+
+AtomIterator AtomIterator::operator++(int)
+{
+ AtomIterator temp = *this;
+ ++(*this);
+ return temp;
+}
+
+bool AtomIterator::operator==(const AtomIterator &o) const
+{
+ return mtop_ == o.mtop_ && globalAtomNumber_ == o.globalAtomNumber_;
+}
+
+bool AtomIterator::operator!=(const AtomIterator &o) const
+{
+ return !(*this == o);
}
-const t_atom &SystemAtomIterator::atom() const
+const t_atom &AtomProxy::atom() const
{
- return atoms_->atom[localAtomNumber_];
+ return it_->atoms_->atom[it_->localAtomNumber_];
}
-int SystemAtomIterator::globalAtomNumber() const
+int AtomProxy::globalAtomNumber() const
{
- return globalAtomNumber_;
+ return it_->globalAtomNumber_;
}
-const char *SystemAtomIterator::atomName() const
+const char *AtomProxy::atomName() const
{
- return *(atoms_->atomname[localAtomNumber_]);
+ return *(it_->atoms_->atomname[it_->localAtomNumber_]);
}
-const char *SystemAtomIterator::residueName() const
+const char *AtomProxy::residueName() const
{
- int residueIndexInMolecule = atoms_->atom[localAtomNumber_].resind;
- return *(atoms_->resinfo[residueIndexInMolecule].name);
+ int residueIndexInMolecule = it_->atoms_->atom[it_->localAtomNumber_].resind;
+ return *(it_->atoms_->resinfo[residueIndexInMolecule].name);
}
-int SystemAtomIterator::residueNumber() const
+int AtomProxy::residueNumber() const
{
- int residueIndexInMolecule = atoms_->atom[localAtomNumber_].resind;
- if (atoms_->nres <= mtop_->maxres_renum)
+ int residueIndexInMolecule = it_->atoms_->atom[it_->localAtomNumber_].resind;
+ if (it_->atoms_->nres <= it_->mtop_->maxres_renum)
{
- return highestResidueNumber_ + 1 + residueIndexInMolecule;
+ return it_->highestResidueNumber_ + 1 + residueIndexInMolecule;
}
else
{
- return atoms_->resinfo[residueIndexInMolecule].nr;
+ return it_->atoms_->resinfo[residueIndexInMolecule].nr;
}
}
-const gmx_moltype_t &SystemAtomIterator::moleculeType() const
+const gmx_moltype_t &AtomProxy::moleculeType() const
{
- return mtop_->moltype[mtop_->molblock[mblock_].type];
+ return it_->mtop_->moltype[it_->mtop_->molblock[it_->mblock_].type];
}
-int SystemAtomIterator::atomNumberInMol() const
+int AtomProxy::atomNumberInMol() const
{
- return localAtomNumber_;
+ return it_->localAtomNumber_;
}
typedef struct gmx_mtop_atomloop_block
{
std::vector<real> qA(mtop.natoms);
std::vector<real> qB(mtop.natoms);
- SystemAtomIterator aloop(mtop);
- while (aloop.nextAtom())
+ for (const AtomProxy &atomP : AtomRange(mtop))
{
- const t_atom &local = aloop.atom();
- int index = aloop.globalAtomNumber();
+ const t_atom &local = atomP.atom();
+ int index = atomP.globalAtomNumber();
qA[index] = local.q;
qB[index] = local.qB;
}
{
std::vector<int> atom_index;
- SystemAtomIterator aloop(*mtop);
- while (aloop.nextAtom())
+ for (const AtomProxy &atomP : AtomRange(*mtop))
{
- const t_atom &local = aloop.atom();
- int index = aloop.globalAtomNumber();
+ const t_atom &local = atomP.atom();
+ int index = atomP.globalAtomNumber();
if (local.ptype == eptAtom)
{
atom_index.push_back(index);
/* Removes the charge groups, i.e. makes single atom charge groups, in mtop */
void gmx_mtop_remove_chargegroups(gmx_mtop_t *mtop);
-/*! \brief
- * Object that allows looping over all atoms in an mtop.
- */
-class SystemAtomIterator
+class AtomIterator;
+
+//! Proxy object returned from AtomIterator
+class AtomProxy
{
public:
//! Default constructor.
- explicit SystemAtomIterator(const gmx_mtop_t &mtop);
-
- /*! \brief
- * Loop to next atoms.
- *
- * When not at the end:
- * returns true.
- * sets internal variables to next atom, and sets new global atom number.
- * When at the end, returns false and resets the loop.
- * Use as:
- * LoopOverAllAtom aloop(mtop);
- * while (aloop.nextAtom()) {
- * ...
- * }
- */
- bool nextAtom();
+ AtomProxy(const AtomIterator* it) : it_(it) {}
//! Access current global atom number.
int globalAtomNumber() const;
//! Access current t_atom struct.
const gmx_moltype_t &moleculeType() const;
//! Access the position of the current atom in the molecule.
int atomNumberInMol() const;
+ private:
+ const AtomIterator* it_;
+};
+
+//! Wrapper around proxy object to implement operator->
+template <typename T>
+class ProxyPtr
+{
+ public:
+ //! Construct with proxy object.
+ ProxyPtr(T t) : t_(t) {}
+ //! Member of pointer operator.
+ T* operator->() { return &t_; }
+ private:
+ T t_;
+};
+
+/*! \brief
+ * Object that allows looping over all atoms in an mtop.
+ */
+class AtomIterator
+{
+ public:
+ //! Construct from topology and optionalally a global atom number.
+ explicit AtomIterator(const gmx_mtop_t &mtop, int globalAtomNumber = 0);
+
+ //! Prefix increment.
+ AtomIterator &operator++();
+ //! Postfix increment.
+ AtomIterator operator++(int);
+
+ //! Equality comparison.
+ bool operator==(const AtomIterator &o) const;
+ //! Non-equal comparison.
+ bool operator!=(const AtomIterator &o) const;
+
+ //! Dereference operator. Returns proxy.
+ AtomProxy operator*() const { return {this}; }
+ //! Member of pointer operator.
+ ProxyPtr<AtomProxy> operator->() const { return {this}; };
private:
//! Global topology.
int localAtomNumber_;
//! Global current atom number.
int globalAtomNumber_;
+
+ friend class AtomProxy;
+};
+
+//! Range over all atoms of topology.
+class AtomRange
+{
+ public:
+ //! Default constructor.
+ explicit AtomRange(const gmx_mtop_t &mtop) :
+ begin_(mtop), end_(mtop, mtop.natoms) {}
+ //! Iterator to begin of range.
+ AtomIterator &begin() { return begin_; }
+ //! Iterator to end of range.
+ AtomIterator &end() { return end_; }
+ private:
+ AtomIterator begin_, end_;
};
/* Abstract type for atom loop over atoms in all molecule blocks */
--- /dev/null
+#
+# This file is part of the GROMACS molecular simulation package.
+#
+# Copyright (c) 2019, 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.
+#
+# GROMACS is free software; you can redistribute it and/or
+# modify it under the terms of the GNU Lesser General Public License
+# as published by the Free Software Foundation; either version 2.1
+# of the License, or (at your option) any later version.
+#
+# GROMACS is distributed in the hope that it will be useful,
+# but WITHOUT ANY WARRANTY; without even the implied warranty of
+# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+# Lesser General Public License for more details.
+#
+# You should have received a copy of the GNU Lesser General Public
+# License along with GROMACS; if not, see
+# http://www.gnu.org/licenses, or write to the Free Software Foundation,
+# Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
+#
+# If you want to redistribute modifications to GROMACS, please
+# consider that scientific software is very special. Version
+# control is crucial - bugs must be traceable. We will be happy to
+# consider code for inclusion in the official distribution, but
+# derived work must not be called official GROMACS. Details are found
+# in the README & COPYING files - if they are missing, get the
+# official version at http://www.gromacs.org.
+#
+# To help us fund GROMACS development, we humbly ask that you cite
+# the research papers on the package. Check out http://www.gromacs.org.
+
+gmx_add_unit_test(TopologyTest topology-test
+ mtop.cpp)
+
--- /dev/null
+/*
+ * This file is part of the GROMACS molecular simulation package.
+ *
+ * Copyright (c) 2019, 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.
+ *
+ * GROMACS is free software; you can redistribute it and/or
+ * modify it under the terms of the GNU Lesser General Public License
+ * as published by the Free Software Foundation; either version 2.1
+ * of the License, or (at your option) any later version.
+ *
+ * GROMACS is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ * Lesser General Public License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public
+ * License along with GROMACS; if not, see
+ * http://www.gnu.org/licenses, or write to the Free Software Foundation,
+ * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
+ *
+ * If you want to redistribute modifications to GROMACS, please
+ * consider that scientific software is very special. Version
+ * control is crucial - bugs must be traceable. We will be happy to
+ * consider code for inclusion in the official distribution, but
+ * derived work must not be called official GROMACS. Details are found
+ * in the README & COPYING files - if they are missing, get the
+ * official version at http://www.gromacs.org.
+ *
+ * To help us fund GROMACS development, we humbly ask that you cite
+ * the research papers on the package. Check out http://www.gromacs.org.
+ */
+/*! \internal \file
+ * \brief
+ * Implements test of mtop routines
+ *
+ * \author Roland Schulz <roland.schulz@intel.com>
+ * \ingroup module_topology
+ */
+#include "gmxpre.h"
+
+#include <gtest/gtest.h>
+
+#include "gromacs/topology/mtop_util.h"
+
+namespace gmx
+{
+namespace
+{
+
+/*! \brief Initializes a basic topology with 9 atoms with settle*/
+void createBasicTop(gmx_mtop_t* mtop)
+{
+ gmx_moltype_t moltype;
+ moltype.atoms.nr = NRAL(F_SETTLE);
+ std::vector<int> &iatoms = moltype.ilist[F_SETTLE].iatoms;
+ const int settleType = 0;
+ iatoms.push_back(settleType);
+ iatoms.push_back(0);
+ iatoms.push_back(1);
+ iatoms.push_back(2);
+ mtop->moltype.push_back(moltype);
+
+ mtop->molblock.resize(1);
+ mtop->molblock[0].type = 0;
+ mtop->molblock[0].nmol = 3;
+ mtop->natoms = moltype.atoms.nr * mtop->molblock[0].nmol;
+ gmx_mtop_finalize(mtop);
+}
+
+TEST(MtopTest, RangeBasedLoop)
+{
+ gmx_mtop_t mtop;
+ createBasicTop(&mtop);
+ int count = 0;
+ for (const AtomProxy &atomP : AtomRange(mtop))
+ {
+ EXPECT_EQ(atomP.globalAtomNumber(), count);
+ ++count;
+ }
+ EXPECT_EQ(count, 9);
+}
+
+TEST(MtopTest, Operators)
+{
+ gmx_mtop_t mtop;
+ createBasicTop(&mtop);
+ AtomIterator it(mtop);
+ AtomIterator otherIt(mtop);
+ EXPECT_EQ((*it).globalAtomNumber(), 0);
+ EXPECT_EQ(it->globalAtomNumber(), 0);
+ EXPECT_TRUE (it == otherIt);
+ EXPECT_FALSE(it != otherIt);
+ ++it;
+ EXPECT_EQ(it->globalAtomNumber(), 1);
+ it++;
+ EXPECT_EQ(it->globalAtomNumber(), 2);
+ EXPECT_TRUE (it != otherIt);
+ EXPECT_FALSE(it == otherIt);
+}
+
+} // namespace
+
+} // namespace gmx