From 62343cded1fa4455de3241d345ec82f88be60d33 Mon Sep 17 00:00:00 2001 From: Berk Hess Date: Fri, 6 Dec 2019 10:11:55 +0100 Subject: [PATCH] Use ListOfLists in gmx_mtop_t and gmx_localtop_t This change uses ListOfLists for exclusions in gmx_mtop_t and gmx_localtop_t. Exclusions in t_topolgy as well as other atom and index lists still use t_blocka. Change-Id: Idfd532ad2a4b19d68e60b05dbd0a3082650ec20a --- src/gromacs/domdec/domdec_topology.cpp | 182 +++++------------- src/gromacs/fileio/tpxio.cpp | 26 ++- src/gromacs/gmxpreprocess/gpp_nextnb.cpp | 21 +- src/gromacs/gmxpreprocess/gpp_nextnb.h | 9 +- src/gromacs/gmxpreprocess/grompp.cpp | 2 +- src/gromacs/gmxpreprocess/grompp_impl.h | 3 +- src/gromacs/gmxpreprocess/topio.cpp | 2 +- src/gromacs/gmxpreprocess/toppush.cpp | 42 ++-- src/gromacs/mdlib/dispersioncorrection.cpp | 11 +- src/gromacs/mdlib/forcerec.cpp | 6 +- src/gromacs/mdlib/perf_est.cpp | 2 +- src/gromacs/mdlib/sim_util.cpp | 4 +- src/gromacs/mdrun/tpi.cpp | 2 +- src/gromacs/nbnxm/benchmark/bench_setup.cpp | 2 +- src/gromacs/nbnxm/benchmark/bench_system.cpp | 9 +- src/gromacs/nbnxm/benchmark/bench_system.h | 4 +- src/gromacs/nbnxm/nbnxm.h | 8 +- src/gromacs/nbnxm/pairlist.cpp | 48 ++--- src/gromacs/nbnxm/pairlistset.h | 9 +- src/gromacs/nbnxm/pairlistsets.h | 18 +- src/gromacs/selection/nbsearch.cpp | 44 ++--- src/gromacs/selection/nbsearch.h | 5 +- src/gromacs/selection/tests/nbsearch.cpp | 44 ++--- src/gromacs/tools/convert_tpr.cpp | 46 ++--- src/gromacs/topology/block.cpp | 51 +++++ src/gromacs/topology/block.h | 4 + src/gromacs/topology/exclusionblocks.cpp | 71 +++++-- src/gromacs/topology/exclusionblocks.h | 4 +- src/gromacs/topology/mtop_util.cpp | 63 ++++-- .../topology/tests/exclusionblocks.cpp | 36 +++- src/gromacs/topology/topology.cpp | 26 +-- src/gromacs/topology/topology.h | 11 +- .../trajectoryanalysis/modules/rdf.cpp | 2 +- 33 files changed, 433 insertions(+), 384 deletions(-) diff --git a/src/gromacs/domdec/domdec_topology.cpp b/src/gromacs/domdec/domdec_topology.cpp index 19acbac64b..23346db87c 100644 --- a/src/gromacs/domdec/domdec_topology.cpp +++ b/src/gromacs/domdec/domdec_topology.cpp @@ -75,6 +75,7 @@ #include "gromacs/utility/exceptions.h" #include "gromacs/utility/fatalerror.h" #include "gromacs/utility/gmxassert.h" +#include "gromacs/utility/listoflists.h" #include "gromacs/utility/logger.h" #include "gromacs/utility/smalloc.h" #include "gromacs/utility/strconvert.h" @@ -87,6 +88,8 @@ #include "domdec_vsite.h" #include "dump.h" +using gmx::ListOfLists; + /*! \brief The number of integer item in the local state, used for broadcasting of the state */ #define NITEM_DD_INIT_LOCAL_STATE 5 @@ -111,7 +114,7 @@ struct thread_work_t t_idef idef; /**< Partial local topology */ std::unique_ptr vsitePbc; /**< vsite PBC structure */ int nbonded; /**< The number of bondeds in this struct */ - t_blocka excl; /**< List of exclusions */ + ListOfLists excl; /**< List of exclusions */ int excl_count; /**< The total exclusion count for \p excl */ }; @@ -119,8 +122,6 @@ struct thread_work_t struct gmx_reverse_top_t { //! @cond Doxygen_Suppress - //! \brief The maximum number of exclusions one atom can have - int n_excl_at_max = 0; //! \brief Are there constraints in this revserse top? bool bConstr = false; //! \brief Are there settles in this revserse top? @@ -453,14 +454,15 @@ static void global_atomnr_to_moltype_ind(const gmx_reverse_top_t* rt, int i_gl, } /*! \brief Returns the maximum number of exclusions per atom */ -static int getMaxNumExclusionsPerAtom(const t_blocka& excls) +static int getMaxNumExclusionsPerAtom(const ListOfLists& excls) { int maxNumExcls = 0; - for (int at = 0; at < excls.nr; at++) + for (gmx::index at = 0; at < excls.ssize(); at++) { - const int numExcls = excls.index[at + 1] - excls.index[at]; + const auto list = excls[at]; + const int numExcls = list.ssize(); - GMX_RELEASE_ASSERT(numExcls != 1 || excls.a[excls.index[at]] == at, + GMX_RELEASE_ASSERT(numExcls != 1 || list[0] == at, "With 1 exclusion we expect a self-exclusion"); maxNumExcls = std::max(maxNumExcls, numExcls); @@ -717,10 +719,7 @@ void dd_make_reverse_top(FILE* fplog, make_reverse_top(mtop, ir->efep != efepNO, !dd->comm->systemInfo.haveSplitConstraints, !dd->comm->systemInfo.haveSplitSettles, bBCheck, &dd->nbonded_global); - gmx_reverse_top_t* rt = dd->reverse_top; - dd->haveExclusions = false; - rt->n_excl_at_max = 0; for (const gmx_molblock_t& molb : mtop->molblock) { const int maxNumExclusionsPerAtom = getMaxNumExclusionsPerAtom(mtop->moltype[molb.type].excls); @@ -729,7 +728,6 @@ void dd_make_reverse_top(FILE* fplog, { dd->haveExclusions = true; } - rt->n_excl_at_max = std::max(rt->n_excl_at_max, maxNumExclusionsPerAtom); } if (vsite && vsite->numInterUpdategroupVsites > 0) @@ -995,37 +993,12 @@ static real dd_dist2(t_pbc* pbc_null, const rvec* x, const int i, int j) return norm2(dx); } -/*! \brief Append t_blocka block structures 1 to nsrc in src to *dest */ -static void combine_blocka(t_blocka* dest, gmx::ArrayRef src) +/*! \brief Append ListOfLists exclusion objects 1 to nsrc in \p src to \p *dest */ +static void combineExclusions(ListOfLists* dest, gmx::ArrayRef src) { - int ni = src.back().excl.nr; - int na = 0; - for (const thread_work_t& th_work : src) - { - na += th_work.excl.nra; - } - if (ni + 1 > dest->nalloc_index) - { - dest->nalloc_index = over_alloc_large(ni + 1); - srenew(dest->index, dest->nalloc_index); - } - if (dest->nra + na > dest->nalloc_a) - { - dest->nalloc_a = over_alloc_large(dest->nra + na); - srenew(dest->a, dest->nalloc_a); - } for (gmx::index s = 1; s < src.ssize(); s++) { - for (int i = dest->nr + 1; i < src[s].excl.nr + 1; i++) - { - dest->index[i] = dest->nra + src[s].excl.index[i]; - } - for (int i = 0; i < src[s].excl.nra; i++) - { - dest->a[dest->nra + i] = src[s].excl.a[i]; - } - dest->nr = src[s].excl.nr; - dest->nra += src[s].excl.nra; + dest->appendListOfLists(src[s].excl); } } @@ -1377,11 +1350,11 @@ static int make_bondeds_zone(gmx_domdec_t* dd, } /*! \brief Set the exclusion data for i-zone \p iz for the case of no exclusions */ -static void set_no_exclusions_zone(const gmx_domdec_zones_t* zones, int iz, t_blocka* lexcls) +static void set_no_exclusions_zone(const gmx_domdec_zones_t* zones, int iz, ListOfLists* lexcls) { for (int a = zones->cg_range[iz]; a < zones->cg_range[iz + 1]; a++) { - lexcls->index[a + 1] = lexcls->nra; + lexcls->pushBack({}); } } @@ -1390,52 +1363,33 @@ static void make_exclusions_zone(gmx_domdec_t* dd, gmx_domdec_zones_t* zones, const std::vector& moltype, const int* cginfo, - t_blocka* lexcls, + ListOfLists* lexcls, int iz, int at_start, int at_end, const gmx::ArrayRef intermolecularExclusionGroup) { - int n_excl_at_max, n, at; - const gmx_ga2la_t& ga2la = *dd->ga2la; const auto& jAtomRange = zones->iZones[iz].jAtomRange; - n_excl_at_max = dd->reverse_top->n_excl_at_max; - - /* We set the end index, but note that we might not start at zero here */ - lexcls->nr = at_end; + const gmx::index oldNumLists = lexcls->ssize(); - n = lexcls->nra; - for (at = at_start; at < at_end; at++) + std::vector exclusionsForAtom; + for (int at = at_start; at < at_end; at++) { - if (n + 1000 > lexcls->nalloc_a) - { - lexcls->nalloc_a = over_alloc_large(n + 1000); - srenew(lexcls->a, lexcls->nalloc_a); - } + exclusionsForAtom.clear(); if (GET_CGINFO_EXCL_INTER(cginfo[at])) { - int a_gl, mb, mt, mol, a_mol, j; - const t_blocka* excls; - - if (n + n_excl_at_max > lexcls->nalloc_a) - { - lexcls->nalloc_a = over_alloc_large(n + n_excl_at_max); - srenew(lexcls->a, lexcls->nalloc_a); - } + int a_gl, mb, mt, mol, a_mol; /* Copy the exclusions from the global top */ - lexcls->index[at] = n; - a_gl = dd->globalAtomIndices[at]; + a_gl = dd->globalAtomIndices[at]; global_atomnr_to_moltype_ind(dd->reverse_top, a_gl, &mb, &mt, &mol, &a_mol); - excls = &moltype[mt].excls; - for (j = excls->index[a_mol]; j < excls->index[a_mol + 1]; j++) + const auto excls = moltype[mt].excls[a_mol]; + for (const int aj_mol : excls) { - const int aj_mol = excls->a[j]; - if (const auto* jEntry = ga2la.find(a_gl + aj_mol - a_mol)) { /* This check is not necessary, but it can reduce @@ -1444,16 +1398,11 @@ static void make_exclusions_zone(gmx_domdec_t* dd, */ if (jAtomRange.isInRange(jEntry->la)) { - lexcls->a[n++] = jEntry->la; + exclusionsForAtom.push_back(jEntry->la); } } } } - else - { - /* We don't need exclusions for this atom */ - lexcls->index[at] = n; - } bool isExcludedAtom = !intermolecularExclusionGroup.empty() && std::find(intermolecularExclusionGroup.begin(), @@ -1462,51 +1411,26 @@ static void make_exclusions_zone(gmx_domdec_t* dd, if (isExcludedAtom) { - if (n + intermolecularExclusionGroup.ssize() > lexcls->nalloc_a) - { - lexcls->nalloc_a = over_alloc_large(n + intermolecularExclusionGroup.size()); - srenew(lexcls->a, lexcls->nalloc_a); - } for (int qmAtomGlobalIndex : intermolecularExclusionGroup) { if (const auto* entry = dd->ga2la->find(qmAtomGlobalIndex)) { - lexcls->a[n++] = entry->la; + exclusionsForAtom.push_back(entry->la); } } } - } - - lexcls->index[lexcls->nr] = n; - lexcls->nra = n; -} - -/*! \brief Ensure we have enough space in \p ba for \p nindex_max indices */ -static void check_alloc_index(t_blocka* ba, int nindex_max) -{ - if (nindex_max + 1 > ba->nalloc_index) - { - ba->nalloc_index = over_alloc_dd(nindex_max + 1); - srenew(ba->index, ba->nalloc_index); + /* Append the exclusions for this atom to the topology */ + lexcls->pushBack(exclusionsForAtom); } -} -/*! \brief Ensure that we have enough space for exclusion storate in \p lexcls */ -static void check_exclusions_alloc(const gmx_domdec_t* dd, const gmx_domdec_zones_t* zones, t_blocka* lexcls) -{ - const int nr = zones->iZones.back().iAtomRange.end(); - - check_alloc_index(lexcls, nr); - - for (size_t thread = 1; thread < dd->reverse_top->th_work.size(); thread++) - { - check_alloc_index(&dd->reverse_top->th_work[thread].excl, nr); - } + GMX_RELEASE_ASSERT( + lexcls->ssize() - oldNumLists == at_end - at_start, + "The number of exclusion list should match the number of atoms in the range"); } /*! \brief Set the total count indexes for the local exclusions, needed by several functions */ -static void finish_local_exclusions(gmx_domdec_t* dd, gmx_domdec_zones_t* zones, t_blocka* lexcls) +static void finish_local_exclusions(gmx_domdec_t* dd, gmx_domdec_zones_t* zones, ListOfLists* lexcls) { const gmx::Range nonhomeIzonesAtomRange(zones->iZones[0].iAtomRange.end(), zones->iZones.back().iAtomRange.end()); @@ -1516,19 +1440,10 @@ static void finish_local_exclusions(gmx_domdec_t* dd, gmx_domdec_zones_t* zones, /* There are no exclusions involving non-home charge groups, * but we need to set the indices for neighborsearching. */ - for (int la : nonhomeIzonesAtomRange) + for (int gmx_unused la : nonhomeIzonesAtomRange) { - lexcls->index[la] = lexcls->nra; + lexcls->pushBack({}); } - - /* nr is only used to loop over the exclusions for Ewald and RF, - * so we can set it to the number of home atoms for efficiency. - */ - lexcls->nr = nonhomeIzonesAtomRange.begin(); - } - else - { - lexcls->nr = nonhomeIzonesAtomRange.end(); } } @@ -1556,7 +1471,7 @@ static int make_local_bondeds_excls(gmx_domdec_t* dd, t_pbc* pbc_null, rvec* cg_cm, t_idef* idef, - t_blocka* lexcls, + ListOfLists* lexcls, int* excl_count) { int nzone_bondeds, nzone_excl; @@ -1588,8 +1503,6 @@ static int make_local_bondeds_excls(gmx_domdec_t* dd, nzone_excl = 1; } - check_exclusions_alloc(dd, zones, lexcls); - rt = dd->reverse_top; rc2 = rc * rc; @@ -1598,8 +1511,7 @@ static int make_local_bondeds_excls(gmx_domdec_t* dd, clear_idef(idef); nbonded_local = 0; - lexcls->nr = 0; - lexcls->nra = 0; + lexcls->clear(); *excl_count = 0; for (int izone = 0; izone < nzone_bondeds; izone++) @@ -1613,9 +1525,8 @@ static int make_local_bondeds_excls(gmx_domdec_t* dd, { try { - int cg0t, cg1t; - t_idef* idef_t; - t_blocka* excl_t; + int cg0t, cg1t; + t_idef* idef_t; cg0t = cg0 + ((cg1 - cg0) * thread) / numThreads; cg1t = cg0 + ((cg1 - cg0) * (thread + 1)) / numThreads; @@ -1636,15 +1547,17 @@ static int make_local_bondeds_excls(gmx_domdec_t* dd, if (izone < nzone_excl) { + ListOfLists* excl_t; if (thread == 0) { + // Thread 0 stores exclusions directly in the final storage excl_t = lexcls; } else { - excl_t = &rt->th_work[thread].excl; - excl_t->nr = 0; - excl_t->nra = 0; + // Threads > 0 store in temporary storage, starting at list index 0 + excl_t = &rt->th_work[thread].excl; + excl_t->clear(); } /* No charge groups and no distance check required */ @@ -1669,7 +1582,7 @@ static int make_local_bondeds_excls(gmx_domdec_t* dd, { if (rt->th_work.size() > 1) { - combine_blocka(lexcls, rt->th_work); + combineExclusions(lexcls, rt->th_work); } for (const thread_work_t& th_work : rt->th_work) @@ -1690,7 +1603,7 @@ static int make_local_bondeds_excls(gmx_domdec_t* dd, finish_local_exclusions(dd, zones, lexcls); if (debug) { - fprintf(debug, "We have %d exclusions, check count %d\n", lexcls->nra, *excl_count); + fprintf(debug, "We have %d exclusions, check count %d\n", lexcls->numElements(), *excl_count); } return nbonded_local; @@ -2062,12 +1975,11 @@ static void bonded_cg_distance_mol(const gmx_moltype_t* molt, } if (bExcl) { - const t_blocka* excls = &molt->excls; - for (int ai = 0; ai < excls->nr; ai++) + const auto& excls = molt->excls; + for (gmx::index ai = 0; ai < excls.ssize(); ai++) { - for (int j = excls->index[ai]; j < excls->index[ai + 1]; j++) + for (const int aj : excls[ai]) { - int aj = excls->a[j]; if (ai != aj) { real rij2 = distance2(cg_cm[ai], cg_cm[aj]); diff --git a/src/gromacs/fileio/tpxio.cpp b/src/gromacs/fileio/tpxio.cpp index 86698e5e82..801174e5e5 100644 --- a/src/gromacs/fileio/tpxio.cpp +++ b/src/gromacs/fileio/tpxio.cpp @@ -2125,19 +2125,25 @@ static void do_block(gmx::ISerializer* serializer, t_block* block) serializer->doIntArray(block->index, block->nr + 1); } -static void do_blocka(gmx::ISerializer* serializer, t_blocka* block) +static void doListOfLists(gmx::ISerializer* serializer, gmx::ListOfLists* listOfLists) { - serializer->doInt(&block->nr); - serializer->doInt(&block->nra); + int numLists = listOfLists->ssize(); + serializer->doInt(&numLists); + int numElements = listOfLists->elementsView().ssize(); + serializer->doInt(&numElements); if (serializer->reading()) { - block->nalloc_index = block->nr + 1; - snew(block->index, block->nalloc_index); - block->nalloc_a = block->nra; - snew(block->a, block->nalloc_a); + std::vector listRanges(numLists + 1); + serializer->doIntArray(listRanges.data(), numLists + 1); + std::vector elements(numElements); + serializer->doIntArray(elements.data(), numElements); + *listOfLists = gmx::ListOfLists(std::move(listRanges), std::move(elements)); + } + else + { + serializer->doIntArray(const_cast(listOfLists->listRangesView().data()), numLists + 1); + serializer->doIntArray(const_cast(listOfLists->elementsView().data()), numElements); } - serializer->doIntArray(block->index, block->nr + 1); - serializer->doIntArray(block->a, block->nra); } /* This is a primitive routine to make it possible to translate atomic numbers @@ -2450,7 +2456,7 @@ static void do_moltype(gmx::ISerializer* serializer, gmx_moltype_t* molt, t_symt sfree(cgs.index); /* This used to be in the atoms struct */ - do_blocka(serializer, &molt->excls); + doListOfLists(serializer, &molt->excls); } static void do_molblock(gmx::ISerializer* serializer, gmx_molblock_t* molb, int numAtomsPerMolecule) diff --git a/src/gromacs/gmxpreprocess/gpp_nextnb.cpp b/src/gromacs/gmxpreprocess/gpp_nextnb.cpp index 19aa460d27..b411d6a602 100644 --- a/src/gromacs/gmxpreprocess/gpp_nextnb.cpp +++ b/src/gromacs/gmxpreprocess/gpp_nextnb.cpp @@ -172,14 +172,15 @@ void __print_nnb(t_nextnb* nnb, char* s) } #endif -static void nnb2excl(t_nextnb* nnb, t_blocka* excl) +static void nnb2excl(t_nextnb* nnb, gmx::ListOfLists* excls) { int i, j, j_index; int nre, nrx, nrs, nr_of_sortables; sortable* s; - srenew(excl->index, nnb->nr + 1); - excl->index[0] = 0; + excls->clear(); + + std::vector exclusionsForAtom; for (i = 0; (i < nnb->nr); i++) { /* calculate the total number of exclusions for atom i */ @@ -230,16 +231,13 @@ static void nnb2excl(t_nextnb* nnb, t_blocka* excl) nr_of_sortables = j_index; prints("after rm-double", j_index, s); - /* make space for arrays */ - srenew(excl->a, excl->nra + nr_of_sortables); - /* put the sorted exclusions in the target list */ + exclusionsForAtom.clear(); for (nrs = 0; (nrs < nr_of_sortables); nrs++) { - excl->a[excl->nra + nrs] = s[nrs].aj; + exclusionsForAtom.push_back(s[nrs].aj); } - excl->nra += nr_of_sortables; - excl->index[i + 1] = excl->nra; + excls->pushBack(exclusionsForAtom); /* cleanup temporary space */ sfree(s); @@ -432,7 +430,7 @@ static void sort_and_purge_nnb(t_nextnb* nnb) } -void generate_excl(int nrexcl, int nratoms, gmx::ArrayRef plist, t_blocka* excl) +void generate_excl(int nrexcl, int nratoms, gmx::ArrayRef plist, gmx::ListOfLists* excls) { t_nextnb nnb; if (nrexcl < 0) @@ -441,8 +439,7 @@ void generate_excl(int nrexcl, int nratoms, gmx::ArrayRef pl } init_nnb(&nnb, nratoms, nrexcl); gen_nnb(&nnb, plist); - excl->nr = nratoms; sort_and_purge_nnb(&nnb); - nnb2excl(&nnb, excl); + nnb2excl(&nnb, excls); done_nnb(&nnb); } diff --git a/src/gromacs/gmxpreprocess/gpp_nextnb.h b/src/gromacs/gmxpreprocess/gpp_nextnb.h index 4d014cdf97..8e640c1bab 100644 --- a/src/gromacs/gmxpreprocess/gpp_nextnb.h +++ b/src/gromacs/gmxpreprocess/gpp_nextnb.h @@ -40,9 +40,14 @@ #include "gromacs/utility/arrayref.h" -struct t_blocka; struct InteractionsOfType; +namespace gmx +{ +template +class ListOfLists; +} + struct t_nextnb { int nr; /* nr atoms (0 <= i < nr) (atoms->nr) */ @@ -75,7 +80,7 @@ void gen_nnb(t_nextnb* nnb, gmx::ArrayRef plist); * initiated using init_nnb. */ -void generate_excl(int nrexcl, int nratoms, gmx::ArrayRef plist, t_blocka* excl); +void generate_excl(int nrexcl, int nratoms, gmx::ArrayRef plist, gmx::ListOfLists* excls); /* Generate an exclusion block from bonds and constraints in * plist. */ diff --git a/src/gromacs/gmxpreprocess/grompp.cpp b/src/gromacs/gmxpreprocess/grompp.cpp index 3922eeb01c..735b5c37d6 100644 --- a/src/gromacs/gmxpreprocess/grompp.cpp +++ b/src/gromacs/gmxpreprocess/grompp.cpp @@ -234,7 +234,7 @@ void InteractionOfType::setForceParameter(int pos, real value) void MoleculeInformation::initMolInfo() { init_block(&mols); - init_blocka(&excls); + excls.clear(); init_t_atoms(&atoms, 0, FALSE); } diff --git a/src/gromacs/gmxpreprocess/grompp_impl.h b/src/gromacs/gmxpreprocess/grompp_impl.h index 1e06a7d6da..7eaa8f1c23 100644 --- a/src/gromacs/gmxpreprocess/grompp_impl.h +++ b/src/gromacs/gmxpreprocess/grompp_impl.h @@ -47,6 +47,7 @@ #include "gromacs/utility/arrayref.h" #include "gromacs/utility/basedefinitions.h" #include "gromacs/utility/exceptions.h" +#include "gromacs/utility/listoflists.h" #include "gromacs/utility/real.h" /*! \libinternal \brief @@ -176,7 +177,7 @@ struct MoleculeInformation //! Molecules separated in datastructure. t_block mols; //! Exclusions in the molecule. - t_blocka excls; + gmx::ListOfLists excls; //! Interactions of a defined type. std::array interactions; diff --git a/src/gromacs/gmxpreprocess/topio.cpp b/src/gromacs/gmxpreprocess/topio.cpp index 7058eb623b..041cb95136 100644 --- a/src/gromacs/gmxpreprocess/topio.cpp +++ b/src/gromacs/gmxpreprocess/topio.cpp @@ -1376,7 +1376,7 @@ void generate_qmexcl(gmx_mtop_t* sys, t_inputrec* ir, warninp* wi, GmxQmmmMode q /* Copy the exclusions to a new array, since this is the only * thing that needs to be modified for QMMM. */ - copy_blocka(&sys->moltype[molb->type].excls, &sys->moltype.back().excls); + sys->moltype.back().excls = sys->moltype[molb->type].excls; /* Set the molecule type for the QMMM molblock */ molb->type = sys->moltype.size() - 1; } diff --git a/src/gromacs/gmxpreprocess/toppush.cpp b/src/gromacs/gmxpreprocess/toppush.cpp index af199cdc79..f78a72af02 100644 --- a/src/gromacs/gmxpreprocess/toppush.cpp +++ b/src/gromacs/gmxpreprocess/toppush.cpp @@ -2615,10 +2615,8 @@ static void convert_pairs_to_pairsQ(gmx::ArrayRef interactio static void generate_LJCpairsNB(MoleculeInformation* mol, int nb_funct, InteractionsOfType* nbp, warninp* wi) { - int n, ntype; - t_atom* atom; - t_blocka* excl; - bool bExcl; + int n, ntype; + t_atom* atom; n = mol->atoms.nr; atom = mol->atoms.atom; @@ -2628,20 +2626,20 @@ static void generate_LJCpairsNB(MoleculeInformation* mol, int nb_funct, Interact "Number of pairs of generated non-bonded parameters should be a perfect square"); /* Add a pair interaction for all non-excluded atom pairs */ - excl = &mol->excls; + const auto& excls = mol->excls; for (int i = 0; i < n; i++) { for (int j = i + 1; j < n; j++) { - bExcl = FALSE; - for (int k = excl->index[i]; k < excl->index[i + 1]; k++) + bool pairIsExcluded = false; + for (const int atomK : excls[i]) { - if (excl->a[k] == j) + if (atomK == j) { - bExcl = TRUE; + pairIsExcluded = true; } } - if (!bExcl) + if (!pairIsExcluded) { if (nb_funct != F_LJ) { @@ -2661,24 +2659,20 @@ static void generate_LJCpairsNB(MoleculeInformation* mol, int nb_funct, Interact } } -static void set_excl_all(t_blocka* excl) +static void set_excl_all(gmx::ListOfLists* excl) { - int nat, i, j, k; - /* Get rid of the current exclusions and exclude all atom pairs */ - nat = excl->nr; - excl->nra = nat * nat; - srenew(excl->a, excl->nra); - k = 0; - for (i = 0; i < nat; i++) + const int numAtoms = excl->ssize(); + std::vector exclusionsForAtom(numAtoms); + for (int i = 0; i < numAtoms; i++) { - excl->index[i] = k; - for (j = 0; j < nat; j++) - { - excl->a[k++] = j; - } + exclusionsForAtom[i] = i; + } + excl->clear(); + for (int i = 0; i < numAtoms; i++) + { + excl->pushBack(exclusionsForAtom); } - excl->index[nat] = k; } static void decouple_atoms(t_atoms* atoms, diff --git a/src/gromacs/mdlib/dispersioncorrection.cpp b/src/gromacs/mdlib/dispersioncorrection.cpp index 7cbd9e41c7..e436bc9d42 100644 --- a/src/gromacs/mdlib/dispersioncorrection.cpp +++ b/src/gromacs/mdlib/dispersioncorrection.cpp @@ -186,17 +186,14 @@ DispersionCorrection::TopologyParams::TopologyParams(const gmx_mtop_t& m */ for (const gmx_molblock_t& molb : mtop.molblock) { - const int nmol = molb.nmol; - const t_atoms* atoms = &mtop.moltype[molb.type].atoms; - const t_blocka* excl = &mtop.moltype[molb.type].excls; + const int nmol = molb.nmol; + const t_atoms* atoms = &mtop.moltype[molb.type].atoms; + const auto& excl = mtop.moltype[molb.type].excls; for (int i = 0; (i < atoms->nr); i++) { const int tpi = atomtypeAOrB(atoms->atom[i], q); - const int j1 = excl->index[i]; - const int j2 = excl->index[i + 1]; - for (int j = j1; j < j2; j++) + for (const int k : excl[i]) { - const int k = excl->a[j]; if (k > i) { const int tpj = atomtypeAOrB(atoms->atom[k], q); diff --git a/src/gromacs/mdlib/forcerec.cpp b/src/gromacs/mdlib/forcerec.cpp index 66c4c33b56..9440cf8f08 100644 --- a/src/gromacs/mdlib/forcerec.cpp +++ b/src/gromacs/mdlib/forcerec.cpp @@ -212,7 +212,7 @@ static std::vector init_cginfo_mb(const gmx_mtop_t* mtop, const t_f { const gmx_molblock_t& molb = mtop->molblock[mb]; const gmx_moltype_t& molt = mtop->moltype[molb.type]; - const t_blocka& excl = molt.excls; + const auto& excl = molt.excls; /* Check if the cginfo is identical for all molecules in this block. * If so, we only need an array of the size of one molecule. @@ -285,9 +285,9 @@ static std::vector init_cginfo_mb(const gmx_mtop_t* mtop, const t_f bool haveExclusions = false; /* Loop over all the exclusions of atom ai */ - for (int j = excl.index[a]; j < excl.index[a + 1]; j++) + for (const int j : excl[a]) { - if (excl.a[j] != a) + if (j != a) { haveExclusions = true; break; diff --git a/src/gromacs/mdlib/perf_est.cpp b/src/gromacs/mdlib/perf_est.cpp index c93de7f105..b2b061061f 100644 --- a/src/gromacs/mdlib/perf_est.cpp +++ b/src/gromacs/mdlib/perf_est.cpp @@ -233,7 +233,7 @@ void count_bonded_distances(const gmx_mtop_t& mtop, const t_inputrec& ir, double } if (bExcl) { - ndtot_c += molb.nmol * (molt->excls.nra - molt->atoms.nr) / 2.; + ndtot_c += molb.nmol * (molt->excls.numElements() - molt->atoms.nr) / 2.; } } diff --git a/src/gromacs/mdlib/sim_util.cpp b/src/gromacs/mdlib/sim_util.cpp index 692c2570b6..8e0541bc86 100644 --- a/src/gromacs/mdlib/sim_util.cpp +++ b/src/gromacs/mdlib/sim_util.cpp @@ -1143,7 +1143,7 @@ void do_force(FILE* fplog, wallcycle_start_nocount(wcycle, ewcNS); wallcycle_sub_start(wcycle, ewcsNBS_SEARCH_LOCAL); /* Note that with a GPU the launch overhead of the list transfer is not timed separately */ - nbv->constructPairlist(InteractionLocality::Local, &top->excls, step, nrnb); + nbv->constructPairlist(InteractionLocality::Local, top->excls, step, nrnb); nbv->setupGpuShortRangeWork(fr->gpuBonded, InteractionLocality::Local); @@ -1242,7 +1242,7 @@ void do_force(FILE* fplog, wallcycle_start_nocount(wcycle, ewcNS); wallcycle_sub_start(wcycle, ewcsNBS_SEARCH_NONLOCAL); /* Note that with a GPU the launch overhead of the list transfer is not timed separately */ - nbv->constructPairlist(InteractionLocality::NonLocal, &top->excls, step, nrnb); + nbv->constructPairlist(InteractionLocality::NonLocal, top->excls, step, nrnb); nbv->setupGpuShortRangeWork(fr->gpuBonded, InteractionLocality::NonLocal); wallcycle_sub_stop(wcycle, ewcsNBS_SEARCH_NONLOCAL); diff --git a/src/gromacs/mdrun/tpi.cpp b/src/gromacs/mdrun/tpi.cpp index 2ed8bd4586..1c9066c7be 100644 --- a/src/gromacs/mdrun/tpi.cpp +++ b/src/gromacs/mdrun/tpi.cpp @@ -666,7 +666,7 @@ void LegacySimulator::do_tpi() /* TODO: Avoid updating all atoms at every bNS step */ fr->nbv->setAtomProperties(*mdatoms, fr->cginfo); - fr->nbv->constructPairlist(InteractionLocality::Local, &top.excls, step, nrnb); + fr->nbv->constructPairlist(InteractionLocality::Local, top.excls, step, nrnb); bNS = FALSE; } diff --git a/src/gromacs/nbnxm/benchmark/bench_setup.cpp b/src/gromacs/nbnxm/benchmark/bench_setup.cpp index f3d6b27526..c2c2328f66 100644 --- a/src/gromacs/nbnxm/benchmark/bench_setup.cpp +++ b/src/gromacs/nbnxm/benchmark/bench_setup.cpp @@ -223,7 +223,7 @@ static std::unique_ptr setupNbnxmForBenchInstance(const Kern { 0, int(system.coordinates.size()) }, atomDensity, atomInfo, system.coordinates, 0, nullptr); - nbv->constructPairlist(gmx::InteractionLocality::Local, &system.excls, 0, &nrnb); + nbv->constructPairlist(gmx::InteractionLocality::Local, system.excls, 0, &nrnb); t_mdatoms mdatoms; // We only use (read) the atom type and charge from mdatoms diff --git a/src/gromacs/nbnxm/benchmark/bench_system.cpp b/src/gromacs/nbnxm/benchmark/bench_system.cpp index fba71d2b9a..8ed07ece85 100644 --- a/src/gromacs/nbnxm/benchmark/bench_system.cpp +++ b/src/gromacs/nbnxm/benchmark/bench_system.cpp @@ -167,10 +167,8 @@ BenchmarkSystem::BenchmarkSystem(const int multiplicationFactor) charges.resize(numAtoms); atomInfoAllVdw.resize(numAtoms); atomInfoOxygenVdw.resize(numAtoms); - snew(excls.index, numAtoms + 1); - snew(excls.a, numAtoms * numAtomsInMolecule); - excls.index[0] = 0; + std::vector exclusionsForAtom; for (int a = 0; a < numAtoms; a++) { if (a % numAtomsInMolecule == 0) @@ -191,12 +189,13 @@ BenchmarkSystem::BenchmarkSystem(const int multiplicationFactor) SET_CGINFO_HAS_Q(atomInfoAllVdw[a]); SET_CGINFO_HAS_Q(atomInfoOxygenVdw[a]); + exclusionsForAtom.clear(); const int firstAtomInMolecule = a - (a % numAtomsInMolecule); for (int aj = 0; aj < numAtomsInMolecule; aj++) { - excls.a[a * numAtomsInMolecule + aj] = firstAtomInMolecule + aj; + exclusionsForAtom.push_back(firstAtomInMolecule + aj); } - excls.index[a + 1] = (a + 1) * numAtomsInMolecule; + excls.pushBack(exclusionsForAtom); } forceRec.ntype = numAtomTypes; diff --git a/src/gromacs/nbnxm/benchmark/bench_system.h b/src/gromacs/nbnxm/benchmark/bench_system.h index 512368ddda..adcc85d4ff 100644 --- a/src/gromacs/nbnxm/benchmark/bench_system.h +++ b/src/gromacs/nbnxm/benchmark/bench_system.h @@ -48,7 +48,7 @@ #include "gromacs/math/vectypes.h" #include "gromacs/mdtypes/forcerec.h" -#include "gromacs/topology/block.h" +#include "gromacs/utility/listoflists.h" #include "gromacs/utility/smalloc.h" namespace gmx @@ -80,7 +80,7 @@ struct BenchmarkSystem //! Atom info where only oxygen atoms are marked to have Van der Waals interactions std::vector atomInfoOxygenVdw; //! Information about exclusions. - t_blocka excls; + ListOfLists excls; //! Storage for atom positions. std::vector coordinates; //! System simulation box. diff --git a/src/gromacs/nbnxm/nbnxm.h b/src/gromacs/nbnxm/nbnxm.h index 1548f3704f..cb5da28eed 100644 --- a/src/gromacs/nbnxm/nbnxm.h +++ b/src/gromacs/nbnxm/nbnxm.h @@ -133,7 +133,6 @@ struct interaction_const_t; struct nonbonded_verlet_t; class PairSearch; class PairlistSets; -struct t_blocka; struct t_commrec; struct t_lambda; struct t_mdatoms; @@ -153,6 +152,8 @@ class GpuEventSynchronizer; namespace gmx { class ForceWithShiftForces; +template +class ListOfLists; class MDLogger; class UpdateGroupsCog; } // namespace gmx @@ -251,7 +252,10 @@ public: gmx::ArrayRef getGridIndices() const; //! Constructs the pairlist for the given locality - void constructPairlist(gmx::InteractionLocality iLocality, const t_blocka* excl, int64_t step, t_nrnb* nrnb); + void constructPairlist(gmx::InteractionLocality iLocality, + const gmx::ListOfLists& exclusions, + int64_t step, + t_nrnb* nrnb); //! Updates all the atom properties in Nbnxm void setAtomProperties(const t_mdatoms& mdatoms, gmx::ArrayRef atomInfo); diff --git a/src/gromacs/nbnxm/pairlist.cpp b/src/gromacs/nbnxm/pairlist.cpp index 6174905d57..f5efbef638 100644 --- a/src/gromacs/nbnxm/pairlist.cpp +++ b/src/gromacs/nbnxm/pairlist.cpp @@ -58,10 +58,10 @@ #include "gromacs/pbcutil/pbc.h" #include "gromacs/simd/simd.h" #include "gromacs/simd/vector_operations.h" -#include "gromacs/topology/block.h" #include "gromacs/utility/exceptions.h" #include "gromacs/utility/fatalerror.h" #include "gromacs/utility/gmxomp.h" +#include "gromacs/utility/listoflists.h" #include "gromacs/utility/smalloc.h" #include "atomdata.h" @@ -1363,12 +1363,12 @@ static nbnxn_sci_t* getOpenIEntry(NbnxnPairlistGpu* nbl) * Set all atom-pair exclusions from the topology stored in exclusions * as masks in the pair-list for simple list entry iEntry. */ -static void setExclusionsForIEntry(const Nbnxm::GridSet& gridSet, - NbnxnPairlistCpu* nbl, - gmx_bool diagRemoved, - int na_cj_2log, - const nbnxn_ci_t& iEntry, - const t_blocka& exclusions) +static void setExclusionsForIEntry(const Nbnxm::GridSet& gridSet, + NbnxnPairlistCpu* nbl, + gmx_bool diagRemoved, + int na_cj_2log, + const nbnxn_ci_t& iEntry, + const ListOfLists& exclusions) { if (iEntry.cj_ind_end == iEntry.cj_ind_start) { @@ -1391,11 +1391,8 @@ static void setExclusionsForIEntry(const Nbnxm::GridSet& gridSet, if (iAtom >= 0) { /* Loop over the topology-based exclusions for this i-atom */ - for (int exclIndex = exclusions.index[iAtom]; exclIndex < exclusions.index[iAtom + 1]; - exclIndex++) + for (const int jAtom : exclusions[iAtom]) { - const int jAtom = exclusions.a[exclIndex]; - if (jAtom == iAtom) { /* The self exclusion are already set, save some time */ @@ -1876,9 +1873,9 @@ static void make_fep_list(gmx::ArrayRef atomIndices, static void setExclusionsForIEntry(const Nbnxm::GridSet& gridSet, NbnxnPairlistGpu* nbl, gmx_bool diagRemoved, - int gmx_unused na_cj_2log, - const nbnxn_sci_t& iEntry, - const t_blocka& exclusions) + int gmx_unused na_cj_2log, + const nbnxn_sci_t& iEntry, + const ListOfLists& exclusions) { if (iEntry.numJClusterGroups() == 0) { @@ -1912,11 +1909,8 @@ static void setExclusionsForIEntry(const Nbnxm::GridSet& gridSet, const int iCluster = i / c_clusterSize; /* Loop over the topology-based exclusions for this i-atom */ - for (int exclIndex = exclusions.index[iAtom]; exclIndex < exclusions.index[iAtom + 1]; - exclIndex++) + for (const int jAtom : exclusions[iAtom]) { - const int jAtom = exclusions.a[exclIndex]; - if (jAtom == iAtom) { /* The self exclusions are already set, save some time */ @@ -3083,7 +3077,7 @@ static void nbnxn_make_pairlist_part(const Nbnxm::GridSet& gridSet, const Grid& jGrid, PairsearchWork* work, const nbnxn_atomdata_t* nbat, - const t_blocka& exclusions, + const ListOfLists& exclusions, real rlist, const PairlistType pairlistType, int ci_block, @@ -3924,7 +3918,7 @@ static void prepareListsForDynamicPruning(gmx::ArrayRef lists) void PairlistSet::constructPairlists(const Nbnxm::GridSet& gridSet, gmx::ArrayRef searchWork, nbnxn_atomdata_t* nbat, - const t_blocka* excl, + const ListOfLists& exclusions, const int minimumIlistCountForGpuBalancing, t_nrnb* nrnb, SearchCycleCounting* searchCycleCounting) @@ -4037,14 +4031,14 @@ void PairlistSet::constructPairlists(const Nbnxm::GridSet& gridSet, /* Divide the i cells equally over the pairlists */ if (isCpuType_) { - nbnxn_make_pairlist_part(gridSet, iGrid, jGrid, &work, nbat, *excl, rlist, + nbnxn_make_pairlist_part(gridSet, iGrid, jGrid, &work, nbat, exclusions, rlist, params_.pairlistType, ci_block, nbat->bUseBufferFlags, nsubpair_target, progBal, nsubpair_tot_est, th, numLists, &cpuLists_[th], fepListPtr); } else { - nbnxn_make_pairlist_part(gridSet, iGrid, jGrid, &work, nbat, *excl, rlist, + nbnxn_make_pairlist_part(gridSet, iGrid, jGrid, &work, nbat, exclusions, rlist, params_.pairlistType, ci_block, nbat->bUseBufferFlags, nsubpair_target, progBal, nsubpair_tot_est, th, numLists, &gpuLists_[th], fepListPtr); @@ -4202,12 +4196,12 @@ void PairlistSet::constructPairlists(const Nbnxm::GridSet& gridSet, void PairlistSets::construct(const InteractionLocality iLocality, PairSearch* pairSearch, nbnxn_atomdata_t* nbat, - const t_blocka* excl, + const ListOfLists& exclusions, const int64_t step, t_nrnb* nrnb) { - pairlistSet(iLocality).constructPairlists(pairSearch->gridSet(), pairSearch->work(), nbat, excl, - minimumIlistCountForGpuBalancing_, nrnb, + pairlistSet(iLocality).constructPairlists(pairSearch->gridSet(), pairSearch->work(), nbat, + exclusions, minimumIlistCountForGpuBalancing_, nrnb, &pairSearch->cycleCounting_); if (iLocality == InteractionLocality::Local) @@ -4234,11 +4228,11 @@ void PairlistSets::construct(const InteractionLocality iLocality, } void nonbonded_verlet_t::constructPairlist(const InteractionLocality iLocality, - const t_blocka* excl, + const ListOfLists& exclusions, int64_t step, t_nrnb* nrnb) { - pairlistSets_->construct(iLocality, pairSearch_.get(), nbat.get(), excl, step, nrnb); + pairlistSets_->construct(iLocality, pairSearch_.get(), nbat.get(), exclusions, step, nrnb); if (useGpu()) { diff --git a/src/gromacs/nbnxm/pairlistset.h b/src/gromacs/nbnxm/pairlistset.h index d0e930330c..1133846a6b 100644 --- a/src/gromacs/nbnxm/pairlistset.h +++ b/src/gromacs/nbnxm/pairlistset.h @@ -62,9 +62,14 @@ struct nbnxn_atomdata_t; struct PairlistParams; struct PairsearchWork; struct SearchCycleCounting; -struct t_blocka; struct t_nrnb; +namespace gmx +{ +template +class ListOfLists; +} + namespace Nbnxm { class GridSet; @@ -85,7 +90,7 @@ public: void constructPairlists(const Nbnxm::GridSet& gridSet, gmx::ArrayRef searchWork, nbnxn_atomdata_t* nbat, - const t_blocka* excl, + const gmx::ListOfLists& exclusions, int minimumIlistCountForGpuBalancing, t_nrnb* nrnb, SearchCycleCounting* searchCycleCounting); diff --git a/src/gromacs/nbnxm/pairlistsets.h b/src/gromacs/nbnxm/pairlistsets.h index 2349159a1d..960f5f9679 100644 --- a/src/gromacs/nbnxm/pairlistsets.h +++ b/src/gromacs/nbnxm/pairlistsets.h @@ -57,9 +57,13 @@ struct nbnxn_atomdata_t; class PairlistSet; enum class PairlistType; class PairSearch; -struct t_blocka; struct t_nrnb; +namespace gmx +{ +template +class ListOfLists; +} class PairlistSets { @@ -69,12 +73,12 @@ public: int minimumIlistCountForGpuBalancing); //! Construct the pairlist set for the given locality - void construct(gmx::InteractionLocality iLocality, - PairSearch* pairSearch, - nbnxn_atomdata_t* nbat, - const t_blocka* excl, - int64_t step, - t_nrnb* nrnb); + void construct(gmx::InteractionLocality iLocality, + PairSearch* pairSearch, + nbnxn_atomdata_t* nbat, + const gmx::ListOfLists& exclusions, + int64_t step, + t_nrnb* nrnb); //! Dispatches the dynamic pruning kernel for the given locality void dispatchPruneKernel(gmx::InteractionLocality iLocality, diff --git a/src/gromacs/selection/nbsearch.cpp b/src/gromacs/selection/nbsearch.cpp index ddc8bd67d3..08d8f6cc7b 100644 --- a/src/gromacs/selection/nbsearch.cpp +++ b/src/gromacs/selection/nbsearch.cpp @@ -62,10 +62,10 @@ #include "gromacs/math/functions.h" #include "gromacs/math/vec.h" #include "gromacs/pbcutil/pbc.h" -#include "gromacs/topology/block.h" #include "gromacs/utility/arrayref.h" #include "gromacs/utility/exceptions.h" #include "gromacs/utility/gmxassert.h" +#include "gromacs/utility/listoflists.h" #include "gromacs/utility/mutex.h" #include "gromacs/utility/stringutil.h" @@ -137,7 +137,7 @@ public: */ void init(AnalysisNeighborhood::SearchMode mode, bool bXY, - const t_blocka* excls, + const ListOfLists* excls, const t_pbc* pbc, const AnalysisNeighborhoodPositions& positions); PairSearchImplPointer getPairSearch(); @@ -280,7 +280,7 @@ private: //! Reference position indices (NULL if no indices). const int* refIndices_; //! Exclusions. - const t_blocka* excls_; + const ListOfLists* excls_; //! PBC data. t_pbc pbc_; @@ -337,8 +337,6 @@ public: testPositions_ = nullptr; testExclusionIds_ = nullptr; testIndices_ = nullptr; - nexcl_ = 0; - excl_ = nullptr; clear_rvec(xtest_); clear_rvec(testcell_); clear_ivec(currCell_); @@ -376,10 +374,8 @@ private: const int* testExclusionIds_; //! Reference to the test position indices. const int* testIndices_; - //! Number of excluded reference positions for current test particle. - int nexcl_; //! Exclusions for current test particle. - const int* excl_; + ArrayRef excl_; //! Index of the currently active test position in \p testPositions_. int testIndex_; //! Stores test position during a pair loop. @@ -862,7 +858,7 @@ int AnalysisNeighborhoodSearchImpl::shiftCell(const ivec cell, rvec shift) const void AnalysisNeighborhoodSearchImpl::init(AnalysisNeighborhood::SearchMode mode, bool bXY, - const t_blocka* excls, + const ListOfLists* excls, const t_pbc* pbc, const AnalysisNeighborhoodPositions& positions) { @@ -988,16 +984,13 @@ void AnalysisNeighborhoodPairSearchImpl::reset(int testIndex) if (search_.excls_ != nullptr) { const int exclIndex = testExclusionIds_[index]; - if (exclIndex < search_.excls_->nr) + if (exclIndex < search_.excls_->ssize()) { - const int startIndex = search_.excls_->index[exclIndex]; - nexcl_ = search_.excls_->index[exclIndex + 1] - startIndex; - excl_ = &search_.excls_->a[startIndex]; + excl_ = (*search_.excls_)[exclIndex]; } else { - nexcl_ = 0; - excl_ = nullptr; + excl_ = ArrayRef(); } } } @@ -1014,15 +1007,16 @@ void AnalysisNeighborhoodPairSearchImpl::nextTestPosition() bool AnalysisNeighborhoodPairSearchImpl::isExcluded(int j) { - if (exclind_ < nexcl_) + const int nexcl = excl_.ssize(); + if (exclind_ < nexcl) { const int index = (search_.refIndices_ != nullptr ? search_.refIndices_[j] : j); const int refId = search_.refExclusionIds_[index]; - while (exclind_ < nexcl_ && excl_[exclind_] < refId) + while (exclind_ < nexcl && excl_[exclind_] < refId) { ++exclind_; } - if (exclind_ < nexcl_ && refId == excl_[exclind_]) + if (exclind_ < nexcl && refId == excl_[exclind_]) { ++exclind_; return true; @@ -1258,12 +1252,12 @@ public: SearchImplPointer getSearch(); - Mutex createSearchMutex_; - SearchList searchList_; - real cutoff_; - const t_blocka* excls_; - SearchMode mode_; - bool bXY_; + Mutex createSearchMutex_; + SearchList searchList_; + real cutoff_; + const ListOfLists* excls_; + SearchMode mode_; + bool bXY_; }; AnalysisNeighborhood::Impl::SearchImplPointer AnalysisNeighborhood::Impl::getSearch() @@ -1304,7 +1298,7 @@ void AnalysisNeighborhood::setXYMode(bool bXY) impl_->bXY_ = bXY; } -void AnalysisNeighborhood::setTopologyExclusions(const t_blocka* excls) +void AnalysisNeighborhood::setTopologyExclusions(const ListOfLists* excls) { GMX_RELEASE_ASSERT(impl_->searchList_.empty(), "Changing the exclusions after initSearch() not currently supported"); diff --git a/src/gromacs/selection/nbsearch.h b/src/gromacs/selection/nbsearch.h index ffe93854b6..7066cafd12 100644 --- a/src/gromacs/selection/nbsearch.h +++ b/src/gromacs/selection/nbsearch.h @@ -59,11 +59,12 @@ #include "gromacs/utility/gmxassert.h" #include "gromacs/utility/real.h" -struct t_blocka; struct t_pbc; namespace gmx { +template +class ListOfLists; namespace internal { @@ -282,7 +283,7 @@ public: * * \see AnalysisNeighborhoodPositions::exclusionIds() */ - void setTopologyExclusions(const t_blocka* excls); + void setTopologyExclusions(const ListOfLists* excls); /*! \brief * Sets the algorithm to use for searching. * diff --git a/src/gromacs/selection/tests/nbsearch.cpp b/src/gromacs/selection/tests/nbsearch.cpp index 7e85995a1d..ec81ab5588 100644 --- a/src/gromacs/selection/tests/nbsearch.cpp +++ b/src/gromacs/selection/tests/nbsearch.cpp @@ -64,6 +64,7 @@ #include "gromacs/random/threefry.h" #include "gromacs/random/uniformrealdistribution.h" #include "gromacs/topology/block.h" +#include "gromacs/utility/listoflists.h" #include "gromacs/utility/smalloc.h" #include "gromacs/utility/stringutil.h" @@ -309,13 +310,13 @@ void NeighborhoodSearchTestData::computeReferencesInternal(t_pbc* pbc, bool bXY) class ExclusionsHelper { public: - static void markExcludedPairs(RefPairList* refPairs, int testIndex, const t_blocka* excls); + static void markExcludedPairs(RefPairList* refPairs, int testIndex, const gmx::ListOfLists* excls); ExclusionsHelper(int refPosCount, int testPosCount); void generateExclusions(); - const t_blocka* exclusions() const { return &excls_; } + const gmx::ListOfLists* exclusions() const { return &excls_; } gmx::ArrayRef refPosIds() const { @@ -327,21 +328,18 @@ public: } private: - int refPosCount_; - int testPosCount_; - std::vector exclusionIds_; - std::vector exclsIndex_; - std::vector exclsAtoms_; - t_blocka excls_; + int refPosCount_; + int testPosCount_; + std::vector exclusionIds_; + gmx::ListOfLists excls_; }; // static -void ExclusionsHelper::markExcludedPairs(RefPairList* refPairs, int testIndex, const t_blocka* excls) +void ExclusionsHelper::markExcludedPairs(RefPairList* refPairs, int testIndex, const gmx::ListOfLists* excls) { int count = 0; - for (int i = excls->index[testIndex]; i < excls->index[testIndex + 1]; ++i) + for (const int excludedIndex : (*excls)[testIndex]) { - const int excludedIndex = excls->a[i]; NeighborhoodSearchTestData::RefPair searchPair(excludedIndex, 0.0); RefPairList::iterator excludedRefPair = std::lower_bound(refPairs->begin(), refPairs->end(), searchPair); @@ -364,13 +362,6 @@ ExclusionsHelper::ExclusionsHelper(int refPosCount, int testPosCount) : exclusionIds_.resize(std::max(refPosCount, testPosCount), 1); exclusionIds_[0] = 0; std::partial_sum(exclusionIds_.begin(), exclusionIds_.end(), exclusionIds_.begin()); - - excls_.nr = 0; - excls_.index = nullptr; - excls_.nra = 0; - excls_.a = nullptr; - excls_.nalloc_index = 0; - excls_.nalloc_a = 0; } void ExclusionsHelper::generateExclusions() @@ -379,21 +370,16 @@ void ExclusionsHelper::generateExclusions() // particles would be higher, or where the exclusions would not be random, // to make a higher percentage of the exclusions to actually be within the // cutoff. - exclsIndex_.reserve(testPosCount_ + 1); - exclsAtoms_.reserve(testPosCount_ * 20); - exclsIndex_.push_back(0); + std::vector exclusionsForAtom; for (int i = 0; i < testPosCount_; ++i) { + exclusionsForAtom.clear(); for (int j = 0; j < 20; ++j) { - exclsAtoms_.push_back(i + j * 3); + exclusionsForAtom.push_back(i + j * 3); } - exclsIndex_.push_back(exclsAtoms_.size()); + excls_.pushBack(exclusionsForAtom); } - excls_.nr = exclsIndex_.size(); - excls_.index = exclsIndex_.data(); - excls_.nra = exclsAtoms_.size(); - excls_.a = exclsAtoms_.data(); } /******************************************************************** @@ -414,7 +400,7 @@ public: void testPairSearchFull(gmx::AnalysisNeighborhoodSearch* search, const NeighborhoodSearchTestData& data, const gmx::AnalysisNeighborhoodPositions& pos, - const t_blocka* excls, + const gmx::ListOfLists* excls, const gmx::ArrayRef& refIndices, const gmx::ArrayRef& testIndices, bool selfPairs); @@ -522,7 +508,7 @@ void NeighborhoodSearchTest::testPairSearchIndexed(gmx::AnalysisNeighborhood* void NeighborhoodSearchTest::testPairSearchFull(gmx::AnalysisNeighborhoodSearch* search, const NeighborhoodSearchTestData& data, const gmx::AnalysisNeighborhoodPositions& pos, - const t_blocka* excls, + const gmx::ListOfLists* excls, const gmx::ArrayRef& refIndices, const gmx::ArrayRef& testIndices, bool selfPairs) diff --git a/src/gromacs/tools/convert_tpr.cpp b/src/gromacs/tools/convert_tpr.cpp index 09adddc795..d4d4685a58 100644 --- a/src/gromacs/tools/convert_tpr.cpp +++ b/src/gromacs/tools/convert_tpr.cpp @@ -129,39 +129,33 @@ static void reduce_block(const gmx_bool bKeep[], t_block* block, const char* nam block->nr = newi; } -static void reduce_blocka(const int invindex[], const gmx_bool bKeep[], t_blocka* block, const char* name) +static gmx::ListOfLists +reduce_blocka(const int invindex[], const gmx_bool bKeep[], const t_blocka& block, const char* name) { - int *index, *a; - int i, j, k, newi, newj; + gmx::ListOfLists lists; - snew(index, block->nr); - snew(a, block->nra); - - newi = newj = 0; - for (i = 0; (i < block->nr); i++) + std::vector exclusionsForAtom; + for (int i = 0; i < block.nr; i++) { - for (j = block->index[i]; (j < block->index[i + 1]); j++) + if (bKeep[i]) { - k = block->a[j]; - if (bKeep[k]) + exclusionsForAtom.clear(); + for (int j = block.index[i]; j < block.index[i + 1]; j++) { - a[newj] = invindex[k]; - newj++; + const int k = block.a[j]; + if (bKeep[k]) + { + exclusionsForAtom.push_back(invindex[k]); + } } - } - if (newj > index[newi]) - { - newi++; - index[newi] = newj; + lists.pushBack(exclusionsForAtom); } } - fprintf(stderr, "Reduced block %8s from %6d to %6d index-, %6d to %6d a-entries\n", name, - block->nr, newi, block->nra, newj); - block->index = index; - block->a = a; - block->nr = newi; - block->nra = newj; + fprintf(stderr, "Reduced block %8s from %6d to %6zu index-, %6d to %6d a-entries\n", name, + block.nr, lists.size(), block.nra, lists.numElements()); + + return lists; } static void reduce_rvec(int gnx, const int index[], rvec vv[]) @@ -271,7 +265,7 @@ static void reduce_topology_x(int gnx, int index[], gmx_mtop_t* mtop, rvec x[], invindex = invind(gnx, top.atoms.nr, index); reduce_block(bKeep, &(top.mols), "mols"); - reduce_blocka(invindex, bKeep, &(top.excls), "excls"); + gmx::ListOfLists excls = reduce_blocka(invindex, bKeep, top.excls, "excls"); reduce_rvec(gnx, index, x); reduce_rvec(gnx, index, v); reduce_atom(gnx, index, top.atoms.atom, top.atoms.atomname, &(top.atoms.nres), top.atoms.resinfo); @@ -297,7 +291,7 @@ static void reduce_topology_x(int gnx, int index[], gmx_mtop_t* mtop, rvec x[], } } mtop->moltype[0].atoms = top.atoms; - mtop->moltype[0].excls = top.excls; + mtop->moltype[0].excls = excls; mtop->molblock.resize(1); mtop->molblock[0].type = 0; diff --git a/src/gromacs/topology/block.cpp b/src/gromacs/topology/block.cpp index 7c1cc84357..f0a91414fe 100644 --- a/src/gromacs/topology/block.cpp +++ b/src/gromacs/topology/block.cpp @@ -42,6 +42,7 @@ #include +#include "gromacs/utility/listoflists.h" #include "gromacs/utility/smalloc.h" #include "gromacs/utility/txtdump.h" @@ -218,6 +219,19 @@ static int pr_blocka_title(FILE* fp, int indent, const char* title, const t_bloc return indent; } +static int pr_listoflists_title(FILE* fp, int indent, const char* title, const gmx::ListOfLists* lists) +{ + if (available(fp, lists, indent, title)) + { + indent = pr_title(fp, indent, title); + pr_indent(fp, indent); + fprintf(fp, "numLists=%zu\n", lists->size()); + pr_indent(fp, indent); + fprintf(fp, "numElements=%d\n", lists->numElements()); + } + return indent; +} + static void low_pr_blocka(FILE* fp, int indent, const char* title, const t_blocka* block, gmx_bool bShowNumbers) { int i; @@ -325,6 +339,43 @@ void pr_blocka(FILE* fp, int indent, const char* title, const t_blocka* block, g } } +void pr_listoflists(FILE* fp, int indent, const char* title, const gmx::ListOfLists* lists, gmx_bool bShowNumbers) +{ + if (available(fp, lists, indent, title)) + { + indent = pr_listoflists_title(fp, indent, title, lists); + for (gmx::index i = 0; i < lists->ssize(); i++) + { + int size = pr_indent(fp, indent); + gmx::ArrayRef list = (*lists)[i]; + if (list.empty()) + { + size += fprintf(fp, "%s[%d]={", title, int(i)); + } + else + { + size += fprintf(fp, "%s[%d][num=%zu]={", title, bShowNumbers ? int(i) : -1, list.size()); + } + bool isFirst = true; + for (const int j : list) + { + if (!isFirst) + { + size += fprintf(fp, ", "); + } + if ((size) > (USE_WIDTH)) + { + fprintf(fp, "\n"); + size = pr_indent(fp, indent + INDENT); + } + size += fprintf(fp, "%d", j); + isFirst = false; + } + fprintf(fp, "}\n"); + } + } +} + void copy_block(const t_block* src, t_block* dst) { dst->nr = src->nr; diff --git a/src/gromacs/topology/block.h b/src/gromacs/topology/block.h index 8db0164b90..fe0a30a56a 100644 --- a/src/gromacs/topology/block.h +++ b/src/gromacs/topology/block.h @@ -48,6 +48,9 @@ namespace gmx { +template +class ListOfLists; + /*! \brief Division of a range of indices into consecutive blocks * * A range of consecutive indices 0 to full.range.end() is divided @@ -217,5 +220,6 @@ void stupid_fill_blocka(t_blocka* grp, int natom); void pr_block(FILE* fp, int indent, const char* title, const t_block* block, gmx_bool bShowNumbers); void pr_blocka(FILE* fp, int indent, const char* title, const t_blocka* block, gmx_bool bShowNumbers); +void pr_listoflists(FILE* fp, int indent, const char* title, const gmx::ListOfLists* block, gmx_bool bShowNumbers); #endif diff --git a/src/gromacs/topology/exclusionblocks.cpp b/src/gromacs/topology/exclusionblocks.cpp index 9472c43781..b73b68097e 100644 --- a/src/gromacs/topology/exclusionblocks.cpp +++ b/src/gromacs/topology/exclusionblocks.cpp @@ -42,12 +42,41 @@ #include "gromacs/topology/block.h" #include "gromacs/utility/fatalerror.h" +#include "gromacs/utility/listoflists.h" #include "gromacs/utility/smalloc.h" #include "gromacs/utility/stringutil.h" namespace gmx { +namespace +{ + +//! Converts ListOfLists to a list of ExclusionBlocks +void listOfListsToExclusionBlocks(const ListOfLists& b, gmx::ArrayRef b2) +{ + for (gmx::index i = 0; i < b.ssize(); i++) + { + for (int jAtom : b[i]) + { + b2[i].atomNumber.push_back(jAtom); + } + } +} + +//! Converts a list of ExclusionBlocks to ListOfLists +void exclusionBlocksToListOfLists(gmx::ArrayRef b2, ListOfLists* b) +{ + b->clear(); + + for (const auto& block : b2) + { + b->pushBack(block.atomNumber); + } +} + +} // namespace + void blockaToExclusionBlocks(const t_blocka* b, gmx::ArrayRef b2) { for (int i = 0; (i < b->nr); i++) @@ -79,21 +108,12 @@ void exclusionBlocksToBlocka(gmx::ArrayRef b2, t_blocka* b b->index[i] = nra; } -void mergeExclusions(t_blocka* excl, gmx::ArrayRef b2) +namespace { - if (b2.empty()) - { - return; - } - GMX_RELEASE_ASSERT(b2.ssize() == excl->nr, - "Cannot merge exclusions for " - "blocks that do not describe the same number " - "of particles"); - /* Convert the t_blocka entries to ExclusionBlock form */ - blockaToExclusionBlocks(excl, b2); - - /* Count and sort the exclusions */ +//! Counts and sorts the exclusions +int countAndSortExclusions(gmx::ArrayRef b2) +{ int nra = 0; for (auto& block : b2) { @@ -118,10 +138,29 @@ void mergeExclusions(t_blocka* excl, gmx::ArrayRef b2) nra += block.nra(); } } - excl->nra = nra; - srenew(excl->a, excl->nra); - exclusionBlocksToBlocka(b2, excl); + return nra; +} + +} // namespace + +void mergeExclusions(ListOfLists* excl, gmx::ArrayRef b2) +{ + if (b2.empty()) + { + return; + } + GMX_RELEASE_ASSERT(b2.ssize() == excl->ssize(), + "Cannot merge exclusions for " + "blocks that do not describe the same number " + "of particles"); + + /* Convert the t_blocka entries to ExclusionBlock form */ + listOfListsToExclusionBlocks(*excl, b2); + + countAndSortExclusions(b2); + + exclusionBlocksToListOfLists(b2, excl); } } // namespace gmx diff --git a/src/gromacs/topology/exclusionblocks.h b/src/gromacs/topology/exclusionblocks.h index b27db36093..69c9c54c21 100644 --- a/src/gromacs/topology/exclusionblocks.h +++ b/src/gromacs/topology/exclusionblocks.h @@ -44,6 +44,8 @@ struct t_blocka; namespace gmx { +template +class ListOfLists; /*! \libinternal \brief * Describes exclusions for a single atom. @@ -61,7 +63,7 @@ struct ExclusionBlock * Requires that \c b2 and \c excl describe the same number of * particles, if \c b2 describes a non-zero number. */ -void mergeExclusions(t_blocka* excl, gmx::ArrayRef b2); +void mergeExclusions(ListOfLists* excl, gmx::ArrayRef b2); /*! \brief * Convert the exclusions. diff --git a/src/gromacs/topology/mtop_util.cpp b/src/gromacs/topology/mtop_util.cpp index 4bc2e2946e..541788aa60 100644 --- a/src/gromacs/topology/mtop_util.cpp +++ b/src/gromacs/topology/mtop_util.cpp @@ -650,39 +650,41 @@ t_atoms gmx_mtop_global_atoms(const gmx_mtop_t* mtop) * The cat routines below are old code from src/kernel/topcat.c */ -static void blockacat(t_blocka* dest, const t_blocka* src, int copies, int dnum, int snum) +static void blockacat(t_blocka* dest, const gmx::ListOfLists& src, int copies, int dnum, int snum) { - int i, j, l, size; + int j, l; int destnr = dest->nr; int destnra = dest->nra; - if (src->nr) + if (!src.empty()) { - size = (dest->nr + copies * src->nr + 1); + size_t size = (dest->nr + copies * src.size() + 1); srenew(dest->index, size); } - if (src->nra) + gmx::ArrayRef srcElements = src.elementsView(); + if (!srcElements.empty()) { - size = (dest->nra + copies * src->nra); + size_t size = (dest->nra + copies * srcElements.size()); srenew(dest->a, size); } + gmx::ArrayRef srcListRanges = src.listRangesView(); for (l = destnr, j = 0; (j < copies); j++) { - for (i = 0; (i < src->nr); i++) + for (gmx::index i = 0; i < src.ssize(); i++) { - dest->index[l++] = dest->nra + src->index[i]; + dest->index[l++] = dest->nra + srcListRanges[i]; } - dest->nra += src->nra; + dest->nra += src.ssize(); } for (l = destnra, j = 0; (j < copies); j++) { - for (i = 0; (i < src->nra); i++) + for (const int srcElement : srcElements) { - dest->a[l++] = dnum + src->a[i]; + dest->a[l++] = dnum + srcElement; } dnum += snum; - dest->nr += src->nr; + dest->nr += src.ssize(); } dest->index[dest->nr] = dest->nra; dest->nalloc_index = dest->nr; @@ -941,6 +943,31 @@ static void copyAtomtypesFromMtop(const gmx_mtop_t& mtop, t_atomtypes* atomtypes * \param[in] mtop Reference to input mtop. * \param[in] excls Pointer to final excls data structure. */ +static void copyExclsFromMtop(const gmx_mtop_t& mtop, gmx::ListOfLists* excls) +{ + excls->clear(); + int atomIndex = 0; + for (const gmx_molblock_t& molb : mtop.molblock) + { + const gmx_moltype_t& molt = mtop.moltype[molb.type]; + + for (int mol = 0; mol < molb.nmol; mol++) + { + excls->appendListOfLists(molt.excls, atomIndex); + + atomIndex += molt.atoms.nr; + } + } +} + +/*! \brief Copy excls from mtop. + * + * Makes a deep copy of excls(t_blocka) from gmx_mtop_t. + * Used to initialize legacy t_topology. + * + * \param[in] mtop Reference to input mtop. + * \param[in] excls Pointer to final excls data structure. + */ static void copyExclsFromMtop(const gmx_mtop_t& mtop, t_blocka* excls) { init_blocka(excls); @@ -952,7 +979,7 @@ static void copyExclsFromMtop(const gmx_mtop_t& mtop, t_blocka* excls) int srcnr = molt.atoms.nr; int destnr = natoms; - blockacat(excls, &molt.excls, molb.nmol, destnr, srcnr); + blockacat(excls, molt.excls, molb.nmol, destnr, srcnr); natoms += molb.nmol * srcnr; } @@ -966,21 +993,21 @@ static void copyExclsFromMtop(const gmx_mtop_t& mtop, t_blocka* excls) * \param[inout] excls existing exclusions in local topology * \param[in] ids list of global IDs of atoms */ -static void addMimicExclusions(t_blocka* excls, const gmx::ArrayRef ids) +static void addMimicExclusions(gmx::ListOfLists* excls, const gmx::ArrayRef ids) { t_blocka inter_excl{}; init_blocka(&inter_excl); size_t n_q = ids.size(); - inter_excl.nr = excls->nr; + inter_excl.nr = excls->ssize(); inter_excl.nra = n_q * n_q; size_t total_nra = n_q * n_q; - snew(inter_excl.index, excls->nr + 1); + snew(inter_excl.index, excls->ssize() + 1); snew(inter_excl.a, total_nra); - for (int i = 0; i < excls->nr; ++i) + for (int i = 0; i < inter_excl.nr; ++i) { inter_excl.index[i] = 0; } @@ -1012,7 +1039,7 @@ static void addMimicExclusions(t_blocka* excls, const gmx::ArrayRef i inter_excl.index[inter_excl.nr] = n_q * n_q; - std::vector qmexcl2(excls->nr); + std::vector qmexcl2(excls->size()); gmx::blockaToExclusionBlocks(&inter_excl, qmexcl2); // Merge the created exclusion list with the existing one diff --git a/src/gromacs/topology/tests/exclusionblocks.cpp b/src/gromacs/topology/tests/exclusionblocks.cpp index ab9f1ad122..53a2a05fbc 100644 --- a/src/gromacs/topology/tests/exclusionblocks.cpp +++ b/src/gromacs/topology/tests/exclusionblocks.cpp @@ -47,6 +47,7 @@ #include "gromacs/topology/block.h" #include "gromacs/utility/arrayref.h" +#include "gromacs/utility/listoflists.h" #include "gromacs/utility/smalloc.h" #include "testutils/cmdlinetest.h" @@ -101,6 +102,21 @@ void makeTestBlockAData(t_blocka* ba) addGroupToBlocka(ba, indices); } +//! Return ListOfLists filled with some datastructures +ListOfLists makeTestListOfLists() +{ + ListOfLists list; + + std::vector indices = { 12, 11, 9, 6, 2 }; + list.pushBack(indices); + indices = { 10, 8, 5, 1 }; + list.pushBack(indices); + indices = { 7, 4, 0 }; + list.pushBack(indices); + + return list; +} + class ExclusionBlockTest : public ::testing::Test { public: @@ -108,6 +124,7 @@ public: { const int natom = 3; makeTestBlockAData(&ba_); + list_ = makeTestListOfLists(); b_.resize(natom); } ~ExclusionBlockTest() override { done_blocka(&ba_); } @@ -126,8 +143,23 @@ public: } } + void compareBlocksAndList() + { + GMX_RELEASE_ASSERT(ssize(b_) == list_.ssize(), "The list counts should match"); + for (index i = 0; i < ssize(b_); i++) + { + gmx::ArrayRef jList = list_[i]; + EXPECT_EQ(b_[i].nra(), jList.ssize()) << "Block size mismatch at " << i << "."; + for (int j = 0; j < b_[i].nra(); j++) + { + EXPECT_EQ(b_[i].atomNumber[j], jList[j]) << "Block mismatch at " << i << " , " << j << "."; + } + } + } + protected: t_blocka ba_; + ListOfLists list_; std::vector b_; }; @@ -148,8 +180,8 @@ TEST_F(ExclusionBlockTest, ConvertExclusionBlockToBlocka) TEST_F(ExclusionBlockTest, MergeExclusions) { - mergeExclusions(&ba_, b_); - compareBlocks(); + mergeExclusions(&list_, b_); + compareBlocksAndList(); } } // namespace diff --git a/src/gromacs/topology/topology.cpp b/src/gromacs/topology/topology.cpp index e01993300b..3a1b82eb50 100644 --- a/src/gromacs/topology/topology.cpp +++ b/src/gromacs/topology/topology.cpp @@ -76,7 +76,7 @@ void init_top(t_topology* top) } -gmx_moltype_t::gmx_moltype_t() : name(nullptr), excls() +gmx_moltype_t::gmx_moltype_t() : name(nullptr) { init_t_atoms(&atoms, 0, FALSE); } @@ -84,7 +84,6 @@ gmx_moltype_t::gmx_moltype_t() : name(nullptr), excls() gmx_moltype_t::~gmx_moltype_t() { done_atom(&atoms); - done_blocka(&excls); } gmx_mtop_t::gmx_mtop_t() @@ -136,7 +135,6 @@ void done_top_mtop(t_topology* top, gmx_mtop_t* mtop) gmx_localtop_t::gmx_localtop_t() { - init_blocka_null(&excls); init_idef(&idef); init_atomtypes(&atomtypes); } @@ -146,7 +144,6 @@ gmx_localtop_t::~gmx_localtop_t() if (!useInDomainDecomp_) { done_idef(&idef); - done_blocka(&excls); done_atomtypes(&atomtypes); } } @@ -287,7 +284,7 @@ static void pr_moltype(FILE* fp, pr_indent(fp, indent); fprintf(fp, "name=\"%s\"\n", *(molt->name)); pr_atoms(fp, indent, "atoms", &(molt->atoms), bShowNumbers); - pr_blocka(fp, indent, "excls", &molt->excls, bShowNumbers); + pr_listoflists(fp, indent, "excls", &molt->excls, bShowNumbers); for (j = 0; (j < F_NRE); j++) { pr_ilist(fp, indent, interaction_function[j].longname, ffparams->functype.data(), @@ -457,15 +454,18 @@ static void cmp_cmap(FILE* fp, const gmx_cmap_t* cmap1, const gmx_cmap_t* cmap2, } } -static void cmp_blocka(FILE* fp, const t_blocka* b1, const t_blocka* b2, const char* s) +static void cmp_listoflists(FILE* fp, + const gmx::ListOfLists& list1, + const gmx::ListOfLists& list2, + const char* s) { char buf[32]; fprintf(fp, "comparing blocka %s\n", s); - sprintf(buf, "%s.nr", s); - cmp_int(fp, buf, -1, b1->nr, b2->nr); - sprintf(buf, "%s.nra", s); - cmp_int(fp, buf, -1, b1->nra, b2->nra); + sprintf(buf, "%s.numLists", s); + cmp_int(fp, buf, -1, list1.ssize(), list2.ssize()); + sprintf(buf, "%s.numElements", s); + cmp_int(fp, buf, -1, list1.numElements(), list2.numElements()); } static void compareFfparams(FILE* fp, @@ -534,7 +534,7 @@ static void compareMoltypes(FILE* fp, compareAtoms(fp, &mt1[i].atoms, &mt2[i].atoms, relativeTolerance, absoluteTolerance); compareInteractionLists(fp, &mt1[i].ilist, &mt2[i].ilist); std::string buf = gmx::formatString("excls[%d]", i); - cmp_blocka(fp, &mt1[i].excls, &mt2[i].excls, buf.c_str()); + cmp_listoflists(fp, mt1[i].excls, mt2[i].excls, buf.c_str()); } } @@ -674,8 +674,8 @@ int getGroupType(const SimulationGroups& group, SimulationAtomGroupType type, in void copy_moltype(const gmx_moltype_t* src, gmx_moltype_t* dst) { - dst->name = src->name; - copy_blocka(&src->excls, &dst->excls); + dst->name = src->name; + dst->excls = src->excls; t_atoms* atomsCopy = copy_t_atoms(&src->atoms); dst->atoms = *atomsCopy; sfree(atomsCopy); diff --git a/src/gromacs/topology/topology.h b/src/gromacs/topology/topology.h index 710bbbd2d8..0e5c0d9e17 100644 --- a/src/gromacs/topology/topology.h +++ b/src/gromacs/topology/topology.h @@ -48,6 +48,7 @@ #include "gromacs/topology/idef.h" #include "gromacs/topology/symtab.h" #include "gromacs/utility/enumerationhelpers.h" +#include "gromacs/utility/listoflists.h" #include "gromacs/utility/unique_cptr.h" enum class SimulationAtomGroupType : int @@ -83,10 +84,10 @@ struct gmx_moltype_t /*! \brief Default copy constructor */ gmx_moltype_t(const gmx_moltype_t&) = default; - char** name; /**< Name of the molecule type */ - t_atoms atoms; /**< The atoms in this molecule */ - InteractionLists ilist; /**< Interaction list with local indices */ - t_blocka excls; /**< The exclusions */ + char** name; /**< Name of the molecule type */ + t_atoms atoms; /**< The atoms in this molecule */ + InteractionLists ilist; /**< Interaction list with local indices */ + gmx::ListOfLists excls; /**< The exclusions */ }; /*! \brief Block of molecules of the same type, used in gmx_mtop_t */ @@ -214,7 +215,7 @@ struct gmx_localtop_t //! Atomtype properties t_atomtypes atomtypes; //! The exclusions - t_blocka excls; + gmx::ListOfLists excls; //! Flag for domain decomposition so we don't free already freed memory. bool useInDomainDecomp_ = false; }; diff --git a/src/gromacs/trajectoryanalysis/modules/rdf.cpp b/src/gromacs/trajectoryanalysis/modules/rdf.cpp index cc45dadebc..f9a58d2f2f 100644 --- a/src/gromacs/trajectoryanalysis/modules/rdf.cpp +++ b/src/gromacs/trajectoryanalysis/modules/rdf.cpp @@ -386,7 +386,7 @@ void Rdf::initAnalysis(const TrajectoryAnalysisSettings& settings, const Topolog } } localTop_ = top.expandedTopology(); - if (localTop_->excls.nr == 0) + if (localTop_->excls.empty()) { GMX_THROW(InconsistentInputError( "-excl is set, but the file provided to -s does not define exclusions")); -- 2.22.0