Speed up mtop_util atom lookup
authorBerk Hess <hess@kth.se>
Tue, 27 Sep 2016 12:50:05 +0000 (14:50 +0200)
committerBerk Hess <hess@kth.se>
Tue, 11 Oct 2016 09:02:53 +0000 (11:02 +0200)
The lookup of atom indices and properties on global atom index have
been sped up by moving functions to a new header file mtop_lookup.h
and by storing start and end global atom indices in gmx_mtop_t.
Another performance improvement is that the previous molblock index is
used as starting value for the next search.
The atom+residue lookup function now also returns the reside index.
This change also simplifies the code, since we no longer need a lookup
data structure.
A large number of files are touched because the t_atom return pointer
is changed to const also in the atomloop functions.

Change-Id: I185b8c2e614604e9561190dd5e447077d88933ca

28 files changed:
src/gromacs/domdec/domdec.cpp
src/gromacs/domdec/domdec_constraints.cpp
src/gromacs/essentialdynamics/edsam.cpp
src/gromacs/fileio/groio.cpp
src/gromacs/gmxana/gmx_clustsize.cpp
src/gromacs/gmxana/gmx_energy.cpp
src/gromacs/gmxana/gmx_pme_error.cpp
src/gromacs/gmxpreprocess/gen_maxwell_velocities.cpp
src/gromacs/gmxpreprocess/grompp.cpp
src/gromacs/gmxpreprocess/readir.cpp
src/gromacs/listed-forces/orires.cpp
src/gromacs/mdlib/broadcaststructs.cpp
src/gromacs/mdlib/constr.cpp
src/gromacs/mdlib/mdatoms.cpp
src/gromacs/mdlib/qmmm.cpp
src/gromacs/mdlib/shellfc.cpp
src/gromacs/mdlib/tgroup.cpp
src/gromacs/pulling/pull.cpp
src/gromacs/pulling/pull_rotation.cpp
src/gromacs/swap/swapcoords.cpp
src/gromacs/swap/swapcoords.h
src/gromacs/tools/check.cpp
src/gromacs/topology/mtop_lookup.h [new file with mode: 0644]
src/gromacs/topology/mtop_util.cpp
src/gromacs/topology/mtop_util.h
src/gromacs/topology/topology.h
src/programs/mdrun/md.cpp
src/programs/mdrun/membed.cpp

index b612a1ebdf88f462a24c5184172c4f1c558e101b..1a3cbe0303bc641dd67e7f126d966bbf3688191f 100644 (file)
@@ -92,6 +92,7 @@
 #include "gromacs/topology/block.h"
 #include "gromacs/topology/idef.h"
 #include "gromacs/topology/ifunc.h"
+#include "gromacs/topology/mtop_lookup.h"
 #include "gromacs/topology/mtop_util.h"
 #include "gromacs/topology/topology.h"
 #include "gromacs/utility/basedefinitions.h"
