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, 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.
37 /* This file is completely threadsafe - keep it that way! */
40 #include "broadcaststructs.h"
44 #include "gromacs/compat/make_unique.h"
45 #include "gromacs/gmxlib/network.h"
46 #include "gromacs/math/vec.h"
47 #include "gromacs/mdlib/mdrun.h"
48 #include "gromacs/mdlib/tgroup.h"
49 #include "gromacs/mdtypes/awh-params.h"
50 #include "gromacs/mdtypes/commrec.h"
51 #include "gromacs/mdtypes/inputrec.h"
52 #include "gromacs/mdtypes/md_enums.h"
53 #include "gromacs/mdtypes/pull-params.h"
54 #include "gromacs/mdtypes/state.h"
55 #include "gromacs/topology/mtop_util.h"
56 #include "gromacs/topology/symtab.h"
57 #include "gromacs/topology/topology.h"
58 #include "gromacs/utility/fatalerror.h"
59 #include "gromacs/utility/inmemoryserializer.h"
60 #include "gromacs/utility/keyvaluetree.h"
61 #include "gromacs/utility/keyvaluetreeserializer.h"
62 #include "gromacs/utility/smalloc.h"
64 static void bc_cstring(const t_commrec *cr, char **s)
68 if (MASTER(cr) && *s != nullptr)
70 /* Size of the char buffer is string length + 1 for '\0' */
71 size = strlen(*s) + 1;
80 nblock_bc(cr, size, *s);
82 else if (!MASTER(cr) && *s != nullptr)
89 static void bc_string(const t_commrec *cr, t_symtab *symtab, char ***s)
95 handle = lookup_symtab(symtab, *s);
100 *s = get_symtab_handle(symtab, handle);
104 static void bc_strings(const t_commrec *cr, t_symtab *symtab, int nr, char ****nm)
112 for (i = 0; (i < nr); i++)
114 handle[i] = lookup_symtab(symtab, (*nm)[i]);
117 nblock_bc(cr, nr, handle);
121 snew_bc(cr, *nm, nr);
122 for (i = 0; (i < nr); i++)
124 (*nm)[i] = get_symtab_handle(symtab, handle[i]);
130 static void bc_strings_resinfo(const t_commrec *cr, t_symtab *symtab,
131 int nr, t_resinfo *resinfo)
139 for (i = 0; (i < nr); i++)
141 handle[i] = lookup_symtab(symtab, resinfo[i].name);
144 nblock_bc(cr, nr, handle);
148 for (i = 0; (i < nr); i++)
150 resinfo[i].name = get_symtab_handle(symtab, handle[i]);
156 static void bc_symtab(const t_commrec *cr, t_symtab *symtab)
161 block_bc(cr, symtab->nr);
163 snew_bc(cr, symtab->symbuf, 1);
164 symbuf = symtab->symbuf;
165 symbuf->bufsize = nr;
166 snew_bc(cr, symbuf->buf, nr);
167 for (i = 0; i < nr; i++)
171 len = strlen(symbuf->buf[i]) + 1;
174 snew_bc(cr, symbuf->buf[i], len);
175 nblock_bc(cr, len, symbuf->buf[i]);
179 static void bc_block(const t_commrec *cr, t_block *block)
181 block_bc(cr, block->nr);
182 snew_bc(cr, block->index, block->nr+1);
183 nblock_bc(cr, block->nr+1, block->index);
186 static void bc_blocka(const t_commrec *cr, t_blocka *block)
188 block_bc(cr, block->nr);
189 snew_bc(cr, block->index, block->nr+1);
190 nblock_bc(cr, block->nr+1, block->index);
191 block_bc(cr, block->nra);
194 snew_bc(cr, block->a, block->nra);
195 nblock_bc(cr, block->nra, block->a);
199 static void bc_grps(const t_commrec *cr, t_grps grps[])
203 for (i = 0; (i < egcNR); i++)
205 block_bc(cr, grps[i].nr);
206 snew_bc(cr, grps[i].nm_ind, grps[i].nr);
207 nblock_bc(cr, grps[i].nr, grps[i].nm_ind);
211 static void bc_atoms(const t_commrec *cr, t_symtab *symtab, t_atoms *atoms)
213 block_bc(cr, atoms->nr);
214 snew_bc(cr, atoms->atom, atoms->nr);
215 nblock_bc(cr, atoms->nr, atoms->atom);
216 bc_strings(cr, symtab, atoms->nr, &atoms->atomname);
217 block_bc(cr, atoms->nres);
218 snew_bc(cr, atoms->resinfo, atoms->nres);
219 nblock_bc(cr, atoms->nres, atoms->resinfo);
220 bc_strings_resinfo(cr, symtab, atoms->nres, atoms->resinfo);
221 /* QMMM requires atomtypes to be known on all nodes as well */
222 bc_strings(cr, symtab, atoms->nr, &atoms->atomtype);
223 bc_strings(cr, symtab, atoms->nr, &atoms->atomtypeB);
226 static void bc_groups(const t_commrec *cr, t_symtab *symtab,
227 int natoms, gmx_groups_t *groups)
231 bc_grps(cr, groups->grps);
232 block_bc(cr, groups->ngrpname);
233 bc_strings(cr, symtab, groups->ngrpname, &groups->grpname);
234 for (g = 0; g < egcNR; g++)
238 if (groups->grpnr[g])
250 groups->grpnr[g] = nullptr;
254 snew_bc(cr, groups->grpnr[g], n);
255 nblock_bc(cr, n, groups->grpnr[g]);
260 fprintf(debug, "after bc_groups\n");
264 template <typename AllocatorType>
265 static void bcastPaddedRVecVector(const t_commrec *cr, std::vector<gmx::RVec, AllocatorType> *v, int numAtoms)
267 v->resize(gmx::paddedRVecVectorSize(numAtoms));
268 nblock_bc(cr, numAtoms, as_rvec_array(v->data()));
271 void broadcastStateWithoutDynamics(const t_commrec *cr, t_state *state)
273 GMX_RELEASE_ASSERT(!DOMAINDECOMP(cr), "broadcastStateWithoutDynamics should only be used for special cases without domain decomposition");
280 /* Broadcasts the state sizes and flags from the master to all ranks
281 * in cr->mpi_comm_mygroup.
283 block_bc(cr, state->natoms);
284 block_bc(cr, state->flags);
286 for (int i = 0; i < estNR; i++)
288 if (state->flags & (1 << i))
293 nblock_bc(cr, efptNR, state->lambda.data());
296 block_bc(cr, state->fep_state);
299 block_bc(cr, state->box);
302 bcastPaddedRVecVector(cr, &state->x, state->natoms);
305 GMX_RELEASE_ASSERT(false, "The state has a dynamic entry, while no dynamic entries should be present");
312 static void bc_ilists(const t_commrec *cr, InteractionList *ilist)
316 /* Here we only communicate the non-zero length ilists */
319 for (ftype = 0; ftype < F_NRE; ftype++)
321 if (ilist[ftype].size() > 0)
324 int nr = ilist[ftype].size();
326 nblock_bc(cr, nr, ilist[ftype].iatoms.data());
334 for (ftype = 0; ftype < F_NRE; ftype++)
336 ilist[ftype].iatoms.clear();
345 ilist[ftype].iatoms.resize(nr);
346 nblock_bc(cr, nr, ilist[ftype].iatoms.data());
354 fprintf(debug, "after bc_ilists\n");
358 static void bc_cmap(const t_commrec *cr, gmx_cmap_t *cmap_grid)
362 block_bc(cr, cmap_grid->ngrid);
363 block_bc(cr, cmap_grid->grid_spacing);
365 ngrid = cmap_grid->ngrid;
366 nelem = cmap_grid->grid_spacing * cmap_grid->grid_spacing;
370 snew_bc(cr, cmap_grid->cmapdata, ngrid);
372 for (i = 0; i < ngrid; i++)
374 snew_bc(cr, cmap_grid->cmapdata[i].cmap, 4*nelem);
375 nblock_bc(cr, 4*nelem, cmap_grid->cmapdata[i].cmap);
380 static void bc_ffparams(const t_commrec *cr, gmx_ffparams_t *ffp)
382 block_bc(cr, ffp->ntypes);
383 block_bc(cr, ffp->atnr);
384 snew_bc(cr, ffp->functype, ffp->ntypes);
385 snew_bc(cr, ffp->iparams, ffp->ntypes);
386 nblock_bc(cr, ffp->ntypes, ffp->functype);
387 nblock_bc(cr, ffp->ntypes, ffp->iparams);
388 block_bc(cr, ffp->reppow);
389 block_bc(cr, ffp->fudgeQQ);
390 bc_cmap(cr, &ffp->cmap_grid);
393 static void bc_grpopts(const t_commrec *cr, t_grpopts *g)
397 block_bc(cr, g->ngtc);
398 block_bc(cr, g->ngacc);
399 block_bc(cr, g->ngfrz);
400 block_bc(cr, g->ngener);
401 snew_bc(cr, g->nrdf, g->ngtc);
402 snew_bc(cr, g->tau_t, g->ngtc);
403 snew_bc(cr, g->ref_t, g->ngtc);
404 snew_bc(cr, g->acc, g->ngacc);
405 snew_bc(cr, g->nFreeze, g->ngfrz);
406 snew_bc(cr, g->egp_flags, g->ngener*g->ngener);
408 nblock_bc(cr, g->ngtc, g->nrdf);
409 nblock_bc(cr, g->ngtc, g->tau_t);
410 nblock_bc(cr, g->ngtc, g->ref_t);
411 nblock_bc(cr, g->ngacc, g->acc);
412 nblock_bc(cr, g->ngfrz, g->nFreeze);
413 nblock_bc(cr, g->ngener*g->ngener, g->egp_flags);
414 snew_bc(cr, g->annealing, g->ngtc);
415 snew_bc(cr, g->anneal_npoints, g->ngtc);
416 snew_bc(cr, g->anneal_time, g->ngtc);
417 snew_bc(cr, g->anneal_temp, g->ngtc);
418 nblock_bc(cr, g->ngtc, g->annealing);
419 nblock_bc(cr, g->ngtc, g->anneal_npoints);
420 for (i = 0; (i < g->ngtc); i++)
422 n = g->anneal_npoints[i];
425 snew_bc(cr, g->anneal_time[i], n);
426 snew_bc(cr, g->anneal_temp[i], n);
427 nblock_bc(cr, n, g->anneal_time[i]);
428 nblock_bc(cr, n, g->anneal_temp[i]);
432 /* QMMM stuff, see inputrec */
433 block_bc(cr, g->ngQM);
434 snew_bc(cr, g->QMmethod, g->ngQM);
435 snew_bc(cr, g->QMbasis, g->ngQM);
436 snew_bc(cr, g->QMcharge, g->ngQM);
437 snew_bc(cr, g->QMmult, g->ngQM);
438 snew_bc(cr, g->bSH, g->ngQM);
439 snew_bc(cr, g->CASorbitals, g->ngQM);
440 snew_bc(cr, g->CASelectrons, g->ngQM);
441 snew_bc(cr, g->SAon, g->ngQM);
442 snew_bc(cr, g->SAoff, g->ngQM);
443 snew_bc(cr, g->SAsteps, g->ngQM);
447 nblock_bc(cr, g->ngQM, g->QMmethod);
448 nblock_bc(cr, g->ngQM, g->QMbasis);
449 nblock_bc(cr, g->ngQM, g->QMcharge);
450 nblock_bc(cr, g->ngQM, g->QMmult);
451 nblock_bc(cr, g->ngQM, g->bSH);
452 nblock_bc(cr, g->ngQM, g->CASorbitals);
453 nblock_bc(cr, g->ngQM, g->CASelectrons);
454 nblock_bc(cr, g->ngQM, g->SAon);
455 nblock_bc(cr, g->ngQM, g->SAoff);
456 nblock_bc(cr, g->ngQM, g->SAsteps);
457 /* end of QMMM stuff */
461 static void bc_awhBias(const t_commrec *cr, gmx::AwhBiasParams *awhBiasParams)
463 block_bc(cr, *awhBiasParams);
465 snew_bc(cr, awhBiasParams->dimParams, awhBiasParams->ndim);
466 nblock_bc(cr, awhBiasParams->ndim, awhBiasParams->dimParams);
469 static void bc_awh(const t_commrec *cr, gmx::AwhParams *awhParams)
473 block_bc(cr, *awhParams);
474 snew_bc(cr, awhParams->awhBiasParams, awhParams->numBias);
475 for (k = 0; k < awhParams->numBias; k++)
477 bc_awhBias(cr, &awhParams->awhBiasParams[k]);
481 static void bc_pull_group(const t_commrec *cr, t_pull_group *pgrp)
486 snew_bc(cr, pgrp->ind, pgrp->nat);
487 nblock_bc(cr, pgrp->nat, pgrp->ind);
489 if (pgrp->nweight > 0)
491 snew_bc(cr, pgrp->weight, pgrp->nweight);
492 nblock_bc(cr, pgrp->nweight, pgrp->weight);
496 static void bc_pull(const t_commrec *cr, pull_params_t *pull)
501 snew_bc(cr, pull->group, pull->ngroup);
502 for (g = 0; g < pull->ngroup; g++)
504 bc_pull_group(cr, &pull->group[g]);
506 snew_bc(cr, pull->coord, pull->ncoord);
507 nblock_bc(cr, pull->ncoord, pull->coord);
508 for (int c = 0; c < pull->ncoord; c++)
512 pull->coord[c].externalPotentialProvider = nullptr;
514 if (pull->coord[c].eType == epullEXTERNAL)
516 bc_cstring(cr, &pull->coord[c].externalPotentialProvider);
521 static void bc_rotgrp(const t_commrec *cr, t_rotgrp *rotg)
526 snew_bc(cr, rotg->ind, rotg->nat);
527 nblock_bc(cr, rotg->nat, rotg->ind);
528 snew_bc(cr, rotg->x_ref, rotg->nat);
529 nblock_bc(cr, rotg->nat, rotg->x_ref);
533 static void bc_rot(const t_commrec *cr, t_rot *rot)
538 snew_bc(cr, rot->grp, rot->ngrp);
539 for (g = 0; g < rot->ngrp; g++)
541 bc_rotgrp(cr, &rot->grp[g]);
545 static void bc_imd(const t_commrec *cr, t_IMD *imd)
548 snew_bc(cr, imd->ind, imd->nat);
549 nblock_bc(cr, imd->nat, imd->ind);
552 static void bc_fepvals(const t_commrec *cr, t_lambda *fep)
556 block_bc(cr, fep->nstdhdl);
557 block_bc(cr, fep->init_lambda);
558 block_bc(cr, fep->init_fep_state);
559 block_bc(cr, fep->delta_lambda);
560 block_bc(cr, fep->edHdLPrintEnergy);
561 block_bc(cr, fep->n_lambda);
562 if (fep->n_lambda > 0)
564 snew_bc(cr, fep->all_lambda, efptNR);
565 nblock_bc(cr, efptNR, fep->all_lambda);
566 for (i = 0; i < efptNR; i++)
568 snew_bc(cr, fep->all_lambda[i], fep->n_lambda);
569 nblock_bc(cr, fep->n_lambda, fep->all_lambda[i]);
572 block_bc(cr, fep->sc_alpha);
573 block_bc(cr, fep->sc_power);
574 block_bc(cr, fep->sc_r_power);
575 block_bc(cr, fep->sc_sigma);
576 block_bc(cr, fep->sc_sigma_min);
577 block_bc(cr, fep->bScCoul);
578 nblock_bc(cr, efptNR, &(fep->separate_dvdl[0]));
579 block_bc(cr, fep->dhdl_derivatives);
580 block_bc(cr, fep->dh_hist_size);
581 block_bc(cr, fep->dh_hist_spacing);
584 fprintf(debug, "after bc_fepvals\n");
588 static void bc_expandedvals(const t_commrec *cr, t_expanded *expand, int n_lambda)
590 block_bc(cr, expand->nstexpanded);
591 block_bc(cr, expand->elamstats);
592 block_bc(cr, expand->elmcmove);
593 block_bc(cr, expand->elmceq);
594 block_bc(cr, expand->equil_n_at_lam);
595 block_bc(cr, expand->equil_wl_delta);
596 block_bc(cr, expand->equil_ratio);
597 block_bc(cr, expand->equil_steps);
598 block_bc(cr, expand->equil_samples);
599 block_bc(cr, expand->lmc_seed);
600 block_bc(cr, expand->minvar);
601 block_bc(cr, expand->minvar_const);
602 block_bc(cr, expand->c_range);
603 block_bc(cr, expand->bSymmetrizedTMatrix);
604 block_bc(cr, expand->nstTij);
605 block_bc(cr, expand->lmc_repeats);
606 block_bc(cr, expand->lmc_forced_nstart);
607 block_bc(cr, expand->gibbsdeltalam);
608 block_bc(cr, expand->wl_scale);
609 block_bc(cr, expand->wl_ratio);
610 block_bc(cr, expand->init_wl_delta);
611 block_bc(cr, expand->bInit_weights);
612 snew_bc(cr, expand->init_lambda_weights, n_lambda);
613 nblock_bc(cr, n_lambda, expand->init_lambda_weights);
614 block_bc(cr, expand->mc_temp);
617 fprintf(debug, "after bc_expandedvals\n");
621 static void bc_simtempvals(const t_commrec *cr, t_simtemp *simtemp, int n_lambda)
623 block_bc(cr, simtemp->simtemp_low);
624 block_bc(cr, simtemp->simtemp_high);
625 block_bc(cr, simtemp->eSimTempScale);
626 snew_bc(cr, simtemp->temperatures, n_lambda);
627 nblock_bc(cr, n_lambda, simtemp->temperatures);
630 fprintf(debug, "after bc_simtempvals\n");
635 static void bc_swapions(const t_commrec *cr, t_swapcoords *swap)
639 /* Broadcast atom indices for split groups, solvent group, and for all user-defined swap groups */
640 snew_bc(cr, swap->grp, swap->ngrp);
641 for (int i = 0; i < swap->ngrp; i++)
643 t_swapGroup *g = &swap->grp[i];
646 snew_bc(cr, g->ind, g->nat);
647 nblock_bc(cr, g->nat, g->ind);
652 len = strlen(g->molname);
655 snew_bc(cr, g->molname, len);
656 nblock_bc(cr, len, g->molname);
661 static void bc_inputrec(const t_commrec *cr, t_inputrec *inputrec)
663 // Note that this overwrites pointers in inputrec, so all pointer fields
664 // Must be initialized separately below.
665 block_bc(cr, *inputrec);
668 gmx::InMemorySerializer serializer;
669 gmx::serializeKeyValueTree(*inputrec->params, &serializer);
670 std::vector<char> buffer = serializer.finishAndGetBuffer();
671 size_t size = buffer.size();
673 nblock_bc(cr, size, buffer.data());
677 // block_bc() above overwrites the old pointer, so set it to a
678 // reasonable value in case code below throws.
679 inputrec->params = nullptr;
680 std::vector<char> buffer;
683 nblock_abc(cr, size, &buffer);
684 gmx::InMemoryDeserializer serializer(buffer);
685 inputrec->params = new gmx::KeyValueTreeObject(
686 gmx::deserializeKeyValueTree(&serializer));
689 bc_grpopts(cr, &(inputrec->opts));
691 /* even if efep is efepNO, we need to initialize to make sure that
692 * n_lambda is set to zero */
694 snew_bc(cr, inputrec->fepvals, 1);
695 if (inputrec->efep != efepNO || inputrec->bSimTemp)
697 bc_fepvals(cr, inputrec->fepvals);
699 /* need to initialize this as well because of data checked for in the logic */
700 snew_bc(cr, inputrec->expandedvals, 1);
701 if (inputrec->bExpanded)
703 bc_expandedvals(cr, inputrec->expandedvals, inputrec->fepvals->n_lambda);
705 snew_bc(cr, inputrec->simtempvals, 1);
706 if (inputrec->bSimTemp)
708 bc_simtempvals(cr, inputrec->simtempvals, inputrec->fepvals->n_lambda);
712 snew_bc(cr, inputrec->pull, 1);
713 bc_pull(cr, inputrec->pull);
715 if (inputrec->bDoAwh)
717 snew_bc(cr, inputrec->awhParams, 1);
718 bc_awh(cr, inputrec->awhParams);
723 snew_bc(cr, inputrec->rot, 1);
724 bc_rot(cr, inputrec->rot);
728 snew_bc(cr, inputrec->imd, 1);
729 bc_imd(cr, inputrec->imd);
731 if (inputrec->eSwapCoords != eswapNO)
733 snew_bc(cr, inputrec->swap, 1);
734 bc_swapions(cr, inputrec->swap);
738 static void bc_moltype(const t_commrec *cr, t_symtab *symtab,
739 gmx_moltype_t *moltype)
741 bc_string(cr, symtab, &moltype->name);
742 bc_atoms(cr, symtab, &moltype->atoms);
745 fprintf(debug, "after bc_atoms\n");
748 bc_ilists(cr, moltype->ilist);
749 bc_block(cr, &moltype->cgs);
750 bc_blocka(cr, &moltype->excls);
753 static void bc_vector_of_rvec(const t_commrec *cr, std::vector<gmx::RVec> *vec)
755 int numElements = vec->size();
756 block_bc(cr, numElements);
759 vec->resize(numElements);
763 nblock_bc(cr, numElements, as_rvec_array(vec->data()));
767 static void bc_molblock(const t_commrec *cr, gmx_molblock_t *molb)
769 block_bc(cr, molb->type);
770 block_bc(cr, molb->nmol);
771 bc_vector_of_rvec(cr, &molb->posres_xA);
772 bc_vector_of_rvec(cr, &molb->posres_xB);
775 fprintf(debug, "after bc_molblock\n");
779 static void bc_atomtypes(const t_commrec *cr, t_atomtypes *atomtypes)
781 block_bc(cr, atomtypes->nr);
784 /*! \brief Broadcasts ir and mtop from the master to all nodes in
785 * cr->mpi_comm_mygroup. */
787 void bcast_ir_mtop(const t_commrec *cr, t_inputrec *inputrec, gmx_mtop_t *mtop)
791 fprintf(debug, "in bc_data\n");
793 bc_inputrec(cr, inputrec);
796 fprintf(debug, "after bc_inputrec\n");
798 bc_symtab(cr, &mtop->symtab);
801 fprintf(debug, "after bc_symtab\n");
803 bc_string(cr, &mtop->symtab, &mtop->name);
806 fprintf(debug, "after bc_name\n");
809 bc_ffparams(cr, &mtop->ffparams);
811 int nmoltype = mtop->moltype.size();
812 block_bc(cr, nmoltype);
813 mtop->moltype.resize(nmoltype);
814 for (gmx_moltype_t &moltype : mtop->moltype)
816 bc_moltype(cr, &mtop->symtab, &moltype);
819 block_bc(cr, mtop->bIntermolecularInteractions);
820 if (mtop->bIntermolecularInteractions)
822 mtop->intermolecular_ilist = gmx::compat::make_unique<InteractionLists>();
823 bc_ilists(cr, mtop->intermolecular_ilist->data());
826 int nmolblock = mtop->molblock.size();
827 block_bc(cr, nmolblock);
828 mtop->molblock.resize(nmolblock);
829 for (gmx_molblock_t &molblock : mtop->molblock)
831 bc_molblock(cr, &molblock);
834 block_bc(cr, mtop->natoms);
836 bc_atomtypes(cr, &mtop->atomtypes);
838 bc_groups(cr, &mtop->symtab, mtop->natoms, &mtop->groups);
840 GMX_RELEASE_ASSERT(!MASTER(cr) || mtop->haveMoleculeIndices, "mtop should have valid molecule indices");
843 mtop->haveMoleculeIndices = true;
845 gmx_mtop_finalize(mtop);
849 void init_parallel(t_commrec *cr, t_inputrec *inputrec,
852 bcast_ir_mtop(cr, inputrec, mtop);