2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
5 * Copyright (c) 2001-2004, The GROMACS development team.
6 * Copyright (c) 2013,2014,2015,2016,2017,2018,2019, by the GROMACS development team, led by
7 * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
8 * and including many others, as listed in the AUTHORS file in the
9 * top-level source directory and at http://www.gromacs.org.
11 * GROMACS is free software; you can redistribute it and/or
12 * modify it under the terms of the GNU Lesser General Public License
13 * as published by the Free Software Foundation; either version 2.1
14 * of the License, or (at your option) any later version.
16 * GROMACS is distributed in the hope that it will be useful,
17 * but WITHOUT ANY WARRANTY; without even the implied warranty of
18 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
19 * Lesser General Public License for more details.
21 * You should have received a copy of the GNU Lesser General Public
22 * License along with GROMACS; if not, see
23 * http://www.gnu.org/licenses, or write to the Free Software Foundation,
24 * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
26 * If you want to redistribute modifications to GROMACS, please
27 * consider that scientific software is very special. Version
28 * control is crucial - bugs must be traceable. We will be happy to
29 * consider code for inclusion in the official distribution, but
30 * derived work must not be called official GROMACS. Details are found
31 * in the README & COPYING files - if they are missing, get the
32 * official version at http://www.gromacs.org.
34 * To help us fund GROMACS development, we humbly ask that you cite
35 * the research papers on the package. Check out http://www.gromacs.org.
45 #include "gromacs/math/vecdump.h"
46 #include "gromacs/topology/atoms.h"
47 #include "gromacs/topology/idef.h"
48 #include "gromacs/topology/ifunc.h"
49 #include "gromacs/topology/symtab.h"
50 #include "gromacs/utility/arrayref.h"
51 #include "gromacs/utility/compare.h"
52 #include "gromacs/utility/gmxassert.h"
53 #include "gromacs/utility/smalloc.h"
54 #include "gromacs/utility/strconvert.h"
55 #include "gromacs/utility/txtdump.h"
57 const char *gtypes[egcNR+1] = {
58 "T-Coupling", "Energy Mon.", "Acceleration", "Freeze",
59 "User1", "User2", "VCM", "Compressed X", "Or. Res. Fit", "QMMM", nullptr
62 static void init_groups(gmx_groups_t *groups)
65 groups->grpname = nullptr;
66 for (int g = 0; g < egcNR; g++)
68 groups->grps[g].nr = 0;
69 groups->grps[g].nm_ind = nullptr;
70 groups->ngrpnr[g] = 0;
71 groups->grpnr[g] = nullptr;
76 void init_top(t_topology *top)
79 init_idef(&top->idef);
80 init_atom(&(top->atoms));
81 init_atomtypes(&(top->atomtypes));
82 init_block(&top->cgs);
83 init_block(&top->mols);
84 init_blocka(&top->excls);
85 open_symtab(&top->symtab);
89 gmx_moltype_t::gmx_moltype_t() :
94 init_t_atoms(&atoms, 0, FALSE);
97 gmx_moltype_t::~gmx_moltype_t()
104 void done_gmx_groups_t(gmx_groups_t *g)
108 for (i = 0; (i < egcNR); i++)
110 if (nullptr != g->grps[i].nm_ind)
112 sfree(g->grps[i].nm_ind);
113 g->grps[i].nm_ind = nullptr;
115 if (nullptr != g->grpnr[i])
118 g->grpnr[i] = nullptr;
121 /* The contents of this array is in symtab, don't free it here */
125 gmx_mtop_t::gmx_mtop_t()
127 init_atomtypes(&atomtypes);
128 init_groups(&groups);
129 open_symtab(&symtab);
132 gmx_mtop_t::~gmx_mtop_t()
134 done_symtab(&symtab);
138 done_atomtypes(&atomtypes);
139 done_gmx_groups_t(&groups);
142 void done_top(t_topology *top)
144 done_idef(&top->idef);
145 done_atom(&(top->atoms));
148 done_atomtypes(&(top->atomtypes));
150 done_symtab(&(top->symtab));
151 done_block(&(top->cgs));
152 done_block(&(top->mols));
153 done_blocka(&(top->excls));
156 void done_top_mtop(t_topology *top, gmx_mtop_t *mtop)
162 done_idef(&top->idef);
163 done_atom(&top->atoms);
164 done_block(&top->cgs);
165 done_blocka(&top->excls);
166 done_block(&top->mols);
167 done_symtab(&top->symtab);
168 open_symtab(&mtop->symtab);
169 done_atomtypes(&(top->atomtypes));
172 // Note that the rest of mtop will be freed by the destructor
176 gmx_localtop_t::gmx_localtop_t()
178 init_block_null(&cgs);
179 init_blocka_null(&excls);
181 init_atomtypes(&atomtypes);
184 gmx_localtop_t::~gmx_localtop_t()
186 if (!useInDomainDecomp_)
191 done_atomtypes(&atomtypes);
195 bool gmx_mtop_has_masses(const gmx_mtop_t *mtop)
201 return mtop->moltype.empty() || mtop->moltype[0].atoms.haveMass;
204 bool gmx_mtop_has_charges(const gmx_mtop_t *mtop)
210 return mtop->moltype.empty() || mtop->moltype[0].atoms.haveCharge;
213 bool gmx_mtop_has_perturbed_charges(const gmx_mtop_t &mtop)
215 for (const gmx_moltype_t &moltype : mtop.moltype)
217 const t_atoms &atoms = moltype.atoms;
218 if (atoms.haveBState)
220 for (int a = 0; a < atoms.nr; a++)
222 if (atoms.atom[a].q != atoms.atom[a].qB)
232 bool gmx_mtop_has_atomtypes(const gmx_mtop_t *mtop)
238 return mtop->moltype.empty() || mtop->moltype[0].atoms.haveType;
241 bool gmx_mtop_has_pdbinfo(const gmx_mtop_t *mtop)
247 return mtop->moltype.empty() || mtop->moltype[0].atoms.havePdbInfo;
250 static void pr_grps(FILE *fp, const char *title, const t_grps grps[], char **grpname[])
254 for (i = 0; (i < egcNR); i++)
256 fprintf(fp, "%s[%-12s] nr=%d, name=[", title, gtypes[i], grps[i].nr);
257 for (j = 0; (j < grps[i].nr); j++)
259 fprintf(fp, " %s", *(grpname[grps[i].nm_ind[j]]));
265 static void pr_groups(FILE *fp, int indent,
266 const gmx_groups_t *groups,
267 gmx_bool bShowNumbers)
271 pr_grps(fp, "grp", groups->grps, groups->grpname);
272 pr_strings(fp, indent, "grpname", groups->grpname, groups->ngrpname, bShowNumbers);
274 pr_indent(fp, indent);
275 fprintf(fp, "groups ");
276 for (g = 0; g < egcNR; g++)
278 printf(" %5.5s", gtypes[g]);
282 pr_indent(fp, indent);
283 fprintf(fp, "allocated ");
285 for (g = 0; g < egcNR; g++)
287 printf(" %5d", groups->ngrpnr[g]);
288 nat_max = std::max(nat_max, groups->ngrpnr[g]);
294 pr_indent(fp, indent);
295 fprintf(fp, "groupnr[%5s] =", "*");
296 for (g = 0; g < egcNR; g++)
298 fprintf(fp, " %3d ", 0);
304 for (i = 0; i < nat_max; i++)
306 pr_indent(fp, indent);
307 fprintf(fp, "groupnr[%5d] =", i);
308 for (g = 0; g < egcNR; g++)
311 groups->grpnr[g] ? groups->grpnr[g][i] : 0);
318 static void pr_moltype(FILE *fp, int indent, const char *title,
319 const gmx_moltype_t *molt, int n,
320 const gmx_ffparams_t *ffparams,
321 gmx_bool bShowNumbers, gmx_bool bShowParameters)
325 indent = pr_title_n(fp, indent, title, n);
326 pr_indent(fp, indent);
327 fprintf(fp, "name=\"%s\"\n", *(molt->name));
328 pr_atoms(fp, indent, "atoms", &(molt->atoms), bShowNumbers);
329 pr_block(fp, indent, "cgs", &molt->cgs, bShowNumbers);
330 pr_blocka(fp, indent, "excls", &molt->excls, bShowNumbers);
331 for (j = 0; (j < F_NRE); j++)
333 pr_ilist(fp, indent, interaction_function[j].longname,
334 ffparams->functype.data(), molt->ilist[j],
335 bShowNumbers, bShowParameters, ffparams->iparams.data());
339 static void pr_molblock(FILE *fp, int indent, const char *title,
340 const gmx_molblock_t *molb, int n,
341 const std::vector<gmx_moltype_t> &molt)
343 indent = pr_title_n(fp, indent, title, n);
344 pr_indent(fp, indent);
345 fprintf(fp, "%-20s = %d \"%s\"\n",
346 "moltype", molb->type, *(molt[molb->type].name));
347 pr_int(fp, indent, "#molecules", molb->nmol);
348 pr_int(fp, indent, "#posres_xA", molb->posres_xA.size());
349 if (!molb->posres_xA.empty())
351 pr_rvecs(fp, indent, "posres_xA", as_rvec_array(molb->posres_xA.data()), molb->posres_xA.size());
353 pr_int(fp, indent, "#posres_xB", molb->posres_xB.size());
354 if (!molb->posres_xB.empty())
356 pr_rvecs(fp, indent, "posres_xB", as_rvec_array(molb->posres_xB.data()), molb->posres_xB.size());
360 void pr_mtop(FILE *fp, int indent, const char *title, const gmx_mtop_t *mtop,
361 gmx_bool bShowNumbers, gmx_bool bShowParameters)
363 if (available(fp, mtop, indent, title))
365 indent = pr_title(fp, indent, title);
366 pr_indent(fp, indent);
367 fprintf(fp, "name=\"%s\"\n", *(mtop->name));
368 pr_int(fp, indent, "#atoms", mtop->natoms);
369 pr_int(fp, indent, "#molblock", mtop->molblock.size());
370 for (size_t mb = 0; mb < mtop->molblock.size(); mb++)
372 pr_molblock(fp, indent, "molblock", &mtop->molblock[mb], mb, mtop->moltype);
374 pr_str(fp, indent, "bIntermolecularInteractions",
375 gmx::boolToString(mtop->bIntermolecularInteractions));
376 if (mtop->bIntermolecularInteractions)
378 for (int j = 0; j < F_NRE; j++)
380 pr_ilist(fp, indent, interaction_function[j].longname,
381 mtop->ffparams.functype.data(),
382 (*mtop->intermolecular_ilist)[j],
383 bShowNumbers, bShowParameters, mtop->ffparams.iparams.data());
386 pr_ffparams(fp, indent, "ffparams", &(mtop->ffparams), bShowNumbers);
387 pr_atomtypes(fp, indent, "atomtypes", &(mtop->atomtypes), bShowNumbers);
388 for (size_t mt = 0; mt < mtop->moltype.size(); mt++)
390 pr_moltype(fp, indent, "moltype", &mtop->moltype[mt], mt,
391 &mtop->ffparams, bShowNumbers, bShowParameters);
393 pr_groups(fp, indent, &mtop->groups, bShowNumbers);
397 void pr_top(FILE *fp, int indent, const char *title, const t_topology *top,
398 gmx_bool bShowNumbers, gmx_bool bShowParameters)
400 if (available(fp, top, indent, title))
402 indent = pr_title(fp, indent, title);
403 pr_indent(fp, indent);
404 fprintf(fp, "name=\"%s\"\n", *(top->name));
405 pr_atoms(fp, indent, "atoms", &(top->atoms), bShowNumbers);
406 pr_atomtypes(fp, indent, "atomtypes", &(top->atomtypes), bShowNumbers);
407 pr_block(fp, indent, "cgs", &top->cgs, bShowNumbers);
408 pr_block(fp, indent, "mols", &top->mols, bShowNumbers);
409 pr_str(fp, indent, "bIntermolecularInteractions",
410 gmx::boolToString(top->bIntermolecularInteractions));
411 pr_blocka(fp, indent, "excls", &top->excls, bShowNumbers);
412 pr_idef(fp, indent, "idef", &top->idef, bShowNumbers, bShowParameters);
416 static void cmp_iparm(FILE *fp, const char *s, t_functype ft,
417 const t_iparams &ip1, const t_iparams &ip2, real relativeTolerance, real absoluteTolerance)
423 for (i = 0; i < MAXFORCEPARAM && !bDiff; i++)
425 bDiff = !equal_real(ip1.generic.buf[i], ip2.generic.buf[i], relativeTolerance, absoluteTolerance);
429 fprintf(fp, "%s1: ", s);
430 pr_iparams(fp, ft, &ip1);
431 fprintf(fp, "%s2: ", s);
432 pr_iparams(fp, ft, &ip2);
436 static void cmp_iparm_AB(FILE *fp, const char *s, t_functype ft,
437 const t_iparams &ip1, real relativeTolerance, real absoluteTolerance)
439 int nrfpA, nrfpB, p0, i;
442 /* Normally the first parameter is perturbable */
444 nrfpA = interaction_function[ft].nrfpA;
445 nrfpB = interaction_function[ft].nrfpB;
450 else if (interaction_function[ft].flags & IF_TABULATED)
452 /* For tabulated interactions only the second parameter is perturbable */
457 for (i = 0; i < nrfpB && !bDiff; i++)
459 bDiff = !equal_real(ip1.generic.buf[p0+i], ip1.generic.buf[nrfpA+i], relativeTolerance, absoluteTolerance);
463 fprintf(fp, "%s: ", s);
464 pr_iparams(fp, ft, &ip1);
468 static void cmp_cmap(FILE *fp, const gmx_cmap_t *cmap1, const gmx_cmap_t *cmap2, real relativeTolerance, real absoluteTolerance)
470 int cmap1_ngrid = (cmap1 ? cmap1->cmapdata.size() : 0);
471 int cmap2_ngrid = (cmap2 ? cmap2->cmapdata.size() : 0);
473 cmp_int(fp, "cmap ngrid", -1, cmap1_ngrid, cmap2_ngrid);
475 if (cmap1 == nullptr || cmap2 == nullptr)
480 cmp_int(fp, "cmap grid_spacing", -1, cmap1->grid_spacing, cmap2->grid_spacing);
481 if (cmap1->cmapdata.size() == cmap2->cmapdata.size() &&
482 cmap1->grid_spacing == cmap2->grid_spacing)
484 for (size_t g = 0; g < cmap1->cmapdata.size(); g++)
488 fprintf(fp, "comparing cmap %zu\n", g);
490 for (i = 0; i < 4*cmap1->grid_spacing*cmap1->grid_spacing; i++)
492 cmp_real(fp, "", i, cmap1->cmapdata[g].cmap[i], cmap2->cmapdata[g].cmap[i], relativeTolerance, absoluteTolerance);
498 static void cmp_block(FILE *fp, const t_block *b1, const t_block *b2, const char *s)
502 fprintf(fp, "comparing block %s\n", s);
503 sprintf(buf, "%s.nr", s);
504 cmp_int(fp, buf, -1, b1->nr, b2->nr);
507 static void cmp_blocka(FILE *fp, const t_blocka *b1, const t_blocka *b2, const char *s)
511 fprintf(fp, "comparing blocka %s\n", s);
512 sprintf(buf, "%s.nr", s);
513 cmp_int(fp, buf, -1, b1->nr, b2->nr);
514 sprintf(buf, "%s.nra", s);
515 cmp_int(fp, buf, -1, b1->nra, b2->nra);
518 static void compareFfparams(FILE *fp, const gmx_ffparams_t &ff1, const gmx_ffparams_t &ff2, real relativeTolerance, real absoluteTolerance)
520 fprintf(fp, "comparing force field parameters\n");
521 cmp_int(fp, "numTypes", -1, ff1.numTypes(), ff2.numTypes());
522 cmp_int(fp, "atnr", -1, ff1.atnr, ff1.atnr);
523 cmp_double(fp, "reppow", -1, ff1.reppow, ff2.reppow, relativeTolerance, absoluteTolerance);
524 cmp_real(fp, "fudgeQQ", -1, ff1.fudgeQQ, ff2.fudgeQQ, relativeTolerance, absoluteTolerance);
525 cmp_cmap(fp, &ff1.cmap_grid, &ff2.cmap_grid, relativeTolerance, absoluteTolerance);
526 for (int i = 0; i < std::min(ff1.numTypes(), ff2.numTypes()); i++)
528 std::string buf = gmx::formatString("ffparams->functype[%d]", i);
529 cmp_int(fp, buf.c_str(), i, ff1.functype[i], ff2.functype[i]);
530 buf = gmx::formatString("ffparams->iparams[%d]", i);
531 cmp_iparm(fp, buf.c_str(), ff1.functype[i], ff1.iparams[i], ff2.iparams[i], relativeTolerance, absoluteTolerance);
536 static void compareFfparamAB(FILE *fp, const gmx_ffparams_t &ff1, real relativeTolerance, real absoluteTolerance)
538 fprintf(fp, "comparing free energy parameters\n");
539 for (int i = 0; i < ff1.numTypes(); i++)
541 std::string buf = gmx::formatString("ffparams->iparams[%d]", i);
542 cmp_iparm_AB(fp, buf.c_str(), ff1.functype[i], ff1.iparams[i], relativeTolerance, absoluteTolerance);
545 static void compareInteractionLists(FILE *fp, const InteractionLists *il1, const InteractionLists *il2)
547 fprintf(fp, "comparing InteractionLists\n");
548 if ((il1 || il2) && (!il1 || !il2))
550 fprintf(fp, "InteractionLists are present in topology %d but not in the other\n", il1 ? 1 : 2);
554 for (int i = 0; i < F_NRE; i++)
556 cmp_int(fp, "InteractionList size", i, il1->at(i).size(), il2->at(i).size());
557 int nr = std::min(il1->at(i).size(), il2->at(i).size());
558 for (int j = 0; j < nr; j++)
560 cmp_int(fp, "InteractionList entry", j, il1->at(i).iatoms.at(j), il2->at(i).iatoms.at(j));
566 static void compareMoltypes(FILE *fp, gmx::ArrayRef<const gmx_moltype_t> mt1, gmx::ArrayRef<const gmx_moltype_t> mt2, real relativeTolerance, real absoluteTolerance)
568 fprintf(fp, "comparing molecule types\n");
569 cmp_int(fp, "moltype size", -1, mt1.size(), mt2.size());
570 for (int i = 0; i < std::min(mt1.ssize(), mt2.ssize()); i++)
572 cmp_str(fp, "Name", i, *mt1[i].name, *mt2[i].name);
573 compareAtoms(fp, &mt1[i].atoms, &mt2[i].atoms, relativeTolerance, absoluteTolerance);
574 compareInteractionLists(fp, &mt1[i].ilist, &mt2[i].ilist);
575 std::string buf = gmx::formatString("cgs[%d]", i);
576 cmp_block(fp, &mt1[i].cgs, &mt2[i].cgs, buf.c_str());
577 buf = gmx::formatString("excls[%d]", i);
578 cmp_blocka(fp, &mt1[i].excls, &mt2[i].excls, buf.c_str());
582 static void compareMoletypeAB(FILE *fp, gmx::ArrayRef<const gmx_moltype_t> mt1, real relativeTolerance, real absoluteTolerance)
584 fprintf(fp, "comparing free energy molecule types\n");
585 for (int i = 0; i < mt1.ssize(); i++)
587 compareAtoms(fp, &mt1[i].atoms, nullptr, relativeTolerance, absoluteTolerance);
590 static void compareMolblocks(FILE *fp, gmx::ArrayRef<const gmx_molblock_t> mb1, gmx::ArrayRef<const gmx_molblock_t> mb2)
592 fprintf(fp, "comparing molecule blocks\n");
593 cmp_int(fp, "molblock size", -1, mb1.size(), mb2.size());
594 int nr = std::min(mb1.size(), mb2.size());
595 for (int i = 0; i < nr; i++)
597 cmp_int(fp, "type", i, mb1[i].type, mb2[i].type);
598 cmp_int(fp, "nmol", i, mb1[i].nmol, mb2[i].nmol);
599 // Only checking size of restraint vectors for now
600 cmp_int(fp, "posres_xA size", i, mb1[i].posres_xA.size(), mb2[i].posres_xA.size());
601 cmp_int(fp, "posres_xB size", i, mb1[i].posres_xB.size(), mb2[i].posres_xB.size());
606 static void compareAtomtypes(FILE *fp, const t_atomtypes &at1, const t_atomtypes &at2)
608 fprintf(fp, "comparing atomtypes\n");
609 cmp_int(fp, "nr", -1, at1.nr, at2.nr);
610 int nr = std::min(at1.nr, at2.nr);
611 for (int i = 0; i < nr; i++)
613 cmp_int(fp, "atomtype", i, at1.atomnumber[i], at2.atomnumber[i]);
617 static void compareIntermolecularExclusions(FILE *fp, gmx::ArrayRef<const int> ime1, gmx::ArrayRef<const int> ime2)
619 fprintf(fp, "comparing intermolecular exclusions\n");
620 cmp_int(fp, "exclusion number", -1, ime1.size(), ime2.size());
621 int nr = std::min(ime1.size(), ime2.size());
622 for (int i = 0; i < nr; i++)
624 cmp_int(fp, "exclusion", i, ime1[i], ime2[i]);
628 static void compareBlockIndices(FILE *fp, gmx::ArrayRef<const MoleculeBlockIndices> mbi1, gmx::ArrayRef<const MoleculeBlockIndices> mbi2)
630 fprintf(fp, "comparing moleculeBlockIndices\n");
631 cmp_int(fp, "size", -1, mbi1.size(), mbi2.size());
632 int nr = std::min(mbi1.size(), mbi2.size());
633 for (int i = 0; i < nr; i++)
635 cmp_int(fp, "numAtomsPerMolecule", i, mbi1[i].numAtomsPerMolecule, mbi2[i].numAtomsPerMolecule);
636 cmp_int(fp, "globalAtomStart", i, mbi1[i].globalAtomStart, mbi2[i].globalAtomStart);
637 cmp_int(fp, "globalAtomEnd", i, mbi1[i].globalAtomEnd, mbi2[i].globalAtomEnd);
638 cmp_int(fp, "globalResidueStart", i, mbi1[i].globalResidueStart, mbi2[i].globalResidueStart);
639 cmp_int(fp, "moleculeIndexStart", i, mbi1[i].moleculeIndexStart, mbi2[i].moleculeIndexStart);
643 void compareMtop(FILE *fp, const gmx_mtop_t &mtop1, const gmx_mtop_t &mtop2, real relativeTolerance, real absoluteTolerance)
645 fprintf(fp, "comparing mtop topology\n");
646 cmp_str(fp, "Name", -1, *mtop1.name, *mtop2.name);
647 cmp_int(fp, "natoms", -1, mtop1.natoms, mtop2.natoms);
648 cmp_int(fp, "maxres_renum", -1, mtop1.maxres_renum, mtop2.maxres_renum);
649 cmp_int(fp, "maxresnr", -1, mtop1.maxresnr, mtop2.maxresnr);
650 cmp_bool(fp, "bIntermolecularInteractions", -1, mtop1.bIntermolecularInteractions, mtop2.bIntermolecularInteractions);
651 cmp_bool(fp, "haveMoleculeIndices", -1, mtop1.haveMoleculeIndices, mtop2.haveMoleculeIndices);
653 compareFfparams(fp, mtop1.ffparams, mtop2.ffparams, relativeTolerance, absoluteTolerance);
654 compareMoltypes(fp, mtop1.moltype, mtop2.moltype, relativeTolerance, absoluteTolerance);
655 compareMolblocks(fp, mtop1.molblock, mtop2.molblock);
656 compareInteractionLists(fp, mtop1.intermolecular_ilist.get(), mtop2.intermolecular_ilist.get());
657 compareAtomtypes(fp, mtop1.atomtypes, mtop2.atomtypes);
658 compareAtomGroups(fp, mtop1.groups, mtop2.groups, mtop1.natoms, mtop2.natoms);
659 compareIntermolecularExclusions(fp, mtop1.intermolecularExclusionGroup, mtop2.intermolecularExclusionGroup);
660 compareBlockIndices(fp, mtop1.moleculeBlockIndices, mtop2.moleculeBlockIndices);
663 void compareMtopAB(FILE *fp, const gmx_mtop_t &mtop1, real relativeTolerance, real absoluteTolerance)
665 fprintf(fp, "comparing topAB\n");
666 compareFfparamAB(fp, mtop1.ffparams, relativeTolerance, absoluteTolerance);
667 compareMoletypeAB(fp, mtop1.moltype, relativeTolerance, absoluteTolerance);
670 void compareAtomGroups(FILE *fp, const gmx_groups_t &g0, const gmx_groups_t &g1,
671 int natoms0, int natoms1)
673 fprintf(fp, "comparing groups\n");
675 for (int i = 0; i < egcNR; i++)
677 std::string buf = gmx::formatString("grps[%d].nr", i);
678 cmp_int(fp, buf.c_str(), -1, g0.grps[i].nr, g1.grps[i].nr);
679 if (g0.grps[i].nr == g1.grps[i].nr)
681 for (int j = 0; j < g0.grps[i].nr; j++)
683 buf = gmx::formatString("grps[%d].name[%d]", i, j);
684 cmp_str(fp, buf.c_str(), -1,
685 *g0.grpname[g0.grps[i].nm_ind[j]],
686 *g1.grpname[g1.grps[i].nm_ind[j]]);
689 cmp_int(fp, "ngrpnr", i, g0.ngrpnr[i], g1.ngrpnr[i]);
690 if (g0.ngrpnr[i] == g1.ngrpnr[i] && natoms0 == natoms1 &&
691 (g0.grpnr[i] != nullptr || g1.grpnr[i] != nullptr))
693 for (int j = 0; j < natoms0; j++)
695 cmp_int(fp, gtypes[i], j, getGroupType(g0, i, j), getGroupType(g1, i, j));
699 /* We have compared the names in the groups lists,
700 * so we can skip the grpname list comparison.
704 int getGroupType(const gmx_groups_t &group, int type, int atom)
706 return (group.grpnr[type] ? group.grpnr[type][atom] : 0);
709 void copy_moltype(const gmx_moltype_t *src, gmx_moltype_t *dst)
711 dst->name = src->name;
712 copy_blocka(&src->excls, &dst->excls);
713 copy_block(&src->cgs, &dst->cgs);
714 t_atoms *atomsCopy = copy_t_atoms(&src->atoms);
715 dst->atoms = *atomsCopy;
718 for (int i = 0; i < F_NRE; ++i)
720 dst->ilist[i] = src->ilist[i];