@@ -1753,7 +1754,7 @@ void write_dd_pdb(const char *fn, gmx_int64_t step, const char *title,
     char          fname[STRLEN], buf[22];
     FILE         *out;
     int           i, ii, resnr, c;
-    char         *atomname, *resname;
+    const char   *atomname, *resname;
     real          b;
     gmx_domdec_t *dd;
 
@@ -1769,10 +1770,11 @@ void write_dd_pdb(const char *fn, gmx_int64_t step, const char *title,
 
     fprintf(out, "TITLE     %s\n", title);
     gmx_write_pdb_box(out, dd->bScrewPBC ? epbcSCREW : epbcXYZ, box);
+    int molb = 0;
     for (i = 0; i < natoms; i++)
     {
         ii = dd->gatindex[i];
-        gmx_mtop_atominfo_global(mtop, ii, &atomname, &resnr, &resname);
+        mtopGetAtomAndResidueName(mtop, ii, &molb, &atomname, &resnr, &resname, NULL);
         if (i < dd->comm->nat[ddnatZONE])
         {
             c = 0;
index 5cece422a6e8e294fee875c2480bbcadc6b158a1..91565d67117875cf2991b7fc20432ed0c4826863 100644 (file)
@@ -1,7 +1,7 @@
 /*
  * This file is part of the GROMACS molecular simulation package.
  *
- * Copyright (c) 2006,2007,2008,2009,2010,2012,2013,2014,2015, by the GROMACS development team, led by
+ * Copyright (c) 2006,2007,2008,2009,2010,2012,2013,2014,2015,2016, 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.
@@ -58,7 +58,7 @@
 #include "gromacs/mdlib/gmx_omp_nthreads.h"
 #include "gromacs/mdtypes/commrec.h"
 #include "gromacs/pbcutil/ishift.h"
-#include "gromacs/topology/mtop_util.h"
+#include "gromacs/topology/mtop_lookup.h"
 #include "gromacs/utility/exceptions.h"
 #include "gromacs/utility/fatalerror.h"
 #include "gromacs/utility/gmxassert.h"
@@ -242,51 +242,37 @@ static void atoms_to_settles(gmx_domdec_t *dd,
                              t_ilist *ils_local,
                              ind_req_t *ireq)
 {
-    gmx_ga2la_t            *ga2la;
-    gmx_mtop_atomlookup_t   alook;
-    int                     settle;
-    int                     nral, sa;
-    int                     cg, a, a_gl, a_glsa, a_gls[3], a_locs[3];
-    int                     mb, molnr, a_mol, offset;
-    const gmx_molblock_t   *molb;
-    const t_iatom          *ia1;
-    gmx_bool                a_home[3];
-    int                     nlocal;
-    gmx_bool                bAssign;
+    gmx_ga2la_t *ga2la = dd->ga2la;
+    int          nral  = NRAL(F_SETTLE);
 
-    ga2la  = dd->ga2la;
-
-    alook = gmx_mtop_atomlookup_settle_init(mtop);
-
-    nral = NRAL(F_SETTLE);
-
-    for (cg = cg_start; cg < cg_end; cg++)
+    int          mb    = 0;
+    for (int cg = cg_start; cg < cg_end; cg++)
     {
         if (GET_CGINFO_SETTLE(cginfo[cg]))
         {
-            for (a = dd->cgindex[cg]; a < dd->cgindex[cg+1]; a++)
+            for (int a = dd->cgindex[cg]; a < dd->cgindex[cg+1]; a++)
             {
-                a_gl = dd->gatindex[a];
-
-                gmx_mtop_atomnr_to_molblock_ind(alook, a_gl, &mb, &molnr, &a_mol);
-                molb = &mtop->molblock[mb];
+                int a_gl = dd->gatindex[a];
+                int a_mol;
+                mtopGetMolblockIndex(mtop, a_gl, &mb, NULL, &a_mol);
 
-                settle = at2settle_mt[molb->type][a_mol];
+                const gmx_molblock_t *molb   = &mtop->molblock[mb];
+                int                   settle = at2settle_mt[molb->type][a_mol];
 
                 if (settle >= 0)
                 {
-                    offset = a_gl - a_mol;
+                    int      offset  = a_gl - a_mol;
 
-                    ia1 = mtop->moltype[molb->type].ilist[F_SETTLE].iatoms;
+                    t_iatom *ia1     = mtop->moltype[molb->type].ilist[F_SETTLE].iatoms;
 
-                    bAssign = FALSE;
-                    nlocal  = 0;
-                    for (sa = 0; sa < nral; sa++)
+                    int      a_gls[3], a_locs[3];
+                    gmx_bool bAssign = FALSE;
+                    int      nlocal  = 0;
+                    for (int sa = 0; sa < nral; sa++)
                     {
-                        a_glsa     = offset + ia1[settle*(1+nral)+1+sa];
+                        int a_glsa = offset + ia1[settle*(1+nral)+1+sa];
                         a_gls[sa]  = a_glsa;
-                        a_home[sa] = ga2la_get_home(ga2la, a_glsa, &a_locs[sa]);
-                        if (a_home[sa])
+                        if (ga2la_get_home(ga2la, a_glsa, &a_locs[sa]))
                         {
                             if (nlocal == 0 && a_gl == a_glsa)
                             {
@@ -306,7 +292,7 @@ static void atoms_to_settles(gmx_domdec_t *dd,
 
                         ils_local->iatoms[ils_local->nr++] = ia1[settle*4];
 
-                        for (sa = 0; sa < nral; sa++)
+                        for (int sa = 0; sa < nral; sa++)
                         {
                             if (ga2la_get_home(ga2la, a_gls[sa], &a_locs[sa]))
                             {
@@ -332,8 +318,6 @@ static void atoms_to_settles(gmx_domdec_t *dd,
             }
         }
     }
-
-    gmx_mtop_atomlookup_destroy(alook);
 }
 
 /*! \brief Looks up constraint for the local atoms */
@@ -345,33 +329,28 @@ static void atoms_to_constraints(gmx_domdec_t *dd,
                                  ind_req_t *ireq)
 {
     const t_blocka             *at2con;
-    gmx_ga2la_t                *ga2la;
-    gmx_mtop_atomlookup_t       alook;
     int                         ncon1;
-    gmx_molblock_t             *molb;
     t_iatom                    *ia1, *ia2, *iap;
-    int                         nhome, cg, a, a_gl, a_mol, a_loc, b_lo, offset, mb, molnr, b_mol, i, con, con_offset;
-    gmx_domdec_constraints_t   *dc;
-    gmx_domdec_specat_comm_t   *dcc;
+    int                         a_loc, b_lo, offset, b_mol, i, con, con_offset;
 
-    dc  = dd->constraints;
-    dcc = dd->constraint_comm;
+    gmx_domdec_constraints_t   *dc     = dd->constraints;
+    gmx_domdec_specat_comm_t   *dcc    = dd->constraint_comm;
 
-    ga2la  = dd->ga2la;
+    gmx_ga2la_t                *ga2la  = dd->ga2la;
 
-    alook = gmx_mtop_atomlookup_init(mtop);
-
-    nhome = 0;
-    for (cg = 0; cg < dd->ncg_home; cg++)
+    int mb    = 0;
+    int nhome = 0;
+    for (int cg = 0; cg < dd->ncg_home; cg++)
     {
         if (GET_CGINFO_CONSTR(cginfo[cg]))
         {
-            for (a = dd->cgindex[cg]; a < dd->cgindex[cg+1]; a++)
+            for (int a = dd->cgindex[cg]; a < dd->cgindex[cg+1]; a++)
             {
-                a_gl = dd->gatindex[a];
+                int a_gl = dd->gatindex[a];
+                int molnr, a_mol;
+                mtopGetMolblockIndex(mtop, a_gl, &mb, &molnr, &a_mol);
 
-                gmx_mtop_atomnr_to_molblock_ind(alook, a_gl, &mb, &molnr, &a_mol);
-                molb = &mtop->molblock[mb];
+                const gmx_molblock_t *molb = &mtop->molblock[mb];
 
                 ncon1 = mtop->moltype[molb->type].ilist[F_CONSTR].nr/NRAL(F_SETTLE);
 
@@ -443,8 +422,6 @@ static void atoms_to_constraints(gmx_domdec_t *dd,
         }
     }
 
-    gmx_mtop_atomlookup_destroy(alook);
-
     if (debug)
     {
         fprintf(debug,
index 6dec771bddba3027fffbfb9edc7701da8164cc50..3586977a7e90e3bddc9f228fd72c1a6786673e4a 100644 (file)
@@ -61,7 +61,7 @@
 #include "gromacs/mdtypes/inputrec.h"
 #include "gromacs/mdtypes/md_enums.h"
 #include "gromacs/pbcutil/pbc.h"
-#include "gromacs/topology/mtop_util.h"
+#include "gromacs/topology/mtop_lookup.h"
 #include "gromacs/topology/topology.h"
 #include "gromacs/utility/cstringutil.h"
 #include "gromacs/utility/fatalerror.h"
@@ -1367,22 +1367,20 @@ static void init_edi(const gmx_mtop_t *mtop, t_edpar *edi)
     int                   i;
     real                  totalmass = 0.0;
     rvec                  com;
-    gmx_mtop_atomlookup_t alook = NULL;
-    t_atom               *atom;
 
     /* NOTE Init_edi is executed on the master process only
      * The initialized data sets are then transmitted to the
      * other nodes in broadcast_ed_data */
 
-    alook = gmx_mtop_atomlookup_init(mtop);
-
     /* evaluate masses (reference structure) */
     snew(edi->sref.m, edi->sref.nr);
+    int molb = 0;
     for (i = 0; i < edi->sref.nr; i++)
     {
         if (edi->fitmas)
         {
-            gmx_mtop_atomnr_to_atom(alook, edi->sref.anrs[i], &atom);
+            const t_atom *atom;
+            mtopGetAtomParameters(mtop, edi->sref.anrs[i], &molb, &atom);
             edi->sref.m[i] = atom->m;
         }
         else
@@ -1410,7 +1408,8 @@ static void init_edi(const gmx_mtop_t *mtop, t_edpar *edi)
     snew(edi->sav.m, edi->sav.nr );
     for (i = 0; i < edi->sav.nr; i++)
     {
-        gmx_mtop_atomnr_to_atom(alook, edi->sav.anrs[i], &atom);
+        const t_atom *atom;
+        mtopGetAtomParameters(mtop, edi->sav.anrs[i], &molb, &atom);
         edi->sav.m[i] = atom->m;
         if (edi->pcamas)
         {
@@ -1432,8 +1431,6 @@ static void init_edi(const gmx_mtop_t *mtop, t_edpar *edi)
         }
     }
 
-    gmx_mtop_atomlookup_destroy(alook);
-
     /* put reference structure in origin */
     get_center(edi->sref.x, edi->sref.m, edi->sref.nr, com);
     com[XX] = -com[XX];
index 46f5cd5b093ef6d0bd614d80c5dc665db886f81a..0aa130216cb5d1c1884c50ebf7b4d84795ca2459 100644 (file)
@@ -507,7 +507,7 @@ void write_hconf_mtop(FILE *out, const char *title, gmx_mtop_t *mtop,
 {
     int                     i, resnr;
     gmx_mtop_atomloop_all_t aloop;
-    t_atom                 *atom;
+    const t_atom           *atom;
     char                   *atomname, *resname;
 
     fprintf(out, "%s\n", (title && title[0]) ? title : gmx::bromacs().c_str());
index 0bf6fb501eecf4d0387853e1c7f19528c708dbbe..8b14533fe059352ee97c15ab9c167a08fea7a85c 100644 (file)
@@ -3,7 +3,7 @@
  *
  * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
  * Copyright (c) 2001-2007, The GROMACS development team.
- * Copyright (c) 2013,2014,2015, by the GROMACS development team, led by
+ * Copyright (c) 2013,2014,2015,2016, 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.
@@ -52,7 +52,7 @@
 #include "gromacs/math/vec.h"
 #include "gromacs/pbcutil/pbc.h"
 #include "gromacs/topology/index.h"
-#include "gromacs/topology/mtop_util.h"
+#include "gromacs/topology/mtop_lookup.h"
 #include "gromacs/topology/topology.h"
 #include "gromacs/trajectory/trajectoryframe.h"
 #include "gromacs/utility/arraysize.h"
@@ -85,8 +85,6 @@ static void clust_size(const char *ndx, const char *trx, const char *xpm,
     gmx_mtop_t           *mtop = NULL;
     int                   ePBC = -1;
     t_block              *mols = NULL;
-    gmx_mtop_atomlookup_t alook;
-    t_atom               *atom;
     int                   ii, jj;
     real                  temp, tfac;
     /* Cluster size distribution (matrix) */
@@ -157,8 +155,6 @@ static void clust_size(const char *ndx, const char *trx, const char *xpm,
         rd_index(ndx, 1, &nindex, &index, &gname);
     }
 
-    alook = gmx_mtop_atomlookup_init(mtop);
-
     snew(clust_index, nindex);
     snew(clust_size, nindex);
     cut2   = cut*cut;
@@ -171,6 +167,7 @@ static void clust_size(const char *ndx, const char *trx, const char *xpm,
     }
     max_clust_size = 1;
     max_clust_ind  = -1;
+    int molb       = 0;
     do
     {
         if ((nskip == 0) || ((nskip > 0) && ((nframe % nskip) == 0)))
@@ -324,7 +321,8 @@ static void clust_size(const char *ndx, const char *trx, const char *xpm,
                         if (clust_index[i] == max_clust_ind)
                         {
                             ai    = index[i];
-                            gmx_mtop_atomnr_to_atom(alook, ai, &atom);
+                            const t_atom       *atom;
+                            mtopGetAtomParameters(mtop, ai, &molb, &atom);
                             ekin += 0.5*atom->m*iprod(v[ai], v[ai]);
                         }
                     }
@@ -342,8 +340,6 @@ static void clust_size(const char *ndx, const char *trx, const char *xpm,
     xvgrclose(hp);
     xvgrclose(tp);
 
-    gmx_mtop_atomlookup_destroy(alook);
-
     if (max_clust_ind >= 0)
     {
         fp = gmx_ffopen(mcn, "w");
index 55e8b0b628886956f67ba93af0b050ee65cf16ba..2ecc85a4b14643ec0407aec4f62433d69c603510 100644 (file)
@@ -60,6 +60,7 @@
 #include "gromacs/mdtypes/inputrec.h"
 #include "gromacs/mdtypes/md_enums.h"
 #include "gromacs/topology/ifunc.h"
+#include "gromacs/topology/mtop_lookup.h"
 #include "gromacs/topology/mtop_util.h"
 #include "gromacs/topology/topology.h"
 #include "gromacs/utility/arraysize.h"
@@ -2013,7 +2014,7 @@ int gmx_energy(int argc, char *argv[])
     gmx_bool          *bIsEner = NULL;
     char             **pairleg, **odtleg, **otenleg;
     char             **leg = NULL;
-    char              *anm_j, *anm_k, *resnm_j, *resnm_k;
+    const char        *anm_j, *anm_k, *resnm_j, *resnm_k;
     int                resnr_j, resnr_k;
     const char        *orinst_sub = "@ subtitle \"instantaneous\"\n";
     char               buf[256];
@@ -2476,13 +2477,14 @@ int gmx_energy(int argc, char *argv[])
                               ndisre, top->idef.il[F_DISRES].nr/3);
                 }
                 snew(pairleg, ndisre);
+                int molb = 0;
                 for (i = 0; i < ndisre; i++)
                 {
                     snew(pairleg[i], 30);
                     j = fa[3*i+1];
                     k = fa[3*i+2];
-                    gmx_mtop_atominfo_global(&mtop, j, &anm_j, &resnr_j, &resnm_j);
-                    gmx_mtop_atominfo_global(&mtop, k, &anm_k, &resnr_k, &resnm_k);
+                    mtopGetAtomAndResidueName(&mtop, j, &molb, &anm_j, &resnr_j, &resnm_j, NULL);
+                    mtopGetAtomAndResidueName(&mtop, k, &molb, &anm_k, &resnr_k, &resnm_k, NULL);
                     sprintf(pairleg[i], "%d %s %d %s (%d)",
                             resnr_j, anm_j, resnr_k, anm_k,
                             ip[fa[3*i]].disres.label);
index 4025d73b7865b80fa7dceaa79b21db7782a8f96d..9dece7a6bdafea0d38be6bf0d04de1ee7d54bbb3 100644 (file)
@@ -818,7 +818,6 @@ static int prepare_x_q(real *q[], rvec *x[], const gmx_mtop_t *mtop, const rvec
     int                     i;
     int                     nq; /* number of charged particles */
     gmx_mtop_atomloop_all_t aloop;
-    t_atom                 *atom;
 
 
     if (MASTER(cr))
@@ -828,7 +827,7 @@ static int prepare_x_q(real *q[], rvec *x[], const gmx_mtop_t *mtop, const rvec
         nq = 0;
 
         aloop = gmx_mtop_atomloop_all_init(mtop);
-
+        const t_atom *atom;
         while (gmx_mtop_atomloop_all_next(aloop, &i, &atom))
         {
             if (is_charge(atom->q))
index 6f25a14da8bc8516f72d0df0d6055c4bdbfdce17..350e1ad95176a82c4f356cf4e4cd1e7e346be393 100644 (file)
@@ -57,7 +57,7 @@ static void low_mspeed(real tempi,
     real                                    boltz, sd;
     real                                    ekin, temp, mass, scal;
     gmx_mtop_atomloop_all_t                 aloop;
-    t_atom                                 *atom;
+    const t_atom                           *atom;
     gmx::TabulatedNormalDistribution<real>  normalDist;
 
     boltz = BOLTZ*tempi;
index ca768804ea8738da69b7b4074fe985482cb3932c..a6ce086301db5ffc109a38cc7a86f0f1bde96d81 100644 (file)
@@ -363,7 +363,7 @@ static void check_bonds_timestep(gmx_mtop_t *mtop, double dt, warninp_t wi)
 static void check_vel(gmx_mtop_t *mtop, rvec v[])
 {
     gmx_mtop_atomloop_all_t aloop;
-    t_atom                 *atom;
+    const t_atom           *atom;
     int                     a;
 
     aloop = gmx_mtop_atomloop_all_init(mtop);
@@ -383,7 +383,7 @@ static void check_shells_inputrec(gmx_mtop_t *mtop,
                                   warninp_t   wi)
 {
     gmx_mtop_atomloop_all_t aloop;
-    t_atom                 *atom;
+    const t_atom           *atom;
     int                     a, nshells = 0;
     char                    warn_buf[STRLEN];
 
@@ -670,7 +670,7 @@ new_status(const char *topfile, const char *topppfile, const char *confin,
     {
         real                   *mass;
         gmx_mtop_atomloop_all_t aloop;
-        t_atom                 *atom;
+        const t_atom           *atom;
 
         snew(mass, state->natoms);
         aloop = gmx_mtop_atomloop_all_init(sys);
@@ -1326,7 +1326,7 @@ static real calc_temp(const gmx_mtop_t *mtop,
                       rvec             *v)
 {
     gmx_mtop_atomloop_all_t aloop;
-    t_atom                 *atom;
+    const t_atom           *atom;
     int                     a;
 
     double                  sum_mv2 = 0;
index 6f1c559f8c2938fb420e5db3c45b7a361b02a017..2fd9c2e3537263d87d8721f3ecc7a13495a5c9c7 100644 (file)
@@ -2735,7 +2735,6 @@ static void calc_nrdf(gmx_mtop_t *mtop, t_inputrec *ir, char **gnames)
     double                 *nrdf_tc, *nrdf_vcm, nrdf_uc, *nrdf_vcm_sub;
     ivec                   *dof_vcm;
     gmx_mtop_atomloop_all_t aloop;
-    t_atom                 *atom;
     int                     mb, mol, ftype, as;
     gmx_molblock_t         *molb;
     gmx_moltype_t          *molt;
@@ -2776,6 +2775,7 @@ static void calc_nrdf(gmx_mtop_t *mtop, t_inputrec *ir, char **gnames)
 
     snew(nrdf2, natoms);
     aloop = gmx_mtop_atomloop_all_init(mtop);
+    const t_atom *atom;
     while (gmx_mtop_atomloop_all_next(aloop, &i, &atom))
     {
         nrdf2[i] = 0;
@@ -4007,7 +4007,6 @@ void triple_check(const char *mdparin, t_inputrec *ir, gmx_mtop_t *sys,
     rvec                      acc;
     gmx_mtop_atomloop_block_t aloopb;
     gmx_mtop_atomloop_all_t   aloop;
-    t_atom                   *atom;
     ivec                      AbsRef;
     char                      warn_buf[STRLEN];
 
@@ -4104,6 +4103,7 @@ void triple_check(const char *mdparin, t_inputrec *ir, gmx_mtop_t *sys,
 
     bCharge = FALSE;
     aloopb  = gmx_mtop_atomloop_block_init(sys);
+    const t_atom *atom;
     while (gmx_mtop_atomloop_block_next(aloopb, &atom, &nmol))
     {
         if (atom->q != 0 || atom->qB != 0)
@@ -4172,6 +4172,7 @@ void triple_check(const char *mdparin, t_inputrec *ir, gmx_mtop_t *sys,
         clear_rvec(acc);
         snew(mgrp, sys->groups.grps[egcACC].nr);
         aloop = gmx_mtop_atomloop_all_init(sys);
+        const t_atom *atom;
         while (gmx_mtop_atomloop_all_next(aloop, &i, &atom))
         {
             mgrp[ggrpnr(&sys->groups, egcACC, i)] += atom->m;
index 44ac6e3f25462fb2f79388d452bc82db2a9caece..93be131e79b92208f5dbdaad910b82511a9a2fd4 100644 (file)
@@ -76,7 +76,6 @@ void init_orires(FILE *fplog, const gmx_mtop_t *mtop,
     gmx_mtop_ilistloop_t    iloop;
     t_ilist                *il;
     gmx_mtop_atomloop_all_t aloop;
-    t_atom                 *atom;
     const gmx_multisim_t   *ms;
 
     od->nr = gmx_mtop_ftype_count(mtop, F_ORIRES);
@@ -209,6 +208,7 @@ void init_orires(FILE *fplog, const gmx_mtop_t *mtop,
     mtot  = 0.0;
     j     = 0;
     aloop = gmx_mtop_atomloop_all_init(mtop);
+    const t_atom *atom;
     while (gmx_mtop_atomloop_all_next(aloop, &i, &atom))
     {
         if (mtop->groups.grpnr[egcORFIT] == NULL ||
index e4cc11daafe4aa9d14c892ef6cbc642dc32e0294..c33b16b3279d455d9855d00c62aaf0647a383aaf 100644 (file)
@@ -762,16 +762,12 @@ static void bc_moltype(const t_commrec *cr, t_symtab *symtab,
 
 static void bc_molblock(const t_commrec *cr, gmx_molblock_t *molb)
 {
-    block_bc(cr, molb->type);
-    block_bc(cr, molb->nmol);
-    block_bc(cr, molb->natoms_mol);
-    block_bc(cr, molb->nposres_xA);
+    block_bc(cr, *molb);
     if (molb->nposres_xA > 0)
     {
         snew_bc(cr, molb->posres_xA, molb->nposres_xA);
         nblock_bc(cr, molb->nposres_xA*DIM, molb->posres_xA[0]);
     }
-    block_bc(cr, molb->nposres_xB);
     if (molb->nposres_xB > 0)
     {
         snew_bc(cr, molb->posres_xB, molb->nposres_xB);
index 5ed9837feaf0adef94304e07783cec665c3f3945..a4b9f1e9a158d7b2c24d93cfdbc74a81485afba2 100644 (file)
@@ -64,6 +64,7 @@
 #include "gromacs/pulling/pull.h"
 #include "gromacs/topology/block.h"
 #include "gromacs/topology/invblock.h"
+#include "gromacs/topology/mtop_lookup.h"
 #include "gromacs/topology/mtop_util.h"
 #include "gromacs/utility/exceptions.h"
 #include "gromacs/utility/fatalerror.h"
@@ -183,7 +184,7 @@ static void write_constr_pdb(const char *fn, const char *title,
     FILE         *out;
     int           dd_ac0 = 0, dd_ac1 = 0, i, ii, resnr;
     gmx_domdec_t *dd;
-    char         *anm, *resnm;
+    const char   *anm, *resnm;
 
     dd = NULL;
     if (DOMAINDECOMP(cr))
@@ -207,6 +208,7 @@ static void write_constr_pdb(const char *fn, const char *title,
 
     fprintf(out, "TITLE     %s\n", title);
     gmx_write_pdb_box(out, -1, box);
+    int molb = 0;
     for (i = start; i < start+homenr; i++)
     {
         if (dd != NULL)
@@ -221,7 +223,7 @@ static void write_constr_pdb(const char *fn, const char *title,
         {
             ii = i;
         }
-        gmx_mtop_atominfo_global(mtop, ii, &anm, &resnr, &resnm);
+        mtopGetAtomAndResidueName(mtop, ii, &molb, &anm, &resnr, &resnm, NULL);
         gmx_fprintf_pdb_atomline(out, epdbATOM, ii+1, anm, ' ', resnm, ' ', resnr, ' ',
                                  10*x[i][XX], 10*x[i][YY], 10*x[i][ZZ], 1.0, 0.0, "");
     }
index a035bc79ff89797cc3ac6444876a7cded6917eb5..40abcfcc51d2f38fbd4bca31f3df0651a145820d 100644 (file)
@@ -45,6 +45,7 @@
 #include "gromacs/mdlib/qmmm.h"
 #include "gromacs/mdtypes/inputrec.h"
 #include "gromacs/mdtypes/md_enums.h"
+#include "gromacs/topology/mtop_lookup.h"
 #include "gromacs/topology/mtop_util.h"
 #include "gromacs/topology/topology.h"
 #include "gromacs/utility/exceptions.h"
@@ -56,7 +57,7 @@ t_mdatoms *init_mdatoms(FILE *fp, const gmx_mtop_t *mtop, gmx_bool bFreeEnergy)
 {
     int                     a;
     double                  tmA, tmB;
-    t_atom                 *atom;
+    const t_atom           *atom;
     t_mdatoms              *md;
     gmx_mtop_atomloop_all_t aloop;
 
@@ -117,8 +118,6 @@ void atoms2md(const gmx_mtop_t *mtop, const t_inputrec *ir,
               t_mdatoms *md)
 {
     gmx_bool              bLJPME;
-    gmx_mtop_atomlookup_t alook;
-    int                   i;
     const t_grpopts      *opts;
     const gmx_groups_t   *groups;
     int                   nthreads gmx_unused;
@@ -219,19 +218,18 @@ void atoms2md(const gmx_mtop_t *mtop, const t_inputrec *ir,
         }
     }
 
-    alook = gmx_mtop_atomlookup_init(mtop);
+    int molb = 0;
 
     // cppcheck-suppress unreadVariable
     nthreads = gmx_omp_nthreads_get(emntDefault);
-#pragma omp parallel for num_threads(nthreads) schedule(static)
-    for (i = 0; i < md->nr; i++)
+#pragma omp parallel for num_threads(nthreads) schedule(static) firstprivate(molb)
+    for (int i = 0; i < md->nr; i++)
     {
         try
         {
             int      g, ag;
             real     mA, mB, fac;
             real     c6, c12;
-            t_atom  *atom;
 
             if (index == NULL)
             {
@@ -239,9 +237,10 @@ void atoms2md(const gmx_mtop_t *mtop, const t_inputrec *ir,
             }
             else
             {
-                ag   = index[i];
+                ag = index[i];
             }
-            gmx_mtop_atomnr_to_atom(alook, ag, &atom);
+            const t_atom *atom;
+            mtopGetAtomParameters(mtop, ag, &molb, &atom);
 
             if (md->cFREEZE)
             {
@@ -399,8 +398,6 @@ void atoms2md(const gmx_mtop_t *mtop, const t_inputrec *ir,
         GMX_CATCH_ALL_AND_EXIT_WITH_FATAL_ERROR;
     }
 
-    gmx_mtop_atomlookup_destroy(alook);
-
     md->homenr = homenr;
     md->lambda = 0;
 }
index 98b3a3324411081f70daec3b438a4815b2a95234..ab4a4c0d0b123f3a27b66b680706590ef521c403 100644 (file)
@@ -65,6 +65,7 @@
 #include "gromacs/mdtypes/nblist.h"
 #include "gromacs/pbcutil/ishift.h"
 #include "gromacs/pbcutil/pbc.h"
+#include "gromacs/topology/mtop_lookup.h"
 #include "gromacs/topology/mtop_util.h"
 #include "gromacs/topology/topology.h"
 #include "gromacs/utility/fatalerror.h"
@@ -315,32 +316,26 @@ static void init_QMrec(int grpnr, t_QMrec *qm, int nr, int *atomarray,
 {
     /* fills the t_QMrec struct of QM group grpnr
      */
-    int                   i;
-    gmx_mtop_atomlookup_t alook;
-    t_atom               *atom;
-
 
     qm->nrQMatoms = nr;
     snew(qm->xQM, nr);
     snew(qm->indexQM, nr);
     snew(qm->shiftQM, nr); /* the shifts */
-    for (i = 0; i < nr; i++)
+    for (int i = 0; i < nr; i++)
     {
         qm->indexQM[i] = atomarray[i];
     }
 
-    alook = gmx_mtop_atomlookup_init(mtop);
-
     snew(qm->atomicnumberQM, nr);
-    for (i = 0; i < qm->nrQMatoms; i++)
+    int molb = 0;
+    for (int i = 0; i < qm->nrQMatoms; i++)
     {
-        gmx_mtop_atomnr_to_atom(alook, qm->indexQM[i], &atom);
+        const t_atom *atom;
+        mtopGetAtomParameters(mtop, qm->indexQM[i], &molb, &atom);
         qm->nelectrons       += mtop->atomtypes.atomnumber[atom->type];
         qm->atomicnumberQM[i] = mtop->atomtypes.atomnumber[atom->type];
     }
 
-    gmx_mtop_atomlookup_destroy(alook);
-
     qm->QMcharge       = ir->opts.QMcharge[grpnr];
     qm->multiplicity   = ir->opts.QMmult[grpnr];
     qm->nelectrons    -= ir->opts.QMcharge[grpnr];
@@ -468,11 +463,9 @@ void init_QMMMrec(t_commrec  *cr,
     t_iatom                 *iatoms;
     real                     c12au, c6au;
     gmx_mtop_atomloop_all_t  aloop;
-    t_atom                  *atom;
     gmx_mtop_ilistloop_all_t iloop;
     int                      a_offset;
     t_ilist                 *ilist_mol;
-    gmx_mtop_atomlookup_t    alook;
 
     if (ir->cutoff_scheme != ecutsGROUP)
     {
@@ -535,6 +528,7 @@ void init_QMMMrec(t_commrec  *cr,
     {
         /* new layer */
         aloop = gmx_mtop_atomloop_all_init(mtop);
+        const t_atom *atom;
         while (gmx_mtop_atomloop_all_next(aloop, &i, &atom))
         {
             if (qm_nr >= qm_max)
@@ -661,15 +655,17 @@ void init_QMMMrec(t_commrec  *cr,
 
         /* standard QMMM, all layers are merged together so there is one QM
          * subsystem and one MM subsystem.
-         * Also we set the charges to zero in the md->charge arrays to prevent
-         * the innerloops from doubly counting the electostatic QM MM interaction
+         * Also we set the charges to zero in mtop to prevent the innerloops
+         * from doubly counting the electostatic QM MM interaction
+         * TODO: Consider doing this in grompp instead.
          */
 
-        alook = gmx_mtop_atomlookup_init(mtop);
-
+        int molb = 0;
         for (k = 0; k < qm_nr; k++)
         {
-            gmx_mtop_atomnr_to_atom(alook, qm_arr[k], &atom);
+            int     indexInMolecule;
+            mtopGetMolblockIndex(mtop, qm_arr[k], &molb, NULL, &indexInMolecule);
+            t_atom *atom = &mtop->moltype[mtop->molblock[molb].type].atoms.atom[indexInMolecule];
             atom->q  = 0.0;
             atom->qB = 0.0;
         }
@@ -681,7 +677,8 @@ void init_QMMMrec(t_commrec  *cr,
         {
             for (i = 0; i < qm_nr; i++)
             {
-                gmx_mtop_atomnr_to_atom(alook, qm_arr[i], &atom);
+                const t_atom *atom;
+                mtopGetAtomParameters(mtop, qm_arr[i], &molb, &atom);
                 /* nbfp now includes the 6.0/12.0 derivative prefactors */
                 qr->qm[0]->c6[i]  =  C6(fr->nbfp, mtop->ffparams.atnr, atom->type, atom->type)/c6au/6.0;
                 qr->qm[0]->c12[i] = C12(fr->nbfp, mtop->ffparams.atnr, atom->type, atom->type)/c12au/12.0;
@@ -692,9 +689,10 @@ void init_QMMMrec(t_commrec  *cr,
          */
         for (i = 0; i < qm_nr; i++)
         {
-            gmx_mtop_atomnr_to_ilist(alook, qm_arr[i], &ilist_mol, &a_offset);
-            nrvsite2 = ilist_mol[F_VSITE2].nr;
-            iatoms   = ilist_mol[F_VSITE2].iatoms;
+            mtopGetMolblockIndex(mtop, qm_arr[i], &molb, NULL, &a_offset);
+            ilist_mol = mtop->moltype[mtop->molblock[molb].type].ilist;
+            nrvsite2  = ilist_mol[F_VSITE2].nr;
+            iatoms    = ilist_mol[F_VSITE2].iatoms;
 
             for (k = 0; k < nrvsite2; k += 4)
             {
@@ -722,8 +720,6 @@ void init_QMMMrec(t_commrec  *cr,
             }
         }
 
-        gmx_mtop_atomlookup_destroy(alook);
-
         /* MM rec creation */
         mm               = mk_MMrec();
         mm->scalefactor  = ir->scalefactor;
index 1a86663027320d84f56083715f1eca22b44001e0..2db79b7e8797eab65d70c8891e43b7cc49ce53b9 100644 (file)
@@ -64,6 +64,7 @@
 #include "gromacs/mdtypes/md_enums.h"
 #include "gromacs/pbcutil/mshift.h"
 #include "gromacs/pbcutil/pbc.h"
+#include "gromacs/topology/mtop_lookup.h"
 #include "gromacs/topology/mtop_util.h"
 #include "gromacs/utility/arraysize.h"
 #include "gromacs/utility/cstringutil.h"
@@ -150,13 +151,7 @@ static void predict_shells(FILE *fplog, rvec x[], rvec v[], real dt,
     int                   i, m, s1, n1, n2, n3;
     real                  dt_1, fudge, tm, m1, m2, m3;
     rvec                 *ptr;
-    gmx_mtop_atomlookup_t alook = NULL;
-    t_atom               *atom;
-
-    if (mass == NULL)
-    {
-        alook = gmx_mtop_atomlookup_init(mtop);
-    }
+    const t_atom         *atom;
 
     /* We introduce a fudge factor for performance reasons: with this choice
      * the initial force on the shells is about a factor of two lower than
@@ -179,6 +174,7 @@ static void predict_shells(FILE *fplog, rvec x[], rvec v[], real dt,
         dt_1 = fudge*dt;
     }
 
+    int molb = 0;
     for (i = 0; (i < ns); i++)
     {
         s1 = s[i].shell;
@@ -206,9 +202,9 @@ static void predict_shells(FILE *fplog, rvec x[], rvec v[], real dt,
                 else
                 {
                     /* Not the correct masses with FE, but it is just a prediction... */
-                    gmx_mtop_atomnr_to_atom(alook, n1, &atom);
+                    mtopGetAtomParameters(mtop, n1, &molb, &atom);
                     m1 = atom->m;
-                    gmx_mtop_atomnr_to_atom(alook, n2, &atom);
+                    mtopGetAtomParameters(mtop, n2, &molb, &atom);
                     m2 = atom->m;
                 }
                 tm = dt_1/(m1+m2);
@@ -230,11 +226,11 @@ static void predict_shells(FILE *fplog, rvec x[], rvec v[], real dt,
                 else
                 {
                     /* Not the correct masses with FE, but it is just a prediction... */
-                    gmx_mtop_atomnr_to_atom(alook, n1, &atom);
+                    mtopGetAtomParameters(mtop, n1, &molb, &atom);
                     m1 = atom->m;
-                    gmx_mtop_atomnr_to_atom(alook, n2, &atom);
+                    mtopGetAtomParameters(mtop, n2, &molb, &atom);
                     m2 = atom->m;
-                    gmx_mtop_atomnr_to_atom(alook, n3, &atom);
+                    mtopGetAtomParameters(mtop, n3, &molb, &atom);
                     m3 = atom->m;
                 }
                 tm = dt_1/(m1+m2+m3);
@@ -247,11 +243,6 @@ static void predict_shells(FILE *fplog, rvec x[], rvec v[], real dt,
                 gmx_fatal(FARGS, "Shell %d has %d nuclei!", i, s[i].nnucl);
         }
     }
-
-    if (mass == NULL)
-    {
-        gmx_mtop_atomlookup_destroy(alook);
-    }
 }
 
 /*! \brief Count the different particle types in a system
@@ -274,7 +265,7 @@ static std::array<int, eptNR> countPtypes(FILE       *fplog,
 
     gmx_mtop_atomloop_block_t  aloopb = gmx_mtop_atomloop_block_init(mtop);
     int                        nmol;
-    t_atom                    *atom;
+    const t_atom              *atom;
     while (gmx_mtop_atomloop_block_next(aloopb, &atom, &nmol))
     {
         switch (atom->ptype)
@@ -313,7 +304,7 @@ gmx_shellfc_t *init_shell_flexcon(FILE *fplog,
     gmx_shellfc_t            *shfc;
     t_shell                  *shell;
     int                      *shell_index = NULL, *at2cg;
-    t_atom                   *atom;
+    const t_atom             *atom;
 
     int                       ns, nshell, nsi;
     int                       i, j, type, mb, a_offset, cg, mol, ftype, nra;
index 78ece843ca9321ea3707a4833b826577dbefeb37..0d7712be6dbee67f3a1eb6fcc7f201bdeb0d3b75 100644 (file)
@@ -74,7 +74,7 @@ static void init_grpstat(gmx_mtop_t *mtop, int ngacc, t_grp_acc gstat[])
     gmx_groups_t           *groups;
     gmx_mtop_atomloop_all_t aloop;
     int                     i, grp;
-    t_atom                 *atom;
+    const t_atom           *atom;
 
     if (ngacc > 0)
     {
index 0b53424b2913d38de4da0c0f3eeed0b93d7ac979..332535531631ac283ad8ebc561b7c24cfff97500 100644 (file)
@@ -66,7 +66,7 @@
 #include "gromacs/mdtypes/mdatom.h"
 #include "gromacs/pbcutil/pbc.h"
 #include "gromacs/pulling/pull_internal.h"
-#include "gromacs/topology/mtop_util.h"
+#include "gromacs/topology/mtop_lookup.h"
 #include "gromacs/topology/topology.h"
 #include "gromacs/utility/cstringutil.h"
 #include "gromacs/utility/exceptions.h"
@@ -1920,8 +1920,6 @@ static void init_pull_group_index(FILE *fplog, t_commrec *cr,
     real                  m, w, mbd;
     double                tmass, wmass, wwmass;
     const gmx_groups_t   *groups;
-    gmx_mtop_atomlookup_t alook;
-    t_atom               *atom;
 
     if (EI_ENERGY_MINIMIZATION(ir->eI) || ir->eI == eiBD)
     {
@@ -1959,16 +1957,16 @@ static void init_pull_group_index(FILE *fplog, t_commrec *cr,
 
     groups = &mtop->groups;
 
-    alook = gmx_mtop_atomlookup_init(mtop);
-
-    nfrozen = 0;
-    tmass   = 0;
-    wmass   = 0;
-    wwmass  = 0;
+    nfrozen  = 0;
+    tmass    = 0;
+    wmass    = 0;
+    wwmass   = 0;
+    int molb = 0;
     for (i = 0; i < pg->params.nat; i++)
     {
         ii = pg->params.ind[i];
-        gmx_mtop_atomnr_to_atom(alook, ii, &atom);
+        const t_atom *atom;
+        mtopGetAtomParameters(mtop, ii, &molb, &atom);
         if (bConstraint && ir->opts.nFreeze)
         {
             for (d = 0; d < DIM; d++)
@@ -2029,8 +2027,6 @@ static void init_pull_group_index(FILE *fplog, t_commrec *cr,
         wwmass += m*w*w;
     }
 
-    gmx_mtop_atomlookup_destroy(alook);
-
     if (wmass == 0)
     {
         /* We can have single atom groups with zero mass with potential pulling
index a9e1bac40dedaa320789d02187251c08f1621a34..209029dc154d91d707a0a0595bd038d058683d4c 100644 (file)
@@ -65,7 +65,7 @@
 #include "gromacs/pbcutil/pbc.h"
 #include "gromacs/timing/cyclecounter.h"
 #include "gromacs/timing/wallcycle.h"
-#include "gromacs/topology/mtop_util.h"
+#include "gromacs/topology/mtop_lookup.h"
 #include "gromacs/utility/fatalerror.h"
 #include "gromacs/utility/pleasecite.h"
 #include "gromacs/utility/qsort_threadsafe.h"
@@ -3390,10 +3390,8 @@ static void init_rot_group(FILE *fplog, t_commrec *cr, int g, t_rotgrp *rotg,
     int                   i, ii;
     rvec                  coord, xref, *xdum;
     gmx_bool              bFlex, bColl;
-    t_atom               *atom;
     gmx_enfrotgrp_t       erg; /* Pointer to enforced rotation group data */
     int                   ref_firstindex, ref_lastindex;
-    gmx_mtop_atomlookup_t alook = NULL;
     real                  mass, totalmass;
     real                  start = 0.0;
     double                t_start;
@@ -3461,10 +3459,6 @@ static void init_rot_group(FILE *fplog, t_commrec *cr, int g, t_rotgrp *rotg,
 
     /* Copy the masses so that the center can be determined. For all types of
      * enforced rotation, we store the masses in the erg->mc array. */
-    if (rotg->bMassW)
-    {
-        alook = gmx_mtop_atomlookup_init(mtop);
-    }
     snew(erg->mc, rotg->nat);
     if (bFlex)
     {
@@ -3475,11 +3469,13 @@ static void init_rot_group(FILE *fplog, t_commrec *cr, int g, t_rotgrp *rotg,
         snew(erg->m_loc, rotg->nat);
     }
     totalmass = 0.0;
+    int molb  = 0;
     for (i = 0; i < rotg->nat; i++)
     {
         if (rotg->bMassW)
         {
-            gmx_mtop_atomnr_to_atom(alook, rotg->ind[i], &atom);
+            const t_atom *atom;
+            mtopGetAtomParameters(mtop, rotg->ind[i], &molb, &atom);
             mass = atom->m;
         }
         else
@@ -3491,11 +3487,6 @@ static void init_rot_group(FILE *fplog, t_commrec *cr, int g, t_rotgrp *rotg,
     }
     erg->invmass = 1.0/totalmass;
 
-    if (rotg->bMassW)
-    {
-        gmx_mtop_atomlookup_destroy(alook);
-    }
-
     /* Set xc_ref_center for any rotation potential */
     if ((rotg->eType == erotgISO) || (rotg->eType == erotgPM) || (rotg->eType == erotgRM) || (rotg->eType == erotgRM2))
     {
index ad1e8372517942b22d69ad1c8790fa57a98046ce..e9e1536d97e99b6d6ed6c781bd624df4ce43e05f 100644 (file)
@@ -62,7 +62,7 @@
 #include "gromacs/mdtypes/md_enums.h"
 #include "gromacs/pbcutil/pbc.h"
 #include "gromacs/timing/wallcycle.h"
-#include "gromacs/topology/mtop_util.h"
+#include "gromacs/topology/mtop_lookup.h"
 #include "gromacs/topology/topology.h"
 #include "gromacs/utility/cstringutil.h"
 #include "gromacs/utility/fatalerror.h"
@@ -975,17 +975,16 @@ static int get_group_apm_check(
         int                         igroup,
         t_swap                     *s,
         gmx_bool                    bVerbose,
-        const gmx_mtop_atomlookup_t alook,
         gmx_mtop_t                 *mtop)
 {
-    int        molb, molnr, atnr_mol;
     t_swapgrp *g   = &s->group[igroup];
     int       *ind = s->group[igroup].ind;
     int        nat = s->group[igroup].nat;
 
     /* Determine the number of solvent atoms per solvent molecule from the
      * first solvent atom: */
-    gmx_mtop_atomnr_to_molblock_ind(alook, ind[0], &molb, &molnr, &atnr_mol);
+    int molb = 0;
+    mtopGetMolblockIndex(mtop, ind[0], &molb, NULL, NULL);
     int apm = mtop->molblock[molb].natoms_mol;
 
     if (bVerbose)
@@ -997,7 +996,7 @@ static int get_group_apm_check(
     /* Check whether this is also true for all other solvent atoms */
     for (int i = 1; i < nat; i++)
     {
-        gmx_mtop_atomnr_to_molblock_ind(alook, ind[i], &molb, &molnr, &atnr_mol);
+        mtopGetMolblockIndex(mtop, ind[i], &molb, NULL, NULL);
         if (apm != mtop->molblock[molb].natoms_mol)
         {
             gmx_fatal(FARGS, "Not all molecules of swap group %d consist of %d atoms.",
@@ -1388,8 +1387,6 @@ void convertOldToNewGroupFormat(
         gmx_bool      bVerbose,
         t_commrec    *cr)
 {
-    t_atom                *atom;
-    gmx_mtop_atomlookup_t  alook = gmx_mtop_atomlookup_init(mtop);
     t_swapGroup           *g     = &sc->grp[3];
 
     /* Loop through the atom indices of group #3 (anions) and put all indices
@@ -1402,9 +1399,11 @@ void convertOldToNewGroupFormat(
     snew(indAnions, g->nat);
     snew(indCations, g->nat);
 
+    int molb = 0;
     for (int i = 0; i < g->nat; i++)
     {
-        gmx_mtop_atomnr_to_atom(alook, g->ind[i], &atom);
+        const t_atom *atom;
+        mtopGetAtomParameters(mtop, g->ind[i], &molb, &atom);
         if (atom->q < 0)
         {
             // This is an anion, add it to the list of anions
@@ -1466,16 +1465,13 @@ void init_swapcoords(
 {
     t_swapcoords          *sc;
     t_swap                *s;
-    t_atom                *atom;
+    const t_atom          *atom;
     t_swapgrp             *g;
     swapstateIons_t       *gs;
     gmx_bool               bAppend, bStartFromCpt, bRerun;
-    gmx_mtop_atomlookup_t  alook = NULL;
     matrix                 boxCopy;
 
 
-    alook = gmx_mtop_atomlookup_init(mtop);
-
     if ( (PAR(cr)) && !DOMAINDECOMP(cr) )
     {
         gmx_fatal(FARGS, "Position swapping is only implemented for domain decomposition!");
@@ -1591,15 +1587,16 @@ void init_swapcoords(
         real charge;
 
         g      = &(s->group[ig]);
-        g->apm = get_group_apm_check(ig, s, MASTER(cr) && bVerbose, alook, mtop);
+        g->apm = get_group_apm_check(ig, s, MASTER(cr) && bVerbose, mtop);
 
         /* Since all molecules of a group are equal, we only need enough space
          * to determine properties of a single molecule at at time */
         snew(g->m, g->apm);  /* For the center of mass */
         charge = 0;          /* To determine the charge imbalance */
+        int molb = 0;
         for (int j = 0; j < g->apm; j++)
         {
-            gmx_mtop_atomnr_to_atom(alook, g->ind[j], &atom);
+            mtopGetAtomParameters(mtop, g->ind[j], &molb, &atom);
             g->m[j] = atom->m;
             charge += atom->q;
         }
@@ -1616,9 +1613,10 @@ void init_swapcoords(
         {
             /* Save the split group masses if mass-weighting is requested */
             snew(g->m, g->nat);
+            int molb = 0;
             for (int i = 0; i < g->nat; i++)
             {
-                gmx_mtop_atomnr_to_atom(alook, g->ind[i], &atom);
+                mtopGetAtomParameters(mtop, g->ind[i], &molb, &atom);
                 g->m[i] = atom->m;
             }
         }
@@ -1997,7 +1995,6 @@ gmx_bool do_swapcoords(
         gmx_wallcycle_t   wcycle,
         rvec              x[],
         matrix            box,
-        gmx_mtop_t       *mtop,
         gmx_bool          bVerbose,
         gmx_bool          bRerun)
 {
@@ -2009,7 +2006,6 @@ gmx_bool do_swapcoords(
     t_swapgrp            *g, *gsol;
     int                   isol, iion;
     rvec                  com_solvent, com_particle; /* solvent and swap molecule's center of mass */
-    gmx_mtop_atomlookup_t alook = NULL;
 
 
     wallcycle_start(wcycle, ewcSWAP);
@@ -2089,7 +2085,6 @@ gmx_bool do_swapcoords(
         }
 
         /* Now actually perform the particle exchanges, one swap group after another */
-        alook  = gmx_mtop_atomlookup_init(mtop);
         gsol   = &s->group[eGrpSolvent];
         for (ig = eSwapFixedGrpNR; ig < s->ngrp; ig++)
         {
@@ -2151,7 +2146,6 @@ gmx_bool do_swapcoords(
                         SwS, nswaps, nswaps > 1 ? "s" : "", step, g->molname);
             }
         }
-        gmx_mtop_atomlookup_destroy(alook);
 
         if (s->fpout != NULL)
         {
index 88b353ffa99b3e98fd53ffff247f628a9490e3fb..9f7d556155d0d735c7352a42fae537ef5547e5e1 100644 (file)
@@ -121,7 +121,6 @@ void dd_make_local_swap_groups(gmx_domdec_t *dd, t_swapcoords *si_pub);
  * \param[in] wcycle   Count wallcycles of swap routines for diagnostic output.
  * \param[in] x        Positions of home particles this node owns.
  * \param[in] box      The simulation box.
- * \param[in] mtop     Molecular topology.
  * \param[in] bVerbose Should we be quiet or verbose?
  * \param[in] bRerun   Are we doing a rerun?
  *
@@ -135,7 +134,6 @@ gmx_bool do_swapcoords(
         gmx_wallcycle_t   wcycle,
         rvec              x[],
         matrix            box,
-        gmx_mtop_t       *mtop,
         gmx_bool          bVerbose,
         gmx_bool          bRerun);
 
index 856d9809f7535679eddff8d6db01a37fabcac63f..2238ae7d21dceadce34e8bbb056032b3111b1a83 100644 (file)
@@ -188,7 +188,7 @@ static void tpx2system(FILE *fp, const gmx_mtop_t *mtop)
 {
     int                       nmol, nvsite = 0;
     gmx_mtop_atomloop_block_t aloop;
-    t_atom                   *atom;
+    const t_atom             *atom;
 
     fprintf(fp, "\\subsection{Simulation system}\n");
     aloop = gmx_mtop_atomloop_block_init(mtop);
diff --git a/src/gromacs/topology/mtop_lookup.h b/src/gromacs/topology/mtop_lookup.h
new file mode 100644 (file)
index 0000000..d6cc715
--- /dev/null
@@ -0,0 +1,199 @@
+/*
+ * This file is part of the GROMACS molecular simulation package.
+ *
+ * Copyright (c) 2016, 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.
+ */
+/*! \libinternal \file
+ *
+ * \brief This file contains inline functions to look up atom information
+ * using the global atom index.
+ *
+ * \author Berk Hess <hess@kth.se>
+ * \inlibraryapi
+ * \ingroup module_mtop
+ */
+
+#ifndef GMX_TOPOLOGY_MTOP_LOOKUP_H
+#define GMX_TOPOLOGY_MTOP_LOOKUP_H
+
+#include "gromacs/topology/topology.h"
+#include "gromacs/utility/basedefinitions.h"
+#include "gromacs/utility/gmxassert.h"
+
+struct t_atom;
+
+/*! \brief Look up the molecule block and other indices of a global atom index
+ *
+ * The atom index has to be in range: 0 <= \p globalAtomIndex < \p mtop->natoms.
+ * The input value of moleculeBlock should be in range. Use 0 as starting value.
+ * For subsequent calls to this function, e.g. in a loop, pass in the previously
+ * returned value for best performance. Atoms in a group tend to be in the same
+ * molecule(block), so this minimizes the search time.
+ *
+ * \param[in]     mtop                 The molecule topology
+ * \param[in]     globalAtomIndex      The global atom index to look up
+ * \param[in,out] moleculeBlock        The molecule block index in \p mtop
+ * \param[out]    moleculeIndex        The index of the molecule in the block, can be NULL
+ * \param[out]    atomIndexInMolecule  The atom index in the molecule, can be NULL
+ */
+static void
+mtopGetMolblockIndex(const gmx_mtop_t *mtop,
+                     int               globalAtomIndex,
+                     int              *moleculeBlock,
+                     int              *moleculeIndex,
+                     int              *atomIndexInMolecule)
+{
+    GMX_ASSERT(globalAtomIndex >= 0 && globalAtomIndex < mtop->natoms, "The atom index to look up should be within range");
+    GMX_ASSERT(moleculeBlock != nullptr, "molBlock can not be NULL");
+    GMX_ASSERT(*moleculeBlock >= 0 && *moleculeBlock < mtop->nmolblock, "The starting molecule block index for the search should be within range");
+
+    /* Search the molecue block index using bisection */
+    int molBlock0 = -1;
+    int molBlock1 = mtop->nmolblock;
+
+    int globalAtomStart;
+    while (TRUE)
+    {
+        globalAtomStart = mtop->molblock[*moleculeBlock].globalAtomStart;
+        if (globalAtomIndex < globalAtomStart)
+        {
+            molBlock1 = *moleculeBlock;
+        }
+        else if (globalAtomIndex >= mtop->molblock[*moleculeBlock].globalAtomEnd)
+        {
+            molBlock0 = *moleculeBlock;
+        }
+        else
+        {
+            break;
+        }
+        *moleculeBlock = ((molBlock0 + molBlock1 + 1) >> 1);
+    }
+
+    int molIndex = (globalAtomIndex - globalAtomStart) / mtop->molblock[*moleculeBlock].natoms_mol;
+    if (moleculeIndex != nullptr)
+    {
+        *moleculeIndex = molIndex;
+    }
+    if (atomIndexInMolecule != nullptr)
+    {
+        *atomIndexInMolecule = globalAtomIndex - globalAtomStart - molIndex*mtop->molblock[*moleculeBlock].natoms_mol;
+    }
+}
+
+/*! \brief Look up the molecule block and atom data of a global atom index
+ *
+ * The atom index has to be in range: 0 <= \p globalAtomIndex < \p mtop->natoms.
+ * The input value of moleculeBlock should be in range. Use 0 as starting value.
+ * For subsequent calls to this function, e.g. in a loop, pass in the previously
+ * returned value for best performance. Atoms in a group tend to be in the same
+ * molecule(block), so this minimizes the search time.
+ *
+ * \param[in]     mtop                 The molecule topology
+ * \param[in]     globalAtomIndex      The global atom index to look up
+ * \param[in,out] moleculeBlock        The molecule block index in \p mtop
+ * \param[out]    atom                 Atom data of the global atom index
+ */
+static void
+mtopGetAtomParameters(const gmx_mtop_t  *mtop,
+                      int                globalAtomIndex,
+                      int               *moleculeBlock,
+                      const t_atom     **atom)
+{
+    int atomIndexInMolecule;
+    mtopGetMolblockIndex(mtop, globalAtomIndex, moleculeBlock,
+                         NULL, &atomIndexInMolecule);
+    *atom = &mtop->moltype[mtop->molblock[*moleculeBlock].type].atoms.atom[atomIndexInMolecule];
+}
+
+/*! \brief Look up the atom and residue name and residue number and index of a global atom index
+ *
+ * The atom index has to be in range: 0 <= \p globalAtomIndex < \p mtop->natoms.
+ * The input value of moleculeBlock should be in range. Use 0 as starting value.
+ * For subsequent calls to this function, e.g. in a loop, pass in the previously
+ * returned value for best performance. Atoms in a group tend to be in the same
+ * molecule(block), so this minimizes the search time.
+ * Note that this function does a (somewhat expensive) lookup. If you want
+ * to look up data sequentially for all atoms in a molecule or the system,
+ * use one of the mtop loop functionalities.
+ *
+ * \param[in]     mtop                The molecule topology
+ * \param[in]     globalAtomIndex     The global atom index to look up
+ * \param[in,out] moleculeBlock       The molecule block index in \p mtop
+ * \param[out]    atomName            The atom name, input can be NULL
+ * \param[out]    residueNumber       The residue number, input can be NULL
+ * \param[out]    residueName         The residue name, input can be NULL
+ * \param[out]    globalResidueIndex  The gobal residue index, input can be NULL
+ */
+static void
+mtopGetAtomAndResidueName(const gmx_mtop_t  *mtop,
+                          int                globalAtomIndex,
+                          int               *moleculeBlock,
+                          const char       **atomName,
+                          int               *residueNumber,
+                          const char       **residueName,
+                          int               *globalResidueIndex)
+{
+    int moleculeIndex;
+    int atomIndexInMolecule;
+    mtopGetMolblockIndex(mtop, globalAtomIndex, moleculeBlock,
+                         &moleculeIndex, &atomIndexInMolecule);
+
+    const gmx_molblock_t *molb  = &mtop->molblock[*moleculeBlock];
+    const t_atoms        *atoms = &mtop->moltype[molb->type].atoms;
+    if (atomName != nullptr)
+    {
+        *atomName = *(atoms->atomname[atomIndexInMolecule]);
+    }
+    if (residueNumber != nullptr)
+    {
+        if (atoms->nres > mtop->maxres_renum)
+        {
+            *residueNumber = atoms->resinfo[atoms->atom[atomIndexInMolecule].resind].nr;
+        }
+        else
+        {
+            /* Single residue molecule, keep counting */
+            *residueNumber = molb->residueNumberStart + moleculeIndex*atoms->nres + atoms->atom[atomIndexInMolecule].resind;
+        }
+    }
+    if (residueName != nullptr)
+    {
+        *residueName = *(atoms->resinfo[atoms->atom[atomIndexInMolecule].resind].name);
+    }
+    if (globalResidueIndex != nullptr)
+    {
+        *globalResidueIndex = molb->globalResidueStart + moleculeIndex*atoms->nres + atoms->atom[atomIndexInMolecule].resind;
+    }
+}
+
+#endif
index d4482b564a324da54921662938c4b32119c6f5d5..0f140660761e9be2188045d02c0b3987c8e8edcd 100644 (file)
@@ -78,6 +78,28 @@ static int gmx_mtop_maxresnr(const gmx_mtop_t *mtop, int maxres_renum)
     return maxresnr;
 }
 
+static void finalizeMolblocks(gmx_mtop_t *mtop)
+{
+    int atomIndex          = 0;
+    int residueIndex       = 0;
+    int residueNumberStart = mtop->maxresnr + 1;
+    for (int mb = 0; mb < mtop->nmolblock; mb++)
+    {
+        gmx_molblock_t *molb          = &mtop->molblock[mb];
+        int             numResPerMol  = mtop->moltype[molb->type].atoms.nres;
+        molb->globalAtomStart         = atomIndex;
+        molb->globalResidueStart      = residueIndex;
+        atomIndex                    += molb->nmol*molb->natoms_mol;
+        residueIndex                 += molb->nmol*numResPerMol;
+        molb->globalAtomEnd           = atomIndex;
+        molb->residueNumberStart      = residueNumberStart;
+        if (numResPerMol <= mtop->maxres_renum)
+        {
+            residueNumberStart       += molb->nmol*numResPerMol;
+        }
+    }
+}
+
 void gmx_mtop_finalize(gmx_mtop_t *mtop)
 {
     char *env;
@@ -110,6 +132,8 @@ void gmx_mtop_finalize(gmx_mtop_t *mtop)
     }
 
     mtop->maxresnr = gmx_mtop_maxresnr(mtop, mtop->maxres_renum);
+
+    finalizeMolblocks(mtop);
 }
 
 void gmx_mtop_count_atomtypes(const gmx_mtop_t *mtop, int state, int typecount[])
@@ -177,279 +201,6 @@ void gmx_mtop_remove_chargegroups(gmx_mtop_t *mtop)
     }
 }
 
-
-typedef struct
-{
-    int a_start;
-    int a_end;
-    int na_mol;
-} mb_at_t;
-
-typedef struct gmx_mtop_atomlookup
-{
-    const gmx_mtop_t *mtop;
-    int               nmb;
-    int               mb_start;
-    mb_at_t          *mba;
-} t_gmx_mtop_atomlookup;
-
-
-gmx_mtop_atomlookup_t
-gmx_mtop_atomlookup_init(const gmx_mtop_t *mtop)
-{
-    t_gmx_mtop_atomlookup *alook;
-    int                    mb;
-    int                    a_start, a_end, na, na_start = -1;
-
-    snew(alook, 1);
-
-    alook->mtop     = mtop;
-    alook->nmb      = mtop->nmolblock;
-    alook->mb_start = 0;
-    snew(alook->mba, alook->nmb);
-
-    a_start = 0;
-    for (mb = 0; mb < mtop->nmolblock; mb++)
-    {
-        na    = mtop->molblock[mb].nmol*mtop->molblock[mb].natoms_mol;
-        a_end = a_start + na;
-
-        alook->mba[mb].a_start = a_start;
-        alook->mba[mb].a_end   = a_end;
-        alook->mba[mb].na_mol  = mtop->molblock[mb].natoms_mol;
-
-        /* We start the binary search with the largest block */
-        if (mb == 0 || na > na_start)
-        {
-            alook->mb_start = mb;
-            na_start        = na;
-        }
-
-        a_start = a_end;
-    }
-
-    return alook;
-}
-
-gmx_mtop_atomlookup_t
-gmx_mtop_atomlookup_settle_init(const gmx_mtop_t *mtop)
-{
-    t_gmx_mtop_atomlookup *alook;
-    int                    mb;
-    int                    na, na_start = -1;
-
-    alook = gmx_mtop_atomlookup_init(mtop);
-
-    /* Check if the starting molblock has settle */
-    if (mtop->moltype[mtop->molblock[alook->mb_start].type].ilist[F_SETTLE].nr  == 0)
-    {
-        /* Search the largest molblock with settle */
-        alook->mb_start = -1;
-        for (mb = 0; mb < mtop->nmolblock; mb++)
-        {
-            if (mtop->moltype[mtop->molblock[mb].type].ilist[F_SETTLE].nr > 0)
-            {
-                na = alook->mba[mb].a_end - alook->mba[mb].a_start;
-                if (alook->mb_start == -1 || na > na_start)
-                {
-                    alook->mb_start = mb;
-                    na_start        = na;
-                }
-            }
-        }
-
-        if (alook->mb_start == -1)
-        {
-            gmx_incons("gmx_mtop_atomlookup_settle_init called without settles");
-        }
-    }
-
-    return alook;
-}
-
-void
-gmx_mtop_atomlookup_destroy(gmx_mtop_atomlookup_t alook)
-{
-    sfree(alook->mba);
-    sfree(alook);
-}
-
-void gmx_mtop_atomnr_to_atom(const gmx_mtop_atomlookup_t alook,
-                             int                         atnr_global,
-                             t_atom                    **atom)
-{
-    int mb0, mb1, mb;
-    int a_start, atnr_mol;
-
-#ifdef DEBUG_MTOP
-    if (atnr_global < 0 || atnr_global >= mtop->natoms)
-    {
-        gmx_fatal(FARGS, "gmx_mtop_atomnr_to_moltype was called with atnr_global=%d which is not in the atom range of this system (%d-%d)",
-                  atnr_global, 0, mtop->natoms-1);
-    }
-#endif
-
-    mb0 = -1;
-    mb1 = alook->nmb;
-    mb  = alook->mb_start;
-
-    while (TRUE)
-    {
-        a_start = alook->mba[mb].a_start;
-        if (atnr_global < a_start)
-        {
-            mb1 = mb;
-        }
-        else if (atnr_global >= alook->mba[mb].a_end)
-        {
-            mb0 = mb;
-        }
-        else
-        {
-            break;
-        }
-        mb = ((mb0 + mb1 + 1)>>1);
-    }
-
-    atnr_mol = (atnr_global - a_start) % alook->mba[mb].na_mol;
-
-    *atom = &alook->mtop->moltype[alook->mtop->molblock[mb].type].atoms.atom[atnr_mol];
-}
-
-void gmx_mtop_atomnr_to_ilist(const gmx_mtop_atomlookup_t alook,
-                              int atnr_global,
-                              t_ilist **ilist_mol, int *atnr_offset)
-{
-    int mb0, mb1, mb;
-    int a_start, atnr_local;
-
-#ifdef DEBUG_MTOP
-    if (atnr_global < 0 || atnr_global >= mtop->natoms)
-    {
-        gmx_fatal(FARGS, "gmx_mtop_atomnr_to_moltype was called with atnr_global=%d which is not in the atom range of this system (%d-%d)",
-                  atnr_global, 0, mtop->natoms-1);
-    }
-#endif
-
-    mb0 = -1;
-    mb1 = alook->nmb;
-    mb  = alook->mb_start;
-
-    while (TRUE)
-    {
-        a_start = alook->mba[mb].a_start;
-        if (atnr_global < a_start)
-        {
-            mb1 = mb;
-        }
-        else if (atnr_global >= alook->mba[mb].a_end)
-        {
-            mb0 = mb;
-        }
-        else
-        {
-            break;
-        }
-        mb = ((mb0 + mb1 + 1)>>1);
-    }
-
-    *ilist_mol = alook->mtop->moltype[alook->mtop->molblock[mb].type].ilist;
-
-    atnr_local = (atnr_global - a_start) % alook->mba[mb].na_mol;
-
-    *atnr_offset = atnr_global - atnr_local;
-}
-
-void gmx_mtop_atomnr_to_molblock_ind(const gmx_mtop_atomlookup_t alook,
-                                     int atnr_global,
-                                     int *molb, int *molnr, int *atnr_mol)
-{
-    int mb0, mb1, mb;
-    int a_start;
-
-#ifdef DEBUG_MTOP
-    if (atnr_global < 0 || atnr_global >= mtop->natoms)
-    {
-        gmx_fatal(FARGS, "gmx_mtop_atomnr_to_moltype was called with atnr_global=%d which is not in the atom range of this system (%d-%d)",
-                  atnr_global, 0, mtop->natoms-1);
-    }
-#endif
-
-    mb0 = -1;
-    mb1 = alook->nmb;
-    mb  = alook->mb_start;
-
-    while (TRUE)
-    {
-        a_start = alook->mba[mb].a_start;
-        if (atnr_global < a_start)
-        {
-            mb1 = mb;
-        }
-        else if (atnr_global >= alook->mba[mb].a_end)
-        {
-            mb0 = mb;
-        }
-        else
-        {
-            break;
-        }
-        mb = ((mb0 + mb1 + 1)>>1);
-    }
-
-    *molb     = mb;
-    *molnr    = (atnr_global - a_start) / alook->mba[mb].na_mol;
-    *atnr_mol = atnr_global - a_start - (*molnr)*alook->mba[mb].na_mol;
-}
-
-void gmx_mtop_atominfo_global(const gmx_mtop_t *mtop, int atnr_global,
-                              char **atomname, int *resnr, char **resname)
-{
-    int             mb, a_start, a_end, maxresnr, at_loc;
-    t_atoms        *atoms = NULL;
-
-    if (atnr_global < 0 || atnr_global >= mtop->natoms)
-    {
-        gmx_fatal(FARGS, "gmx_mtop_atominfo_global was called with atnr_global=%d which is not in the atom range of this system (%d-%d)",
-                  atnr_global, 0, mtop->natoms-1);
-    }
-
-    mb       = -1;
-    a_end    = 0;
-    maxresnr = mtop->maxresnr;
-    do
-    {
-        if (mb >= 0)
-        {
-            /* cppcheck-suppress nullPointer #6330 will be fixed in cppcheck 1.73 */
-            if (atoms->nres <= mtop->maxres_renum)
-            {
-                /* Single residue molecule, keep counting */
-                /* cppcheck-suppress nullPointer #6330 will be fixed in cppcheck 1.73 */
-                maxresnr += mtop->molblock[mb].nmol*atoms->nres;
-            }
-        }
-        mb++;
-        atoms   = &mtop->moltype[mtop->molblock[mb].type].atoms;
-        a_start = a_end;
-        a_end   = a_start + mtop->molblock[mb].nmol*atoms->nr;
-    }
-    while (atnr_global >= a_end);
-
-    at_loc    = (atnr_global - a_start) % atoms->nr;
-    *atomname = *(atoms->atomname[at_loc]);
-    if (atoms->nres > mtop->maxres_renum)
-    {
-        *resnr = atoms->resinfo[atoms->atom[at_loc].resind].nr;
-    }
-    else
-    {
-        /* Single residue molecule, keep counting */
-        *resnr = maxresnr + 1 + (atnr_global - a_start)/atoms->nr*atoms->nres + atoms->atom[at_loc].resind;
-    }
-    *resname  = *(atoms->resinfo[atoms->atom[at_loc].resind].name);
-}
-
 typedef struct gmx_mtop_atomloop_all
 {
     const gmx_mtop_t *mtop;
@@ -486,7 +237,7 @@ static void gmx_mtop_atomloop_all_destroy(gmx_mtop_atomloop_all_t aloop)
 }
 
 gmx_bool gmx_mtop_atomloop_all_next(gmx_mtop_atomloop_all_t aloop,
-                                    int *at_global, t_atom **atom)
+                                    int *at_global, const t_atom **atom)
 {
     if (aloop == NULL)
     {
@@ -575,7 +326,7 @@ static void gmx_mtop_atomloop_block_destroy(gmx_mtop_atomloop_block_t aloop)
 }
 
 gmx_bool gmx_mtop_atomloop_block_next(gmx_mtop_atomloop_block_t aloop,
-                                      t_atom **atom, int *nmol)
+                                      const t_atom **atom, int *nmol)
 {
     if (aloop == NULL)
     {
@@ -1084,7 +835,6 @@ static void gen_local_top(const gmx_mtop_t *mtop,
     real                   *qA, *qB;
     gmx_mtop_atomloop_all_t aloop;
     int                     ag;
-    t_atom                 *atom;
 
     top->atomtypes = mtop->atomtypes;
 
@@ -1177,6 +927,7 @@ static void gen_local_top(const gmx_mtop_t *mtop,
         snew(qA, mtop->natoms);
         snew(qB, mtop->natoms);
         aloop = gmx_mtop_atomloop_all_init(mtop);
+        const t_atom *atom;
         while (gmx_mtop_atomloop_all_next(aloop, &ag, &atom))
         {
             qA[ag] = atom->q;
@@ -1250,7 +1001,7 @@ std::vector<size_t> get_atom_index(const gmx_mtop_t *mtop)
 
     std::vector<size_t>       atom_index;
     gmx_mtop_atomloop_block_t aloopb = gmx_mtop_atomloop_block_init(mtop);
-    t_atom                   *atom;
+    const t_atom             *atom;
     int                       nmol, j = 0;
     while (gmx_mtop_atomloop_block_next(aloopb, &atom, &nmol))
     {
index 77b1569d81a37263c47230d02643c799250474b6..622f4cd2c88b7f6adaf8037e163d0a255f6dcd19 100644 (file)
 
 #include <vector>
 
+#include "gromacs/topology/topology.h"
 #include "gromacs/utility/basedefinitions.h"
+#include "gromacs/utility/gmxassert.h"
 
 struct gmx_localtop_t;
-struct gmx_moltype_t;
-struct gmx_mtop_t;
 struct t_atom;
 struct t_atoms;
 struct t_block;
 struct t_ilist;
 struct t_symtab;
-struct t_topology;
 
 /* Should be called after generating or reading mtop,
  * to set some compute intesive variables to avoid
@@ -74,59 +73,6 @@ ncg_mtop(const gmx_mtop_t *mtop);
 /* Removes the charge groups, i.e. makes single atom charge groups, in mtop */
 void gmx_mtop_remove_chargegroups(gmx_mtop_t *mtop);
 
-
-/* Abstract data type for looking up atoms by global atom number */
-typedef struct gmx_mtop_atomlookup *gmx_mtop_atomlookup_t;
-
-/* Initialize atom lookup by global atom number */
-gmx_mtop_atomlookup_t
-gmx_mtop_atomlookup_init(const gmx_mtop_t *mtop);
-
-/* As gmx_mtop_atomlookup_init, but optimized for atoms involved in settle */
-gmx_mtop_atomlookup_t
-gmx_mtop_atomlookup_settle_init(const gmx_mtop_t *mtop);
-
-/* Destroy a gmx_mtop_atomlookup_t data structure */
-void
-gmx_mtop_atomlookup_destroy(gmx_mtop_atomlookup_t alook);
-
-
-/* Returns a pointer to the t_atom struct belonging to atnr_global.
- * This can be an expensive operation, so if possible use
- * one of the atom loop constructs below.
- */
-void
-gmx_mtop_atomnr_to_atom(const gmx_mtop_atomlookup_t alook,
-                        int                         atnr_global,
-                        t_atom                    **atom);
-
-
-/* Returns a pointer to the molecule interaction array ilist_mol[F_NRE]
- * and the local atom number in the molecule belonging to atnr_global.
- */
-void
-gmx_mtop_atomnr_to_ilist(const gmx_mtop_atomlookup_t alook,
-                         int atnr_global,
-                         t_ilist **ilist_mol, int *atnr_offset);
-
-
-/* Returns the molecule block index
- * and the molecule number in the block
- * and the atom number offset for the atom indices in moltype
- * belonging to atnr_global.
- */
-void
-gmx_mtop_atomnr_to_molblock_ind(const gmx_mtop_atomlookup_t alook,
-                                int atnr_global,
-                                int *molb, int *molnr, int *atnr_mol);
-
-
-/* Returns atom name, global resnr and residue name  of atom atnr_global */
-void
-gmx_mtop_atominfo_global(const gmx_mtop_t *mtop, int atnr_global,
-                         char **atomname, int *resnr, char **resname);
-
-
 /* Abstract type for atom loop over all atoms */
 typedef struct gmx_mtop_atomloop_all *gmx_mtop_atomloop_all_t;
 
@@ -154,7 +100,7 @@ gmx_mtop_atomloop_all_init(const gmx_mtop_t *mtop);
  */
 gmx_bool
 gmx_mtop_atomloop_all_next(gmx_mtop_atomloop_all_t aloop,
-                           int *at_global, t_atom **atom);
+                           int *at_global, const t_atom **atom);
 
 /* Return the atomname, the residue number and residue name
  * of the current atom in the loop.
@@ -194,7 +140,7 @@ gmx_mtop_atomloop_block_init(const gmx_mtop_t *mtop);
  */
 gmx_bool
 gmx_mtop_atomloop_block_next(gmx_mtop_atomloop_block_t aloop,
-                             t_atom **atom, int *nmol);
+                             const t_atom **atom, int *nmol);
 
 
 /* Abstract type for ilist loop over all ilists */
index ec16e7cda7a999393ab3763b0d9f6cc928249818..3599ac893a90921dfc01f231825b2efbd964aac7 100644 (file)
@@ -63,15 +63,22 @@ typedef struct gmx_moltype_t
     t_blocka        excls;        /* The exclusions                       */
 } gmx_moltype_t;
 
+/*! \brief Block of molecules of the same type, used in gmx_mtop_t */
 typedef struct gmx_molblock_t
 {
-    int            type;        /* The molcule type index in mtop.moltype */
-    int            nmol;        /* The number of molecules in this block  */
-    int            natoms_mol;  /* The number of atoms in one molecule    */
-    int            nposres_xA;  /* The number of posres coords for top A  */
-    rvec          *posres_xA;   /* The posres coords for top A            */
-    int            nposres_xB;  /* The number of posres coords for top B  */
-    rvec          *posres_xB;   /* The posres coords for top B            */
+    int     type;               /**< The molecule type index in mtop.moltype  */
+    int     nmol;               /**< The number of molecules in this block    */
+    int     nposres_xA;         /**< The number of posres coords for top A    */
+    rvec   *posres_xA;          /**< Position restraint coordinates for top A */
+    int     nposres_xB;         /**< The number of posres coords for top B    */
+    rvec   *posres_xB;          /**< Position restraint coordinates for top B */
+
+    /* Convenience information, derived from other gmx_mtop_t contents     */
+    int     natoms_mol;         /**< The number of atoms in one molecule      */
+    int     globalAtomStart;    /**< Global atom index of the first atom in the block */
+    int     globalAtomEnd;      /**< Global atom index + 1 of the last atom in the block */
+    int     globalResidueStart; /**< Global residue index of the first residue in the block */
+    int     residueNumberStart; /**< Residue numbers start from this value if the number of residues per molecule is <= maxres_renum */
 } gmx_molblock_t;
 
 typedef struct gmx_groups_t
index b7d924fa0273f8cd8b45329716513bbc9718c7f1..d92c53074f6ca9bdd234e9443c184ad2eedc3810 100644 (file)
@@ -1651,7 +1651,7 @@ double gmx::do_md(FILE *fplog, t_commrec *cr, const gmx::MDLogger &mdlog,
             bNeedRepartition = do_swapcoords(cr, step, t, ir, wcycle,
                                              bRerunMD ? rerun_fr.x   : as_rvec_array(state->x.data()),
                                              bRerunMD ? rerun_fr.box : state->box,
-                                             top_global, MASTER(cr) && bVerbose, bRerunMD);
+                                             MASTER(cr) && bVerbose, bRerunMD);
 
             if (bNeedRepartition && DOMAINDECOMP(cr))
             {
index 03b35b200a22f01c5bd376ac575598576552d097..c57dc73c3234dcc05ec1900d4d604081c50381bc 100644 (file)
@@ -50,6 +50,7 @@
 #include "gromacs/mdtypes/state.h"
 #include "gromacs/pbcutil/pbc.h"
 #include "gromacs/topology/index.h"
+#include "gromacs/topology/mtop_lookup.h"
 #include "gromacs/topology/mtop_util.h"
 #include "gromacs/topology/topology.h"
 #include "gromacs/utility/cstringutil.h"
@@ -106,18 +107,15 @@ static int get_mol_id(int at, gmx_mtop_t  *mtop, int *type, int *block)
     int                   mol_id = 0;
     int                   i;
     int                   atnr_mol;
-    gmx_mtop_atomlookup_t alook;
 
-    alook = gmx_mtop_atomlookup_settle_init(mtop);
-    gmx_mtop_atomnr_to_molblock_ind(alook, at, block, &mol_id, &atnr_mol);
+    *block = 0;
+    mtopGetMolblockIndex(mtop, at, block, &mol_id, &atnr_mol);
     for (i = 0; i < *block; i++)
     {
         mol_id += mtop->molblock[i].nmol;
     }
     *type = mtop->molblock[*block].type;
 
-    gmx_mtop_atomlookup_destroy(alook);
-
     return mol_id;
 }