From 5018b5a0f165076593a535f6160749298de5de6d Mon Sep 17 00:00:00 2001 From: Berk Hess Date: Tue, 13 Feb 2018 14:30:22 +0100 Subject: [PATCH] Convert gmx_mtop_t to C++ gmx_mtop_t now uses std vectors for moltype and molblock. gmx_mtop_t, gmx_moltype_t and gmx_molblock_t now have constructors and destructors. Changed declarations of gmx_mtop_t* to not use a pointer or to use std::unique_ptr. Several sub-structs still use C style pointers and will be converted later. Change-Id: Iee802e4382a0a389496bb9395389e5926918d0f8 --- src/gromacs/domdec/domdec.cpp | 2 +- src/gromacs/domdec/domdec_constraints.cpp | 11 +- src/gromacs/domdec/domdec_topology.cpp | 207 +++++++-------- src/gromacs/fileio/confio.cpp | 11 +- src/gromacs/fileio/tngio.cpp | 37 ++- src/gromacs/fileio/tpxio.cpp | 26 +- src/gromacs/gmxana/gmx_clustsize.cpp | 5 +- src/gromacs/gmxana/gmx_disre.cpp | 17 +- src/gromacs/gmxana/gmx_energy.cpp | 1 - src/gromacs/gmxana/gmx_nmeig.cpp | 5 +- src/gromacs/gmxana/gmx_pme_error.cpp | 28 +- src/gromacs/gmxana/gmx_trjconv.cpp | 23 +- src/gromacs/gmxlib/chargegroup.cpp | 21 +- src/gromacs/gmxpreprocess/convparm.cpp | 4 +- src/gromacs/gmxpreprocess/gpp_atomtype.cpp | 6 +- src/gromacs/gmxpreprocess/grompp.cpp | 290 ++++++++++----------- src/gromacs/gmxpreprocess/insert-molecules.cpp | 31 +-- src/gromacs/gmxpreprocess/readir.cpp | 29 +-- src/gromacs/gmxpreprocess/topio.cpp | 120 ++++----- src/gromacs/gmxpreprocess/topio.h | 39 +-- src/gromacs/gmxpreprocess/toputil.cpp | 2 +- src/gromacs/listed-forces/disre.cpp | 2 +- src/gromacs/listed-forces/orires.cpp | 2 +- src/gromacs/mdlib/broadcaststructs.cpp | 20 +- src/gromacs/mdlib/calc_verletbuf.cpp | 19 +- src/gromacs/mdlib/constr.cpp | 20 +- src/gromacs/mdlib/constraintrange.cpp | 9 +- src/gromacs/mdlib/forcerec.cpp | 80 +++--- src/gromacs/mdlib/lincs.cpp | 20 +- src/gromacs/mdlib/ns.cpp | 7 +- src/gromacs/mdlib/perf_est.cpp | 50 ++-- src/gromacs/mdlib/qmmm.cpp | 2 +- src/gromacs/mdlib/rf_util.cpp | 10 +- src/gromacs/mdlib/settle.cpp | 2 +- src/gromacs/mdlib/shellfc.cpp | 21 +- src/gromacs/mdlib/shellfc.h | 2 +- src/gromacs/mdlib/sim_util.cpp | 20 +- src/gromacs/mdlib/tests/settle.cpp | 46 ++-- src/gromacs/mdlib/vsite.cpp | 28 +- src/gromacs/selection/indexutil.cpp | 4 +- src/gromacs/selection/selectionoptionbehavior.cpp | 5 +- src/gromacs/selection/tests/toputils.cpp | 26 +- src/gromacs/selection/tests/toputils.h | 11 +- src/gromacs/tools/convert_tpr.cpp | 18 +- src/gromacs/topology/block.cpp | 5 +- src/gromacs/topology/mtop_lookup.h | 4 +- src/gromacs/topology/mtop_util.cpp | 255 ++++++++---------- src/gromacs/topology/mtop_util.h | 8 +- src/gromacs/topology/topology.cpp | 139 ++++++---- src/gromacs/topology/topology.h | 91 ++++--- src/gromacs/topology/topsort.cpp | 28 +- .../trajectoryanalysis/analysissettings.cpp | 7 +- src/gromacs/trajectoryanalysis/analysissettings.h | 6 +- src/gromacs/trajectoryanalysis/runnercommon.cpp | 12 +- src/programs/mdrun/membed.cpp | 39 ++- src/programs/mdrun/runner.cpp | 48 ++-- 56 files changed, 925 insertions(+), 1056 deletions(-) diff --git a/src/gromacs/domdec/domdec.cpp b/src/gromacs/domdec/domdec.cpp index 801fe558eb..6f6c9f6d62 100644 --- a/src/gromacs/domdec/domdec.cpp +++ b/src/gromacs/domdec/domdec.cpp @@ -6069,7 +6069,7 @@ static int multi_body_bondeds_count(const gmx_mtop_t *mtop) { int n, nmol, ftype; gmx_mtop_ilistloop_t iloop; - t_ilist *il; + const t_ilist *il; n = 0; iloop = gmx_mtop_ilistloop_init(mtop); diff --git a/src/gromacs/domdec/domdec_constraints.cpp b/src/gromacs/domdec/domdec_constraints.cpp index 4174a0bbd0..1b2d141b09 100644 --- a/src/gromacs/domdec/domdec_constraints.cpp +++ b/src/gromacs/domdec/domdec_constraints.cpp @@ -651,7 +651,6 @@ void init_domdec_constraints(gmx_domdec_t *dd, { gmx_domdec_constraints_t *dc; const gmx_molblock_t *molb; - int mb, ncon, c; if (debug) { @@ -661,11 +660,11 @@ void init_domdec_constraints(gmx_domdec_t *dd, snew(dd->constraints, 1); dc = dd->constraints; - snew(dc->molb_con_offset, mtop->nmolblock); - snew(dc->molb_ncon_mol, mtop->nmolblock); + snew(dc->molb_con_offset, mtop->molblock.size()); + snew(dc->molb_ncon_mol, mtop->molblock.size()); - ncon = 0; - for (mb = 0; mb < mtop->nmolblock; mb++) + int ncon = 0; + for (size_t mb = 0; mb < mtop->molblock.size(); mb++) { molb = &mtop->molblock[mb]; dc->molb_con_offset[mb] = ncon; @@ -678,7 +677,7 @@ void init_domdec_constraints(gmx_domdec_t *dd, if (ncon > 0) { snew(dc->gc_req, ncon); - for (c = 0; c < ncon; c++) + for (int c = 0; c < ncon; c++) { dc->gc_req[c] = 0; } diff --git a/src/gromacs/domdec/domdec_topology.cpp b/src/gromacs/domdec/domdec_topology.cpp index b45f2127b9..e6fdd8277f 100644 --- a/src/gromacs/domdec/domdec_topology.cpp +++ b/src/gromacs/domdec/domdec_topology.cpp @@ -314,25 +314,20 @@ static void print_missing_interactions_atoms(FILE *fplog, t_commrec *cr, const gmx_mtop_t *mtop, const t_idef *idef) { - int mb, a_start, a_end; - const gmx_molblock_t *molb; - const gmx_reverse_top_t *rt; - - rt = cr->dd->reverse_top; + const gmx_reverse_top_t *rt = cr->dd->reverse_top; /* Print the atoms in the missing interactions per molblock */ - a_end = 0; - for (mb = 0; mb < mtop->nmolblock; mb++) + int a_end = 0; + for (const gmx_molblock_t &molb : mtop->molblock) { - molb = &mtop->molblock[mb]; - a_start = a_end; - a_end = a_start + molb->nmol*molb->natoms_mol; + int a_start = a_end; + a_end = a_start + molb.nmol*molb.natoms_mol; print_missing_interactions_mb(fplog, cr, rt, - *(mtop->moltype[molb->type].name), - &rt->ril_mt[molb->type], - a_start, a_end, molb->natoms_mol, - molb->nmol, + *(mtop->moltype[molb.type].name), + &rt->ril_mt[molb.type], + a_start, a_end, molb.natoms_mol, + molb.nmol, idef); } } @@ -666,10 +661,8 @@ static gmx_reverse_top_t *make_reverse_top(const gmx_mtop_t *mtop, gmx_bool bFE, gmx_bool bConstr, gmx_bool bSettle, gmx_bool bBCheck, int *nint) { - int mt, i, mb; gmx_reverse_top_t *rt; int *nint_mt; - gmx_moltype_t *molt; int thread; snew(rt, 1); @@ -680,25 +673,25 @@ static gmx_reverse_top_t *make_reverse_top(const gmx_mtop_t *mtop, gmx_bool bFE, rt->bBCheck = bBCheck; rt->bInterCGInteractions = mtop->bIntermolecularInteractions; - snew(nint_mt, mtop->nmoltype); - snew(rt->ril_mt, mtop->nmoltype); + snew(nint_mt, mtop->moltype.size()); + snew(rt->ril_mt, mtop->moltype.size()); rt->ril_mt_tot_size = 0; - for (mt = 0; mt < mtop->nmoltype; mt++) + for (size_t mt = 0; mt < mtop->moltype.size(); mt++) { - molt = &mtop->moltype[mt]; - if (molt->cgs.nr > 1) + const gmx_moltype_t &molt = mtop->moltype[mt]; + if (molt.cgs.nr > 1) { rt->bInterCGInteractions = TRUE; } /* Make the atom to interaction list for this molecule type */ nint_mt[mt] = - make_reverse_ilist(molt->ilist, &molt->atoms, + make_reverse_ilist(molt.ilist, &molt.atoms, vsite_pbc_molt ? vsite_pbc_molt[mt] : nullptr, rt->bConstr, rt->bSettle, rt->bBCheck, FALSE, &rt->ril_mt[mt]); - rt->ril_mt_tot_size += rt->ril_mt[mt].index[molt->atoms.nr]; + rt->ril_mt_tot_size += rt->ril_mt[mt].index[molt.atoms.nr]; } if (debug) { @@ -706,9 +699,9 @@ static gmx_reverse_top_t *make_reverse_top(const gmx_mtop_t *mtop, gmx_bool bFE, } *nint = 0; - for (mb = 0; mb < mtop->nmolblock; mb++) + for (const gmx_molblock_t &molblock : mtop->molblock) { - *nint += mtop->molblock[mb].nmol*nint_mt[mtop->molblock[mb].type]; + *nint += molblock.nmol*nint_mt[molblock.type]; } sfree(nint_mt); @@ -741,10 +734,10 @@ static gmx_reverse_top_t *make_reverse_top(const gmx_mtop_t *mtop, gmx_bool bFE, } /* Make a molblock index for fast searching */ - snew(rt->mbi, mtop->nmolblock); - rt->nmolblock = mtop->nmolblock; - i = 0; - for (mb = 0; mb < mtop->nmolblock; mb++) + snew(rt->mbi, mtop->molblock.size()); + rt->nmolblock = mtop->molblock.size(); + int i = 0; + for (size_t mb = 0; mb < mtop->molblock.size(); mb++) { rt->mbi[mb].a_start = i; i += mtop->molblock[mb].nmol*mtop->molblock[mb].natoms_mol; @@ -801,23 +794,18 @@ void dd_make_reverse_top(FILE *fplog, rt->bExclRequired = (ir->cutoff_scheme == ecutsGROUP && inputrecExclForces(ir)); - int nexcl, mb; - - nexcl = 0; + int nexcl = 0; dd->n_intercg_excl = 0; rt->n_excl_at_max = 0; - for (mb = 0; mb < mtop->nmolblock; mb++) + for (const gmx_molblock_t &molb : mtop->molblock) { - gmx_molblock_t *molb; - gmx_moltype_t *molt; - int n_excl_mol, n_excl_icg, n_excl_at_max; + int n_excl_mol, n_excl_icg, n_excl_at_max; - molb = &mtop->molblock[mb]; - molt = &mtop->moltype[molb->type]; - count_excls(&molt->cgs, &molt->excls, + const gmx_moltype_t &molt = mtop->moltype[molb.type]; + count_excls(&molt.cgs, &molt.excls, &n_excl_mol, &n_excl_icg, &n_excl_at_max); - nexcl += molb->nmol*n_excl_mol; - dd->n_intercg_excl += molb->nmol*n_excl_icg; + nexcl += molb.nmol*n_excl_mol; + dd->n_intercg_excl += molb.nmol*n_excl_icg; rt->n_excl_at_max = std::max(rt->n_excl_at_max, n_excl_at_max); } if (rt->bExclRequired) @@ -1536,7 +1524,7 @@ check_assign_interactions_atom(int i, int i_gl, */ static int make_bondeds_zone(gmx_domdec_t *dd, const gmx_domdec_zones_t *zones, - const gmx_molblock_t *molb, + const std::vector &molb, gmx_bool bRCheckMB, ivec rcheck, gmx_bool bRCheck2B, real rc2, int *la2lc, t_pbc *pbc_null, rvec *cg_cm, @@ -1632,7 +1620,7 @@ static void set_no_exclusions_zone(gmx_domdec_t *dd, gmx_domdec_zones_t *zones, * are supported. */ static int make_exclusions_zone_cg(gmx_domdec_t *dd, gmx_domdec_zones_t *zones, - const gmx_moltype_t *moltype, + const std::vector &moltype, gmx_bool bRCheck, real rc2, int *la2lc, t_pbc *pbc_null, rvec *cg_cm, const int *cginfo, @@ -1771,7 +1759,7 @@ static int make_exclusions_zone_cg(gmx_domdec_t *dd, gmx_domdec_zones_t *zones, /*! \brief Set the exclusion data for i-zone \p iz */ static void make_exclusions_zone(gmx_domdec_t *dd, gmx_domdec_zones_t *zones, - const gmx_moltype_t *moltype, + const std::vector &moltype, const int *cginfo, t_blocka *lexcls, int iz, @@ -2294,15 +2282,13 @@ static void check_link(t_blocka *link, int cg_gl, int cg_gl_j) } } -/*! \brief Allocate and return an array of the charge group index for all atoms */ -static int *make_at2cg(t_block *cgs) +/*! \brief Return a vector of the charge group index for all atoms */ +static std::vector make_at2cg(const t_block &cgs) { - int *at2cg, cg, a; - - snew(at2cg, cgs->index[cgs->nr]); - for (cg = 0; cg < cgs->nr; cg++) + std::vector at2cg(cgs.index[cgs.nr]); + for (int cg = 0; cg < cgs.nr; cg++) { - for (a = cgs->index[cg]; a < cgs->index[cg+1]; a++) + for (int a = cgs.index[cg]; a < cgs.index[cg + 1]; a++) { at2cg[a] = cg; } @@ -2315,12 +2301,6 @@ t_blocka *make_charge_group_links(const gmx_mtop_t *mtop, gmx_domdec_t *dd, cginfo_mb_t *cginfo_mb) { gmx_bool bExclRequired; - int mb, cg_offset, cg, cg_gl, a, aj, i, j, ftype, nral, nlink_mol, mol, ncgi; - gmx_molblock_t *molb; - gmx_moltype_t *molt; - t_block *cgs; - t_blocka *excls; - int *a2c; reverse_ilist_t ril, ril_intermol; t_blocka *link; cginfo_mb_t *cgi_mb; @@ -2354,46 +2334,47 @@ t_blocka *make_charge_group_links(const gmx_mtop_t *mtop, gmx_domdec_t *dd, link->a = nullptr; link->index[0] = 0; - cg_offset = 0; - ncgi = 0; - for (mb = 0; mb < mtop->nmolblock; mb++) + int cg_offset = 0; + int ncgi = 0; + for (size_t mb = 0; mb < mtop->molblock.size(); mb++) { - molb = &mtop->molblock[mb]; - if (molb->nmol == 0) + const gmx_molblock_t &molb = mtop->molblock[mb]; + if (molb.nmol == 0) { continue; } - molt = &mtop->moltype[molb->type]; - cgs = &molt->cgs; - excls = &molt->excls; - a2c = make_at2cg(cgs); + const gmx_moltype_t &molt = mtop->moltype[molb.type]; + const t_block &cgs = molt.cgs; + const t_blocka &excls = molt.excls; + std::vector a2c = make_at2cg(cgs); /* Make a reverse ilist in which the interactions are linked * to all atoms, not only the first atom as in gmx_reverse_top. * The constraints are discarded here. */ - make_reverse_ilist(molt->ilist, &molt->atoms, + make_reverse_ilist(molt.ilist, &molt.atoms, nullptr, FALSE, FALSE, FALSE, TRUE, &ril); cgi_mb = &cginfo_mb[mb]; - for (mol = 0; mol < (mtop->bIntermolecularInteractions ? molb->nmol : 1); mol++) + int mol; + for (mol = 0; mol < (mtop->bIntermolecularInteractions ? molb.nmol : 1); mol++) { - for (cg = 0; cg < cgs->nr; cg++) + for (int cg = 0; cg < cgs.nr; cg++) { - cg_gl = cg_offset + cg; + int cg_gl = cg_offset + cg; link->index[cg_gl+1] = link->index[cg_gl]; - for (a = cgs->index[cg]; a < cgs->index[cg+1]; a++) + for (int a = cgs.index[cg]; a < cgs.index[cg + 1]; a++) { - i = ril.index[a]; + int i = ril.index[a]; while (i < ril.index[a+1]) { - ftype = ril.il[i++]; - nral = NRAL(ftype); + int ftype = ril.il[i++]; + int nral = NRAL(ftype); /* Skip the ifunc index */ i++; - for (j = 0; j < nral; j++) + for (int j = 0; j < nral; j++) { - aj = ril.il[i+j]; + int aj = ril.il[i + j]; if (a2c[aj] != cg) { check_link(link, cg_gl, cg_offset+a2c[aj]); @@ -2404,9 +2385,9 @@ t_blocka *make_charge_group_links(const gmx_mtop_t *mtop, gmx_domdec_t *dd, if (bExclRequired) { /* Exclusions always go both ways */ - for (j = excls->index[a]; j < excls->index[a+1]; j++) + for (int j = excls.index[a]; j < excls.index[a + 1]; j++) { - aj = excls->a[j]; + int aj = excls.a[j]; if (a2c[aj] != cg) { check_link(link, cg_gl, cg_offset+a2c[aj]); @@ -2416,19 +2397,19 @@ t_blocka *make_charge_group_links(const gmx_mtop_t *mtop, gmx_domdec_t *dd, if (mtop->bIntermolecularInteractions) { - i = ril_intermol.index[a]; + int i = ril_intermol.index[a]; while (i < ril_intermol.index[a+1]) { - ftype = ril_intermol.il[i++]; - nral = NRAL(ftype); + int ftype = ril_intermol.il[i++]; + int nral = NRAL(ftype); /* Skip the ifunc index */ i++; - for (j = 0; j < nral; j++) + for (int j = 0; j < nral; j++) { /* Here we assume we have no charge groups; * this has been checked above. */ - aj = ril_intermol.il[i+j]; + int aj = ril_intermol.il[i + j]; check_link(link, cg_gl, aj); } i += nral_rt(ftype); @@ -2442,33 +2423,32 @@ t_blocka *make_charge_group_links(const gmx_mtop_t *mtop, gmx_domdec_t *dd, } } - cg_offset += cgs->nr; + cg_offset += cgs.nr; } - nlink_mol = link->index[cg_offset] - link->index[cg_offset-cgs->nr]; + int nlink_mol = link->index[cg_offset] - link->index[cg_offset - cgs.nr]; destroy_reverse_ilist(&ril); - sfree(a2c); if (debug) { - fprintf(debug, "molecule type '%s' %d cgs has %d cg links through bonded interac.\n", *molt->name, cgs->nr, nlink_mol); + fprintf(debug, "molecule type '%s' %d cgs has %d cg links through bonded interac.\n", *molt.name, cgs.nr, nlink_mol); } - if (molb->nmol > mol) + if (molb.nmol > mol) { /* Copy the data for the rest of the molecules in this block */ - link->nalloc_a += (molb->nmol - mol)*nlink_mol; + link->nalloc_a += (molb.nmol - mol)*nlink_mol; srenew(link->a, link->nalloc_a); - for (; mol < molb->nmol; mol++) + for (; mol < molb.nmol; mol++) { - for (cg = 0; cg < cgs->nr; cg++) + for (int cg = 0; cg < cgs.nr; cg++) { - cg_gl = cg_offset + cg; - link->index[cg_gl+1] = - link->index[cg_gl+1-cgs->nr] + nlink_mol; - for (j = link->index[cg_gl]; j < link->index[cg_gl+1]; j++) + int cg_gl = cg_offset + cg; + link->index[cg_gl + 1] = + link->index[cg_gl + 1 - cgs.nr] + nlink_mol; + for (int j = link->index[cg_gl]; j < link->index[cg_gl+1]; j++) { - link->a[j] = link->a[j-nlink_mol] + cgs->nr; + link->a[j] = link->a[j - nlink_mol] + cgs.nr; } if (link->index[cg_gl+1] - link->index[cg_gl] > 0 && cg_gl - cgi_mb->cg_start < cgi_mb->cg_mod) @@ -2477,7 +2457,7 @@ t_blocka *make_charge_group_links(const gmx_mtop_t *mtop, gmx_domdec_t *dd, ncgi++; } } - cg_offset += cgs->nr; + cg_offset += cgs.nr; } } } @@ -2516,7 +2496,8 @@ static void update_max_bonded_distance(real r2, int ftype, int a1, int a2, } /*! \brief Set the distance, function type and atom indices for the longest distance between charge-groups of molecule type \p molt for two-body and multi-body bonded interactions */ -static void bonded_cg_distance_mol(gmx_moltype_t *molt, int *at2cg, +static void bonded_cg_distance_mol(const gmx_moltype_t *molt, + const std::vector &at2cg, gmx_bool bBCheck, gmx_bool bExcl, rvec *cg_cm, bonded_distance_t *bd_2b, bonded_distance_t *bd_mb) @@ -2688,10 +2669,8 @@ void dd_bonded_cg_distance(FILE *fplog, real *r_2b, real *r_mb) { gmx_bool bExclRequired; - int mb, at_offset, *at2cg, mol; + int at_offset; t_graph graph; - gmx_molblock_t *molb; - gmx_moltype_t *molt; rvec *xs, *cg_cm; bonded_distance_t bd_2b = { 0, -1, -1, -1 }; bonded_distance_t bd_mb = { 0, -1, -1, -1 }; @@ -2701,34 +2680,33 @@ void dd_bonded_cg_distance(FILE *fplog, *r_2b = 0; *r_mb = 0; at_offset = 0; - for (mb = 0; mb < mtop->nmolblock; mb++) + for (const gmx_molblock_t &molb : mtop->molblock) { - molb = &mtop->molblock[mb]; - molt = &mtop->moltype[molb->type]; - if (molt->cgs.nr == 1 || molb->nmol == 0) + const gmx_moltype_t &molt = mtop->moltype[molb.type]; + if (molt.cgs.nr == 1 || molb.nmol == 0) { - at_offset += molb->nmol*molt->atoms.nr; + at_offset += molb.nmol*molt.atoms.nr; } else { if (ir->ePBC != epbcNONE) { - mk_graph_ilist(nullptr, molt->ilist, 0, molt->atoms.nr, FALSE, FALSE, + mk_graph_ilist(nullptr, molt.ilist, 0, molt.atoms.nr, FALSE, FALSE, &graph); } - at2cg = make_at2cg(&molt->cgs); - snew(xs, molt->atoms.nr); - snew(cg_cm, molt->cgs.nr); - for (mol = 0; mol < molb->nmol; mol++) + std::vector at2cg = make_at2cg(molt.cgs); + snew(xs, molt.atoms.nr); + snew(cg_cm, molt.cgs.nr); + for (int mol = 0; mol < molb.nmol; mol++) { - get_cgcm_mol(molt, &mtop->ffparams, ir->ePBC, &graph, box, + get_cgcm_mol(&molt, &mtop->ffparams, ir->ePBC, &graph, box, x+at_offset, xs, cg_cm); bonded_distance_t bd_mol_2b = { 0, -1, -1, -1 }; bonded_distance_t bd_mol_mb = { 0, -1, -1, -1 }; - bonded_cg_distance_mol(molt, at2cg, bBCheck, bExclRequired, cg_cm, + bonded_cg_distance_mol(&molt, at2cg, bBCheck, bExclRequired, cg_cm, &bd_mol_2b, &bd_mol_mb); /* Process the mol data adding the atom index offset */ @@ -2741,11 +2719,10 @@ void dd_bonded_cg_distance(FILE *fplog, at_offset + bd_mol_mb.a2, &bd_mb); - at_offset += molt->atoms.nr; + at_offset += molt.atoms.nr; } sfree(cg_cm); sfree(xs); - sfree(at2cg); if (ir->ePBC != epbcNONE) { done_graph(&graph); diff --git a/src/gromacs/fileio/confio.cpp b/src/gromacs/fileio/confio.cpp index b9015383b4..0eea8ba210 100644 --- a/src/gromacs/fileio/confio.cpp +++ b/src/gromacs/fileio/confio.cpp @@ -3,7 +3,7 @@ * * Copyright (c) 1991-2000, University of Groningen, The Netherlands. * Copyright (c) 2001-2004, The GROMACS development team. - * Copyright (c) 2013,2014,2015,2016,2017, by the GROMACS development team, led by + * Copyright (c) 2013,2014,2015,2016,2017,2018, 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. @@ -436,14 +436,11 @@ gmx_bool read_tps_conf(const char *infile, t_topology *top, int *ePBC, rvec **x, rvec **v, matrix box, gmx_bool requireMasses) { bool haveTopology; - gmx_mtop_t *mtop; + gmx_mtop_t mtop; - // Note: We should have an initializer instead of relying on snew - snew(mtop, 1); - readConfAndTopology(infile, &haveTopology, mtop, ePBC, x, v, box); + readConfAndTopology(infile, &haveTopology, &mtop, ePBC, x, v, box); - *top = gmx_mtop_t_to_t_topology(mtop, true); - sfree(mtop); + *top = gmx_mtop_t_to_t_topology(&mtop, true); tpx_make_chain_identifiers(&top->atoms, &top->mols); diff --git a/src/gromacs/fileio/tngio.cpp b/src/gromacs/fileio/tngio.cpp index 67c0225eae..12ef42ba29 100644 --- a/src/gromacs/fileio/tngio.cpp +++ b/src/gromacs/fileio/tngio.cpp @@ -308,17 +308,17 @@ void gmx_tng_add_mtop(gmx_tng_trajectory_t gmx_tng, atomCharges.reserve(mtop->natoms); atomMasses.reserve(mtop->natoms); - for (int molIndex = 0; molIndex < mtop->nmolblock; molIndex++) + for (const gmx_molblock_t &molBlock : mtop->molblock) { tng_molecule_t tngMol = nullptr; - const gmx_moltype_t *molType = &mtop->moltype[mtop->molblock[molIndex].type]; + const gmx_moltype_t *molType = &mtop->moltype[molBlock.type]; /* Add a molecule to the TNG trajectory with the same name as the * current molecule. */ addTngMoleculeFromTopology(gmx_tng, *(molType->name), &molType->atoms, - mtop->molblock[molIndex].nmol, + molBlock.nmol, &tngMol); /* Bonds have to be deduced from interactions (constraints etc). Different @@ -359,7 +359,7 @@ void gmx_tng_add_mtop(gmx_tng_trajectory_t gmx_tng, atomCharges.push_back(molType->atoms.atom[atomCounter].q); atomMasses.push_back(molType->atoms.atom[atomCounter].m); } - for (int molCounter = 1; molCounter < mtop->molblock[molIndex].nmol; molCounter++) + for (int molCounter = 1; molCounter < molBlock.nmol; molCounter++) { std::copy_n(atomCharges.end() - molType->atoms.nr, molType->atoms.nr, std::back_inserter(atomCharges)); std::copy_n(atomMasses.end() - molType->atoms.nr, molType->atoms.nr, std::back_inserter(atomMasses)); @@ -574,21 +574,18 @@ void gmx_tng_prepare_md_writing(gmx_tng_trajectory_t gmx_tng, static gmx_bool all_atoms_selected(const gmx_mtop_t *mtop, const int gtype) { - const gmx_moltype_t *molType; - const t_atoms *atoms; - /* Iterate through all atoms in the molecule system and * check if they belong to a selection group of the * requested type. */ - for (int molIndex = 0, i = 0; molIndex < mtop->nmoltype; molIndex++) + int i = 0; + for (const gmx_molblock_t &molBlock : mtop->molblock) { - molType = &mtop->moltype[mtop->molblock[molIndex].type]; + const gmx_moltype_t &molType = mtop->moltype[molBlock.type]; + const t_atoms &atoms = molType.atoms; - atoms = &molType->atoms; - - for (int j = 0; j < mtop->molblock[molIndex].nmol; j++) + for (int j = 0; j < molBlock.nmol; j++) { - for (int atomIndex = 0; atomIndex < atoms->nr; atomIndex++, i++) + for (int atomIndex = 0; atomIndex < atoms.nr; atomIndex++, i++) { if (ggrpnr(&mtop->groups, gtype, i) != 0) { @@ -610,7 +607,6 @@ static gmx_bool all_atoms_selected(const gmx_mtop_t *mtop, static void add_selection_groups(gmx_tng_trajectory_t gmx_tng, const gmx_mtop_t *mtop) { - const gmx_moltype_t *molType; const t_atoms *atoms; const t_atom *at; const t_resinfo *resInfo; @@ -653,13 +649,14 @@ static void add_selection_groups(gmx_tng_trajectory_t gmx_tng, tng_molecule_alloc(tng, &mol); tng_molecule_name_set(tng, mol, groupName); tng_molecule_chain_add(tng, mol, "", &chain); - for (int molIndex = 0, i = 0; molIndex < mtop->nmolblock; molIndex++) + int i = 0; + for (const gmx_molblock_t &molBlock : mtop->molblock) { - molType = &mtop->moltype[mtop->molblock[molIndex].type]; + const gmx_moltype_t &molType = mtop->moltype[molBlock.type]; - atoms = &molType->atoms; + atoms = &molType.atoms; - for (int j = 0; j < mtop->molblock[molIndex].nmol; j++) + for (int j = 0; j < molBlock.nmol; j++) { bool bAtomsAdded = FALSE; for (int atomIndex = 0; atomIndex < atoms->nr; atomIndex++, i++) @@ -704,7 +701,7 @@ static void add_selection_groups(gmx_tng_trajectory_t gmx_tng, { if (IS_CHEMBOND(k)) { - ilist = &molType->ilist[k]; + ilist = &molType.ilist[k]; if (ilist) { int l = 1; @@ -725,7 +722,7 @@ static void add_selection_groups(gmx_tng_trajectory_t gmx_tng, } } /* Settle is described using three atoms */ - ilist = &molType->ilist[F_SETTLE]; + ilist = &molType.ilist[F_SETTLE]; if (ilist) { int l = 1; diff --git a/src/gromacs/fileio/tpxio.cpp b/src/gromacs/fileio/tpxio.cpp index 7a841ec380..9615c1b5dc 100644 --- a/src/gromacs/fileio/tpxio.cpp +++ b/src/gromacs/fileio/tpxio.cpp @@ -2449,7 +2449,7 @@ static void set_disres_npair(gmx_mtop_t *mtop) { t_iparams *ip; gmx_mtop_ilistloop_t iloop; - t_ilist *ilist, *il; + const t_ilist *ilist, *il; int nmol, i, npair; t_iatom *a; @@ -2480,8 +2480,6 @@ static void set_disres_npair(gmx_mtop_t *mtop) static void do_mtop(t_fileio *fio, gmx_mtop_t *mtop, gmx_bool bRead, int file_version) { - int mt, mb; - if (bRead) { init_mtop(mtop); @@ -2492,25 +2490,26 @@ static void do_mtop(t_fileio *fio, gmx_mtop_t *mtop, gmx_bool bRead, do_ffparams(fio, &mtop->ffparams, bRead, file_version); - gmx_fio_do_int(fio, mtop->nmoltype); - + int nmoltype = mtop->moltype.size(); + gmx_fio_do_int(fio, nmoltype); if (bRead) { - snew(mtop->moltype, mtop->nmoltype); + mtop->moltype.resize(nmoltype); } - for (mt = 0; mt < mtop->nmoltype; mt++) + for (gmx_moltype_t &moltype : mtop->moltype) { - do_moltype(fio, &mtop->moltype[mt], bRead, &mtop->symtab, file_version); + do_moltype(fio, &moltype, bRead, &mtop->symtab, file_version); } - gmx_fio_do_int(fio, mtop->nmolblock); + int nmolblock = mtop->molblock.size(); + gmx_fio_do_int(fio, nmolblock); if (bRead) { - snew(mtop->molblock, mtop->nmolblock); + mtop->molblock.resize(nmolblock); } - for (mb = 0; mb < mtop->nmolblock; mb++) + for (gmx_molblock_t &molblock : mtop->molblock) { - do_molblock(fio, &mtop->molblock[mb], bRead); + do_molblock(fio, &molblock, bRead); } gmx_fio_do_int(fio, mtop->natoms); @@ -2706,7 +2705,6 @@ static int do_tpx(t_fileio *fio, gmx_bool bRead, gmx_mtop_t *mtop) { t_tpxheader tpx; - gmx_mtop_t dum_top; gmx_bool TopOnlyOK; int ePBC; gmx_bool bPeriodicMols; @@ -2810,8 +2808,8 @@ static int do_tpx(t_fileio *fio, gmx_bool bRead, } else { + gmx_mtop_t dum_top; do_mtop(fio, &dum_top, bRead, fileVersion); - done_mtop(&dum_top); } } do_test(fio, tpx.bX, x); diff --git a/src/gromacs/gmxana/gmx_clustsize.cpp b/src/gromacs/gmxana/gmx_clustsize.cpp index 5ea9d9fa79..3a81b9b736 100644 --- a/src/gromacs/gmxana/gmx_clustsize.cpp +++ b/src/gromacs/gmxana/gmx_clustsize.cpp @@ -114,7 +114,7 @@ static void clust_size(const char *ndx, const char *trx, const char *xpm, if (tpr) { - snew(mtop, 1); + mtop = new gmx_mtop_t; read_tpxheader(tpr, &tpxh, TRUE); if (tpxh.natoms != natoms) { @@ -444,8 +444,7 @@ static void clust_size(const char *ndx, const char *trx, const char *xpm, gmx_ffclose(fp); if (mtop) { - done_mtop(mtop); - sfree(mtop); + delete mtop; } sfree(t_x); sfree(t_y); diff --git a/src/gromacs/gmxana/gmx_disre.cpp b/src/gromacs/gmxana/gmx_disre.cpp index 7498dbf6e1..de2731492e 100644 --- a/src/gromacs/gmxana/gmx_disre.cpp +++ b/src/gromacs/gmxana/gmx_disre.cpp @@ -532,8 +532,7 @@ static void dump_disre_matrix(const char *fn, t_dr_result *dr, int ndr, { FILE *fp; int *resnr; - int n_res, a_offset, mb, mol, a; - t_atoms *atoms; + int n_res, a_offset, mol, a; int i, j, nra, nratoms, tp, ri, rj, index, nlabel, label; int ai, aj, *ptr; real **matrix, *t_res, hi, *w_dr, rav, rviol; @@ -547,17 +546,17 @@ static void dump_disre_matrix(const char *fn, t_dr_result *dr, int ndr, snew(resnr, mtop->natoms); n_res = 0; a_offset = 0; - for (mb = 0; mb < mtop->nmolblock; mb++) + for (const gmx_molblock_t &molb : mtop->molblock) { - atoms = &mtop->moltype[mtop->molblock[mb].type].atoms; - for (mol = 0; mol < mtop->molblock[mb].nmol; mol++) + const t_atoms &atoms = mtop->moltype[molb.type].atoms; + for (mol = 0; mol < molb.nmol; mol++) { - for (a = 0; a < atoms->nr; a++) + for (a = 0; a < atoms.nr; a++) { - resnr[a_offset+a] = n_res + atoms->atom[a].resind; + resnr[a_offset + a] = n_res + atoms.atom[a].resind; } - n_res += atoms->nres; - a_offset += atoms->nr; + n_res += atoms.nres; + a_offset += atoms.nr; } } diff --git a/src/gromacs/gmxana/gmx_energy.cpp b/src/gromacs/gmxana/gmx_energy.cpp index 3132b67322..77de9a3222 100644 --- a/src/gromacs/gmxana/gmx_energy.cpp +++ b/src/gromacs/gmxana/gmx_energy.cpp @@ -312,7 +312,6 @@ static void get_dhdl_parms(const char *topnm, t_inputrec *ir) /* all we need is the ir to be able to write the label */ read_tpx(topnm, ir, box, &natoms, nullptr, nullptr, &mtop); - done_mtop(&mtop); } static void einstein_visco(const char *fn, const char *fni, int nsets, diff --git a/src/gromacs/gmxana/gmx_nmeig.cpp b/src/gromacs/gmxana/gmx_nmeig.cpp index dbddf790ba..1f5fbad25d 100644 --- a/src/gromacs/gmxana/gmx_nmeig.cpp +++ b/src/gromacs/gmxana/gmx_nmeig.cpp @@ -114,10 +114,9 @@ static int get_nharm(const gmx_mtop_t *mtop) { int nh = 0; - for (int j = 0; (j < mtop->nmolblock); j++) + for (const gmx_molblock_t &molb : mtop->molblock) { - int mt = mtop->molblock[j].type; - nh += mtop->molblock[j].nmol * get_nharm_mt(&(mtop->moltype[mt])); + nh += molb.nmol * get_nharm_mt(&(mtop->moltype[molb.type])); } return nh; } diff --git a/src/gromacs/gmxana/gmx_pme_error.cpp b/src/gromacs/gmxana/gmx_pme_error.cpp index e854b1023b..55d35e0836 100644 --- a/src/gromacs/gmxana/gmx_pme_error.cpp +++ b/src/gromacs/gmxana/gmx_pme_error.cpp @@ -1,7 +1,7 @@ /* * This file is part of the GROMACS molecular simulation package. * - * Copyright (c) 2010,2011,2012,2013,2014,2015,2016,2017, by the GROMACS development team, led by + * Copyright (c) 2010,2011,2012,2013,2014,2015,2016,2017,2018, 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. @@ -123,28 +123,24 @@ static gmx_bool is_charge(real charge) static void calc_q2all(const gmx_mtop_t *mtop, /* molecular topology */ real *q2all, real *q2allnr) { - int imol, iatom; /* indices for loops */ real q2_all = 0; /* Sum of squared charges */ int nrq_mol; /* Number of charges in a single molecule */ int nrq_all; /* Total number of charges in the MD system */ real qi, q2_mol; - gmx_moltype_t *molecule; - gmx_molblock_t *molblock; #ifdef DEBUG fprintf(stderr, "\nCharge density:\n"); #endif - q2_all = 0.0; /* total q squared */ - nrq_all = 0; /* total number of charges in the system */ - for (imol = 0; imol < mtop->nmolblock; imol++) /* Loop over molecule types */ + q2_all = 0.0; /* total q squared */ + nrq_all = 0; /* total number of charges in the system */ + for (const gmx_molblock_t &molblock : mtop->molblock) /* Loop over molecule types */ { - q2_mol = 0.0; /* q squared value of this molecule */ - nrq_mol = 0; /* number of charges this molecule carries */ - molecule = &(mtop->moltype[imol]); - molblock = &(mtop->molblock[imol]); - for (iatom = 0; iatom < molblock->natoms_mol; iatom++) /* Loop over atoms in this molecule */ + q2_mol = 0.0; /* q squared value of this molecule */ + nrq_mol = 0; /* number of charges this molecule carries */ + const gmx_moltype_t &molecule = mtop->moltype[molblock.type]; + for (int i = 0; i < molecule.atoms.nr; i++) { - qi = molecule->atoms.atom[iatom].q; /* Charge of this atom */ + qi = molecule.atoms.atom[i].q; /* Is this charge worth to be considered? */ if (is_charge(qi)) { @@ -153,11 +149,11 @@ static void calc_q2all(const gmx_mtop_t *mtop, /* molecular topology */ } } /* Multiply with the number of molecules present of this type and add */ - q2_all += q2_mol*molblock->nmol; - nrq_all += nrq_mol*molblock->nmol; + q2_all += q2_mol*molblock.nmol; + nrq_all += nrq_mol*molblock.nmol; #ifdef DEBUG fprintf(stderr, "Molecule %2d (%5d atoms) q2_mol=%10.3e nr.mol.charges=%5d (%6dx) q2_all=%10.3e tot.charges=%d\n", - imol, molblock->natoms_mol, q2_mol, nrq_mol, molblock->nmol, q2_all, nrq_all); + imol, molecule.atoms.nr, q2_mol, nrq_mol, molblock.nmol, q2_all, nrq_all); #endif } diff --git a/src/gromacs/gmxana/gmx_trjconv.cpp b/src/gromacs/gmxana/gmx_trjconv.cpp index 462e3455b9..f3b1bfdacc 100644 --- a/src/gromacs/gmxana/gmx_trjconv.cpp +++ b/src/gromacs/gmxana/gmx_trjconv.cpp @@ -44,6 +44,7 @@ #include "gromacs/commandline/pargs.h" #include "gromacs/commandline/viewit.h" +#include "gromacs/compat/make_unique.h" #include "gromacs/fileio/confio.h" #include "gromacs/fileio/g96io.h" #include "gromacs/fileio/gmxfio.h" @@ -569,20 +570,21 @@ static void do_trunc(const char *fn, real t0) * molecule information will generally be present if the input TNG * file was written by a GROMACS tool, this seems like reasonable * behaviour. */ -static gmx_mtop_t *read_mtop_for_tng(const char *tps_file, - const char *input_file, - const char *output_file) +static std::unique_ptr +read_mtop_for_tng(const char *tps_file, + const char *input_file, + const char *output_file) { - gmx_mtop_t *mtop = nullptr; + std::unique_ptr mtop; if (fn2bTPX(tps_file) && efTNG != fn2ftp(input_file) && efTNG == fn2ftp(output_file)) { int temp_natoms = -1; - snew(mtop, 1); + mtop = gmx::compat::make_unique(); read_tpx(tps_file, nullptr, nullptr, &temp_natoms, - nullptr, nullptr, mtop); + nullptr, nullptr, mtop.get()); } return mtop; @@ -871,7 +873,6 @@ int gmx_trjconv(int argc, char *argv[]) int m, i, d, frame, outframe, natoms, nout, ncent, newstep = 0, model_nr; #define SKIP 10 t_topology top; - gmx_mtop_t *mtop = nullptr; gmx_conect gc = nullptr; int ePBC = -1; t_atoms *atoms = nullptr, useatoms; @@ -1081,7 +1082,7 @@ int gmx_trjconv(int argc, char *argv[]) { } - mtop = read_mtop_for_tng(top_file, in_file, out_file); + std::unique_ptr mtop = read_mtop_for_tng(top_file, in_file, out_file); /* Determine whether to read a topology */ bTPS = (ftp2bSet(efTPS, NFILE, fnm) || @@ -1365,7 +1366,7 @@ int gmx_trjconv(int argc, char *argv[]) trxin, nullptr, nout, - mtop, + mtop.get(), index, grpnm); break; @@ -1976,10 +1977,6 @@ int gmx_trjconv(int argc, char *argv[]) } } - if (mtop) - { - sfree(mtop); - } if (bTPS) { done_top(&top); diff --git a/src/gromacs/gmxlib/chargegroup.cpp b/src/gromacs/gmxlib/chargegroup.cpp index 424c847944..19f70e9eba 100644 --- a/src/gromacs/gmxlib/chargegroup.cpp +++ b/src/gromacs/gmxlib/chargegroup.cpp @@ -3,7 +3,7 @@ * * Copyright (c) 1991-2000, University of Groningen, The Netherlands. * Copyright (c) 2001-2004, The GROMACS development team. - * Copyright (c) 2013,2014,2015, by the GROMACS development team, led by + * Copyright (c) 2013,2014,2015,2018, 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,12 +52,8 @@ void calc_chargegroup_radii(const gmx_mtop_t *mtop, rvec *x, real *rcoul1, real *rcoul2) { real r2v1, r2v2, r2c1, r2c2, r2; - int ntype, i, j, mb, m, cg, a_mol, a0, a1, a; + int ntype, i, j, m, cg, a_mol, a0, a1, a; gmx_bool *bLJ; - gmx_molblock_t *molb; - gmx_moltype_t *molt; - t_block *cgs; - t_atom *atom; rvec cen; r2v1 = 0; @@ -81,13 +77,12 @@ void calc_chargegroup_radii(const gmx_mtop_t *mtop, rvec *x, } a_mol = 0; - for (mb = 0; mb < mtop->nmolblock; mb++) + for (const gmx_molblock_t &molb : mtop->molblock) { - molb = &mtop->molblock[mb]; - molt = &mtop->moltype[molb->type]; - cgs = &molt->cgs; - atom = molt->atoms.atom; - for (m = 0; m < molb->nmol; m++) + const gmx_moltype_t *molt = &mtop->moltype[molb.type]; + const t_block *cgs = &molt->cgs; + const t_atom *atom = molt->atoms.atom; + for (m = 0; m < molb.nmol; m++) { for (cg = 0; cg < cgs->nr; cg++) { @@ -133,7 +128,7 @@ void calc_chargegroup_radii(const gmx_mtop_t *mtop, rvec *x, } } } - a_mol += molb->natoms_mol; + a_mol += molb.natoms_mol; } } diff --git a/src/gromacs/gmxpreprocess/convparm.cpp b/src/gromacs/gmxpreprocess/convparm.cpp index 619e55b826..7dc03d1e6d 100644 --- a/src/gromacs/gmxpreprocess/convparm.cpp +++ b/src/gromacs/gmxpreprocess/convparm.cpp @@ -553,7 +553,7 @@ void convert_params(int atnr, t_params nbtypes[], int comb, double reppow, real fudgeQQ, gmx_mtop_t *mtop) { - int i, maxtypes, mt; + int i, maxtypes; unsigned long flags; gmx_ffparams_t *ffp; gmx_moltype_t *molt; @@ -573,7 +573,7 @@ void convert_params(int atnr, t_params nbtypes[], enter_function(&(nbtypes[F_BHAM]), (t_functype)F_BHAM, comb, reppow, ffp, nullptr, &maxtypes, TRUE, TRUE); - for (mt = 0; mt < mtop->nmoltype; mt++) + for (size_t mt = 0; mt < mtop->moltype.size(); mt++) { molt = &mtop->moltype[mt]; for (i = 0; (i < F_NRE); i++) diff --git a/src/gromacs/gmxpreprocess/gpp_atomtype.cpp b/src/gromacs/gmxpreprocess/gpp_atomtype.cpp index 50f7d7f124..f27615c981 100644 --- a/src/gromacs/gmxpreprocess/gpp_atomtype.cpp +++ b/src/gromacs/gmxpreprocess/gpp_atomtype.cpp @@ -337,7 +337,7 @@ void renum_atype(t_params plist[], gmx_mtop_t *mtop, int *wall_atomtype, gpp_atomtype_t ga, gmx_bool bVerbose) { - int i, j, k, l, molt, mi, mj, nat, nrfp, ftype, ntype; + int i, j, k, l, mi, mj, nat, nrfp, ftype, ntype; t_atoms *atoms; t_param *nbsnew; int *typelist; @@ -375,9 +375,9 @@ void renum_atype(t_params plist[], gmx_mtop_t *mtop, * can determine if two types should be merged. */ nat = 0; - for (molt = 0; molt < mtop->nmoltype; molt++) + for (gmx_moltype_t &moltype : mtop->moltype) { - atoms = &mtop->moltype[molt].atoms; + atoms = &moltype.atoms; for (i = 0; (i < atoms->nr); i++) { atoms->atom[i].type = diff --git a/src/gromacs/gmxpreprocess/grompp.cpp b/src/gromacs/gmxpreprocess/grompp.cpp index c38b6825a3..8280fcd7a0 100644 --- a/src/gromacs/gmxpreprocess/grompp.cpp +++ b/src/gromacs/gmxpreprocess/grompp.cpp @@ -50,6 +50,7 @@ #include "gromacs/awh/read-params.h" #include "gromacs/commandline/pargs.h" +#include "gromacs/compat/make_unique.h" #include "gromacs/ewald/ewald-utils.h" #include "gromacs/ewald/pme.h" #include "gromacs/fft/calcgrid.h" @@ -119,7 +120,7 @@ static int rm_interactions(int ifunc, int nrmols, t_molinfo mols[]) static int check_atom_names(const char *fn1, const char *fn2, gmx_mtop_t *mtop, const t_atoms *at) { - int mb, m, i, j, nmismatch; + int m, i, j, nmismatch; t_atoms *tat; #define MAXMISMATCH 20 @@ -130,10 +131,10 @@ static int check_atom_names(const char *fn1, const char *fn2, nmismatch = 0; i = 0; - for (mb = 0; mb < mtop->nmolblock; mb++) + for (const gmx_molblock_t &molb : mtop->molblock) { - tat = &mtop->moltype[mtop->molblock[mb].type].atoms; - for (m = 0; m < mtop->molblock[mb].nmol; m++) + tat = &mtop->moltype[molb.type].atoms; + for (m = 0; m < molb.nmol; m++) { for (j = 0; j < tat->nr; j++) { @@ -161,7 +162,7 @@ static int check_atom_names(const char *fn1, const char *fn2, static void check_eg_vs_cg(gmx_mtop_t *mtop) { - int astart, mb, m, cg, j, firstj; + int astart, m, cg, j, firstj; unsigned char firsteg, eg; gmx_moltype_t *molt; @@ -170,10 +171,10 @@ static void check_eg_vs_cg(gmx_mtop_t *mtop) */ astart = 0; - for (mb = 0; mb < mtop->nmolblock; mb++) + for (const gmx_molblock_t &molb : mtop->molblock) { - molt = &mtop->moltype[mtop->molblock[mb].type]; - for (m = 0; m < mtop->molblock[mb].nmol; m++) + molt = &mtop->moltype[molb.type]; + for (m = 0; m < molb.nmol; m++) { for (cg = 0; cg < molt->cgs.nr; cg++) { @@ -195,7 +196,7 @@ static void check_eg_vs_cg(gmx_mtop_t *mtop) } } -static void check_cg_sizes(const char *topfn, t_block *cgs, warninp_t wi) +static void check_cg_sizes(const char *topfn, const t_block *cgs, warninp_t wi) { int maxsize, cg; char warn_buf[STRLEN]; @@ -223,7 +224,7 @@ static void check_cg_sizes(const char *topfn, t_block *cgs, warninp_t wi) } } -static void check_bonds_timestep(gmx_mtop_t *mtop, double dt, warninp_t wi) +static void check_bonds_timestep(const gmx_mtop_t *mtop, double dt, warninp_t wi) { /* This check is not intended to ensure accurate integration, * rather it is to signal mistakes in the mdp settings. @@ -246,10 +247,6 @@ static void check_bonds_timestep(gmx_mtop_t *mtop, double dt, warninp_t wi) int min_steps_warn = 5; int min_steps_note = 10; t_iparams *ip; - int molt; - gmx_moltype_t *moltype, *w_moltype; - t_atom *atom; - t_ilist *ilist, *ilb, *ilc, *ils; int ftype; int i, a1, a2, w_a1, w_a2, j; real twopi2, limit2, fc, re, m1, m2, period2, w_period2; @@ -265,14 +262,13 @@ static void check_bonds_timestep(gmx_mtop_t *mtop, double dt, warninp_t wi) w_a1 = w_a2 = -1; w_period2 = -1.0; - w_moltype = nullptr; - for (molt = 0; molt < mtop->nmoltype; molt++) + const gmx_moltype_t *w_moltype = nullptr; + for (const gmx_moltype_t &moltype : mtop->moltype) { - moltype = &mtop->moltype[molt]; - atom = moltype->atoms.atom; - ilist = moltype->ilist; - ilc = &ilist[F_CONSTR]; - ils = &ilist[F_SETTLE]; + const t_atom *atom = moltype.atoms.atom; + const t_ilist *ilist = moltype.ilist; + const t_ilist *ilc = &ilist[F_CONSTR]; + const t_ilist *ils = &ilist[F_SETTLE]; for (ftype = 0; ftype < F_NRE; ftype++) { if (!(ftype == F_BONDS || ftype == F_G96BONDS || ftype == F_HARMONIC)) @@ -280,7 +276,7 @@ static void check_bonds_timestep(gmx_mtop_t *mtop, double dt, warninp_t wi) continue; } - ilb = &ilist[ftype]; + const t_ilist *ilb = &ilist[ftype]; for (i = 0; i < ilb->nr; i += 3) { fc = ip[ilb->iatoms[i]].harmonic.krA; @@ -329,7 +325,7 @@ static void check_bonds_timestep(gmx_mtop_t *mtop, double dt, warninp_t wi) if (!bFound && (w_moltype == nullptr || period2 < w_period2)) { - w_moltype = moltype; + w_moltype = &moltype; w_a1 = a1; w_a2 = a2; w_period2 = period2; @@ -416,12 +412,10 @@ static void check_shells_inputrec(gmx_mtop_t *mtop, * gmx_mtop_ftype_count */ static gmx_bool nint_ftype(gmx_mtop_t *mtop, t_molinfo *mi, int ftype) { - int nint, mb; - - nint = 0; - for (mb = 0; mb < mtop->nmolblock; mb++) + int nint = 0; + for (const gmx_molblock_t &molb : mtop->molblock) { - nint += mtop->molblock[mb].nmol*mi[mtop->molblock[mb].type].plist[ftype].nr; + nint += molb.nmol*mi[molb.type].plist[ftype].nr; } return nint; @@ -434,17 +428,17 @@ static gmx_bool nint_ftype(gmx_mtop_t *mtop, t_molinfo *mi, int ftype) static void renumber_moltypes(gmx_mtop_t *sys, int *nmolinfo, t_molinfo **molinfo) { - int *order, norder, i; - int mb, mi; + int *order, norder; t_molinfo *minew; snew(order, *nmolinfo); norder = 0; - for (mb = 0; mb < sys->nmolblock; mb++) + for (gmx_molblock_t &molblock : sys->molblock) { + int i; for (i = 0; i < norder; i++) { - if (order[i] == sys->molblock[mb].type) + if (order[i] == molblock.type) { break; } @@ -452,17 +446,18 @@ static void renumber_moltypes(gmx_mtop_t *sys, if (i == norder) { /* This type did not occur yet, add it */ - order[norder] = sys->molblock[mb].type; + order[norder] = molblock.type; /* Renumber the moltype in the topology */ norder++; } - sys->molblock[mb].type = i; + molblock.type = i; } /* We still need to reorder the molinfo structs */ snew(minew, norder); - for (mi = 0; mi < *nmolinfo; mi++) + for (int mi = 0; mi < *nmolinfo; mi++) { + int i; for (i = 0; i < norder; i++) { if (order[i] == mi) @@ -488,19 +483,15 @@ static void renumber_moltypes(gmx_mtop_t *sys, static void molinfo2mtop(int nmi, t_molinfo *mi, gmx_mtop_t *mtop) { - int m; - gmx_moltype_t *molt; - - mtop->nmoltype = nmi; - snew(mtop->moltype, nmi); - for (m = 0; m < nmi; m++) + mtop->moltype.resize(nmi); + for (int m = 0; m < nmi; m++) { - molt = &mtop->moltype[m]; - molt->name = mi[m].name; - molt->atoms = mi[m].atoms; + gmx_moltype_t &molt = mtop->moltype[m]; + molt.name = mi[m].name; + molt.atoms = mi[m].atoms; /* ilists are copied later */ - molt->cgs = mi[m].cgs; - molt->excls = mi[m].excls; + molt.cgs = mi[m].cgs; + molt.excls = mi[m].excls; } } @@ -515,12 +506,11 @@ new_status(const char *topfile, const char *topppfile, const char *confin, gmx_bool bMorse, warninp_t wi) { - t_molinfo *molinfo = nullptr; - int nmolblock; - gmx_molblock_t *molblock, *molbs; - int mb, i, nrmols, nmismatch; - char buf[STRLEN]; - char warn_buf[STRLEN]; + t_molinfo *molinfo = nullptr; + std::vector molblock; + int i, nrmols, nmismatch; + char buf[STRLEN]; + char warn_buf[STRLEN]; init_mtop(sys); @@ -529,35 +519,33 @@ new_status(const char *topfile, const char *topppfile, const char *confin, plist, comb, reppow, fudgeQQ, atype, &nrmols, &molinfo, intermolecular_interactions, ir, - &nmolblock, &molblock, + &molblock, wi); - sys->nmolblock = 0; - snew(sys->molblock, nmolblock); + sys->molblock.clear(); sys->natoms = 0; - for (mb = 0; mb < nmolblock; mb++) + for (const gmx_molblock_t &molb : molblock) { - if (sys->nmolblock > 0 && - molblock[mb].type == sys->molblock[sys->nmolblock-1].type) + if (!sys->molblock.empty() && + molb.type == sys->molblock.back().type) { /* Merge consecutive blocks with the same molecule type */ - sys->molblock[sys->nmolblock-1].nmol += molblock[mb].nmol; - sys->natoms += molblock[mb].nmol*sys->molblock[sys->nmolblock-1].natoms_mol; + sys->molblock.back().nmol += molb.nmol; + sys->natoms += molb.nmol*sys->molblock.back().natoms_mol; } - else if (molblock[mb].nmol > 0) + else if (molb.nmol > 0) { /* Add a new molblock to the topology */ - molbs = &sys->molblock[sys->nmolblock]; - *molbs = molblock[mb]; - molbs->natoms_mol = molinfo[molbs->type].atoms.nr; - molbs->nposres_xA = 0; - molbs->nposres_xB = 0; - sys->natoms += molbs->nmol*molbs->natoms_mol; - sys->nmolblock++; + sys->molblock.push_back(molb); + gmx_molblock_t &molbs = sys->molblock.back(); + molbs.natoms_mol = molinfo[molbs.type].atoms.nr; + molbs.nposres_xA = 0; + molbs.nposres_xB = 0; + sys->natoms += molbs.nmol*molbs.natoms_mol; } } - if (sys->nmolblock == 0) + if (sys->molblock.empty()) { gmx_fatal(FARGS, "No molecules were defined in the system"); } @@ -847,8 +835,7 @@ static void read_posres(gmx_mtop_t *mtop, t_molinfo *molinfo, gmx_bool bTopB, matrix box, invbox; int natoms, npbcdim = 0; char warn_buf[STRLEN]; - int a, i, ai, j, k, mb, nat_molb; - gmx_molblock_t *molb; + int a, i, ai, j, k, nat_molb; t_params *pr, *prfb; t_atom *atom; @@ -881,22 +868,21 @@ static void read_posres(gmx_mtop_t *mtop, t_molinfo *molinfo, gmx_bool bTopB, totmass = 0; a = 0; snew(hadAtom, natoms); - for (mb = 0; mb < mtop->nmolblock; mb++) + for (gmx_molblock_t &molb : mtop->molblock) { - molb = &mtop->molblock[mb]; - nat_molb = molb->nmol*mtop->moltype[molb->type].atoms.nr; - pr = &(molinfo[molb->type].plist[F_POSRES]); - prfb = &(molinfo[molb->type].plist[F_FBPOSRES]); + nat_molb = molb.nmol*mtop->moltype[molb.type].atoms.nr; + pr = &(molinfo[molb.type].plist[F_POSRES]); + prfb = &(molinfo[molb.type].plist[F_FBPOSRES]); if (pr->nr > 0 || prfb->nr > 0) { - atom = mtop->moltype[molb->type].atoms.atom; + atom = mtop->moltype[molb.type].atoms.atom; for (i = 0; (i < pr->nr); i++) { ai = pr->param[i].ai(); if (ai >= natoms) { gmx_fatal(FARGS, "Position restraint atom index (%d) in moltype '%s' is larger than number of atoms in %s (%d).\n", - ai+1, *molinfo[molb->type].name, fn, natoms); + ai+1, *molinfo[molb.type].name, fn, natoms); } hadAtom[ai] = TRUE; if (rc_scaling == erscCOM) @@ -916,7 +902,7 @@ static void read_posres(gmx_mtop_t *mtop, t_molinfo *molinfo, gmx_bool bTopB, if (ai >= natoms) { gmx_fatal(FARGS, "Position restraint atom index (%d) in moltype '%s' is larger than number of atoms in %s (%d).\n", - ai+1, *molinfo[molb->type].name, fn, natoms); + ai+1, *molinfo[molb.type].name, fn, natoms); } if (rc_scaling == erscCOM && hadAtom[ai] == FALSE) { @@ -930,20 +916,20 @@ static void read_posres(gmx_mtop_t *mtop, t_molinfo *molinfo, gmx_bool bTopB, } if (!bTopB) { - molb->nposres_xA = nat_molb; - snew(molb->posres_xA, molb->nposres_xA); + molb.nposres_xA = nat_molb; + snew(molb.posres_xA, molb.nposres_xA); for (i = 0; i < nat_molb; i++) { - copy_rvec(x[a+i], molb->posres_xA[i]); + copy_rvec(x[a+i], molb.posres_xA[i]); } } else { - molb->nposres_xB = nat_molb; - snew(molb->posres_xB, molb->nposres_xB); + molb.nposres_xB = nat_molb; + snew(molb.posres_xB, molb.nposres_xB); for (i = 0; i < nat_molb; i++) { - copy_rvec(x[a+i], molb->posres_xB[i]); + copy_rvec(x[a+i], molb.posres_xB[i]); } } } @@ -966,13 +952,12 @@ static void read_posres(gmx_mtop_t *mtop, t_molinfo *molinfo, gmx_bool bTopB, { GMX_ASSERT(npbcdim <= DIM, "Only DIM dimensions can have PBC"); - for (mb = 0; mb < mtop->nmolblock; mb++) + for (gmx_molblock_t &molb : mtop->molblock) { - molb = &mtop->molblock[mb]; - nat_molb = molb->nmol*mtop->moltype[molb->type].atoms.nr; - if (molb->nposres_xA > 0 || molb->nposres_xB > 0) + nat_molb = molb.nmol*mtop->moltype[molb.type].atoms.nr; + if (molb.nposres_xA > 0 || molb.nposres_xB > 0) { - xp = (!bTopB ? molb->posres_xA : molb->posres_xB); + xp = (!bTopB ? molb.posres_xA : molb.posres_xB); for (i = 0; i < nat_molb; i++) { for (j = 0; j < npbcdim; j++) @@ -1222,19 +1207,17 @@ static void init_cmap_grid(gmx_cmap_t *cmap_grid, int ngrid, int grid_spacing) } -static int count_constraints(gmx_mtop_t *mtop, t_molinfo *mi, warninp_t wi) +static int count_constraints(const gmx_mtop_t *mtop, t_molinfo *mi, warninp_t wi) { - int count, count_mol, i, mb; - gmx_molblock_t *molb; + int count, count_mol, i; t_params *plist; char buf[STRLEN]; count = 0; - for (mb = 0; mb < mtop->nmolblock; mb++) + for (const gmx_molblock_t &molb : mtop->molblock) { count_mol = 0; - molb = &mtop->molblock[mb]; - plist = mi[molb->type].plist; + plist = mi[molb.type].plist; for (i = 0; i < F_NRE; i++) { @@ -1248,16 +1231,16 @@ static int count_constraints(gmx_mtop_t *mtop, t_molinfo *mi, warninp_t wi) } } - if (count_mol > nrdf_internal(&mi[molb->type].atoms)) + if (count_mol > nrdf_internal(&mi[molb.type].atoms)) { sprintf(buf, "Molecule type '%s' has %d constraints.\n" "For stability and efficiency there should not be more constraints than internal number of degrees of freedom: %d.\n", - *mi[molb->type].name, count_mol, - nrdf_internal(&mi[molb->type].atoms)); + *mi[molb.type].name, count_mol, + nrdf_internal(&mi[molb.type].atoms)); warning(wi, buf); } - count += molb->nmol*count_mol; + count += molb.nmol*count_mol; } return count; @@ -1385,9 +1368,9 @@ static void checkForUnboundAtoms(const gmx_mtop_t *mtop, gmx_bool bVerbose, warninp_t wi) { - for (int mt = 0; mt < mtop->nmoltype; mt++) + for (const gmx_moltype_t &molt : mtop->moltype) { - checkForUnboundAtoms(&mtop->moltype[mt], bVerbose, wi); + checkForUnboundAtoms(&molt, bVerbose, wi); } } @@ -1525,9 +1508,9 @@ static void checkDecoupledModeAccuracy(const gmx_mtop_t *mtop, } bool haveDecoupledMode = false; - for (int mt = 0; mt < mtop->nmoltype; mt++) + for (const gmx_moltype_t &molt : mtop->moltype) { - if (haveDecoupledModeInMol(&mtop->moltype[mt], mtop->ffparams.iparams, + if (haveDecoupledModeInMol(&molt, mtop->ffparams.iparams, massFactorThreshold)) { haveDecoupledMode = true; @@ -1708,11 +1691,10 @@ int gmx_grompp(int argc, char *argv[]) "this option." }; t_gromppopts *opts; - gmx_mtop_t *sys; int nmi; t_molinfo *mi, *intermolecular_interactions; gpp_atomtype_t atype; - int nvsite, comb, mt; + int nvsite, comb; t_params *plist; real fudgeQQ; double reppow; @@ -1820,11 +1802,11 @@ int gmx_grompp(int argc, char *argv[]) snew(plist, F_NRE); init_plist(plist); - snew(sys, 1); + gmx_mtop_t sys; atype = init_atomtype(); if (debug) { - pr_symtab(debug, 0, "Just opened", &sys->symtab); + pr_symtab(debug, 0, "Just opened", &sys.symtab); } const char *fn = ftp2fn(efTOP, NFILE, fnm); @@ -1836,30 +1818,30 @@ int gmx_grompp(int argc, char *argv[]) t_state state; new_status(fn, opt2fn_null("-pp", NFILE, fnm), opt2fn("-c", NFILE, fnm), opts, ir, bZero, bGenVel, bVerbose, &state, - atype, sys, &nmi, &mi, &intermolecular_interactions, + atype, &sys, &nmi, &mi, &intermolecular_interactions, plist, &comb, &reppow, &fudgeQQ, opts->bMorse, wi); if (debug) { - pr_symtab(debug, 0, "After new_status", &sys->symtab); + pr_symtab(debug, 0, "After new_status", &sys.symtab); } nvsite = 0; /* set parameters for virtual site construction (not for vsiten) */ - for (mt = 0; mt < sys->nmoltype; mt++) + for (size_t mt = 0; mt < sys.moltype.size(); mt++) { nvsite += - set_vsites(bVerbose, &sys->moltype[mt].atoms, atype, mi[mt].plist); + set_vsites(bVerbose, &sys.moltype[mt].atoms, atype, mi[mt].plist); } /* now throw away all obsolete bonds, angles and dihedrals: */ /* note: constraints are ALWAYS removed */ if (nvsite) { - for (mt = 0; mt < sys->nmoltype; mt++) + for (size_t mt = 0; mt < sys.moltype.size(); mt++) { - clean_vsite_bondeds(mi[mt].plist, sys->moltype[mt].atoms.nr, bRmVSBds); + clean_vsite_bondeds(mi[mt].plist, sys.moltype[mt].atoms.nr, bRmVSBds); } } @@ -1869,10 +1851,10 @@ int gmx_grompp(int argc, char *argv[]) ecutscheme_names[ir->cutoff_scheme]); /* Remove all charge groups */ - gmx_mtop_remove_chargegroups(sys); + gmx_mtop_remove_chargegroups(&sys); } - if (count_constraints(sys, mi, wi) && (ir->eConstrAlg == econtSHAKE)) + if (count_constraints(&sys, mi, wi) && (ir->eConstrAlg == econtSHAKE)) { if (ir->eI == eiCG || ir->eI == eiLBFGS) { @@ -1916,8 +1898,8 @@ int gmx_grompp(int argc, char *argv[]) */ check_warning_error(wi, FARGS); - if (nint_ftype(sys, mi, F_POSRES) > 0 || - nint_ftype(sys, mi, F_FBPOSRES) > 0) + if (nint_ftype(&sys, mi, F_POSRES) > 0 || + nint_ftype(&sys, mi, F_FBPOSRES) > 0) { if (ir->epc == epcPARRINELLORAHMAN || ir->epc == epcMTTK) { @@ -1967,7 +1949,7 @@ int gmx_grompp(int argc, char *argv[]) fprintf(stderr, " and %s\n", fnB); } } - gen_posres(sys, mi, fn, fnB, + gen_posres(&sys, mi, fn, fnB, ir->refcoord_scaling, ir->ePBC, ir->posres_com, ir->posres_comB, wi); @@ -1976,14 +1958,14 @@ int gmx_grompp(int argc, char *argv[]) /* If we are using CMAP, setup the pre-interpolation grid */ if (plist[F_CMAP].ncmap > 0) { - init_cmap_grid(&sys->ffparams.cmap_grid, plist[F_CMAP].nc, plist[F_CMAP].grid_spacing); - setup_cmap(plist[F_CMAP].grid_spacing, plist[F_CMAP].nc, plist[F_CMAP].cmap, &sys->ffparams.cmap_grid); + init_cmap_grid(&sys.ffparams.cmap_grid, plist[F_CMAP].nc, plist[F_CMAP].grid_spacing); + setup_cmap(plist[F_CMAP].grid_spacing, plist[F_CMAP].nc, plist[F_CMAP].cmap, &sys.ffparams.cmap_grid); } set_wall_atomtype(atype, opts, ir, wi); if (bRenum) { - renum_atype(plist, sys, ir->wall_atomtype, atype, bVerbose); + renum_atype(plist, &sys, ir->wall_atomtype, atype, bVerbose); get_atomtype_ntypes(atype); } @@ -1993,11 +1975,11 @@ int gmx_grompp(int argc, char *argv[]) } /* PELA: Copy the atomtype data to the topology atomtype list */ - copy_atomtype_atomtypes(atype, &(sys->atomtypes)); + copy_atomtype_atomtypes(atype, &(sys.atomtypes)); if (debug) { - pr_symtab(debug, 0, "After renum_atype", &sys->symtab); + pr_symtab(debug, 0, "After renum_atype", &sys.symtab); } if (bVerbose) @@ -2007,47 +1989,47 @@ int gmx_grompp(int argc, char *argv[]) ntype = get_atomtype_ntypes(atype); convert_params(ntype, plist, mi, intermolecular_interactions, - comb, reppow, fudgeQQ, sys); + comb, reppow, fudgeQQ, &sys); if (debug) { - pr_symtab(debug, 0, "After convert_params", &sys->symtab); + pr_symtab(debug, 0, "After convert_params", &sys.symtab); } /* set ptype to VSite for virtual sites */ - for (mt = 0; mt < sys->nmoltype; mt++) + for (gmx_moltype_t &moltype : sys.moltype) { - set_vsites_ptype(FALSE, &sys->moltype[mt]); + set_vsites_ptype(FALSE, &moltype); } if (debug) { - pr_symtab(debug, 0, "After virtual sites", &sys->symtab); + pr_symtab(debug, 0, "After virtual sites", &sys.symtab); } /* Check velocity for virtual sites and shells */ if (bGenVel) { - check_vel(sys, as_rvec_array(state.v.data())); + check_vel(&sys, as_rvec_array(state.v.data())); } /* check for shells and inpurecs */ - check_shells_inputrec(sys, ir, wi); + check_shells_inputrec(&sys, ir, wi); /* check masses */ - check_mol(sys, wi); + check_mol(&sys, wi); - checkForUnboundAtoms(sys, bVerbose, wi); + checkForUnboundAtoms(&sys, bVerbose, wi); - for (i = 0; i < sys->nmoltype; i++) + for (const gmx_moltype_t &moltype : sys.moltype) { - check_cg_sizes(ftp2fn(efTOP, NFILE, fnm), &sys->moltype[i].cgs, wi); + check_cg_sizes(ftp2fn(efTOP, NFILE, fnm), &moltype.cgs, wi); } if (EI_DYNAMICS(ir->eI) && ir->eI != eiBD) { - check_bonds_timestep(sys, ir->delta_t, wi); + check_bonds_timestep(&sys, ir->delta_t, wi); } - checkDecoupledModeAccuracy(sys, ir, wi); + checkDecoupledModeAccuracy(&sys, ir, wi); if (EI_ENERGY_MINIMIZATION(ir->eI) && 0 == ir->nsteps) { @@ -2061,7 +2043,7 @@ int gmx_grompp(int argc, char *argv[]) fprintf(stderr, "initialising group options...\n"); } do_index(mdparin, ftp2fn_null(efNDX, NFILE, fnm), - sys, bVerbose, ir, + &sys, bVerbose, ir, wi); if (ir->cutoff_scheme == ecutsVERLET && ir->verletbuf_tol > 0) @@ -2078,7 +2060,7 @@ int gmx_grompp(int argc, char *argv[]) } else { - buffer_temp = calc_temp(sys, ir, as_rvec_array(state.v.data())); + buffer_temp = calc_temp(&sys, ir, as_rvec_array(state.v.data())); } if (buffer_temp > 0) { @@ -2133,7 +2115,7 @@ int gmx_grompp(int argc, char *argv[]) } } - set_verlet_buffer(sys, ir, buffer_temp, state.box, wi); + set_verlet_buffer(&sys, ir, buffer_temp, state.box, wi); } } } @@ -2145,18 +2127,18 @@ int gmx_grompp(int argc, char *argv[]) { fprintf(stderr, "Checking consistency between energy and charge groups...\n"); } - check_eg_vs_cg(sys); + check_eg_vs_cg(&sys); if (debug) { - pr_symtab(debug, 0, "After index", &sys->symtab); + pr_symtab(debug, 0, "After index", &sys.symtab); } - triple_check(mdparin, ir, sys, wi); - close_symtab(&sys->symtab); + triple_check(mdparin, ir, &sys, wi); + close_symtab(&sys.symtab); if (debug) { - pr_symtab(debug, 0, "After close", &sys->symtab); + pr_symtab(debug, 0, "After close", &sys.symtab); } /* make exclusions between QM atoms */ @@ -2168,7 +2150,7 @@ int gmx_grompp(int argc, char *argv[]) } else { - generate_qmexcl(sys, ir, wi); + generate_qmexcl(&sys, ir, wi); } } @@ -2179,7 +2161,7 @@ int gmx_grompp(int argc, char *argv[]) fprintf(stderr, "getting data from old trajectory ...\n"); } cont_status(ftp2fn(efTRN, NFILE, fnm), ftp2fn_null(efEDR, NFILE, fnm), - bNeedVel, bGenVel, fr_time, ir, &state, sys, oenv); + bNeedVel, bGenVel, fr_time, ir, &state, &sys, oenv); } if (ir->ePBC == epbcXY && ir->nwall != 2) @@ -2190,7 +2172,7 @@ int gmx_grompp(int argc, char *argv[]) if (ir->cutoff_scheme != ecutsVERLET && ir->rlist > 0) { set_warning_line(wi, mdparin, -1); - check_chargegroup_radii(sys, ir, as_rvec_array(state.x.data()), wi); + check_chargegroup_radii(&sys, ir, as_rvec_array(state.x.data()), wi); } if (EEL_FULL(ir->coulombtype) || EVDW_PME(ir->vdwtype)) @@ -2252,7 +2234,7 @@ int gmx_grompp(int argc, char *argv[]) if (ir->bPull) { - pull = set_pull_init(ir, sys, as_rvec_array(state.x.data()), state.box, state.lambda[efptMASS]); + pull = set_pull_init(ir, &sys, as_rvec_array(state.x.data()), state.box, state.lambda[efptMASS]); } /* Modules that supply external potential for pull coordinates @@ -2282,7 +2264,7 @@ int gmx_grompp(int argc, char *argv[]) if (EEL_PME(ir->coulombtype)) { - float ratio = pme_load_estimate(sys, ir, state.box); + float ratio = pme_load_estimate(&sys, ir, state.box); fprintf(stderr, "Estimate for the relative computational load of the PME mesh part: %.2f\n", ratio); /* With free energy we might need to do PME both for the A and B state * charges. This will double the cost, but the optimal performance will @@ -2305,7 +2287,7 @@ int gmx_grompp(int argc, char *argv[]) { char warn_buf[STRLEN]; - double cio = compute_io(ir, sys->natoms, &sys->groups, F_NRE, 1); + double cio = compute_io(ir, sys.natoms, &sys.groups, F_NRE, 1); sprintf(warn_buf, "This run will generate roughly %.0f Mb of data", cio); if (cio > 2000) { @@ -2324,10 +2306,10 @@ int gmx_grompp(int argc, char *argv[]) } done_warning(wi, FARGS); - write_tpx_state(ftp2fn(efTPR, NFILE, fnm), ir, &state, sys); + write_tpx_state(ftp2fn(efTPR, NFILE, fnm), ir, &state, &sys); /* Output IMD group, if bIMD is TRUE */ - write_IMDgroup_to_file(ir->bIMD, ir, &state, sys, NFILE, fnm); + write_IMDgroup_to_file(ir->bIMD, ir, &state, &sys, NFILE, fnm); sfree(opts->define); sfree(opts->include); @@ -2336,8 +2318,6 @@ int gmx_grompp(int argc, char *argv[]) sfree(plist); sfree(mi); done_atomtype(atype); - done_mtop(sys); - sfree(sys); done_inputrec_strings(); return 0; diff --git a/src/gromacs/gmxpreprocess/insert-molecules.cpp b/src/gromacs/gmxpreprocess/insert-molecules.cpp index 90de8972be..86b6920762 100644 --- a/src/gromacs/gmxpreprocess/insert-molecules.cpp +++ b/src/gromacs/gmxpreprocess/insert-molecules.cpp @@ -3,7 +3,7 @@ * * Copyright (c) 1991-2000, University of Groningen, The Netherlands. * Copyright (c) 2001-2004, The GROMACS development team. - * Copyright (c) 2013,2014,2015,2016,2017, by the GROMACS development team, led by + * Copyright (c) 2013,2014,2015,2016,2017,2018, 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. @@ -329,23 +329,15 @@ class InsertMolecules : public ICommandLineOptionsModule, public ITopologyProvid InsertMolecules() : bBox_(false), nmolIns_(0), nmolTry_(10), seed_(0), defaultDistance_(0.105), scaleFactor_(0.57), enumRot_(en_rotXYZ), - top_(nullptr), ePBC_(-1) + top_(), ePBC_(-1) { clear_rvec(newBox_); clear_rvec(deltaR_); clear_mat(box_); } - virtual ~InsertMolecules() - { - if (top_ != nullptr) - { - done_mtop(top_); - sfree(top_); - } - } // From ITopologyProvider - virtual gmx_mtop_t *getTopology(bool /*required*/) { return top_; } + virtual gmx_mtop_t *getTopology(bool /*required*/) { return &top_; } virtual int getAtomCount() { return 0; } // From ICommandLineOptionsModule @@ -377,7 +369,7 @@ class InsertMolecules : public ICommandLineOptionsModule, public ITopologyProvid RotationType enumRot_; Selection replaceSel_; - gmx_mtop_t *top_; + gmx_mtop_t top_; std::vector x_; matrix box_; int ePBC_; @@ -508,12 +500,11 @@ void InsertMolecules::optionsFinished() "together with an existing configuration (-f).")); } - snew(top_, 1); if (!inputConfFile_.empty()) { - readConformation(inputConfFile_.c_str(), top_, &x_, nullptr, + readConformation(inputConfFile_.c_str(), &top_, &x_, nullptr, &ePBC_, box_, "solute"); - if (top_->natoms == 0) + if (top_.natoms == 0) { fprintf(stderr, "Note: no atoms in %s\n", inputConfFile_.c_str()); } @@ -568,9 +559,9 @@ int InsertMolecules::run() gmx_fatal(FARGS, "No molecule in %s, please check your input", insertConfFile_.c_str()); } - if (top_->name == nullptr) + if (top_.name == nullptr) { - top_->name = top_insrt->name; + top_.name = top_insrt->name; } if (positionFile_.empty()) { @@ -579,18 +570,18 @@ int InsertMolecules::run() } // TODO: Adapt to use mtop throughout. - t_atoms atoms = gmx_mtop_global_atoms(top_); + t_atoms atoms = gmx_mtop_global_atoms(&top_); /* add nmol_ins molecules of atoms_ins in random orientation at random place */ insert_mols(nmolIns_, nmolTry_, seed_, defaultDistance_, scaleFactor_, - &atoms, &top_->symtab, &x_, removableAtoms, top_insrt->atoms, x_insrt, + &atoms, &top_.symtab, &x_, removableAtoms, top_insrt->atoms, x_insrt, ePBC_, box_, positionFile_, deltaR_, enumRot_); /* write new configuration to file confout */ fprintf(stderr, "Writing generated configuration to %s\n", outputConfFile_.c_str()); - write_sto_conf(outputConfFile_.c_str(), *top_->name, &atoms, + write_sto_conf(outputConfFile_.c_str(), *top_.name, &atoms, as_rvec_array(x_.data()), nullptr, ePBC_, box_); /* print size of generated configuration */ diff --git a/src/gromacs/gmxpreprocess/readir.cpp b/src/gromacs/gmxpreprocess/readir.cpp index 2ac98cf3ef..0d1cd8ae46 100644 --- a/src/gromacs/gmxpreprocess/readir.cpp +++ b/src/gromacs/gmxpreprocess/readir.cpp @@ -2766,10 +2766,10 @@ static gmx_bool do_numbering(int natoms, gmx_groups_t *groups, int ng, char *ptr return (bRest && grptp == egrptpPART); } -static void calc_nrdf(gmx_mtop_t *mtop, t_inputrec *ir, char **gnames) +static void calc_nrdf(const gmx_mtop_t *mtop, t_inputrec *ir, char **gnames) { t_grpopts *opts; - gmx_groups_t *groups; + const gmx_groups_t *groups; pull_params_t *pull; int natoms, ai, aj, i, j, d, g, imin, jmin; t_iatom *ia; @@ -2777,9 +2777,7 @@ 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; - int mb, mol, ftype, as; - gmx_molblock_t *molb; - gmx_moltype_t *molt; + int mol, ftype, as; /* Calculate nrdf. * First calc 3xnr-atoms for each group @@ -2840,17 +2838,16 @@ static void calc_nrdf(gmx_mtop_t *mtop, t_inputrec *ir, char **gnames) } as = 0; - for (mb = 0; mb < mtop->nmolblock; mb++) + for (const gmx_molblock_t &molb : mtop->molblock) { - molb = &mtop->molblock[mb]; - molt = &mtop->moltype[molb->type]; - atom = molt->atoms.atom; - for (mol = 0; mol < molb->nmol; mol++) + const gmx_moltype_t &molt = mtop->moltype[molb.type]; + atom = molt.atoms.atom; + for (mol = 0; mol < molb.nmol; mol++) { for (ftype = F_CONSTR; ftype <= F_CONSTRNC; ftype++) { - ia = molt->ilist[ftype].iatoms; - for (i = 0; i < molt->ilist[ftype].nr; ) + ia = molt.ilist[ftype].iatoms; + for (i = 0; i < molt.ilist[ftype].nr; ) { /* Subtract degrees of freedom for the constraints, * if the particles still have degrees of freedom left. @@ -2895,8 +2892,8 @@ static void calc_nrdf(gmx_mtop_t *mtop, t_inputrec *ir, char **gnames) i += interaction_function[ftype].nratoms+1; } } - ia = molt->ilist[F_SETTLE].iatoms; - for (i = 0; i < molt->ilist[F_SETTLE].nr; ) + ia = molt.ilist[F_SETTLE].iatoms; + for (i = 0; i < molt.ilist[F_SETTLE].nr; ) { /* Subtract 1 dof from every atom in the SETTLE */ for (j = 0; j < 3; j++) @@ -2910,7 +2907,7 @@ static void calc_nrdf(gmx_mtop_t *mtop, t_inputrec *ir, char **gnames) ia += 4; i += 4; } - as += molt->atoms.nr; + as += molt.atoms.nr; } } @@ -3817,7 +3814,7 @@ static gmx_bool absolute_reference(t_inputrec *ir, gmx_mtop_t *sys, { int d, g, i; gmx_mtop_ilistloop_t iloop; - t_ilist *ilist; + const t_ilist *ilist; int nmol; t_iparams *pr; diff --git a/src/gromacs/gmxpreprocess/topio.cpp b/src/gromacs/gmxpreprocess/topio.cpp index 9c466cc29a..32106323da 100644 --- a/src/gromacs/gmxpreprocess/topio.cpp +++ b/src/gromacs/gmxpreprocess/topio.cpp @@ -181,24 +181,22 @@ static void gen_pairs(t_params *nbs, t_params *pairs, real fudge, int comb) } } -double check_mol(gmx_mtop_t *mtop, warninp_t wi) +double check_mol(const gmx_mtop_t *mtop, warninp_t wi) { char buf[256]; - int i, mb, nmol, ri, pt; + int i, ri, pt; double q; real m, mB; - t_atoms *atoms; /* Check mass and charge */ q = 0.0; - for (mb = 0; mb < mtop->nmolblock; mb++) + for (const gmx_molblock_t &molb : mtop->molblock) { - atoms = &mtop->moltype[mtop->molblock[mb].type].atoms; - nmol = mtop->molblock[mb].nmol; + const t_atoms *atoms = &mtop->moltype[molb.type].atoms; for (i = 0; (i < atoms->nr); i++) { - q += nmol*atoms->atom[i].q; + q += molb.nmol*atoms->atom[i].q; m = atoms->atom[i].m; mB = atoms->atom[i].mB; pt = atoms->atom[i].ptype; @@ -395,28 +393,24 @@ static char ** cpp_opts(const char *define, const char *include, } -static void make_atoms_sys(int nmolb, const gmx_molblock_t *molb, - const t_molinfo *molinfo, - t_atoms *atoms) +static void make_atoms_sys(const std::vector &molblock, + const t_molinfo *molinfo, + t_atoms *atoms) { - int mb, m, a; - const t_atoms *mol_atoms; - atoms->nr = 0; atoms->atom = nullptr; - for (mb = 0; mb < nmolb; mb++) + for (const gmx_molblock_t &molb : molblock) { - assert(molb); - mol_atoms = &molinfo[molb[mb].type].atoms; + const t_atoms &mol_atoms = molinfo[molb.type].atoms; - srenew(atoms->atom, atoms->nr + molb[mb].nmol*mol_atoms->nr); + srenew(atoms->atom, atoms->nr + molb.nmol*mol_atoms.nr); - for (m = 0; m < molb[mb].nmol; m++) + for (int m = 0; m < molb.nmol; m++) { - for (a = 0; a < mol_atoms->nr; a++) + for (int a = 0; a < mol_atoms.nr; a++) { - atoms->atom[atoms->nr++] = mol_atoms->atom[a]; + atoms->atom[atoms->nr++] = mol_atoms.atom[a]; } } } @@ -435,8 +429,7 @@ static char **read_topol(const char *infile, const char *outfile, double *reppow, t_gromppopts *opts, real *fudgeQQ, - int *nmolblock, - gmx_molblock_t **molblock, + std::vector *molblock, gmx_bool bFEP, gmx_bool bZero, gmx_bool usingFullRangeElectrostatics, @@ -448,9 +441,8 @@ static char **read_topol(const char *infile, const char *outfile, char line[STRLEN], errbuf[256], comb_str[256], nb_str[256]; char genpairs[32]; char *dirstr, *dummy2; - int nrcopies, nmol, nmolb = 0, nscan, ncombs, ncopy; + int nrcopies, nmol, nscan, ncombs, ncopy; double fLJ, fQQ, fPOW; - gmx_molblock_t *molb = nullptr; t_molinfo *mi0 = nullptr; DirStack *DS; directive d, newd; @@ -634,7 +626,7 @@ static char **read_topol(const char *infile, const char *outfile, snew(*intermolecular_interactions, 1); init_molinfo(*intermolecular_interactions); mi0 = *intermolecular_interactions; - make_atoms_sys(nmolb, molb, *molinfo, + make_atoms_sys(*molblock, *molinfo, &mi0->atoms); } } @@ -858,10 +850,9 @@ static char **read_topol(const char *infile, const char *outfile, push_mol(nmol, *molinfo, pline, &whichmol, &nrcopies, wi); mi0 = &((*molinfo)[whichmol]); - srenew(molb, nmolb+1); - molb[nmolb].type = whichmol; - molb[nmolb].nmol = nrcopies; - nmolb++; + molblock->resize(molblock->size() + 1); + molblock->back().type = whichmol; + molblock->back().nmol = nrcopies; bCouple = (opts->couple_moltype != nullptr && (gmx_strcasecmp("system", opts->couple_moltype) == 0 || @@ -979,30 +970,26 @@ static char **read_topol(const char *infile, const char *outfile, *nrmols = nmol; - *nmolblock = nmolb; - *molblock = molb; - return title; } -char **do_top(gmx_bool bVerbose, - const char *topfile, - const char *topppfile, - t_gromppopts *opts, - gmx_bool bZero, - t_symtab *symtab, - t_params plist[], - int *combination_rule, - double *repulsion_power, - real *fudgeQQ, - gpp_atomtype_t atype, - int *nrmols, - t_molinfo **molinfo, - t_molinfo **intermolecular_interactions, - const t_inputrec *ir, - int *nmolblock, - gmx_molblock_t **molblock, - warninp_t wi) +char **do_top(gmx_bool bVerbose, + const char *topfile, + const char *topppfile, + t_gromppopts *opts, + gmx_bool bZero, + t_symtab *symtab, + t_params plist[], + int *combination_rule, + double *repulsion_power, + real *fudgeQQ, + gpp_atomtype_t atype, + int *nrmols, + t_molinfo **molinfo, + t_molinfo **intermolecular_interactions, + const t_inputrec *ir, + std::vector *molblock, + warninp_t wi) { /* Tmpfile might contain a long path */ const char *tmpfile; @@ -1025,7 +1012,7 @@ char **do_top(gmx_bool bVerbose, symtab, atype, nrmols, molinfo, intermolecular_interactions, plist, combination_rule, repulsion_power, - opts, fudgeQQ, nmolblock, molblock, + opts, fudgeQQ, molblock, ir->efep != efepNO, bZero, EEL_FULL(ir->coulombtype), wi); @@ -1331,20 +1318,20 @@ void generate_qmexcl(gmx_mtop_t *sys, t_inputrec *ir, warninp_t wi) */ unsigned char *grpnr; - int mb, mol, nat_mol, i, nr_mol_with_qm_atoms = 0; + int mol, nat_mol, nr_mol_with_qm_atoms = 0; gmx_molblock_t *molb; gmx_bool bQMMM; grpnr = sys->groups.grpnr[egcQMMM]; - for (mb = 0; mb < sys->nmolblock; mb++) + for (size_t mb = 0; mb < sys->molblock.size(); mb++) { molb = &sys->molblock[mb]; nat_mol = sys->moltype[molb->type].atoms.nr; for (mol = 0; mol < molb->nmol; mol++) { bQMMM = FALSE; - for (i = 0; i < nat_mol; i++) + for (int i = 0; i < nat_mol; i++) { if ((grpnr ? grpnr[i] : 0) < ir->opts.ngQM) { @@ -1360,12 +1347,8 @@ void generate_qmexcl(gmx_mtop_t *sys, t_inputrec *ir, warninp_t wi) if (mol > 0) { /* Split the molblock at this molecule */ - sys->nmolblock++; - srenew(sys->molblock, sys->nmolblock); - for (i = sys->nmolblock-2; i >= mb; i--) - { - sys->molblock[i+1] = sys->molblock[i]; - } + auto pos = sys->molblock.begin() + mb + 1; + sys->molblock.insert(pos, sys->molblock[mb]); sys->molblock[mb ].nmol = mol; sys->molblock[mb+1].nmol -= mol; mb++; @@ -1374,29 +1357,22 @@ void generate_qmexcl(gmx_mtop_t *sys, t_inputrec *ir, warninp_t wi) if (molb->nmol > 1) { /* Split the molblock after this molecule */ - sys->nmolblock++; - srenew(sys->molblock, sys->nmolblock); + auto pos = sys->molblock.begin() + mb + 1; + sys->molblock.insert(pos, sys->molblock[mb]); molb = &sys->molblock[mb]; - for (i = sys->nmolblock-2; i >= mb; i--) - { - sys->molblock[i+1] = sys->molblock[i]; - } sys->molblock[mb ].nmol = 1; sys->molblock[mb+1].nmol -= 1; } /* Add a moltype for the QMMM molecule */ - sys->nmoltype++; - srenew(sys->moltype, sys->nmoltype); - /* Copy the moltype struct */ - sys->moltype[sys->nmoltype-1] = sys->moltype[molb->type]; + sys->moltype.push_back(sys->moltype[molb->type]); /* 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[sys->nmoltype-1].excls); + &sys->moltype.back().excls); /* Set the molecule type for the QMMM molblock */ - molb->type = sys->nmoltype - 1; + molb->type = sys->moltype.size() - 1; } generate_qmexcl_moltype(&sys->moltype[molb->type], grpnr, ir, wi); } diff --git a/src/gromacs/gmxpreprocess/topio.h b/src/gromacs/gmxpreprocess/topio.h index bb79551e62..3d1fac960a 100644 --- a/src/gromacs/gmxpreprocess/topio.h +++ b/src/gromacs/gmxpreprocess/topio.h @@ -38,6 +38,8 @@ #ifndef GMX_GMXPREPROCESS_TOPIO_H #define GMX_GMXPREPROCESS_TOPIO_H +#include + #include "gromacs/gmxpreprocess/gpp_atomtype.h" #include "gromacs/gmxpreprocess/grompp-impl.h" @@ -48,27 +50,26 @@ struct t_inputrec; struct warninp; typedef warninp *warninp_t; -double check_mol(gmx_mtop_t *mtop, warninp_t wi); +double check_mol(const gmx_mtop_t *mtop, warninp_t wi); /* Check mass and charge */ -char **do_top(gmx_bool bVerbose, - const char *topfile, - const char *topppfile, - t_gromppopts *opts, - gmx_bool bZero, - struct t_symtab *symtab, - t_params plist[], - int *combination_rule, - double *repulsion_power, - real *fudgeQQ, - gpp_atomtype_t atype, - int *nrmols, - t_molinfo **molinfo, - t_molinfo **intermolecular_interactions, - const t_inputrec *ir, - int *nmolblock, - gmx_molblock_t **molblock, - warninp_t wi); +char **do_top(gmx_bool bVerbose, + const char *topfile, + const char *topppfile, + t_gromppopts *opts, + gmx_bool bZero, + struct t_symtab *symtab, + t_params plist[], + int *combination_rule, + double *repulsion_power, + real *fudgeQQ, + gpp_atomtype_t atype, + int *nrmols, + t_molinfo **molinfo, + t_molinfo **intermolecular_interactions, + const t_inputrec *ir, + std::vector *molblock, + warninp_t wi); /* This routine expects sys->molt[m].ilist to be of size F_NRE and ordered. */ void generate_qmexcl(gmx_mtop_t *sys, t_inputrec *ir, warninp_t wi); diff --git a/src/gromacs/gmxpreprocess/toputil.cpp b/src/gromacs/gmxpreprocess/toputil.cpp index b35187d69d..330ffea8fb 100644 --- a/src/gromacs/gmxpreprocess/toputil.cpp +++ b/src/gromacs/gmxpreprocess/toputil.cpp @@ -188,7 +188,7 @@ void init_molinfo(t_molinfo *mol) init_block(&mol->cgs); init_block(&mol->mols); init_blocka(&mol->excls); - init_atom(&mol->atoms); + init_t_atoms(&mol->atoms, 0, FALSE); } /* FREEING MEMORY */ diff --git a/src/gromacs/listed-forces/disre.cpp b/src/gromacs/listed-forces/disre.cpp index 913ea1f24b..959babb668 100644 --- a/src/gromacs/listed-forces/disre.cpp +++ b/src/gromacs/listed-forces/disre.cpp @@ -76,7 +76,7 @@ void init_disres(FILE *fplog, const gmx_mtop_t *mtop, t_disresdata *dd; history_t *hist; gmx_mtop_ilistloop_t iloop; - t_ilist *il; + const t_ilist *il; char *ptr; int type_min, type_max; diff --git a/src/gromacs/listed-forces/orires.cpp b/src/gromacs/listed-forces/orires.cpp index 8f437a31bf..b91d5d001a 100644 --- a/src/gromacs/listed-forces/orires.cpp +++ b/src/gromacs/listed-forces/orires.cpp @@ -113,7 +113,7 @@ void init_orires(FILE *fplog, int typeMin = INT_MAX; int typeMax = 0; gmx_mtop_ilistloop_t iloop = gmx_mtop_ilistloop_init(mtop); - t_ilist *il; + const t_ilist *il; int nmol; while (gmx_mtop_ilistloop_next(iloop, &il, &nmol)) { diff --git a/src/gromacs/mdlib/broadcaststructs.cpp b/src/gromacs/mdlib/broadcaststructs.cpp index e47995f19d..f6c14bbe17 100644 --- a/src/gromacs/mdlib/broadcaststructs.cpp +++ b/src/gromacs/mdlib/broadcaststructs.cpp @@ -51,6 +51,7 @@ #include "gromacs/mdtypes/md_enums.h" #include "gromacs/mdtypes/pull-params.h" #include "gromacs/mdtypes/state.h" +#include "gromacs/topology/mtop_util.h" #include "gromacs/topology/symtab.h" #include "gromacs/topology/topology.h" #include "gromacs/utility/fatalerror.h" @@ -777,7 +778,6 @@ static void bc_atomtypes(const t_commrec *cr, t_atomtypes *atomtypes) static void bcast_ir_mtop(const t_commrec *cr, t_inputrec *inputrec, gmx_mtop_t *mtop) { - int i; if (debug) { fprintf(debug, "in bc_data\n"); @@ -800,11 +800,12 @@ void bcast_ir_mtop(const t_commrec *cr, t_inputrec *inputrec, gmx_mtop_t *mtop) bc_ffparams(cr, &mtop->ffparams); - block_bc(cr, mtop->nmoltype); - snew_bc(cr, mtop->moltype, mtop->nmoltype); - for (i = 0; i < mtop->nmoltype; i++) + int nmoltype = mtop->moltype.size(); + block_bc(cr, nmoltype); + mtop->moltype.resize(nmoltype); + for (gmx_moltype_t &moltype : mtop->moltype) { - bc_moltype(cr, &mtop->symtab, &mtop->moltype[i]); + bc_moltype(cr, &mtop->symtab, &moltype); } block_bc(cr, mtop->bIntermolecularInteractions); @@ -814,11 +815,12 @@ void bcast_ir_mtop(const t_commrec *cr, t_inputrec *inputrec, gmx_mtop_t *mtop) bc_ilists(cr, mtop->intermolecular_ilist); } - block_bc(cr, mtop->nmolblock); - snew_bc(cr, mtop->molblock, mtop->nmolblock); - for (i = 0; i < mtop->nmolblock; i++) + int nmolblock = mtop->molblock.size(); + block_bc(cr, nmolblock); + mtop->molblock.resize(nmolblock); + for (gmx_molblock_t &molblock : mtop->molblock) { - bc_molblock(cr, &mtop->molblock[i]); + bc_molblock(cr, &molblock); } block_bc(cr, mtop->natoms); diff --git a/src/gromacs/mdlib/calc_verletbuf.cpp b/src/gromacs/mdlib/calc_verletbuf.cpp index f730510eea..d9339eef09 100644 --- a/src/gromacs/mdlib/calc_verletbuf.cpp +++ b/src/gromacs/mdlib/calc_verletbuf.cpp @@ -1,7 +1,7 @@ /* * This file is part of the GROMACS molecular simulation package. * - * Copyright (c) 2012,2013,2014,2015,2016,2017, by the GROMACS development team, led by + * Copyright (c) 2012,2013,2014,2015,2016,2017,2018, 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. @@ -338,8 +338,7 @@ static void get_verlet_buffer_atomtypes(const gmx_mtop_t *mtop, { verletbuf_atomtype_t *att; int natt; - int mb, nmol, ft, i, a1, a2, a3, a; - const t_atoms *atoms; + int ft, i, a1, a2, a3, a; const t_ilist *il; const t_iparams *ip; atom_nonbonded_kinetic_prop_t *prop; @@ -354,11 +353,11 @@ static void get_verlet_buffer_atomtypes(const gmx_mtop_t *mtop, *n_nonlin_vsite = 0; } - for (mb = 0; mb < mtop->nmolblock; mb++) + for (const gmx_molblock_t &molblock : mtop->molblock) { - nmol = mtop->molblock[mb].nmol; - - atoms = &mtop->moltype[mtop->molblock[mb].type].atoms; + int nmol = molblock.nmol; + const gmx_moltype_t &moltype = mtop->moltype[molblock.type]; + const t_atoms *atoms = &moltype.atoms; /* Check for constraints, as they affect the kinetic energy. * For virtual sites we need the masses and geometry of @@ -369,7 +368,7 @@ static void get_verlet_buffer_atomtypes(const gmx_mtop_t *mtop, for (ft = F_CONSTR; ft <= F_CONSTRNC; ft++) { - il = &mtop->moltype[mtop->molblock[mb].type].ilist[ft]; + il = &moltype.ilist[ft]; for (i = 0; i < il->nr; i += 1+NRAL(ft)) { @@ -389,7 +388,7 @@ static void get_verlet_buffer_atomtypes(const gmx_mtop_t *mtop, } } - il = &mtop->moltype[mtop->molblock[mb].type].ilist[F_SETTLE]; + il = &moltype.ilist[F_SETTLE]; for (i = 0; i < il->nr; i += 1+NRAL(F_SETTLE)) { @@ -411,7 +410,7 @@ static void get_verlet_buffer_atomtypes(const gmx_mtop_t *mtop, prop[a3].con_len = ip->settle.doh; } - get_vsite_masses(&mtop->moltype[mtop->molblock[mb].type], + get_vsite_masses(&moltype, &mtop->ffparams, vsite_m, &n_nonlin_vsite_mol); diff --git a/src/gromacs/mdlib/constr.cpp b/src/gromacs/mdlib/constr.cpp index ab98e6613c..122a54db07 100644 --- a/src/gromacs/mdlib/constr.cpp +++ b/src/gromacs/mdlib/constr.cpp @@ -771,20 +771,20 @@ gmx_constr_t init_constraints(FILE *fplog, constr->nflexcon = 0; if (nconstraints > 0) { - constr->n_at2con_mt = mtop->nmoltype; + constr->n_at2con_mt = mtop->moltype.size(); snew(constr->at2con_mt, constr->n_at2con_mt); - for (int mt = 0; mt < mtop->nmoltype; mt++) + for (int mt = 0; mt < static_cast(mtop->moltype.size()); mt++) { int nflexcon; constr->at2con_mt[mt] = make_at2con(0, mtop->moltype[mt].atoms.nr, mtop->moltype[mt].ilist, mtop->ffparams.iparams, EI_DYNAMICS(ir->eI), &nflexcon); - for (int i = 0; i < mtop->nmolblock; i++) + for (const gmx_molblock_t &molblock : mtop->molblock) { - if (mtop->molblock[i].type == mt) + if (molblock.type == mt) { - constr->nflexcon += mtop->molblock[i].nmol*nflexcon; + constr->nflexcon += molblock.nmol*nflexcon; } } } @@ -848,9 +848,9 @@ gmx_constr_t init_constraints(FILE *fplog, constr->settled = settle_init(mtop); /* Make an atom to settle index for use in domain decomposition */ - constr->n_at2settle_mt = mtop->nmoltype; + constr->n_at2settle_mt = mtop->moltype.size(); snew(constr->at2settle_mt, constr->n_at2settle_mt); - for (int mt = 0; mt < mtop->nmoltype; mt++) + for (size_t mt = 0; mt < mtop->moltype.size(); mt++) { constr->at2settle_mt[mt] = make_at2settle(mtop->moltype[mt].atoms.nr, @@ -924,12 +924,11 @@ gmx_bool inter_charge_group_constraints(const gmx_mtop_t *mtop) const gmx_moltype_t *molt; const t_block *cgs; const t_ilist *il; - int mb; int *at2cg, cg, a, ftype, i; gmx_bool bInterCG; bInterCG = FALSE; - for (mb = 0; mb < mtop->nmolblock && !bInterCG; mb++) + for (size_t mb = 0; mb < mtop->molblock.size() && !bInterCG; mb++) { molt = &mtop->moltype[mtop->molblock[mb].type]; @@ -971,12 +970,11 @@ gmx_bool inter_charge_group_settles(const gmx_mtop_t *mtop) const gmx_moltype_t *molt; const t_block *cgs; const t_ilist *il; - int mb; int *at2cg, cg, a, ftype, i; gmx_bool bInterCG; bInterCG = FALSE; - for (mb = 0; mb < mtop->nmolblock && !bInterCG; mb++) + for (size_t mb = 0; mb < mtop->molblock.size() && !bInterCG; mb++) { molt = &mtop->moltype[mtop->molblock[mb].type]; diff --git a/src/gromacs/mdlib/constraintrange.cpp b/src/gromacs/mdlib/constraintrange.cpp index 0a46e7079d..abe685a04d 100644 --- a/src/gromacs/mdlib/constraintrange.cpp +++ b/src/gromacs/mdlib/constraintrange.cpp @@ -219,14 +219,11 @@ static real constr_r_max_moltype(const gmx_moltype_t *molt, real constr_r_max(FILE *fplog, const gmx_mtop_t *mtop, const t_inputrec *ir) { - int mt; - real rmax; - - rmax = 0; - for (mt = 0; mt < mtop->nmoltype; mt++) + real rmax = 0; + for (const gmx_moltype_t &molt : mtop->moltype) { rmax = std::max(rmax, - constr_r_max_moltype(&mtop->moltype[mt], + constr_r_max_moltype(&molt, mtop->ffparams.iparams, ir)); } diff --git a/src/gromacs/mdlib/forcerec.cpp b/src/gromacs/mdlib/forcerec.cpp index 5446b5151d..70b7530c3d 100644 --- a/src/gromacs/mdlib/forcerec.cpp +++ b/src/gromacs/mdlib/forcerec.cpp @@ -534,7 +534,7 @@ check_solvent(FILE * fp, { const t_block * cgs; const gmx_moltype_t *molt; - int mb, mol, cg_mol, at_offset, am, cgm, i, nmol_ch, nmol; + int mol, cg_mol, at_offset, am, cgm, i, nmol_ch, nmol; int n_solvent_parameters; solvent_parameters_t *solvent_parameters; int **cg_sp; @@ -548,10 +548,10 @@ check_solvent(FILE * fp, n_solvent_parameters = 0; solvent_parameters = nullptr; /* Allocate temporary array for solvent type */ - snew(cg_sp, mtop->nmolblock); + snew(cg_sp, mtop->molblock.size()); at_offset = 0; - for (mb = 0; mb < mtop->nmolblock; mb++) + for (size_t mb = 0; mb < mtop->molblock.size(); mb++) { molt = &mtop->moltype[mtop->molblock[mb].type]; cgs = &molt->cgs; @@ -605,7 +605,7 @@ check_solvent(FILE * fp, } fr->nWatMol = 0; - for (mb = 0; mb < mtop->nmolblock; mb++) + for (size_t mb = 0; mb < mtop->molblock.size(); mb++) { cgs = &mtop->moltype[mtop->molblock[mb].type].cgs; nmol = (mtop->molblock[mb].nmol*cgs->nr)/cginfo_mb[mb].cg_mod; @@ -653,13 +653,13 @@ static cginfo_mb_t *init_cginfo_mb(FILE *fplog, const gmx_mtop_t *mtop, gmx_bool *type_VDW; int *cginfo; int cg_offset, a_offset; - int mb, m, cg, a0, a1, gid, ai, j, aj, excl_nalloc; + int m, cg, a0, a1, gid, ai, j, aj, excl_nalloc; int *a_con; int ftype; int ia; gmx_bool bId, *bExcl, bExclIntraAll, bExclInter, bHaveVDW, bHaveQ, bHavePerturbedAtoms; - snew(cginfo_mb, mtop->nmolblock); + snew(cginfo_mb, mtop->molblock.size()); snew(type_VDW, fr->ntype); for (ai = 0; ai < fr->ntype; ai++) @@ -681,7 +681,7 @@ static cginfo_mb_t *init_cginfo_mb(FILE *fplog, const gmx_mtop_t *mtop, snew(bExcl, excl_nalloc); cg_offset = 0; a_offset = 0; - for (mb = 0; mb < mtop->nmolblock; mb++) + for (size_t mb = 0; mb < mtop->molblock.size(); mb++) { molb = &mtop->molblock[mb]; molt = &mtop->moltype[molb->type]; @@ -889,7 +889,7 @@ static cginfo_mb_t *init_cginfo_mb(FILE *fplog, const gmx_mtop_t *mtop, } if (!fr->solvent_opt) { - for (mb = 0; mb < mtop->nmolblock; mb++) + for (size_t mb = 0; mb < mtop->molblock.size(); mb++) { for (cg = 0; cg < cginfo_mb[mb].cg_mod; cg++) { @@ -942,17 +942,15 @@ static bool set_chargesum(FILE *log, t_forcerec *fr, const gmx_mtop_t *mtop) { /*This now calculates sum for q and c6*/ double qsum, q2sum, q, c6sum, c6; - int mb, nmol, i; - const t_atoms *atoms; qsum = 0; q2sum = 0; c6sum = 0; - for (mb = 0; mb < mtop->nmolblock; mb++) + for (const gmx_molblock_t &molb : mtop->molblock) { - nmol = mtop->molblock[mb].nmol; - atoms = &mtop->moltype[mtop->molblock[mb].type].atoms; - for (i = 0; i < atoms->nr; i++) + int nmol = molb.nmol; + const t_atoms *atoms = &mtop->moltype[molb.type].atoms; + for (int i = 0; i < atoms->nr; i++) { q = atoms->atom[i].q; qsum += nmol*q; @@ -970,11 +968,11 @@ static bool set_chargesum(FILE *log, t_forcerec *fr, const gmx_mtop_t *mtop) qsum = 0; q2sum = 0; c6sum = 0; - for (mb = 0; mb < mtop->nmolblock; mb++) + for (const gmx_molblock_t &molb : mtop->molblock) { - nmol = mtop->molblock[mb].nmol; - atoms = &mtop->moltype[mtop->molblock[mb].type].atoms; - for (i = 0; i < atoms->nr; i++) + int nmol = molb.nmol; + const t_atoms *atoms = &mtop->moltype[molb.type].atoms; + for (int i = 0; i < atoms->nr; i++) { q = atoms->atom[i].qB; qsum += nmol*q; @@ -1025,7 +1023,7 @@ void set_avcsixtwelve(FILE *fplog, t_forcerec *fr, const gmx_mtop_t *mtop) { const t_atoms *atoms, *atoms_tpi; const t_blocka *excl; - int mb, nmol, nmolc, i, j, tpi, tpj, j1, j2, k, nexcl, q; + int nmolc, i, j, tpi, tpj, j1, j2, k, nexcl, q; gmx_int64_t npair, npair_ij, tmpi, tmpj; double csix, ctwelve; int ntp, *typecount; @@ -1103,12 +1101,12 @@ void set_avcsixtwelve(FILE *fplog, t_forcerec *fr, const gmx_mtop_t *mtop) * any value. These unused values should not influence the dispersion * correction. */ - for (mb = 0; mb < mtop->nmolblock; mb++) + for (const gmx_molblock_t &molb : mtop->molblock) { - nmol = mtop->molblock[mb].nmol; - atoms = &mtop->moltype[mtop->molblock[mb].type].atoms; - excl = &mtop->moltype[mtop->molblock[mb].type].excls; - for (i = 0; (i < atoms->nr); i++) + int nmol = molb.nmol; + atoms = &mtop->moltype[molb.type].atoms; + excl = &mtop->moltype[molb.type].excls; + for (int i = 0; (i < atoms->nr); i++) { if (q == 0) { @@ -1144,7 +1142,7 @@ void set_avcsixtwelve(FILE *fplog, t_forcerec *fr, const gmx_mtop_t *mtop) csix -= nmol*C6 (nbfp, ntp, tpi, tpj)/6.0; ctwelve -= nmol*C12(nbfp, ntp, tpi, tpj)/12.0; } - nexcl += nmol; + nexcl += molb.nmol; } } } @@ -1156,24 +1154,24 @@ void set_avcsixtwelve(FILE *fplog, t_forcerec *fr, const gmx_mtop_t *mtop) * with the rest of the system. */ atoms_tpi = - &mtop->moltype[mtop->molblock[mtop->nmolblock-1].type].atoms; + &mtop->moltype[mtop->molblock.back().type].atoms; npair = 0; - for (mb = 0; mb < mtop->nmolblock; mb++) + for (size_t mb = 0; mb < mtop->molblock.size(); mb++) { - nmol = mtop->molblock[mb].nmol; - atoms = &mtop->moltype[mtop->molblock[mb].type].atoms; + const gmx_molblock_t &molb = mtop->molblock[mb]; + atoms = &mtop->moltype[molb.type].atoms; for (j = 0; j < atoms->nr; j++) { - nmolc = nmol; + nmolc = molb.nmol; /* Remove the interaction of the test charge group * with itself. */ - if (mb == mtop->nmolblock-1) + if (mb == mtop->molblock.size() - 1) { nmolc--; - if (mb == 0 && nmol == 1) + if (mb == 0 && molb.nmol == 1) { gmx_fatal(FARGS, "Old format tpr with TPI, please generate a new tpr file"); } @@ -1257,7 +1255,7 @@ void set_avcsixtwelve(FILE *fplog, t_forcerec *fr, const gmx_mtop_t *mtop) static real calcBuckinghamBMax(FILE *fplog, const gmx_mtop_t *mtop) { const t_atoms *at1, *at2; - int mt1, mt2, i, j, tpi, tpj, ntypes; + int i, j, tpi, tpj, ntypes; real b, bmin; if (fplog) @@ -1268,7 +1266,7 @@ static real calcBuckinghamBMax(FILE *fplog, const gmx_mtop_t *mtop) bmin = -1; real bham_b_max = 0; - for (mt1 = 0; mt1 < mtop->nmoltype; mt1++) + for (size_t mt1 = 0; mt1 < mtop->moltype.size(); mt1++) { at1 = &mtop->moltype[mt1].atoms; for (i = 0; (i < at1->nr); i++) @@ -1279,7 +1277,7 @@ static real calcBuckinghamBMax(FILE *fplog, const gmx_mtop_t *mtop) gmx_fatal(FARGS, "Atomtype[%d] = %d, maximum = %d", i, tpi, ntypes); } - for (mt2 = mt1; mt2 < mtop->nmoltype; mt2++) + for (size_t mt2 = mt1; mt2 < mtop->moltype.size(); mt2++) { at2 = &mtop->moltype[mt2].atoms; for (j = 0; (j < at2->nr); j++) @@ -1400,21 +1398,19 @@ static void make_nbf_tables(FILE *fp, static void count_tables(int ftype1, int ftype2, const gmx_mtop_t *mtop, int *ncount, int **count) { - const gmx_moltype_t *molt; const t_ilist *il; - int mt, ftype, stride, i, j, tabnr; + int ftype, stride, i, j, tabnr; // Loop over all moleculetypes - for (mt = 0; mt < mtop->nmoltype; mt++) + for (const gmx_moltype_t &molt : mtop->moltype) { - molt = &mtop->moltype[mt]; // Loop over all interaction types for (ftype = 0; ftype < F_NRE; ftype++) { // If the current interaction type is one of the types whose tables we're trying to count... if (ftype == ftype1 || ftype == ftype2) { - il = &molt->ilist[ftype]; + il = &molt.ilist[ftype]; stride = 1 + NRAL(ftype); // ... and there are actually some interactions for this type for (i = 0; i < il->nr; i += stride) @@ -2379,7 +2375,7 @@ void init_forcerec(FILE *fp, /* Because of old style topologies, we have to use the last cg * instead of the last molecule type. */ - cgs = &mtop->moltype[mtop->molblock[mtop->nmolblock-1].type].cgs; + cgs = &mtop->moltype[mtop->molblock.back().type].cgs; fr->n_tpi = cgs->index[cgs->nr] - cgs->index[cgs->nr-1]; gmx::BlockRanges molecules = gmx_mtop_molecules(*mtop); if (fr->n_tpi != molecules.index[molecules.numBlocks()] - molecules.index[molecules.numBlocks() - 1]) @@ -3059,7 +3055,7 @@ void init_forcerec(FILE *fp, } else { - fr->cginfo = cginfo_expand(mtop->nmolblock, fr->cginfo_mb); + fr->cginfo = cginfo_expand(mtop->molblock.size(), fr->cginfo_mb); } if (!DOMAINDECOMP(cr)) diff --git a/src/gromacs/mdlib/lincs.cpp b/src/gromacs/mdlib/lincs.cpp index 5adaf1bd2d..c9ec014f3c 100644 --- a/src/gromacs/mdlib/lincs.cpp +++ b/src/gromacs/mdlib/lincs.cpp @@ -1412,8 +1412,6 @@ gmx_lincsdata_t init_lincs(FILE *fplog, const gmx_mtop_t *mtop, gmx_bool bPLINCS, int nIter, int nProjOrder) { gmx_lincsdata *li; - int mt, mb; - gmx_moltype_t *molt; gmx_bool bMoreThanTwoSeq; if (fplog) @@ -1433,12 +1431,9 @@ gmx_lincsdata_t init_lincs(FILE *fplog, const gmx_mtop_t *mtop, li->nOrder = nProjOrder; li->max_connect = 0; - for (mt = 0; mt < mtop->nmoltype; mt++) + for (size_t mt = 0; mt < mtop->moltype.size(); mt++) { - int a; - - molt = &mtop->moltype[mt]; - for (a = 0; a < molt->atoms.nr; a++) + for (int a = 0; a < mtop->moltype[mt].atoms.nr; a++) { li->max_connect = std::max(li->max_connect, at2con[mt].index[a + 1] - at2con[mt].index[a]); @@ -1447,17 +1442,16 @@ gmx_lincsdata_t init_lincs(FILE *fplog, const gmx_mtop_t *mtop, li->ncg_triangle = 0; bMoreThanTwoSeq = FALSE; - for (mb = 0; mb < mtop->nmolblock; mb++) + for (const gmx_molblock_t &molb : mtop->molblock) { - molt = &mtop->moltype[mtop->molblock[mb].type]; + const gmx_moltype_t &molt = mtop->moltype[molb.type]; li->ncg_triangle += - mtop->molblock[mb].nmol* - count_triangle_constraints(molt->ilist, - &at2con[mtop->molblock[mb].type]); + molb.nmol* + count_triangle_constraints(molt.ilist, &at2con[molb.type]); if (!bMoreThanTwoSeq && - more_than_two_sequential_constraints(molt->ilist, &at2con[mtop->molblock[mb].type])) + more_than_two_sequential_constraints(molt.ilist, &at2con[molb.type])) { bMoreThanTwoSeq = TRUE; } diff --git a/src/gromacs/mdlib/ns.cpp b/src/gromacs/mdlib/ns.cpp index d61f0c5f60..e3e792f43b 100644 --- a/src/gromacs/mdlib/ns.cpp +++ b/src/gromacs/mdlib/ns.cpp @@ -2047,14 +2047,13 @@ void init_ns(FILE *fplog, const t_commrec *cr, gmx_ns_t *ns, t_forcerec *fr, const gmx_mtop_t *mtop) { - int mt, icg, nr_in_cg, maxcg, i, j, jcg, ngid, ncg; - t_block *cgs; + int icg, nr_in_cg, maxcg, i, j, jcg, ngid, ncg; /* Compute largest charge groups size (# atoms) */ nr_in_cg = 1; - for (mt = 0; mt < mtop->nmoltype; mt++) + for (const gmx_moltype_t &molt : mtop->moltype) { - cgs = &mtop->moltype[mt].cgs; + const t_block *cgs = &molt.cgs; for (icg = 0; (icg < cgs->nr); icg++) { nr_in_cg = std::max(nr_in_cg, (int)(cgs->index[icg+1]-cgs->index[icg])); diff --git a/src/gromacs/mdlib/perf_est.cpp b/src/gromacs/mdlib/perf_est.cpp index 6888f3b0bc..d638869306 100644 --- a/src/gromacs/mdlib/perf_est.cpp +++ b/src/gromacs/mdlib/perf_est.cpp @@ -3,7 +3,7 @@ * * Copyright (c) 1991-2000, University of Groningen, The Netherlands. * Copyright (c) 2001-2008, The GROMACS development team. - * Copyright (c) 2012,2014,2015,2016,2017, by the GROMACS development team, led by + * Copyright (c) 2012,2014,2015,2016,2017,2018, 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. @@ -171,8 +171,7 @@ void count_bonded_distances(const gmx_mtop_t *mtop, const t_inputrec *ir, { gmx_bool bExcl; double nonsimd_step_frac; - int mb, nmol, ftype; - gmx_moltype_t *molt; + int ftype; double ndtot_c, ndtot_simd; #if GMX_SIMD_HAVE_REAL gmx_bool bSimdBondeds = TRUE; @@ -211,10 +210,9 @@ void count_bonded_distances(const gmx_mtop_t *mtop, const t_inputrec *ir, */ ndtot_c = 0; ndtot_simd = 0; - for (mb = 0; mb < mtop->nmolblock; mb++) + for (const gmx_molblock_t &molb : mtop->molblock) { - molt = &mtop->moltype[mtop->molblock[mb].type]; - nmol = mtop->molblock[mb].nmol; + const gmx_moltype_t *molt = &mtop->moltype[molb.type]; for (ftype = 0; ftype < F_NRE; ftype++) { int nbonds; @@ -248,14 +246,14 @@ void count_bonded_distances(const gmx_mtop_t *mtop, const t_inputrec *ir, nd_c = NRAL(ftype) - 1; break; } - nbonds = nmol*molt->ilist[ftype].nr/(1 + NRAL(ftype)); + nbonds = molb.nmol*molt->ilist[ftype].nr/(1 + NRAL(ftype)); ndtot_c += nbonds*nd_c; ndtot_simd += nbonds*nd_simd; } } if (bExcl) { - ndtot_c += nmol*(molt->excls.nra - molt->atoms.nr)/2; + ndtot_c += molb.nmol*(molt->excls.nra - molt->atoms.nr)/2; } } @@ -280,13 +278,11 @@ static void pp_group_load(const gmx_mtop_t *mtop, const t_inputrec *ir, double *cost_pp, gmx_bool *bChargePerturbed, gmx_bool *bTypePerturbed) { - t_atom *atom; - int mb, nmol, atnr, cg, a, a0, ncqlj, ncq, nclj; + int atnr, cg, a0, ncqlj, ncq, nclj; gmx_bool bBHAM, bLJcut, bWater, bQ, bLJ; int nw, nqlj, nq, nlj; double fq, fqlj, flj, fqljw, fqw; t_iparams *iparams; - gmx_moltype_t *molt; bBHAM = (mtop->ffparams.functype[0] == F_BHAM); @@ -314,12 +310,11 @@ static void pp_group_load(const gmx_mtop_t *mtop, const t_inputrec *ir, nlj = 0; *bChargePerturbed = FALSE; *bTypePerturbed = FALSE; - for (mb = 0; mb < mtop->nmolblock; mb++) + for (const gmx_molblock_t &molb : mtop->molblock) { - molt = &mtop->moltype[mtop->molblock[mb].type]; - atom = molt->atoms.atom; - nmol = mtop->molblock[mb].nmol; - a = 0; + const gmx_moltype_t *molt = &mtop->moltype[molb.type]; + const t_atom *atom = molt->atoms.atom; + int a = 0; for (cg = 0; cg < molt->cgs.nr; cg++) { bWater = !bBHAM; @@ -367,13 +362,13 @@ static void pp_group_load(const gmx_mtop_t *mtop, const t_inputrec *ir, } if (bWater) { - nw += nmol; + nw += molb.nmol; } else { - nqlj += nmol*ncqlj; - nq += nmol*ncq; - nlj += nmol*nclj; + nqlj += molb.nmol*ncqlj; + nq += molb.nmol*ncq; + nlj += molb.nmol*nclj; } } } @@ -407,11 +402,9 @@ static void pp_verlet_load(const gmx_mtop_t *mtop, const t_inputrec *ir, double *cost_pp, gmx_bool *bChargePerturbed, gmx_bool *bTypePerturbed) { - t_atom *atom; - int mb, nmol, atnr, a, nqlj, nq, nlj; + int atnr, a, nqlj, nq, nlj; gmx_bool bQRF; t_iparams *iparams; - gmx_moltype_t *molt; real r_eff; double c_qlj, c_q, c_lj; double nppa; @@ -433,11 +426,10 @@ static void pp_verlet_load(const gmx_mtop_t *mtop, const t_inputrec *ir, nq = 0; *bChargePerturbed = FALSE; *bTypePerturbed = FALSE; - for (mb = 0; mb < mtop->nmolblock; mb++) + for (const gmx_molblock_t &molb : mtop->molblock) { - molt = &mtop->moltype[mtop->molblock[mb].type]; - atom = molt->atoms.atom; - nmol = mtop->molblock[mb].nmol; + const gmx_moltype_t *molt = &mtop->moltype[molb.type]; + const t_atom *atom = molt->atoms.atom; for (a = 0; a < molt->atoms.nr; a++) { if (atom[a].q != 0 || atom[a].qB != 0) @@ -445,11 +437,11 @@ static void pp_verlet_load(const gmx_mtop_t *mtop, const t_inputrec *ir, if (iparams[(atnr+1)*atom[a].type].lj.c6 != 0 || iparams[(atnr+1)*atom[a].type].lj.c12 != 0) { - nqlj += nmol; + nqlj += molb.nmol; } else { - nq += nmol; + nq += molb.nmol; } } if (atom[a].q != atom[a].qB) diff --git a/src/gromacs/mdlib/qmmm.cpp b/src/gromacs/mdlib/qmmm.cpp index 7e834b8dfe..22b5543975 100644 --- a/src/gromacs/mdlib/qmmm.cpp +++ b/src/gromacs/mdlib/qmmm.cpp @@ -380,7 +380,7 @@ void init_QMMMrec(const t_commrec *cr, gmx_mtop_atomloop_all_t aloop; gmx_mtop_ilistloop_all_t iloop; int a_offset; - t_ilist *ilist_mol; + const t_ilist *ilist_mol; if (ir->cutoff_scheme != ecutsGROUP) { diff --git a/src/gromacs/mdlib/rf_util.cpp b/src/gromacs/mdlib/rf_util.cpp index 298d6e768a..0b1cb15831 100644 --- a/src/gromacs/mdlib/rf_util.cpp +++ b/src/gromacs/mdlib/rf_util.cpp @@ -3,7 +3,7 @@ * * Copyright (c) 1991-2000, University of Groningen, The Netherlands. * Copyright (c) 2001-2004, The GROMACS development team. - * Copyright (c) 2013,2014,2015,2017, by the GROMACS development team, led by + * Copyright (c) 2013,2014,2015,2017,2018, 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. @@ -280,7 +280,7 @@ void init_generalized_rf(FILE *fplog, const gmx_mtop_t *mtop, const t_inputrec *ir, t_forcerec *fr) { - int mb, i, j; + int i, j; real q, zsq, nrdf, T; const gmx_moltype_t *molt; const t_block *cgs; @@ -290,9 +290,9 @@ void init_generalized_rf(FILE *fplog, fprintf(fplog, "\nWARNING: the generalized reaction field constants are determined from topology A only\n\n"); } zsq = 0.0; - for (mb = 0; mb < mtop->nmolblock; mb++) + for (const gmx_molblock_t &molb : mtop->molblock) { - molt = &mtop->moltype[mtop->molblock[mb].type]; + molt = &mtop->moltype[molb.type]; cgs = &molt->cgs; for (i = 0; (i < cgs->nr); i++) { @@ -301,7 +301,7 @@ void init_generalized_rf(FILE *fplog, { q += molt->atoms.atom[j].q; } - zsq += mtop->molblock[mb].nmol*q*q; + zsq += molb.nmol*q*q; } fr->zsquare = zsq; } diff --git a/src/gromacs/mdlib/settle.cpp b/src/gromacs/mdlib/settle.cpp index 71deaff160..bca8fe38f8 100644 --- a/src/gromacs/mdlib/settle.cpp +++ b/src/gromacs/mdlib/settle.cpp @@ -174,7 +174,7 @@ gmx_settledata_t settle_init(const gmx_mtop_t *mtop) /* Check that we have only one settle type */ int settle_type = -1; gmx_mtop_ilistloop_t iloop = gmx_mtop_ilistloop_init(mtop); - t_ilist *ilist; + const t_ilist *ilist; int nmol; const int nral1 = 1 + NRAL(F_SETTLE); while (gmx_mtop_ilistloop_next(iloop, &ilist, &nmol)) diff --git a/src/gromacs/mdlib/shellfc.cpp b/src/gromacs/mdlib/shellfc.cpp index c40022a11d..b57c3eed98 100644 --- a/src/gromacs/mdlib/shellfc.cpp +++ b/src/gromacs/mdlib/shellfc.cpp @@ -251,8 +251,8 @@ static void predict_shells(FILE *fplog, rvec x[], rvec v[], real dt, * \param[in] mtop Molecular topology. * \returns Array holding the number of particles of a type */ -static std::array countPtypes(FILE *fplog, - gmx_mtop_t *mtop) +static std::array countPtypes(FILE *fplog, + const gmx_mtop_t *mtop) { std::array nptype = { { 0 } }; /* Count number of shells, and find their indices */ @@ -295,7 +295,7 @@ static std::array countPtypes(FILE *fplog, } gmx_shellfc_t *init_shell_flexcon(FILE *fplog, - gmx_mtop_t *mtop, int nflexcon, + const gmx_mtop_t *mtop, int nflexcon, int nstcalcenergy, bool usingDomainDecomposition) { @@ -305,17 +305,14 @@ gmx_shellfc_t *init_shell_flexcon(FILE *fplog, const t_atom *atom; int ns, nshell, nsi; - int i, j, type, mb, a_offset, cg, mol, ftype, nra; + int i, j, type, a_offset, cg, mol, ftype, nra; real qS, alpha; int aS, aN = 0; /* Shell and nucleus */ int bondtypes[] = { F_BONDS, F_HARMONIC, F_CUBICBONDS, F_POLARIZATION, F_ANHARM_POL, F_WATER_POL }; #define NBT asize(bondtypes) t_iatom *ia; gmx_mtop_atomloop_all_t aloop; - gmx_ffparams_t *ffparams; - gmx_molblock_t *molb; - gmx_moltype_t *molt; - t_block *cgs; + const gmx_ffparams_t *ffparams; std::array n = countPtypes(fplog, mtop); nshell = n[eptShell]; @@ -387,12 +384,12 @@ gmx_shellfc_t *init_shell_flexcon(FILE *fplog, shfc->bInterCG = FALSE; ns = 0; a_offset = 0; - for (mb = 0; mb < mtop->nmolblock; mb++) + for (size_t mb = 0; mb < mtop->molblock.size(); mb++) { - molb = &mtop->molblock[mb]; - molt = &mtop->moltype[molb->type]; + const gmx_molblock_t *molb = &mtop->molblock[mb]; + const gmx_moltype_t *molt = &mtop->moltype[molb->type]; + const t_block *cgs = &molt->cgs; - cgs = &molt->cgs; snew(at2cg, molt->atoms.nr); for (cg = 0; cg < cgs->nr; cg++) { diff --git a/src/gromacs/mdlib/shellfc.h b/src/gromacs/mdlib/shellfc.h index b0b95cd9e1..b15c3e17a8 100644 --- a/src/gromacs/mdlib/shellfc.h +++ b/src/gromacs/mdlib/shellfc.h @@ -59,7 +59,7 @@ class t_state; /* Initialization function, also predicts the initial shell postions. */ gmx_shellfc_t *init_shell_flexcon(FILE *fplog, - gmx_mtop_t *mtop, int nflexcon, + const gmx_mtop_t *mtop, int nflexcon, int nstcalcenergy, bool usingDomainDecomposition); diff --git a/src/gromacs/mdlib/sim_util.cpp b/src/gromacs/mdlib/sim_util.cpp index 876243c6bc..184db980eb 100644 --- a/src/gromacs/mdlib/sim_util.cpp +++ b/src/gromacs/mdlib/sim_util.cpp @@ -2582,8 +2582,7 @@ static void low_do_pbc_mtop(FILE *fplog, int ePBC, const matrix box, gmx_bool bFirst) { t_graph *graph; - int mb, as, mol; - gmx_molblock_t *molb; + int as, mol; if (bFirst && fplog) { @@ -2592,22 +2591,21 @@ static void low_do_pbc_mtop(FILE *fplog, int ePBC, const matrix box, snew(graph, 1); as = 0; - for (mb = 0; mb < mtop->nmolblock; mb++) + for (const gmx_molblock_t &molb : mtop->molblock) { - molb = &mtop->molblock[mb]; - if (molb->natoms_mol == 1 || - (!bFirst && mtop->moltype[molb->type].cgs.nr == 1)) + if (molb.natoms_mol == 1 || + (!bFirst && mtop->moltype[molb.type].cgs.nr == 1)) { /* Just one atom or charge group in the molecule, no PBC required */ - as += molb->nmol*molb->natoms_mol; + as += molb.nmol*molb.natoms_mol; } else { /* Pass NULL iso fplog to avoid graph prints for each molecule type */ - mk_graph_ilist(nullptr, mtop->moltype[molb->type].ilist, - 0, molb->natoms_mol, FALSE, FALSE, graph); + mk_graph_ilist(nullptr, mtop->moltype[molb.type].ilist, + 0, molb.natoms_mol, FALSE, FALSE, graph); - for (mol = 0; mol < molb->nmol; mol++) + for (mol = 0; mol < molb.nmol; mol++) { mk_mshift(fplog, graph, ePBC, box, x+as); @@ -2617,7 +2615,7 @@ static void low_do_pbc_mtop(FILE *fplog, int ePBC, const matrix box, * since we no longer need this graph. */ - as += molb->natoms_mol; + as += molb.natoms_mol; } done_graph(graph); } diff --git a/src/gromacs/mdlib/tests/settle.cpp b/src/gromacs/mdlib/tests/settle.cpp index f717fc0207..c10479a12a 100644 --- a/src/gromacs/mdlib/tests/settle.cpp +++ b/src/gromacs/mdlib/tests/settle.cpp @@ -196,37 +196,31 @@ TEST_P(SettleTest, SatisfiesConstraints) const int atomsPerSettle = NRAL(F_SETTLE); ASSERT_LE(numSettles, updatedPositions_.size() / (atomsPerSettle * DIM)) << "cannot test that many SETTLEs " << testDescription; - // Set up the topology. We still have to make some raw pointers, - // but they are put into scope guards for automatic cleanup. - gmx_mtop_t *mtop; - snew(mtop, 1); - const unique_cptr mtopGuard(mtop); - mtop->nmoltype = 1; - snew(mtop->moltype, mtop->nmoltype); - const unique_cptr moltypeGuard(mtop->moltype); - mtop->nmolblock = 1; - snew(mtop->molblock, mtop->nmolblock); - const unique_cptr molblockGuard(mtop->molblock); - mtop->molblock[0].type = 0; - std::vector iatoms; + // Set up the topology. + gmx_mtop_t mtop; + mtop.moltype.resize(1); + mtop.molblock.resize(1); + mtop.molblock[0].type = 0; + const int settleStride = 1 + atomsPerSettle; + int *iatoms; + snew(iatoms, numSettles*settleStride); for (int i = 0; i < numSettles; ++i) { - iatoms.push_back(settleType); - iatoms.push_back(i*atomsPerSettle+0); - iatoms.push_back(i*atomsPerSettle+1); - iatoms.push_back(i*atomsPerSettle+2); + iatoms[i*settleStride + 0] = settleType; + iatoms[i*settleStride + 1] = i*atomsPerSettle + 0; + iatoms[i*settleStride + 2] = i*atomsPerSettle + 1; + iatoms[i*settleStride + 3] = i*atomsPerSettle + 2; } - mtop->moltype[0].ilist[F_SETTLE].iatoms = iatoms.data(); - mtop->moltype[0].ilist[F_SETTLE].nr = iatoms.size(); + mtop.moltype[0].ilist[F_SETTLE].iatoms = iatoms; + mtop.moltype[0].ilist[F_SETTLE].nr = numSettles*settleStride; // Set up the SETTLE parameters. - mtop->ffparams.ntypes = 1; - snew(mtop->ffparams.iparams, mtop->ffparams.ntypes); - const unique_cptr iparamsGuard(mtop->ffparams.iparams); + mtop.ffparams.ntypes = 1; + snew(mtop.ffparams.iparams, mtop.ffparams.ntypes); const real dOH = 0.09572; const real dHH = 0.15139; - mtop->ffparams.iparams[settleType].settle.doh = dOH; - mtop->ffparams.iparams[settleType].settle.dhh = dHH; + mtop.ffparams.iparams[settleType].settle.doh = dOH; + mtop.ffparams.iparams[settleType].settle.dhh = dHH; // Set up the masses. t_mdatoms mdatoms; @@ -246,8 +240,8 @@ TEST_P(SettleTest, SatisfiesConstraints) mdatoms.homenr = numSettles * atomsPerSettle; // Finally make the settle data structures - gmx_settledata_t settled = settle_init(mtop); - settle_set_constraints(settled, &mtop->moltype[0].ilist[F_SETTLE], &mdatoms); + gmx_settledata_t settled = settle_init(&mtop); + settle_set_constraints(settled, &mtop.moltype[0].ilist[F_SETTLE], &mdatoms); // Copy the original positions from the array of doubles to a vector of reals std::vector startingPositions(std::begin(g_positions), std::end(g_positions)); diff --git a/src/gromacs/mdlib/vsite.cpp b/src/gromacs/mdlib/vsite.cpp index 45903d0172..706b07ed80 100644 --- a/src/gromacs/mdlib/vsite.cpp +++ b/src/gromacs/mdlib/vsite.cpp @@ -769,9 +769,8 @@ void constructVsitesGlobal(const gmx_mtop_t &mtop, { GMX_ASSERT(x.size() >= static_cast(mtop.natoms), "x should contain the whole system"); - for (int mb = 0; mb < mtop.nmolblock; mb++) + for (const gmx_molblock_t &molb : mtop.molblock) { - const gmx_molblock_t &molb = mtop.molblock[mb]; const gmx_moltype_t &molt = mtop.moltype[molb.type]; if (vsiteIlistNrCount(molt.ilist) > 0) { @@ -1844,30 +1843,25 @@ static std::vector atom2cg(const t_block &chargeGroups) int count_intercg_vsites(const gmx_mtop_t *mtop) { - gmx_molblock_t *molb; - gmx_moltype_t *molt; - int n_intercg_vsite; - - n_intercg_vsite = 0; - for (int mb = 0; mb < mtop->nmolblock; mb++) + int n_intercg_vsite = 0; + for (const gmx_molblock_t &molb : mtop->molblock) { - molb = &mtop->molblock[mb]; - molt = &mtop->moltype[molb->type]; + const gmx_moltype_t &molt = mtop->moltype[molb.type]; - std::vector a2cg = atom2cg(molt->cgs); + std::vector a2cg = atom2cg(molt.cgs); for (int ftype = c_ftypeVsiteStart; ftype < c_ftypeVsiteEnd; ftype++) { int nral = NRAL(ftype); - t_ilist *il = &molt->ilist[ftype]; - const t_iatom *ia = il->iatoms; - for (int i = 0; i < il->nr; i += 1 + nral) + const t_ilist &il = molt.ilist[ftype]; + const t_iatom *ia = il.iatoms; + for (int i = 0; i < il.nr; i += 1 + nral) { int cg = a2cg[ia[1+i]]; for (int a = 1; a < nral; a++) { if (a2cg[ia[1+a]] != cg) { - n_intercg_vsite += molb->nmol; + n_intercg_vsite += molb.nmol; break; } } @@ -2060,9 +2054,9 @@ gmx_vsite_t *initVsite(const gmx_mtop_t &mtop, vsite->n_intercg_vsite > 0 && DOMAINDECOMP(cr)) { - vsite->nvsite_pbc_molt = mtop.nmoltype; + vsite->nvsite_pbc_molt = mtop.moltype.size(); snew(vsite->vsite_pbc_molt, vsite->nvsite_pbc_molt); - for (int mt = 0; mt < mtop.nmoltype; mt++) + for (size_t mt = 0; mt < mtop.moltype.size(); mt++) { const gmx_moltype_t &molt = mtop.moltype[mt]; vsite->vsite_pbc_molt[mt] = get_vsite_pbc(mtop.ffparams.iparams, diff --git a/src/gromacs/selection/indexutil.cpp b/src/gromacs/selection/indexutil.cpp index b2c1b5b353..ed0e5fe36a 100644 --- a/src/gromacs/selection/indexutil.cpp +++ b/src/gromacs/selection/indexutil.cpp @@ -953,8 +953,8 @@ gmx_ana_index_make_block(t_blocka *t, const gmx_mtop_t *top, gmx_ana_index_t *g, } case INDEX_MOL: { - int molb = 0; - while (molb + 1 < top->nmolblock && id >= top->molblock[molb].moleculeIndexStart) + size_t molb = 0; + while (molb + 1 < top->molblock.size() && id >= top->molblock[molb].moleculeIndexStart) { ++molb; } diff --git a/src/gromacs/selection/selectionoptionbehavior.cpp b/src/gromacs/selection/selectionoptionbehavior.cpp index f188638646..5cc4a4abcd 100644 --- a/src/gromacs/selection/selectionoptionbehavior.cpp +++ b/src/gromacs/selection/selectionoptionbehavior.cpp @@ -1,7 +1,7 @@ /* * This file is part of the GROMACS molecular simulation package. * - * Copyright (c) 2015,2016,2017, by the GROMACS development team, led by + * Copyright (c) 2015,2016,2017,2018, 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. @@ -160,9 +160,8 @@ class SelectionOptionBehavior::Impl // when the user has not provided the topology. GMX_RELEASE_ASSERT(top != nullptr, "Masses are required, but no topology is loaded"); - for (int i = 0; i < top->nmoltype; ++i) + for (gmx_moltype_t &moltype : top->moltype) { - gmx_moltype_t &moltype = top->moltype[i]; if (!moltype.atoms.haveMass) { atomsSetMassesBasedOnNames(&moltype.atoms, TRUE); diff --git a/src/gromacs/selection/tests/toputils.cpp b/src/gromacs/selection/tests/toputils.cpp index 234f466309..c33343e66a 100644 --- a/src/gromacs/selection/tests/toputils.cpp +++ b/src/gromacs/selection/tests/toputils.cpp @@ -47,6 +47,7 @@ #include +#include "gromacs/compat/make_unique.h" #include "gromacs/fileio/confio.h" #include "gromacs/fileio/trxio.h" #include "gromacs/math/vec.h" @@ -67,18 +68,12 @@ namespace test { TopologyManager::TopologyManager() - : mtop_(nullptr), frame_(nullptr) + : mtop_(), frame_(nullptr) { } TopologyManager::~TopologyManager() { - if (mtop_ != nullptr) - { - done_mtop(mtop_); - sfree(mtop_); - } - if (frame_ != nullptr) { sfree(frame_->x); @@ -134,14 +129,15 @@ void TopologyManager::loadTopology(const char *filename) matrix box; GMX_RELEASE_ASSERT(mtop_ == nullptr, "Topology initialized more than once"); - snew(mtop_, 1); + mtop_ = gmx::compat::make_unique(); readConfAndTopology( gmx::test::TestFileManager::getInputFilePath(filename).c_str(), - &fullTopology, mtop_, &ePBC, frame_ != nullptr ? &xtop : nullptr, + &fullTopology, mtop_.get(), &ePBC, frame_ != nullptr ? &xtop : nullptr, nullptr, box); if (frame_ != nullptr) { + GMX_ASSERT(xtop != nullptr, "Keep the static analyzer happy"); frame_->natoms = mtop_->natoms; frame_->bX = TRUE; snew(frame_->x, frame_->natoms); @@ -156,18 +152,16 @@ void TopologyManager::loadTopology(const char *filename) void TopologyManager::initAtoms(int count) { GMX_RELEASE_ASSERT(mtop_ == nullptr, "Topology initialized more than once"); - snew(mtop_, 1); - mtop_->nmoltype = 1; - snew(mtop_->moltype, 1); + mtop_ = gmx::compat::make_unique(); + mtop_->moltype.resize(1); init_t_atoms(&mtop_->moltype[0].atoms, count, FALSE); - mtop_->nmolblock = 1; - snew(mtop_->molblock, 1); + mtop_->molblock.resize(1); mtop_->molblock[0].type = 0; mtop_->molblock[0].nmol = 1; mtop_->molblock[0].natoms_mol = count; mtop_->natoms = count; mtop_->maxres_renum = 0; - gmx_mtop_finalize(mtop_); + gmx_mtop_finalize(mtop_.get()); GMX_RELEASE_ASSERT(mtop_->maxres_renum == 0, "maxres_renum in mtop can be modified by an env.var., that is not supported in this test"); t_atoms &atoms = this->atoms(); for (int i = 0; i < count; ++i) @@ -232,7 +226,7 @@ void TopologyManager::initUniformResidues(int residueSize) void TopologyManager::initUniformMolecules(int moleculeSize) { GMX_RELEASE_ASSERT(mtop_ != nullptr, "Topology not initialized"); - GMX_RELEASE_ASSERT(mtop_->nmolblock == 1, "initUniformMolecules only implemented for a single molblock"); + GMX_RELEASE_ASSERT(mtop_->molblock.size() == 1, "initUniformMolecules only implemented for a single molblock"); gmx_molblock_t &molblock = mtop_->molblock[0]; t_atoms &atoms = mtop_->moltype[molblock.type].atoms; GMX_RELEASE_ASSERT(atoms.nr % moleculeSize == 0, diff --git a/src/gromacs/selection/tests/toputils.h b/src/gromacs/selection/tests/toputils.h index ac2faf4159..33707c53dd 100644 --- a/src/gromacs/selection/tests/toputils.h +++ b/src/gromacs/selection/tests/toputils.h @@ -1,7 +1,7 @@ /* * This file is part of the GROMACS molecular simulation package. * - * Copyright (c) 2013,2014,2015,2016,2017, by the GROMACS development team, led by + * Copyright (c) 2013,2014,2015,2016,2017,2018, 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. @@ -42,6 +42,7 @@ #ifndef GMX_SELECTION_TESTS_TOPUTILS_H #define GMX_SELECTION_TESTS_TOPUTILS_H +#include #include struct gmx_mtop_t; @@ -75,14 +76,14 @@ class TopologyManager void initFrameIndices(const ArrayRef &index); - gmx_mtop_t *topology() { return mtop_; } + gmx_mtop_t *topology() { return mtop_.get(); } t_atoms &atoms(); t_trxframe *frame() { return frame_; } private: - gmx_mtop_t *mtop_; - t_trxframe *frame_; - std::vector atomtypes_; + std::unique_ptr mtop_; + t_trxframe *frame_; + std::vector atomtypes_; }; } // namespace test diff --git a/src/gromacs/tools/convert_tpr.cpp b/src/gromacs/tools/convert_tpr.cpp index a2b32eb56a..7e574d9bc8 100644 --- a/src/gromacs/tools/convert_tpr.cpp +++ b/src/gromacs/tools/convert_tpr.cpp @@ -3,7 +3,7 @@ * * Copyright (c) 1991-2000, University of Groningen, The Netherlands. * Copyright (c) 2001-2004, The GROMACS development team. - * Copyright (c) 2013,2014,2015,2016,2017, by the GROMACS development team, led by + * Copyright (c) 2013,2014,2015,2016,2017,2018, 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. @@ -288,8 +288,7 @@ static void reduce_topology_x(int gnx, int index[], top.atoms.nr = gnx; - mtop->nmoltype = 1; - snew(mtop->moltype, mtop->nmoltype); + mtop->moltype.resize(1); mtop->moltype[0].name = mtop->name; mtop->moltype[0].atoms = top.atoms; for (i = 0; i < F_NRE; i++) @@ -300,8 +299,7 @@ static void reduce_topology_x(int gnx, int index[], mtop->moltype[0].cgs = top.cgs; mtop->moltype[0].excls = top.excls; - mtop->nmolblock = 1; - snew(mtop->molblock, mtop->nmolblock); + mtop->molblock.resize(1); mtop->molblock[0].type = 0; mtop->molblock[0].nmol = 1; mtop->molblock[0].natoms_mol = top.atoms.nr; @@ -313,14 +311,12 @@ static void reduce_topology_x(int gnx, int index[], static void zeroq(int index[], gmx_mtop_t *mtop) { - int mt, i; - - for (mt = 0; mt < mtop->nmoltype; mt++) + for (gmx_moltype_t &moltype : mtop->moltype) { - for (i = 0; (i < mtop->moltype[mt].atoms.nr); i++) + for (int i = 0; i < moltype.atoms.nr; i++) { - mtop->moltype[mt].atoms.atom[index[i]].q = 0; - mtop->moltype[mt].atoms.atom[index[i]].qB = 0; + moltype.atoms.atom[index[i]].q = 0; + moltype.atoms.atom[index[i]].qB = 0; } } } diff --git a/src/gromacs/topology/block.cpp b/src/gromacs/topology/block.cpp index 73c7e47cbb..937d065abb 100644 --- a/src/gromacs/topology/block.cpp +++ b/src/gromacs/topology/block.cpp @@ -3,7 +3,7 @@ * * Copyright (c) 1991-2000, University of Groningen, The Netherlands. * Copyright (c) 2001-2004, The GROMACS development team. - * Copyright (c) 2013,2014,2015,2017, by the GROMACS development team, led by + * Copyright (c) 2013,2014,2015,2017,2018, 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. @@ -76,8 +76,9 @@ t_blocka *new_blocka(void) void done_block(t_block *block) { - block->nr = 0; + block->nr = 0; sfree(block->index); + block->index = nullptr; block->nalloc_index = 0; } diff --git a/src/gromacs/topology/mtop_lookup.h b/src/gromacs/topology/mtop_lookup.h index 70853e0e70..fd4208c9f5 100644 --- a/src/gromacs/topology/mtop_lookup.h +++ b/src/gromacs/topology/mtop_lookup.h @@ -74,11 +74,11 @@ mtopGetMolblockIndex(const gmx_mtop_t *mtop, { 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"); + GMX_ASSERT(*moleculeBlock >= 0 && *moleculeBlock < static_cast(mtop->molblock.size()), "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 molBlock1 = mtop->molblock.size(); int globalAtomStart; while (TRUE) diff --git a/src/gromacs/topology/mtop_util.cpp b/src/gromacs/topology/mtop_util.cpp index 96dfa454b2..4ab201b766 100644 --- a/src/gromacs/topology/mtop_util.cpp +++ b/src/gromacs/topology/mtop_util.cpp @@ -56,21 +56,18 @@ static int gmx_mtop_maxresnr(const gmx_mtop_t *mtop, int maxres_renum) { - int maxresnr, mt, r; - const t_atoms *atoms; + int maxresnr = 0; - maxresnr = 0; - - for (mt = 0; mt < mtop->nmoltype; mt++) + for (const gmx_moltype_t &moltype : mtop->moltype) { - atoms = &mtop->moltype[mt].atoms; - if (atoms->nres > maxres_renum) + const t_atoms &atoms = moltype.atoms; + if (atoms.nres > maxres_renum) { - for (r = 0; r < atoms->nres; r++) + for (int r = 0; r < atoms.nres; r++) { - if (atoms->resinfo[r].nr > maxresnr) + if (atoms.resinfo[r].nr > maxresnr) { - maxresnr = atoms->resinfo[r].nr; + maxresnr = atoms.resinfo[r].nr; } } } @@ -85,22 +82,21 @@ static void finalizeMolblocks(gmx_mtop_t *mtop) int residueIndex = 0; int residueNumberStart = mtop->maxresnr + 1; int moleculeIndexStart = 0; - for (int mb = 0; mb < mtop->nmolblock; mb++) + for (gmx_molblock_t &molb : mtop->molblock) { - 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; + const 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; + residueNumberStart += molb.nmol*numResPerMol; } - molb->moleculeIndexStart = moleculeIndexStart; - moleculeIndexStart += molb->nmol; + molb.moleculeIndexStart = moleculeIndexStart; + moleculeIndexStart += molb.nmol; } } @@ -108,7 +104,7 @@ void gmx_mtop_finalize(gmx_mtop_t *mtop) { char *env; - if (mtop->nmolblock == 1 && mtop->molblock[0].nmol == 1) + if (mtop->molblock.size() == 1 && mtop->molblock[0].nmol == 1) { /* We have a single molecule only, no renumbering needed. * This case also covers an mtop converted from pdb/gro/... input, @@ -142,43 +138,35 @@ void gmx_mtop_finalize(gmx_mtop_t *mtop) void gmx_mtop_count_atomtypes(const gmx_mtop_t *mtop, int state, int typecount[]) { - int i, mb, nmol, tpi; - t_atoms *atoms; - - for (i = 0; i < mtop->ffparams.atnr; ++i) + for (int i = 0; i < mtop->ffparams.atnr; ++i) { typecount[i] = 0; } - for (mb = 0; mb < mtop->nmolblock; ++mb) + for (const gmx_molblock_t &molb : mtop->molblock) { - nmol = mtop->molblock[mb].nmol; - atoms = &mtop->moltype[mtop->molblock[mb].type].atoms; - for (i = 0; i < atoms->nr; ++i) + const t_atoms &atoms = mtop->moltype[molb.type].atoms; + for (int i = 0; i < atoms.nr; ++i) { + int tpi; if (state == 0) { - tpi = atoms->atom[i].type; + tpi = atoms.atom[i].type; } else { - tpi = atoms->atom[i].typeB; + tpi = atoms.atom[i].typeB; } - typecount[tpi] += nmol; + typecount[tpi] += molb.nmol; } } } int ncg_mtop(const gmx_mtop_t *mtop) { - int ncg; - int mb; - - ncg = 0; - for (mb = 0; mb < mtop->nmolblock; mb++) + int ncg = 0; + for (const gmx_molblock_t &molb : mtop->molblock) { - ncg += - mtop->molblock[mb].nmol* - mtop->moltype[mtop->molblock[mb].type].cgs.nr; + ncg += molb.nmol*mtop->moltype[molb.type].cgs.nr; } return ncg; @@ -187,9 +175,9 @@ int ncg_mtop(const gmx_mtop_t *mtop) int gmx_mtop_num_molecules(const gmx_mtop_t &mtop) { int numMolecules = 0; - for (int mb = 0; mb < mtop.nmolblock; mb++) + for (const gmx_molblock_t &molb : mtop.molblock) { - numMolecules += mtop.molblock[mb].nmol; + numMolecules += molb.nmol; } return numMolecules; } @@ -197,31 +185,25 @@ int gmx_mtop_num_molecules(const gmx_mtop_t &mtop) int gmx_mtop_nres(const gmx_mtop_t *mtop) { int nres = 0; - for (int mb = 0; mb < mtop->nmolblock; ++mb) + for (const gmx_molblock_t &molb : mtop->molblock) { - nres += - mtop->molblock[mb].nmol* - mtop->moltype[mtop->molblock[mb].type].atoms.nres; + nres += molb.nmol*mtop->moltype[molb.type].atoms.nres; } return nres; } void gmx_mtop_remove_chargegroups(gmx_mtop_t *mtop) { - int mt; - t_block *cgs; - int i; - - for (mt = 0; mt < mtop->nmoltype; mt++) + for (gmx_moltype_t &molt : mtop->moltype) { - cgs = &mtop->moltype[mt].cgs; - if (cgs->nr < mtop->moltype[mt].atoms.nr) + t_block &cgs = molt.cgs; + if (cgs.nr < molt.atoms.nr) { - cgs->nr = mtop->moltype[mt].atoms.nr; - srenew(cgs->index, cgs->nr+1); - for (i = 0; i < cgs->nr+1; i++) + cgs.nr = molt.atoms.nr; + srenew(cgs.index, cgs.nr + 1); + for (int i = 0; i < cgs.nr + 1; i++) { - cgs->index[i] = i; + cgs.index[i] = i; } } } @@ -230,8 +212,8 @@ void gmx_mtop_remove_chargegroups(gmx_mtop_t *mtop) typedef struct gmx_mtop_atomloop_all { const gmx_mtop_t *mtop; - int mblock; - t_atoms *atoms; + size_t mblock; + const t_atoms *atoms; int mol; int maxresnr; int at_local; @@ -285,7 +267,7 @@ gmx_bool gmx_mtop_atomloop_all_next(gmx_mtop_atomloop_all_t aloop, if (aloop->mol >= aloop->mtop->molblock[aloop->mblock].nmol) { aloop->mblock++; - if (aloop->mblock >= aloop->mtop->nmolblock) + if (aloop->mblock >= aloop->mtop->molblock.size()) { gmx_mtop_atomloop_all_destroy(aloop); return FALSE; @@ -316,8 +298,9 @@ void gmx_mtop_atomloop_all_names(gmx_mtop_atomloop_all_t aloop, *resname = *(aloop->atoms->resinfo[resind_mol].name); } -void gmx_mtop_atomloop_all_moltype(gmx_mtop_atomloop_all_t aloop, - gmx_moltype_t **moltype, int *at_mol) +void gmx_mtop_atomloop_all_moltype(gmx_mtop_atomloop_all_t aloop, + const gmx_moltype_t **moltype, + int *at_mol) { *moltype = &aloop->mtop->moltype[aloop->mtop->molblock[aloop->mblock].type]; *at_mol = aloop->at_local; @@ -326,8 +309,8 @@ void gmx_mtop_atomloop_all_moltype(gmx_mtop_atomloop_all_t aloop, typedef struct gmx_mtop_atomloop_block { const gmx_mtop_t *mtop; - int mblock; - t_atoms *atoms; + size_t mblock; + const t_atoms *atoms; int at_local; } t_gmx_mtop_atomloop_block; @@ -364,7 +347,7 @@ gmx_bool gmx_mtop_atomloop_block_next(gmx_mtop_atomloop_block_t aloop, if (aloop->at_local >= aloop->atoms->nr) { aloop->mblock++; - if (aloop->mblock >= aloop->mtop->nmolblock) + if (aloop->mblock >= aloop->mtop->molblock.size()) { gmx_mtop_atomloop_block_destroy(aloop); return FALSE; @@ -404,7 +387,7 @@ static void gmx_mtop_ilistloop_destroy(gmx_mtop_ilistloop_t iloop) } gmx_bool gmx_mtop_ilistloop_next(gmx_mtop_ilistloop_t iloop, - t_ilist **ilist_mol, int *nmol) + const t_ilist **ilist_mol, int *nmol) { if (iloop == nullptr) { @@ -412,9 +395,9 @@ gmx_bool gmx_mtop_ilistloop_next(gmx_mtop_ilistloop_t iloop, } iloop->mblock++; - if (iloop->mblock >= iloop->mtop->nmolblock) + if (iloop->mblock >= static_cast(iloop->mtop->molblock.size())) { - if (iloop->mblock == iloop->mtop->nmolblock && + if (iloop->mblock == static_cast(iloop->mtop->molblock.size()) && iloop->mtop->bIntermolecularInteractions) { *ilist_mol = iloop->mtop->intermolecular_ilist; @@ -436,7 +419,7 @@ gmx_bool gmx_mtop_ilistloop_next(gmx_mtop_ilistloop_t iloop, typedef struct gmx_mtop_ilistloop_all { const gmx_mtop_t *mtop; - int mblock; + size_t mblock; int mol; int a_offset; } t_gmx_mtop_ilist_all; @@ -462,7 +445,7 @@ static void gmx_mtop_ilistloop_all_destroy(gmx_mtop_ilistloop_all_t iloop) } gmx_bool gmx_mtop_ilistloop_all_next(gmx_mtop_ilistloop_all_t iloop, - t_ilist **ilist_mol, int *atnr_offset) + const t_ilist **ilist_mol, int *atnr_offset) { if (iloop == nullptr) @@ -481,14 +464,14 @@ gmx_bool gmx_mtop_ilistloop_all_next(gmx_mtop_ilistloop_all_t iloop, * iloop->mblock == iloop->mtop->nmolblock, thus we should separately * check for this value in this conditional. */ - if (iloop->mblock == iloop->mtop->nmolblock || + if (iloop->mblock == iloop->mtop->molblock.size() || iloop->mol >= iloop->mtop->molblock[iloop->mblock].nmol) { iloop->mblock++; iloop->mol = 0; - if (iloop->mblock >= iloop->mtop->nmolblock) + if (iloop->mblock >= iloop->mtop->molblock.size()) { - if (iloop->mblock == iloop->mtop->nmolblock && + if (iloop->mblock == iloop->mtop->molblock.size() && iloop->mtop->bIntermolecularInteractions) { *ilist_mol = iloop->mtop->intermolecular_ilist; @@ -512,7 +495,7 @@ gmx_bool gmx_mtop_ilistloop_all_next(gmx_mtop_ilistloop_all_t iloop, int gmx_mtop_ftype_count(const gmx_mtop_t *mtop, int ftype) { gmx_mtop_ilistloop_t iloop; - t_ilist *il; + const t_ilist *il; int n, nmol; n = 0; @@ -533,26 +516,22 @@ int gmx_mtop_ftype_count(const gmx_mtop_t *mtop, int ftype) t_block gmx_mtop_global_cgs(const gmx_mtop_t *mtop) { - t_block cgs_gl, *cgs_mol; - int mb, mol, cg; - gmx_molblock_t *molb; - + t_block cgs_gl; /* In most cases this is too much, but we realloc at the end */ snew(cgs_gl.index, mtop->natoms+1); cgs_gl.nr = 0; cgs_gl.index[0] = 0; - for (mb = 0; mb < mtop->nmolblock; mb++) + for (const gmx_molblock_t &molb : mtop->molblock) { - molb = &mtop->molblock[mb]; - cgs_mol = &mtop->moltype[molb->type].cgs; - for (mol = 0; mol < molb->nmol; mol++) + const t_block &cgs_mol = mtop->moltype[molb.type].cgs; + for (int mol = 0; mol < molb.nmol; mol++) { - for (cg = 0; cg < cgs_mol->nr; cg++) + for (int cg = 0; cg < cgs_mol.nr; cg++) { - cgs_gl.index[cgs_gl.nr+1] = + cgs_gl.index[cgs_gl.nr + 1] = cgs_gl.index[cgs_gl.nr] + - cgs_mol->index[cg+1] - cgs_mol->index[cg]; + cgs_mol.index[cg + 1] - cgs_mol.index[cg]; cgs_gl.nr++; } } @@ -563,7 +542,7 @@ t_block gmx_mtop_global_cgs(const gmx_mtop_t *mtop) return cgs_gl; } -static void atomcat(t_atoms *dest, t_atoms *src, int copies, +static void atomcat(t_atoms *dest, const t_atoms *src, int copies, int maxres_renum, int *maxresnr) { int i, j, l, size; @@ -670,16 +649,13 @@ static void atomcat(t_atoms *dest, t_atoms *src, int copies, t_atoms gmx_mtop_global_atoms(const gmx_mtop_t *mtop) { t_atoms atoms; - int maxresnr, mb; - gmx_molblock_t *molb; init_t_atoms(&atoms, 0, FALSE); - maxresnr = mtop->maxresnr; - for (mb = 0; mb < mtop->nmolblock; mb++) + int maxresnr = mtop->maxresnr; + for (const gmx_molblock_t &molb : mtop->molblock) { - molb = &mtop->molblock[mb]; - atomcat(&atoms, &mtop->moltype[molb->type].atoms, molb->nmol, + atomcat(&atoms, &mtop->moltype[molb.type].atoms, molb.nmol, mtop->maxres_renum, &maxresnr); } @@ -690,7 +666,7 @@ 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 blockcat(t_block *dest, t_block *src, int copies) +static void blockcat(t_block *dest, const t_block *src, int copies) { int i, j, l, nra, size; @@ -707,13 +683,16 @@ static void blockcat(t_block *dest, t_block *src, int copies) { dest->index[l++] = nra + src->index[i]; } - nra += src->index[src->nr]; + if (src->nr > 0) + { + nra += src->index[src->nr]; + } } dest->nr += copies*src->nr; dest->index[dest->nr] = nra; } -static void blockacat(t_blocka *dest, t_blocka *src, int copies, +static void blockacat(t_blocka *dest, const t_blocka *src, int copies, int dnum, int snum) { int i, j, l, size; @@ -751,7 +730,7 @@ static void blockacat(t_blocka *dest, t_blocka *src, int copies, dest->index[dest->nr] = dest->nra; } -static void ilistcat(int ftype, t_ilist *dest, t_ilist *src, int copies, +static void ilistcat(int ftype, t_ilist *dest, const t_ilist *src, int copies, int dnum, int snum) { int nral, c, i, a; @@ -775,7 +754,7 @@ static void ilistcat(int ftype, t_ilist *dest, t_ilist *src, int copies, } } -static void set_posres_params(t_idef *idef, gmx_molblock_t *molb, +static void set_posres_params(t_idef *idef, const gmx_molblock_t *molb, int i0, int a_offset) { t_ilist *il; @@ -816,7 +795,7 @@ static void set_posres_params(t_idef *idef, gmx_molblock_t *molb, } } -static void set_fbposres_params(t_idef *idef, gmx_molblock_t *molb, +static void set_fbposres_params(t_idef *idef, const gmx_molblock_t *molb, int i0, int a_offset) { t_ilist *il; @@ -853,9 +832,7 @@ static void gen_local_top(const gmx_mtop_t *mtop, bool bMergeConstr, gmx_localtop_t *top) { - int mb, srcnr, destnr, ftype, natoms, mol, nposre_old, nfbposre_old; - gmx_molblock_t *molb; - gmx_moltype_t *molt; + int srcnr, destnr, ftype, natoms, nposre_old, nfbposre_old; const gmx_ffparams_t *ffp; t_idef *idef; real *qA, *qB; @@ -889,54 +866,53 @@ static void gen_local_top(const gmx_mtop_t *mtop, } natoms = 0; - for (mb = 0; mb < mtop->nmolblock; mb++) + for (const gmx_molblock_t &molb : mtop->molblock) { - molb = &mtop->molblock[mb]; - molt = &mtop->moltype[molb->type]; + const gmx_moltype_t &molt = mtop->moltype[molb.type]; - srcnr = molt->atoms.nr; + srcnr = molt.atoms.nr; destnr = natoms; - blockcat(&top->cgs, &molt->cgs, molb->nmol); + blockcat(&top->cgs, &molt.cgs, molb.nmol); - blockacat(&top->excls, &molt->excls, molb->nmol, destnr, srcnr); + blockacat(&top->excls, &molt.excls, molb.nmol, destnr, srcnr); nposre_old = idef->il[F_POSRES].nr; nfbposre_old = idef->il[F_FBPOSRES].nr; for (ftype = 0; ftype < F_NRE; ftype++) { if (bMergeConstr && - ftype == F_CONSTR && molt->ilist[F_CONSTRNC].nr > 0) + ftype == F_CONSTR && molt.ilist[F_CONSTRNC].nr > 0) { /* Merge all constrains into one ilist. * This simplifies the constraint code. */ - for (mol = 0; mol < molb->nmol; mol++) + for (int mol = 0; mol < molb.nmol; mol++) { - ilistcat(ftype, &idef->il[F_CONSTR], &molt->ilist[F_CONSTR], - 1, destnr+mol*srcnr, srcnr); - ilistcat(ftype, &idef->il[F_CONSTR], &molt->ilist[F_CONSTRNC], - 1, destnr+mol*srcnr, srcnr); + ilistcat(ftype, &idef->il[F_CONSTR], &molt.ilist[F_CONSTR], + 1, destnr + mol*srcnr, srcnr); + ilistcat(ftype, &idef->il[F_CONSTR], &molt.ilist[F_CONSTRNC], + 1, destnr + mol*srcnr, srcnr); } } else if (!(bMergeConstr && ftype == F_CONSTRNC)) { - ilistcat(ftype, &idef->il[ftype], &molt->ilist[ftype], - molb->nmol, destnr, srcnr); + ilistcat(ftype, &idef->il[ftype], &molt.ilist[ftype], + molb.nmol, destnr, srcnr); } } if (idef->il[F_POSRES].nr > nposre_old) { /* Executing this line line stops gmxdump -sys working * correctly. I'm not aware there's an elegant fix. */ - set_posres_params(idef, molb, nposre_old/2, natoms); + set_posres_params(idef, &molb, nposre_old/2, natoms); } if (idef->il[F_FBPOSRES].nr > nfbposre_old) { - set_fbposres_params(idef, molb, nfbposre_old/2, natoms); + set_fbposres_params(idef, &molb, nfbposre_old/2, natoms); } - natoms += molb->nmol*srcnr; + natoms += molb.nmol*srcnr; } if (mtop->bIntermolecularInteractions) @@ -993,9 +969,8 @@ static void fillMoleculeIndices(const gmx_mtop_t &mtop, int globalAtomIndex = 0; int globalMolIndex = 0; index[globalMolIndex] = globalAtomIndex; - for (int mb = 0; mb < mtop.nmolblock; mb++) + for (const gmx_molblock_t &molb : mtop.molblock) { - const gmx_molblock_t &molb = mtop.molblock[mb]; for (int mol = 0; mol < molb.nmol; mol++) { globalAtomIndex += molb.natoms_mol; @@ -1035,7 +1010,6 @@ static t_block gmx_mtop_molecules_t_block(const gmx_mtop_t &mtop) t_topology gmx_mtop_t_to_t_topology(gmx_mtop_t *mtop, bool freeMTop) { - int mt, mb; gmx_localtop_t ltop; t_topology top; @@ -1054,20 +1028,10 @@ t_topology gmx_mtop_t_to_t_topology(gmx_mtop_t *mtop, bool freeMTop) if (freeMTop) { - // Free pointers that have not been copied to top. - for (mt = 0; mt < mtop->nmoltype; mt++) - { - done_moltype(&mtop->moltype[mt]); - } - sfree(mtop->moltype); - - for (mb = 0; mb < mtop->nmolblock; mb++) - { - done_molblock(&mtop->molblock[mb]); - } - sfree(mtop->molblock); - - done_gmx_groups_t(&mtop->groups); + // Clear pointers and counts, such that the pointers copied to top + // keep pointing to valid data after destroying mtop. + mtop->symtab.symbuf = nullptr; + mtop->symtab.nr = 0; } return top; @@ -1100,23 +1064,18 @@ void convertAtomsToMtop(t_symtab *symtab, mtop->name = name; - mtop->nmoltype = 1; - // This snew clears all entries, we should replace it by an initializer - snew(mtop->moltype, mtop->nmoltype); - mtop->moltype[0].atoms = *atoms; - init_block(&mtop->moltype[0].cgs); - init_blocka(&mtop->moltype[0].excls); + mtop->moltype.clear(); + mtop->moltype.resize(1); + mtop->moltype.back().atoms = *atoms; - mtop->nmolblock = 1; - // This snew clears all entries, we should replace it by an initializer - snew(mtop->molblock, mtop->nmolblock); + mtop->molblock.resize(1); mtop->molblock[0].type = 0; mtop->molblock[0].nmol = 1; - mtop->molblock[0].natoms_mol = atoms->nr; + mtop->molblock[0].natoms_mol = mtop->moltype[mtop->molblock[0].type].atoms.nr; mtop->bIntermolecularInteractions = FALSE; - mtop->natoms = atoms->nr; + mtop->natoms = mtop->molblock[0].natoms_mol; mtop->haveMoleculeIndices = false; diff --git a/src/gromacs/topology/mtop_util.h b/src/gromacs/topology/mtop_util.h index 918bb3dcfd..61f4b1a904 100644 --- a/src/gromacs/topology/mtop_util.h +++ b/src/gromacs/topology/mtop_util.h @@ -123,7 +123,7 @@ gmx_mtop_atomloop_all_names(gmx_mtop_atomloop_all_t aloop, */ void gmx_mtop_atomloop_all_moltype(gmx_mtop_atomloop_all_t aloop, - gmx_moltype_t **moltype, int *at_mol); + const gmx_moltype_t **moltype, int *at_mol); /* Abstract type for atom loop over atoms in all molecule blocks */ @@ -168,7 +168,7 @@ gmx_mtop_ilistloop_init(const gmx_mtop_t *mtop); */ gmx_bool gmx_mtop_ilistloop_next(gmx_mtop_ilistloop_t iloop, - t_ilist **ilist_mol, int *nmol); + const t_ilist **ilist_mol, int *nmol); /* Abstract type for ilist loop over all ilists of all molecules */ @@ -190,7 +190,7 @@ gmx_mtop_ilistloop_all_init(const gmx_mtop_t *mtop); */ gmx_bool gmx_mtop_ilistloop_all_next(gmx_mtop_ilistloop_all_t iloop, - t_ilist **ilist_mol, int *atnr_offset); + const t_ilist **ilist_mol, int *atnr_offset); /* Returns the total number of interactions in the system of type ftype */ @@ -228,6 +228,8 @@ gmx::BlockRanges gmx_mtop_molecules(const gmx_mtop_t &mtop); /* Converts a gmx_mtop_t struct to t_topology. * + * If the lifetime of the returned topology should be longer than that + * of mtop, your need to pass freeMtop==true. * If freeMTop == true, memory related to mtop will be freed so that done_top() * on the result value will free all memory. * If freeMTop == false, mtop and the return value will share some of their diff --git a/src/gromacs/topology/topology.cpp b/src/gromacs/topology/topology.cpp index a9b6b312eb..198e1a5af2 100644 --- a/src/gromacs/topology/topology.cpp +++ b/src/gromacs/topology/topology.cpp @@ -63,6 +63,7 @@ static void init_groups(gmx_groups_t *groups) groups->grpname = nullptr; for (int g = 0; g < egcNR; g++) { + groups->grps[g].nr = 0; groups->grps[g].nm_ind = nullptr; groups->ngrpnr[g] = 0; groups->grpnr[g] = nullptr; @@ -72,13 +73,24 @@ static void init_groups(gmx_groups_t *groups) void init_mtop(gmx_mtop_t *mtop) { - mtop->name = nullptr; - mtop->nmoltype = 0; - mtop->moltype = nullptr; - mtop->nmolblock = 0; - mtop->molblock = nullptr; + mtop->name = nullptr; + + // TODO: Move to ffparams when that is converted to C++ + mtop->ffparams.functype = nullptr; + mtop->ffparams.iparams = nullptr; + mtop->ffparams.cmap_grid.ngrid = 0; + mtop->ffparams.cmap_grid.grid_spacing = 0; + mtop->ffparams.cmap_grid.cmapdata = nullptr; + + mtop->moltype.clear(); + mtop->molblock.clear(); + mtop->bIntermolecularInteractions = FALSE; + mtop->intermolecular_ilist = nullptr; + + mtop->natoms = 0; mtop->maxres_renum = 0; mtop->maxresnr = -1; + init_atomtypes(&mtop->atomtypes); init_groups(&mtop->groups); open_symtab(&mtop->symtab); } @@ -95,30 +107,63 @@ void init_top(t_topology *top) } -void done_moltype(gmx_moltype_t *molt) +gmx_moltype_t::gmx_moltype_t() : + name(nullptr), + cgs(), + excls() +{ + init_t_atoms(&atoms, 0, FALSE); + + for (int ftype = 0; ftype < F_NRE; ftype++) + { + ilist[ftype].nr = 0; + ilist[ftype].nr_nonperturbed = 0; + ilist[ftype].iatoms = nullptr; + ilist[ftype].nalloc = 0; + } +} + +gmx_moltype_t::~gmx_moltype_t() { - done_atom(&molt->atoms); - done_block(&molt->cgs); - done_blocka(&molt->excls); + done_atom(&atoms); + done_block(&cgs); + done_blocka(&excls); for (int f = 0; f < F_NRE; f++) { - sfree(molt->ilist[f].iatoms); - molt->ilist[f].nalloc = 0; + sfree(ilist[f].iatoms); + ilist[f].nalloc = 0; } } -void done_molblock(gmx_molblock_t *molb) +gmx_molblock_t::gmx_molblock_t() { - if (molb->nposres_xA > 0) + type = -1; + nmol = 0; + nposres_xA = 0; + posres_xA = nullptr; + nposres_xB = 0; + posres_xB = nullptr; + + natoms_mol = 0; + globalAtomStart = -1; + globalAtomEnd = -1; + globalResidueStart = -1; + residueNumberStart = -1; + moleculeIndexStart = -1; +} + +gmx_molblock_t::~gmx_molblock_t() +{ + if (nposres_xA > 0) { - molb->nposres_xA = 0; - sfree(molb->posres_xA); + nposres_xA = 0; + sfree(posres_xA); } - if (molb->nposres_xB > 0) + if (nposres_xB > 0) { - molb->nposres_xB = 0; - sfree(molb->posres_xB); + nposres_xB = 0; + sfree(posres_xB); } } @@ -143,30 +188,27 @@ void done_gmx_groups_t(gmx_groups_t *g) sfree(g->grpname); } -void done_mtop(gmx_mtop_t *mtop) +gmx_mtop_t::gmx_mtop_t() { - done_symtab(&mtop->symtab); + init_mtop(this); +} - sfree(mtop->ffparams.functype); - sfree(mtop->ffparams.iparams); - for (int i = 0; i < mtop->ffparams.cmap_grid.ngrid; i++) - { - sfree(mtop->ffparams.cmap_grid.cmapdata[i].cmap); - } - sfree(mtop->ffparams.cmap_grid.cmapdata); +gmx_mtop_t::~gmx_mtop_t() +{ + done_symtab(&symtab); - for (int i = 0; i < mtop->nmoltype; i++) + sfree(ffparams.functype); + sfree(ffparams.iparams); + for (int i = 0; i < ffparams.cmap_grid.ngrid; i++) { - done_moltype(&mtop->moltype[i]); + sfree(ffparams.cmap_grid.cmapdata[i].cmap); } - sfree(mtop->moltype); - for (int i = 0; i < mtop->nmolblock; i++) - { - done_molblock(&mtop->molblock[i]); - } - sfree(mtop->molblock); - done_atomtypes(&mtop->atomtypes); - done_gmx_groups_t(&mtop->groups); + sfree(ffparams.cmap_grid.cmapdata); + + moltype.clear(); + molblock.clear(); + done_atomtypes(&atomtypes); + done_gmx_groups_t(&groups); } void done_top(t_topology *top) @@ -210,7 +252,8 @@ void done_top_mtop(t_topology *top, gmx_mtop_t *mtop) done_symtab(&top->symtab); open_symtab(&mtop->symtab); } - done_mtop(mtop); + + // Note that the rest of mtop will be freed by the destructor } } @@ -220,7 +263,7 @@ bool gmx_mtop_has_masses(const gmx_mtop_t *mtop) { return false; } - return mtop->nmoltype == 0 || mtop->moltype[0].atoms.haveMass; + return mtop->moltype.empty() || mtop->moltype[0].atoms.haveMass; } bool gmx_mtop_has_charges(const gmx_mtop_t *mtop) @@ -229,7 +272,7 @@ bool gmx_mtop_has_charges(const gmx_mtop_t *mtop) { return false; } - return mtop->nmoltype == 0 || mtop->moltype[0].atoms.haveCharge; + return mtop->moltype.empty() || mtop->moltype[0].atoms.haveCharge; } bool gmx_mtop_has_atomtypes(const gmx_mtop_t *mtop) @@ -238,7 +281,7 @@ bool gmx_mtop_has_atomtypes(const gmx_mtop_t *mtop) { return false; } - return mtop->nmoltype == 0 || mtop->moltype[0].atoms.haveType; + return mtop->moltype.empty() || mtop->moltype[0].atoms.haveType; } bool gmx_mtop_has_pdbinfo(const gmx_mtop_t *mtop) @@ -247,7 +290,7 @@ bool gmx_mtop_has_pdbinfo(const gmx_mtop_t *mtop) { return false; } - return mtop->nmoltype == 0 || mtop->moltype[0].atoms.havePdbInfo; + return mtop->moltype.empty() || mtop->moltype[0].atoms.havePdbInfo; } static void pr_grps(FILE *fp, const char *title, const t_grps grps[], char **grpname[]) @@ -341,7 +384,7 @@ static void pr_moltype(FILE *fp, int indent, const char *title, static void pr_molblock(FILE *fp, int indent, const char *title, const gmx_molblock_t *molb, int n, - const gmx_moltype_t *molt) + const std::vector &molt) { indent = pr_title_n(fp, indent, title, n); pr_indent(fp, indent); @@ -364,16 +407,14 @@ static void pr_molblock(FILE *fp, int indent, const char *title, void pr_mtop(FILE *fp, int indent, const char *title, const gmx_mtop_t *mtop, gmx_bool bShowNumbers, gmx_bool bShowParameters) { - int mt, mb, j; - if (available(fp, mtop, indent, title)) { indent = pr_title(fp, indent, title); pr_indent(fp, indent); fprintf(fp, "name=\"%s\"\n", *(mtop->name)); pr_int(fp, indent, "#atoms", mtop->natoms); - pr_int(fp, indent, "#molblock", mtop->nmolblock); - for (mb = 0; mb < mtop->nmolblock; mb++) + pr_int(fp, indent, "#molblock", mtop->molblock.size()); + for (size_t mb = 0; mb < mtop->molblock.size(); mb++) { pr_molblock(fp, indent, "molblock", &mtop->molblock[mb], mb, mtop->moltype); } @@ -381,7 +422,7 @@ void pr_mtop(FILE *fp, int indent, const char *title, const gmx_mtop_t *mtop, gmx::boolToString(mtop->bIntermolecularInteractions)); if (mtop->bIntermolecularInteractions) { - for (j = 0; (j < F_NRE); j++) + for (int j = 0; j < F_NRE; j++) { pr_ilist(fp, indent, interaction_function[j].longname, mtop->ffparams.functype, @@ -391,7 +432,7 @@ void pr_mtop(FILE *fp, int indent, const char *title, const gmx_mtop_t *mtop, } pr_ffparams(fp, indent, "ffparams", &(mtop->ffparams), bShowNumbers); pr_atomtypes(fp, indent, "atomtypes", &(mtop->atomtypes), bShowNumbers); - for (mt = 0; mt < mtop->nmoltype; mt++) + for (size_t mt = 0; mt < mtop->moltype.size(); mt++) { pr_moltype(fp, indent, "moltype", &mtop->moltype[mt], mt, &mtop->ffparams, bShowNumbers, bShowParameters); diff --git a/src/gromacs/topology/topology.h b/src/gromacs/topology/topology.h index c9a226bd3e..04de807ffc 100644 --- a/src/gromacs/topology/topology.h +++ b/src/gromacs/topology/topology.h @@ -39,6 +39,8 @@ #include +#include + #include "gromacs/math/vectypes.h" #include "gromacs/topology/atoms.h" #include "gromacs/topology/block.h" @@ -54,18 +56,46 @@ enum { /* Names corresponding to groups */ extern const char *gtypes[egcNR+1]; -typedef struct gmx_moltype_t +/*! \brief Molecules type data: atoms, interactions and exclusions */ +struct gmx_moltype_t { - char **name; /* Name of the molecule type */ - t_atoms atoms; /* The atoms in this molecule */ - t_ilist ilist[F_NRE]; /* Interaction list with local indices */ - t_block cgs; /* The charge groups */ - t_blocka excls; /* The exclusions */ -} gmx_moltype_t; + /*! \brief Constructor */ + gmx_moltype_t(); + + /*! \brief Destructor */ + ~gmx_moltype_t(); + + /*! \brief Deleted copy assignment operator to avoid (not) freeing pointers */ + gmx_moltype_t &operator=(const gmx_moltype_t &) = delete; + + /*! \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 */ + t_ilist ilist[F_NRE]; /**< Interaction list with local indices */ + t_block cgs; /**< The charge groups */ + t_blocka excls; /**< The exclusions */ +}; /*! \brief Block of molecules of the same type, used in gmx_mtop_t */ -typedef struct gmx_molblock_t +struct gmx_molblock_t { + /*! \brief Constructor */ + gmx_molblock_t(); + + /*! \brief Destructor */ + ~gmx_molblock_t(); + + /*! \brief Default copy assignment operator. + * + * NOTE: Does not free the old pointers. + */ + gmx_molblock_t &operator=(const gmx_molblock_t &) = default; + + /*! \brief Default copy constructor */ + gmx_molblock_t(const gmx_molblock_t &) = default; + 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 */ @@ -80,7 +110,7 @@ typedef struct gmx_molblock_t 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 */ int moleculeIndexStart; /**< Global molecule indexing starts from this value */ -} gmx_molblock_t; +}; typedef struct gmx_groups_t { @@ -99,27 +129,31 @@ typedef struct gmx_groups_t /* The global, complete system topology struct, based on molecule types. This structure should contain no data that is O(natoms) in memory. */ -typedef struct gmx_mtop_t +struct gmx_mtop_t { - char **name; /* Name of the topology */ - gmx_ffparams_t ffparams; - int nmoltype; - gmx_moltype_t *moltype; - int nmolblock; - gmx_molblock_t *molblock; - gmx_bool bIntermolecularInteractions; /* Are there intermolecular - * interactions? */ - t_ilist *intermolecular_ilist; /* List of intermolecular interactions - * using system wide atom indices, - * either NULL or size F_NRE */ + /* Constructor */ + gmx_mtop_t(); + + /* Destructor */ + ~gmx_mtop_t(); + + char **name; /* Name of the topology */ + gmx_ffparams_t ffparams; + std::vector moltype; + std::vector molblock; + gmx_bool bIntermolecularInteractions; /* Are there intermolecular + * interactions? */ + t_ilist *intermolecular_ilist; /* List of intermolecular interactions + * using system wide atom indices, + * either NULL or size F_NRE */ int natoms; - int maxres_renum; /* Parameter for residue numbering */ - int maxresnr; /* The maximum residue number in moltype */ - t_atomtypes atomtypes; /* Atomtype properties */ + int maxres_renum; /* Parameter for residue numbering */ + int maxresnr; /* The maximum residue number in moltype */ + t_atomtypes atomtypes; /* Atomtype properties */ gmx_groups_t groups; - t_symtab symtab; /* The symbol table */ - bool haveMoleculeIndices; /* Tells whether we have valid molecule indices */ -} gmx_mtop_t; + t_symtab symtab; /* The symbol table */ + bool haveMoleculeIndices; /* Tells whether we have valid molecule indices */ +}; /* The mdrun node-local topology struct, completely written out */ typedef struct gmx_localtop_t @@ -146,10 +180,7 @@ typedef struct t_topology void init_mtop(gmx_mtop_t *mtop); void init_top(t_topology *top); -void done_moltype(gmx_moltype_t *molt); -void done_molblock(gmx_molblock_t *molb); void done_gmx_groups_t(gmx_groups_t *g); -void done_mtop(gmx_mtop_t *mtop); void done_top(t_topology *top); // Frees both t_topology and gmx_mtop_t when the former has been created from // the latter. diff --git a/src/gromacs/topology/topsort.cpp b/src/gromacs/topology/topsort.cpp index 4c3a7e80ce..d893701070 100644 --- a/src/gromacs/topology/topsort.cpp +++ b/src/gromacs/topology/topsort.cpp @@ -3,7 +3,7 @@ * * Copyright (c) 1991-2000, University of Groningen, The Netherlands. * Copyright (c) 2001-2008, The GROMACS development team. - * Copyright (c) 2013,2014,2015,2017, by the GROMACS development team, led by + * Copyright (c) 2013,2014,2015,2017,2018, 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. @@ -157,21 +157,13 @@ static gmx_bool ip_q_pert(int ftype, const t_iatom *ia, gmx_bool gmx_mtop_bondeds_free_energy(const gmx_mtop_t *mtop) { - const gmx_ffparams_t *ffparams; - int i, ftype; - int mb; - t_atom *atom; - t_ilist *il; - t_iatom *ia; - gmx_bool bPert; - - ffparams = &mtop->ffparams; + const gmx_ffparams_t *ffparams = &mtop->ffparams; /* Loop over all the function types and compare the A/B parameters */ - bPert = FALSE; - for (i = 0; i < ffparams->ntypes; i++) + gmx_bool bPert = FALSE; + for (int i = 0; i < ffparams->ntypes; i++) { - ftype = ffparams->functype[i]; + int ftype = ffparams->functype[i]; if (interaction_function[ftype].flags & IF_BOND) { if (ip_pert(ftype, &ffparams->iparams[i])) @@ -182,12 +174,12 @@ gmx_bool gmx_mtop_bondeds_free_energy(const gmx_mtop_t *mtop) } /* Check perturbed charges for 1-4 interactions */ - for (mb = 0; mb < mtop->nmolblock; mb++) + for (const gmx_molblock_t &molb : mtop->molblock) { - atom = mtop->moltype[mtop->molblock[mb].type].atoms.atom; - il = &mtop->moltype[mtop->molblock[mb].type].ilist[F_LJ14]; - ia = il->iatoms; - for (i = 0; i < il->nr; i += 3) + const t_atom *atom = mtop->moltype[molb.type].atoms.atom; + const t_ilist &il = mtop->moltype[molb.type].ilist[F_LJ14]; + const int *ia = il.iatoms; + for (int i = 0; i < il.nr; i += 3) { if (atom[ia[i+1]].q != atom[ia[i+1]].qB || atom[ia[i+2]].q != atom[ia[i+2]].qB) diff --git a/src/gromacs/trajectoryanalysis/analysissettings.cpp b/src/gromacs/trajectoryanalysis/analysissettings.cpp index 56df3ab97e..b6bb0d64db 100644 --- a/src/gromacs/trajectoryanalysis/analysissettings.cpp +++ b/src/gromacs/trajectoryanalysis/analysissettings.cpp @@ -1,7 +1,7 @@ /* * This file is part of the GROMACS molecular simulation package. * - * Copyright (c) 2010,2011,2012,2013,2014,2015,2016,2017, by the GROMACS development team, led by + * Copyright (c) 2010,2011,2012,2013,2014,2015,2016,2017,2018, 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. @@ -194,8 +194,7 @@ TopologyInformation::TopologyInformation() TopologyInformation::~TopologyInformation() { - done_top_mtop(top_, mtop_); - sfree(mtop_); + done_top_mtop(top_, mtop_.get()); sfree(top_); sfree(xtop_); } @@ -206,7 +205,7 @@ t_topology *TopologyInformation::topology() const if (top_ == nullptr && mtop_ != nullptr) { snew(top_, 1); - *top_ = gmx_mtop_t_to_t_topology(mtop_, false); + *top_ = gmx_mtop_t_to_t_topology(mtop_.get(), false); } return top_; } diff --git a/src/gromacs/trajectoryanalysis/analysissettings.h b/src/gromacs/trajectoryanalysis/analysissettings.h index 1e8c879374..cc7033c7c3 100644 --- a/src/gromacs/trajectoryanalysis/analysissettings.h +++ b/src/gromacs/trajectoryanalysis/analysissettings.h @@ -1,7 +1,7 @@ /* * This file is part of the GROMACS molecular simulation package. * - * Copyright (c) 2010,2011,2012,2014,2015,2016,2017, by the GROMACS development team, led by + * Copyright (c) 2010,2011,2012,2014,2015,2016,2017,2018, 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. @@ -258,7 +258,7 @@ class TopologyInformation //! Returns true if a full topology file was loaded. bool hasFullTopology() const { return bTop_; } //! Returns the loaded topology, or NULL if not loaded. - const gmx_mtop_t *mtop() const { return mtop_; } + const gmx_mtop_t *mtop() const { return mtop_.get(); } //! Returns the loaded topology, or NULL if not loaded. t_topology *topology() const; //! Returns the ePBC field from the topology. @@ -284,7 +284,7 @@ class TopologyInformation TopologyInformation(); ~TopologyInformation(); - gmx_mtop_t *mtop_; + std::unique_ptr mtop_; //! The topology structure, or NULL if no topology loaded. // TODO: Replace fully with mtop. mutable t_topology *top_; diff --git a/src/gromacs/trajectoryanalysis/runnercommon.cpp b/src/gromacs/trajectoryanalysis/runnercommon.cpp index 38c74f1f96..80e70bcd32 100644 --- a/src/gromacs/trajectoryanalysis/runnercommon.cpp +++ b/src/gromacs/trajectoryanalysis/runnercommon.cpp @@ -1,7 +1,7 @@ /* * This file is part of the GROMACS molecular simulation package. * - * Copyright (c) 2010,2011,2012,2013,2014,2015,2016,2017, by the GROMACS development team, led by + * Copyright (c) 2010,2011,2012,2013,2014,2015,2016,2017,2018, 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. @@ -48,6 +48,7 @@ #include #include +#include "gromacs/compat/make_unique.h" #include "gromacs/fileio/confio.h" #include "gromacs/fileio/oenv.h" #include "gromacs/fileio/timecontrol.h" @@ -93,7 +94,7 @@ class TrajectoryAnalysisRunnerCommon::Impl : public ITopologyProvider virtual gmx_mtop_t *getTopology(bool required) { initTopology(required); - return topInfo_.mtop_; + return topInfo_.mtop_.get(); } virtual int getAtomCount() { @@ -180,15 +181,14 @@ TrajectoryAnalysisRunnerCommon::Impl::initTopology(bool required) // Load the topology if requested. if (!topfile_.empty()) { - snew(topInfo_.mtop_, 1); - readConfAndTopology(topfile_.c_str(), &topInfo_.bTop_, topInfo_.mtop_, + topInfo_.mtop_ = gmx::compat::make_unique(); + readConfAndTopology(topfile_.c_str(), &topInfo_.bTop_, topInfo_.mtop_.get(), &topInfo_.ePBC_, &topInfo_.xtop_, nullptr, topInfo_.boxtop_); // TODO: Only load this here if the tool actually needs it; selections // take care of themselves. - for (int i = 0; i < topInfo_.mtop_->nmoltype; ++i) + for (gmx_moltype_t &moltype : topInfo_.mtop_->moltype) { - gmx_moltype_t &moltype = topInfo_.mtop_->moltype[i]; if (!moltype.atoms.haveMass) { // Try to read masses from database, be silent about missing masses diff --git a/src/programs/mdrun/membed.cpp b/src/programs/mdrun/membed.cpp index e6f09b7380..3779c96ea8 100644 --- a/src/programs/mdrun/membed.cpp +++ b/src/programs/mdrun/membed.cpp @@ -122,12 +122,11 @@ static int get_mol_id(int at, gmx_mtop_t *mtop, int *type, int *block) } /* Get the id of the molblock from a global molecule id */ -static int get_molblock(int mol_id, int nmblock, gmx_molblock_t *mblock) +static int get_molblock(int mol_id, const std::vector &mblock) { - int i; int nmol = 0; - for (i = 0; i < nmblock; i++) + for (size_t i = 0; i < mblock.size(); i++) { nmol += mblock[i].nmol; if (mol_id < nmol) @@ -637,7 +636,7 @@ static int gen_rm_list(rm_t *rm_p, t_block *ins_at, t_block *rest_at, t_pbc *pbc while (nupper != nlower) { mol_id = mem_p->mol_id[order[i]]; - block = get_molblock(mol_id, mtop->nmolblock, mtop->molblock); + block = get_molblock(mol_id, mtop->molblock); bRM = TRUE; for (l = 0; l < nrm; l++) { @@ -692,7 +691,7 @@ static int gen_rm_list(rm_t *rm_p, t_block *ins_at, t_block *rest_at, t_pbc *pbc static void rm_group(gmx_groups_t *groups, gmx_mtop_t *mtop, rm_t *rm_p, t_state *state, t_block *ins_at, pos_ins_t *pos_ins) { - int i, j, k, n, rm, mol_id, at, block; + int j, k, n, rm, mol_id, at, block; rvec *x_tmp, *v_tmp; int *list; unsigned char *new_egrp[egcNR]; @@ -704,7 +703,7 @@ static void rm_group(gmx_groups_t *groups, gmx_mtop_t *mtop, rm_t *rm_p, t_state snew(list, state->natoms); n = 0; - for (i = 0; i < rm_p->nr; i++) + for (int i = 0; i < rm_p->nr; i++) { mol_id = rm_p->mol[i]; at = molecules.index[mol_id]; @@ -722,7 +721,7 @@ static void rm_group(gmx_groups_t *groups, gmx_mtop_t *mtop, rm_t *rm_p, t_state snew(x_tmp, state->natoms); snew(v_tmp, state->natoms); - for (i = 0; i < egcNR; i++) + for (int i = 0; i < egcNR; i++) { if (groups->grpnr[i] != nullptr) { @@ -732,7 +731,7 @@ static void rm_group(gmx_groups_t *groups, gmx_mtop_t *mtop, rm_t *rm_p, t_state } rm = 0; - for (i = 0; i < state->natoms+n; i++) + for (int i = 0; i < state->natoms+n; i++) { bRM = FALSE; for (j = 0; j < n; j++) @@ -786,7 +785,7 @@ static void rm_group(gmx_groups_t *groups, gmx_mtop_t *mtop, rm_t *rm_p, t_state } sfree(v_tmp); - for (i = 0; i < egcNR; i++) + for (int i = 0; i < egcNR; i++) { if (groups->grpnr[i] != nullptr) { @@ -797,7 +796,7 @@ static void rm_group(gmx_groups_t *groups, gmx_mtop_t *mtop, rm_t *rm_p, t_state /* remove empty molblocks */ RMmolblock = 0; - for (i = 0; i < mtop->nmolblock; i++) + for (size_t i = 0; i < mtop->molblock.size(); i++) { if (mtop->molblock[i].nmol == 0) { @@ -808,13 +807,13 @@ static void rm_group(gmx_groups_t *groups, gmx_mtop_t *mtop, rm_t *rm_p, t_state mtop->molblock[i-RMmolblock] = mtop->molblock[i]; } } - mtop->nmolblock -= RMmolblock; + mtop->molblock.resize(mtop->molblock.size() - RMmolblock); } /* remove al bonded interactions from mtop for the molecule to be embedded */ static int rm_bonded(t_block *ins_at, gmx_mtop_t *mtop) { - int i, j, m; + int j, m; int type, natom, nmol, at, atom1 = 0, rm_at = 0; gmx_bool *bRM, bINS; /*this routine lives dangerously by assuming that all molecules of a given type are in order in the structure*/ @@ -823,13 +822,13 @@ static int rm_bonded(t_block *ins_at, gmx_mtop_t *mtop) * ins_at index group. MGWolf 050710 */ - snew(bRM, mtop->nmoltype); - for (i = 0; i < mtop->nmoltype; i++) + snew(bRM, mtop->moltype.size()); + for (size_t i = 0; i < mtop->moltype.size(); i++) { bRM[i] = TRUE; } - for (i = 0; i < mtop->nmolblock; i++) + for (size_t i = 0; i < mtop->molblock.size(); i++) { /*loop over molecule blocks*/ type = mtop->molblock[i].type; @@ -859,7 +858,7 @@ static int rm_bonded(t_block *ins_at, gmx_mtop_t *mtop) } } - for (i = 0; i < mtop->nmoltype; i++) + for (size_t i = 0; i < mtop->moltype.size(); i++) { if (bRM[i]) { @@ -893,7 +892,7 @@ static void top_update(const char *topfile, rm_t *rm_p, gmx_mtop_t *mtop) gmx_tmpnam(temporary_filename); fpout = gmx_ffopen(temporary_filename, "w"); - snew(nmol_rm, mtop->nmoltype); + snew(nmol_rm, mtop->moltype.size()); for (i = 0; i < rm_p->nr; i++) { nmol_rm[rm_p->block[i]]++; @@ -933,7 +932,7 @@ static void top_update(const char *topfile, rm_t *rm_p, gmx_mtop_t *mtop) } else if (bMolecules == 1) { - for (i = 0; i < mtop->nmolblock; i++) + for (size_t i = 0; i < mtop->molblock.size(); i++) { nmol = mtop->molblock[i].nmol; sprintf(buf, "%-15s %5d\n", *(mtop->moltype[mtop->molblock[i].type].name), nmol); @@ -1268,12 +1267,12 @@ gmx_membed_t *init_membed(FILE *fplog, int nfile, const t_filenm fnm[], gmx_mtop } } - for (i = 0; i < mtop->nmolblock; i++) + for (size_t i = 0; i < mtop->molblock.size(); i++) { ntype = 0; for (j = 0; j < rm_p->nr; j++) { - if (rm_p->block[j] == i) + if (rm_p->block[j] == static_cast(i)) { ntype++; } diff --git a/src/programs/mdrun/runner.cpp b/src/programs/mdrun/runner.cpp index 9e8fe71c58..61410e91e3 100644 --- a/src/programs/mdrun/runner.cpp +++ b/src/programs/mdrun/runner.cpp @@ -425,7 +425,6 @@ int Mdrunner::mdrunner() gmx_ddbox_t ddbox = {0}; int npme_major, npme_minor; t_nrnb *nrnb; - gmx_mtop_t *mtop = nullptr; t_forcerec *fr = nullptr; t_fcdata *fcd = nullptr; real ewaldcoeff_q = 0; @@ -446,7 +445,7 @@ int Mdrunner::mdrunner() std::unique_ptr mdModules(new gmx::MDModules); t_inputrec inputrecInstance; t_inputrec *inputrec = &inputrecInstance; - snew(mtop, 1); + gmx_mtop_t mtop; if (mdrunOptions.continuationOptions.appendFiles) { @@ -555,7 +554,7 @@ int Mdrunner::mdrunner() globalState = std::unique_ptr(new t_state); /* Read (nearly) all data required for the simulation */ - read_tpx_state(ftp2fn(efTPR, nfile, fnm), inputrec, globalState.get(), mtop); + read_tpx_state(ftp2fn(efTPR, nfile, fnm), inputrec, globalState.get(), &mtop); if (inputrec->cutoff_scheme != ecutsVERLET) { @@ -625,7 +624,7 @@ int Mdrunner::mdrunner() gpuIdsToUse, useGpuForNonbonded, useGpuForPme, - inputrec, mtop, + inputrec, &mtop, mdlog, doMembed); @@ -643,7 +642,7 @@ int Mdrunner::mdrunner() if (PAR(cr)) { /* now broadcast everything to the non-master nodes/threads: */ - init_parallel(cr, inputrec, mtop); + init_parallel(cr, inputrec, &mtop); } // Now each rank knows the inputrec that SIMMASTER read and used, @@ -799,9 +798,9 @@ int Mdrunner::mdrunner() snew(fcd, 1); /* This needs to be called before read_checkpoint to extend the state */ - init_disres(fplog, mtop, inputrec, cr, ms, fcd, globalState.get(), replExParams.exchangeInterval > 0); + init_disres(fplog, &mtop, inputrec, cr, ms, fcd, globalState.get(), replExParams.exchangeInterval > 0); - init_orires(fplog, mtop, inputrec, cr, ms, globalState.get(), &(fcd->orires)); + init_orires(fplog, &mtop, inputrec, cr, ms, globalState.get(), &(fcd->orires)); if (inputrecDeform(inputrec)) { @@ -881,7 +880,7 @@ int Mdrunner::mdrunner() /* Update rlist and nstlist. */ if (inputrec->cutoff_scheme == ecutsVERLET) { - prepare_verlet_scheme(fplog, cr, inputrec, nstlist_cmdline, mtop, box, + prepare_verlet_scheme(fplog, cr, inputrec, nstlist_cmdline, &mtop, box, useGpuForNonbonded || (emulateGpuNonbonded == EmulateGpuNonbonded::Yes), *hwinfo->cpuInfo); } @@ -891,7 +890,7 @@ int Mdrunner::mdrunner() const rvec *xOnMaster = (SIMMASTER(cr) ? as_rvec_array(globalState->x.data()) : nullptr); cr->dd = init_domain_decomposition(fplog, cr, domdecOptions, mdrunOptions, - mtop, inputrec, + &mtop, inputrec, box, xOnMaster, &ddbox, &npme_major, &npme_minor); // Note that local state still does not exist yet. @@ -945,7 +944,7 @@ int Mdrunner::mdrunner() /* Check and update the number of OpenMP threads requested */ checkAndUpdateRequestedNumOpenmpThreads(&hw_opt, *hwinfo, cr, ms, physicalNodeComm.size_, - pmeRunMode, *mtop); + pmeRunMode, mtop); gmx_omp_nthreads_init(mdlog, cr, hwinfo->nthreads_hw_avail, @@ -1157,7 +1156,7 @@ int Mdrunner::mdrunner() /* Note that membed cannot work in parallel because mtop is * changed here. Fix this if we ever want to make it run with * multiple ranks. */ - membed = init_membed(fplog, nfile, fnm, mtop, inputrec, globalState.get(), cr, &mdrunOptions.checkpointOptions.period); + membed = init_membed(fplog, nfile, fnm, &mtop, inputrec, globalState.get(), cr, &mdrunOptions.checkpointOptions.period); } std::unique_ptr mdAtoms; @@ -1169,7 +1168,7 @@ int Mdrunner::mdrunner() fr = mk_forcerec(); fr->forceProviders = mdModules->initForceProviders(); init_forcerec(fplog, mdlog, fr, fcd, - inputrec, mtop, cr, box, + inputrec, &mtop, cr, box, opt2fn("-table", nfile, fnm), opt2fn("-tablep", nfile, fnm), opt2fns("-tableb", nfile, fnm), @@ -1184,7 +1183,7 @@ int Mdrunner::mdrunner() appendText("Large parts of the QM/MM support is deprecated, and may be removed in a future " "version. Please get in touch with the developers if you find the support useful, " "as help is needed if the functionality is to continue to be available."); - init_QMMMrec(cr, mtop, inputrec, fr); + init_QMMMrec(cr, &mtop, inputrec, fr); } /* Initialize the mdAtoms structure. @@ -1192,7 +1191,7 @@ int Mdrunner::mdrunner() * as this can not be done now with domain decomposition. */ const bool useGpuForPme = (pmeRunMode == PmeRunMode::GPU) || (pmeRunMode == PmeRunMode::Mixed); - mdAtoms = makeMDAtoms(fplog, *mtop, *inputrec, useGpuForPme && thisRankHasDuty(cr, DUTY_PME)); + mdAtoms = makeMDAtoms(fplog, mtop, *inputrec, useGpuForPme && thisRankHasDuty(cr, DUTY_PME)); if (globalState) { // The pinning of coordinates in the global state object works, because we only use @@ -1204,7 +1203,7 @@ int Mdrunner::mdrunner() } /* Initialize the virtual site communication */ - vsite = initVsite(*mtop, cr); + vsite = initVsite(mtop, cr); calc_shifts(box, fr->shift_vec); @@ -1218,7 +1217,7 @@ int Mdrunner::mdrunner() if (fr->ePBC != epbcNONE) { rvec *xGlobal = as_rvec_array(globalState->x.data()); - do_pbc_first_mtop(fplog, inputrec->ePBC, box, mtop, xGlobal); + do_pbc_first_mtop(fplog, inputrec->ePBC, box, &mtop, xGlobal); } if (vsite) { @@ -1226,7 +1225,7 @@ int Mdrunner::mdrunner() * for the initial distribution in the domain decomposition * and for the initial shell prediction. */ - constructVsitesGlobal(*mtop, globalState->x); + constructVsitesGlobal(mtop, globalState->x); } } @@ -1275,7 +1274,7 @@ int Mdrunner::mdrunner() try { pmedata = gmx_pme_init(cr, npme_major, npme_minor, inputrec, - mtop ? mtop->natoms : 0, nChargePerturbed, nTypePerturbed, + mtop.natoms, nChargePerturbed, nTypePerturbed, mdrunOptions.reproducible, ewaldcoeff_q, ewaldcoeff_lj, nthreads_pme, @@ -1306,7 +1305,7 @@ int Mdrunner::mdrunner() /* Initialize pull code */ inputrec->pull_work = init_pull(fplog, inputrec->pull, inputrec, - mtop, cr, inputrec->fepvals->init_lambda); + &mtop, cr, inputrec->fepvals->init_lambda); if (EI_DYNAMICS(inputrec->eI) && MASTER(cr)) { init_pull_output_files(inputrec->pull_work, @@ -1318,7 +1317,7 @@ int Mdrunner::mdrunner() if (inputrec->bRot) { /* Initialize enforced rotation code */ - init_rot(fplog, inputrec, nfile, fnm, cr, globalState.get(), mtop, oenv, mdrunOptions); + init_rot(fplog, inputrec, nfile, fnm, cr, globalState.get(), &mtop, oenv, mdrunOptions); } /* Let init_constraints know whether we have essential dynamics constraints. @@ -1326,7 +1325,7 @@ int Mdrunner::mdrunner() */ bool doEdsam = (opt2fn_null("-ei", nfile, fnm) != nullptr || observablesHistory.edsamHistory); - constr = init_constraints(fplog, mtop, inputrec, doEdsam, cr); + constr = init_constraints(fplog, &mtop, inputrec, doEdsam, cr); if (DOMAINDECOMP(cr)) { @@ -1334,7 +1333,7 @@ int Mdrunner::mdrunner() /* This call is not included in init_domain_decomposition mainly * because fr->cginfo_mb is set later. */ - dd_init_bondeds(fplog, cr->dd, mtop, vsite, inputrec, + dd_init_bondeds(fplog, cr->dd, &mtop, vsite, inputrec, domdecOptions.checkBondedInteractions, fr->cginfo_mb); } @@ -1345,7 +1344,7 @@ int Mdrunner::mdrunner() mdrunOptions, vsite, constr, mdModules->outputProvider(), - inputrec, mtop, + inputrec, &mtop, fcd, globalState.get(), &observablesHistory, @@ -1403,8 +1402,7 @@ int Mdrunner::mdrunner() free_gpu_resources(fr, physicalNodeComm); free_gpu(nonbondedDeviceInfo); free_gpu(pmeDeviceInfo); - done_forcerec(fr, mtop->nmolblock, mtop->groups.grps[egcENER].nr); - done_mtop(mtop); + done_forcerec(fr, mtop.molblock.size(), mtop.groups.grps[egcENER].nr); sfree(fcd); if (doMembed) -- 2.11.4.GIT