Tidy: modernize-use-nullptr
[gromacs.git] / src / gromacs / mdlib / broadcaststructs.cpp
blobc2a09715f5867a5b325a2e849d49625f406654ac
1 /*
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, 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! */
38 #include "gmxpre.h"
40 #include "broadcaststructs.h"
42 #include <string.h>
44 #include "gromacs/gmxlib/network.h"
45 #include "gromacs/math/vec.h"
46 #include "gromacs/mdlib/mdrun.h"
47 #include "gromacs/mdlib/tgroup.h"
48 #include "gromacs/mdtypes/commrec.h"
49 #include "gromacs/mdtypes/inputrec.h"
50 #include "gromacs/mdtypes/md_enums.h"
51 #include "gromacs/mdtypes/pull-params.h"
52 #include "gromacs/mdtypes/state.h"
53 #include "gromacs/topology/symtab.h"
54 #include "gromacs/topology/topology.h"
55 #include "gromacs/utility/fatalerror.h"
56 #include "gromacs/utility/smalloc.h"
58 static void bc_cstring(const t_commrec *cr, char **s)
60 int size = 0;
62 if (MASTER(cr) && *s != nullptr)
64 /* Size of the char buffer is string length + 1 for '\0' */
65 size = strlen(*s) + 1;
67 block_bc(cr, size);
68 if (size > 0)
70 if (!MASTER(cr))
72 srenew(*s, size);
74 nblock_bc(cr, size, *s);
76 else if (!MASTER(cr) && *s != nullptr)
78 sfree(*s);
79 *s = nullptr;
83 static void bc_string(const t_commrec *cr, t_symtab *symtab, char ***s)
85 int handle;
87 if (MASTER(cr))
89 handle = lookup_symtab(symtab, *s);
91 block_bc(cr, handle);
92 if (!MASTER(cr))
94 *s = get_symtab_handle(symtab, handle);
98 static void bc_strings(const t_commrec *cr, t_symtab *symtab, int nr, char ****nm)
100 int i;
101 int *handle;
103 snew(handle, nr);
104 if (MASTER(cr))
106 for (i = 0; (i < nr); i++)
108 handle[i] = lookup_symtab(symtab, (*nm)[i]);
111 nblock_bc(cr, nr, handle);
113 if (!MASTER(cr))
115 snew_bc(cr, *nm, nr);
116 for (i = 0; (i < nr); i++)
118 (*nm)[i] = get_symtab_handle(symtab, handle[i]);
121 sfree(handle);
124 static void bc_strings_resinfo(const t_commrec *cr, t_symtab *symtab,
125 int nr, t_resinfo *resinfo)
127 int i;
128 int *handle;
130 snew(handle, nr);
131 if (MASTER(cr))
133 for (i = 0; (i < nr); i++)
135 handle[i] = lookup_symtab(symtab, resinfo[i].name);
138 nblock_bc(cr, nr, handle);
140 if (!MASTER(cr))
142 for (i = 0; (i < nr); i++)
144 resinfo[i].name = get_symtab_handle(symtab, handle[i]);
147 sfree(handle);
150 static void bc_symtab(const t_commrec *cr, t_symtab *symtab)
152 int i, nr, len;
153 t_symbuf *symbuf;
155 block_bc(cr, symtab->nr);
156 nr = symtab->nr;
157 snew_bc(cr, symtab->symbuf, 1);
158 symbuf = symtab->symbuf;
159 symbuf->bufsize = nr;
160 snew_bc(cr, symbuf->buf, nr);
161 for (i = 0; i < nr; i++)
163 if (MASTER(cr))
165 len = strlen(symbuf->buf[i]) + 1;
167 block_bc(cr, len);
168 snew_bc(cr, symbuf->buf[i], len);
169 nblock_bc(cr, len, symbuf->buf[i]);
173 static void bc_block(const t_commrec *cr, t_block *block)
175 block_bc(cr, block->nr);
176 snew_bc(cr, block->index, block->nr+1);
177 nblock_bc(cr, block->nr+1, block->index);
180 static void bc_blocka(const t_commrec *cr, t_blocka *block)
182 block_bc(cr, block->nr);
183 snew_bc(cr, block->index, block->nr+1);
184 nblock_bc(cr, block->nr+1, block->index);
185 block_bc(cr, block->nra);
186 if (block->nra)
188 snew_bc(cr, block->a, block->nra);
189 nblock_bc(cr, block->nra, block->a);
193 static void bc_grps(const t_commrec *cr, t_grps grps[])
195 int i;
197 for (i = 0; (i < egcNR); i++)
199 block_bc(cr, grps[i].nr);
200 snew_bc(cr, grps[i].nm_ind, grps[i].nr);
201 nblock_bc(cr, grps[i].nr, grps[i].nm_ind);
205 static void bc_atoms(const t_commrec *cr, t_symtab *symtab, t_atoms *atoms)
207 block_bc(cr, atoms->nr);
208 snew_bc(cr, atoms->atom, atoms->nr);
209 nblock_bc(cr, atoms->nr, atoms->atom);
210 bc_strings(cr, symtab, atoms->nr, &atoms->atomname);
211 block_bc(cr, atoms->nres);
212 snew_bc(cr, atoms->resinfo, atoms->nres);
213 nblock_bc(cr, atoms->nres, atoms->resinfo);
214 bc_strings_resinfo(cr, symtab, atoms->nres, atoms->resinfo);
215 /* QMMM requires atomtypes to be known on all nodes as well */
216 bc_strings(cr, symtab, atoms->nr, &atoms->atomtype);
217 bc_strings(cr, symtab, atoms->nr, &atoms->atomtypeB);
220 static void bc_groups(const t_commrec *cr, t_symtab *symtab,
221 int natoms, gmx_groups_t *groups)
223 int g, n;
225 bc_grps(cr, groups->grps);
226 block_bc(cr, groups->ngrpname);
227 bc_strings(cr, symtab, groups->ngrpname, &groups->grpname);
228 for (g = 0; g < egcNR; g++)
230 if (MASTER(cr))
232 if (groups->grpnr[g])
234 n = natoms;
236 else
238 n = 0;
241 block_bc(cr, n);
242 if (n == 0)
244 groups->grpnr[g] = nullptr;
246 else
248 snew_bc(cr, groups->grpnr[g], n);
249 nblock_bc(cr, n, groups->grpnr[g]);
252 if (debug)
254 fprintf(debug, "after bc_groups\n");
258 static void bcastPaddedRVecVector(const t_commrec *cr, PaddedRVecVector *v, unsigned int n)
260 (*v).resize(n + 1);
261 nblock_bc(cr, n, as_rvec_array(v->data()));
264 void bcast_state(const t_commrec *cr, t_state *state)
266 int i, nnht, nnhtp;
268 if (!PAR(cr) || (cr->nnodes - cr->npmenodes <= 1))
270 return;
273 /* Broadcasts the state sizes and flags from the master to all nodes
274 * in cr->mpi_comm_mygroup. The arrays are not broadcasted. */
275 block_bc(cr, state->natoms);
276 block_bc(cr, state->ngtc);
277 block_bc(cr, state->nnhpres);
278 block_bc(cr, state->nhchainlength);
279 block_bc(cr, state->flags);
280 state->lambda.resize(efptNR);
282 if (cr->dd)
284 /* We allocate dynamically in dd_partition_system. */
285 return;
287 /* The code below is reachable only by TPI and NM, so it is not
288 tested by anything. */
290 nnht = (state->ngtc)*(state->nhchainlength);
291 nnhtp = (state->nnhpres)*(state->nhchainlength);
293 for (i = 0; i < estNR; i++)
295 if (state->flags & (1<<i))
297 switch (i)
299 case estLAMBDA: nblock_bc(cr, efptNR, state->lambda.data()); break;
300 case estFEPSTATE: block_bc(cr, state->fep_state); break;
301 case estBOX: block_bc(cr, state->box); break;
302 case estBOX_REL: block_bc(cr, state->box_rel); break;
303 case estBOXV: block_bc(cr, state->boxv); break;
304 case estPRES_PREV: block_bc(cr, state->pres_prev); break;
305 case estSVIR_PREV: block_bc(cr, state->svir_prev); break;
306 case estFVIR_PREV: block_bc(cr, state->fvir_prev); break;
307 case estNH_XI: nblock_abc(cr, nnht, &state->nosehoover_xi); break;
308 case estNH_VXI: nblock_abc(cr, nnht, &state->nosehoover_vxi); break;
309 case estNHPRES_XI: nblock_abc(cr, nnhtp, &state->nhpres_xi); break;
310 case estNHPRES_VXI: nblock_abc(cr, nnhtp, &state->nhpres_vxi); break;
311 case estTC_INT: nblock_abc(cr, state->ngtc, &state->therm_integral); break;
312 case estVETA: block_bc(cr, state->veta); break;
313 case estVOL0: block_bc(cr, state->vol0); break;
314 case estX: bcastPaddedRVecVector(cr, &state->x, state->natoms);
315 case estV: bcastPaddedRVecVector(cr, &state->v, state->natoms);
316 case estCGP: bcastPaddedRVecVector(cr, &state->cg_p, state->natoms);
317 case estDISRE_INITF: block_bc(cr, state->hist.disre_initf); break;
318 case estDISRE_RM3TAV:
319 block_bc(cr, state->hist.ndisrepairs);
320 nblock_abc(cr, state->hist.ndisrepairs, &state->hist.disre_rm3tav);
321 break;
322 case estORIRE_INITF: block_bc(cr, state->hist.orire_initf); break;
323 case estORIRE_DTAV:
324 block_bc(cr, state->hist.norire_Dtav);
325 nblock_abc(cr, state->hist.norire_Dtav, &state->hist.orire_Dtav);
326 break;
327 default:
328 gmx_fatal(FARGS,
329 "Communication is not implemented for %s in bcast_state",
330 est_names[i]);
336 static void bc_ilists(const t_commrec *cr, t_ilist *ilist)
338 int ftype;
340 /* Here we only communicate the non-zero length ilists */
341 if (MASTER(cr))
343 for (ftype = 0; ftype < F_NRE; ftype++)
345 if (ilist[ftype].nr > 0)
347 block_bc(cr, ftype);
348 block_bc(cr, ilist[ftype].nr);
349 nblock_bc(cr, ilist[ftype].nr, ilist[ftype].iatoms);
352 ftype = -1;
353 block_bc(cr, ftype);
355 else
357 for (ftype = 0; ftype < F_NRE; ftype++)
359 ilist[ftype].nr = 0;
363 block_bc(cr, ftype);
364 if (ftype >= 0)
366 block_bc(cr, ilist[ftype].nr);
367 snew_bc(cr, ilist[ftype].iatoms, ilist[ftype].nr);
368 nblock_bc(cr, ilist[ftype].nr, ilist[ftype].iatoms);
371 while (ftype >= 0);
374 if (debug)
376 fprintf(debug, "after bc_ilists\n");
380 static void bc_cmap(const t_commrec *cr, gmx_cmap_t *cmap_grid)
382 int i, nelem, ngrid;
384 block_bc(cr, cmap_grid->ngrid);
385 block_bc(cr, cmap_grid->grid_spacing);
387 ngrid = cmap_grid->ngrid;
388 nelem = cmap_grid->grid_spacing * cmap_grid->grid_spacing;
390 if (ngrid > 0)
392 snew_bc(cr, cmap_grid->cmapdata, ngrid);
394 for (i = 0; i < ngrid; i++)
396 snew_bc(cr, cmap_grid->cmapdata[i].cmap, 4*nelem);
397 nblock_bc(cr, 4*nelem, cmap_grid->cmapdata[i].cmap);
402 static void bc_ffparams(const t_commrec *cr, gmx_ffparams_t *ffp)
404 block_bc(cr, ffp->ntypes);
405 block_bc(cr, ffp->atnr);
406 snew_bc(cr, ffp->functype, ffp->ntypes);
407 snew_bc(cr, ffp->iparams, ffp->ntypes);
408 nblock_bc(cr, ffp->ntypes, ffp->functype);
409 nblock_bc(cr, ffp->ntypes, ffp->iparams);
410 block_bc(cr, ffp->reppow);
411 block_bc(cr, ffp->fudgeQQ);
412 bc_cmap(cr, &ffp->cmap_grid);
415 static void bc_grpopts(const t_commrec *cr, t_grpopts *g)
417 int i, n;
419 block_bc(cr, g->ngtc);
420 block_bc(cr, g->ngacc);
421 block_bc(cr, g->ngfrz);
422 block_bc(cr, g->ngener);
423 snew_bc(cr, g->nrdf, g->ngtc);
424 snew_bc(cr, g->tau_t, g->ngtc);
425 snew_bc(cr, g->ref_t, g->ngtc);
426 snew_bc(cr, g->acc, g->ngacc);
427 snew_bc(cr, g->nFreeze, g->ngfrz);
428 snew_bc(cr, g->egp_flags, g->ngener*g->ngener);
430 nblock_bc(cr, g->ngtc, g->nrdf);
431 nblock_bc(cr, g->ngtc, g->tau_t);
432 nblock_bc(cr, g->ngtc, g->ref_t);
433 nblock_bc(cr, g->ngacc, g->acc);
434 nblock_bc(cr, g->ngfrz, g->nFreeze);
435 nblock_bc(cr, g->ngener*g->ngener, g->egp_flags);
436 snew_bc(cr, g->annealing, g->ngtc);
437 snew_bc(cr, g->anneal_npoints, g->ngtc);
438 snew_bc(cr, g->anneal_time, g->ngtc);
439 snew_bc(cr, g->anneal_temp, g->ngtc);
440 nblock_bc(cr, g->ngtc, g->annealing);
441 nblock_bc(cr, g->ngtc, g->anneal_npoints);
442 for (i = 0; (i < g->ngtc); i++)
444 n = g->anneal_npoints[i];
445 if (n > 0)
447 snew_bc(cr, g->anneal_time[i], n);
448 snew_bc(cr, g->anneal_temp[i], n);
449 nblock_bc(cr, n, g->anneal_time[i]);
450 nblock_bc(cr, n, g->anneal_temp[i]);
454 /* QMMM stuff, see inputrec */
455 block_bc(cr, g->ngQM);
456 snew_bc(cr, g->QMmethod, g->ngQM);
457 snew_bc(cr, g->QMbasis, g->ngQM);
458 snew_bc(cr, g->QMcharge, g->ngQM);
459 snew_bc(cr, g->QMmult, g->ngQM);
460 snew_bc(cr, g->bSH, g->ngQM);
461 snew_bc(cr, g->CASorbitals, g->ngQM);
462 snew_bc(cr, g->CASelectrons, g->ngQM);
463 snew_bc(cr, g->SAon, g->ngQM);
464 snew_bc(cr, g->SAoff, g->ngQM);
465 snew_bc(cr, g->SAsteps, g->ngQM);
467 if (g->ngQM)
469 nblock_bc(cr, g->ngQM, g->QMmethod);
470 nblock_bc(cr, g->ngQM, g->QMbasis);
471 nblock_bc(cr, g->ngQM, g->QMcharge);
472 nblock_bc(cr, g->ngQM, g->QMmult);
473 nblock_bc(cr, g->ngQM, g->bSH);
474 nblock_bc(cr, g->ngQM, g->CASorbitals);
475 nblock_bc(cr, g->ngQM, g->CASelectrons);
476 nblock_bc(cr, g->ngQM, g->SAon);
477 nblock_bc(cr, g->ngQM, g->SAoff);
478 nblock_bc(cr, g->ngQM, g->SAsteps);
479 /* end of QMMM stuff */
483 static void bc_pull_group(const t_commrec *cr, t_pull_group *pgrp)
485 block_bc(cr, *pgrp);
486 if (pgrp->nat > 0)
488 snew_bc(cr, pgrp->ind, pgrp->nat);
489 nblock_bc(cr, pgrp->nat, pgrp->ind);
491 if (pgrp->nweight > 0)
493 snew_bc(cr, pgrp->weight, pgrp->nweight);
494 nblock_bc(cr, pgrp->nweight, pgrp->weight);
498 static void bc_pull(const t_commrec *cr, pull_params_t *pull)
500 int g;
502 block_bc(cr, *pull);
503 snew_bc(cr, pull->group, pull->ngroup);
504 for (g = 0; g < pull->ngroup; g++)
506 bc_pull_group(cr, &pull->group[g]);
508 snew_bc(cr, pull->coord, pull->ncoord);
509 nblock_bc(cr, pull->ncoord, pull->coord);
510 for (int c = 0; c < pull->ncoord; c++)
512 if (!MASTER(cr))
514 pull->coord[c].externalPotentialProvider = nullptr;
516 if (pull->coord[c].eType == epullEXTERNAL)
518 bc_cstring(cr, &pull->coord[c].externalPotentialProvider);
523 static void bc_rotgrp(const t_commrec *cr, t_rotgrp *rotg)
525 block_bc(cr, *rotg);
526 if (rotg->nat > 0)
528 snew_bc(cr, rotg->ind, rotg->nat);
529 nblock_bc(cr, rotg->nat, rotg->ind);
530 snew_bc(cr, rotg->x_ref, rotg->nat);
531 nblock_bc(cr, rotg->nat, rotg->x_ref);
535 static void bc_rot(const t_commrec *cr, t_rot *rot)
537 int g;
539 block_bc(cr, *rot);
540 snew_bc(cr, rot->grp, rot->ngrp);
541 for (g = 0; g < rot->ngrp; g++)
543 bc_rotgrp(cr, &rot->grp[g]);
547 static void bc_imd(const t_commrec *cr, t_IMD *imd)
549 block_bc(cr, *imd);
550 snew_bc(cr, imd->ind, imd->nat);
551 nblock_bc(cr, imd->nat, imd->ind);
554 static void bc_fepvals(const t_commrec *cr, t_lambda *fep)
556 int i;
558 block_bc(cr, fep->nstdhdl);
559 block_bc(cr, fep->init_lambda);
560 block_bc(cr, fep->init_fep_state);
561 block_bc(cr, fep->delta_lambda);
562 block_bc(cr, fep->edHdLPrintEnergy);
563 block_bc(cr, fep->n_lambda);
564 if (fep->n_lambda > 0)
566 snew_bc(cr, fep->all_lambda, efptNR);
567 nblock_bc(cr, efptNR, fep->all_lambda);
568 for (i = 0; i < efptNR; i++)
570 snew_bc(cr, fep->all_lambda[i], fep->n_lambda);
571 nblock_bc(cr, fep->n_lambda, fep->all_lambda[i]);
574 block_bc(cr, fep->sc_alpha);
575 block_bc(cr, fep->sc_power);
576 block_bc(cr, fep->sc_r_power);
577 block_bc(cr, fep->sc_sigma);
578 block_bc(cr, fep->sc_sigma_min);
579 block_bc(cr, fep->bScCoul);
580 nblock_bc(cr, efptNR, &(fep->separate_dvdl[0]));
581 block_bc(cr, fep->dhdl_derivatives);
582 block_bc(cr, fep->dh_hist_size);
583 block_bc(cr, fep->dh_hist_spacing);
584 if (debug)
586 fprintf(debug, "after bc_fepvals\n");
590 static void bc_expandedvals(const t_commrec *cr, t_expanded *expand, int n_lambda)
592 block_bc(cr, expand->nstexpanded);
593 block_bc(cr, expand->elamstats);
594 block_bc(cr, expand->elmcmove);
595 block_bc(cr, expand->elmceq);
596 block_bc(cr, expand->equil_n_at_lam);
597 block_bc(cr, expand->equil_wl_delta);
598 block_bc(cr, expand->equil_ratio);
599 block_bc(cr, expand->equil_steps);
600 block_bc(cr, expand->equil_samples);
601 block_bc(cr, expand->lmc_seed);
602 block_bc(cr, expand->minvar);
603 block_bc(cr, expand->minvar_const);
604 block_bc(cr, expand->c_range);
605 block_bc(cr, expand->bSymmetrizedTMatrix);
606 block_bc(cr, expand->nstTij);
607 block_bc(cr, expand->lmc_repeats);
608 block_bc(cr, expand->lmc_forced_nstart);
609 block_bc(cr, expand->gibbsdeltalam);
610 block_bc(cr, expand->wl_scale);
611 block_bc(cr, expand->wl_ratio);
612 block_bc(cr, expand->init_wl_delta);
613 block_bc(cr, expand->bInit_weights);
614 snew_bc(cr, expand->init_lambda_weights, n_lambda);
615 nblock_bc(cr, n_lambda, expand->init_lambda_weights);
616 block_bc(cr, expand->mc_temp);
617 if (debug)
619 fprintf(debug, "after bc_expandedvals\n");
623 static void bc_simtempvals(const t_commrec *cr, t_simtemp *simtemp, int n_lambda)
625 block_bc(cr, simtemp->simtemp_low);
626 block_bc(cr, simtemp->simtemp_high);
627 block_bc(cr, simtemp->eSimTempScale);
628 snew_bc(cr, simtemp->temperatures, n_lambda);
629 nblock_bc(cr, n_lambda, simtemp->temperatures);
630 if (debug)
632 fprintf(debug, "after bc_simtempvals\n");
637 static void bc_swapions(const t_commrec *cr, t_swapcoords *swap)
639 block_bc(cr, *swap);
641 /* Broadcast atom indices for split groups, solvent group, and for all user-defined swap groups */
642 snew_bc(cr, swap->grp, swap->ngrp);
643 for (int i = 0; i < swap->ngrp; i++)
645 t_swapGroup *g = &swap->grp[i];
647 block_bc(cr, *g);
648 snew_bc(cr, g->ind, g->nat);
649 nblock_bc(cr, g->nat, g->ind);
651 int len = 0;
652 if (MASTER(cr))
654 len = strlen(g->molname);
656 block_bc(cr, len);
657 snew_bc(cr, g->molname, len);
658 nblock_bc(cr, len, g->molname);
663 static void bc_inputrec(const t_commrec *cr, t_inputrec *inputrec)
665 /* The statement below is dangerous. It overwrites all structures in inputrec.
666 * If something is added to inputrec, like efield it will need to be
667 * treated here.
669 gmx::IInputRecExtension *eptr = inputrec->efield;
670 block_bc(cr, *inputrec);
671 inputrec->efield = eptr;
673 bc_grpopts(cr, &(inputrec->opts));
675 /* even if efep is efepNO, we need to initialize to make sure that
676 * n_lambda is set to zero */
678 snew_bc(cr, inputrec->fepvals, 1);
679 if (inputrec->efep != efepNO || inputrec->bSimTemp)
681 bc_fepvals(cr, inputrec->fepvals);
683 /* need to initialize this as well because of data checked for in the logic */
684 snew_bc(cr, inputrec->expandedvals, 1);
685 if (inputrec->bExpanded)
687 bc_expandedvals(cr, inputrec->expandedvals, inputrec->fepvals->n_lambda);
689 snew_bc(cr, inputrec->simtempvals, 1);
690 if (inputrec->bSimTemp)
692 bc_simtempvals(cr, inputrec->simtempvals, inputrec->fepvals->n_lambda);
694 if (inputrec->bPull)
696 snew_bc(cr, inputrec->pull, 1);
697 bc_pull(cr, inputrec->pull);
699 if (inputrec->bRot)
701 snew_bc(cr, inputrec->rot, 1);
702 bc_rot(cr, inputrec->rot);
704 if (inputrec->bIMD)
706 snew_bc(cr, inputrec->imd, 1);
707 bc_imd(cr, inputrec->imd);
709 inputrec->efield->broadCast(cr);
710 if (inputrec->eSwapCoords != eswapNO)
712 snew_bc(cr, inputrec->swap, 1);
713 bc_swapions(cr, inputrec->swap);
717 static void bc_moltype(const t_commrec *cr, t_symtab *symtab,
718 gmx_moltype_t *moltype)
720 bc_string(cr, symtab, &moltype->name);
721 bc_atoms(cr, symtab, &moltype->atoms);
722 if (debug)
724 fprintf(debug, "after bc_atoms\n");
727 bc_ilists(cr, moltype->ilist);
728 bc_block(cr, &moltype->cgs);
729 bc_blocka(cr, &moltype->excls);
732 static void bc_molblock(const t_commrec *cr, gmx_molblock_t *molb)
734 block_bc(cr, *molb);
735 if (molb->nposres_xA > 0)
737 snew_bc(cr, molb->posres_xA, molb->nposres_xA);
738 nblock_bc(cr, molb->nposres_xA*DIM, molb->posres_xA[0]);
740 if (molb->nposres_xB > 0)
742 snew_bc(cr, molb->posres_xB, molb->nposres_xB);
743 nblock_bc(cr, molb->nposres_xB*DIM, molb->posres_xB[0]);
745 if (debug)
747 fprintf(debug, "after bc_molblock\n");
751 static void bc_atomtypes(const t_commrec *cr, t_atomtypes *atomtypes)
753 int nr;
755 block_bc(cr, atomtypes->nr);
757 nr = atomtypes->nr;
759 snew_bc(cr, atomtypes->radius, nr);
760 snew_bc(cr, atomtypes->vol, nr);
761 snew_bc(cr, atomtypes->surftens, nr);
762 snew_bc(cr, atomtypes->gb_radius, nr);
763 snew_bc(cr, atomtypes->S_hct, nr);
765 nblock_bc(cr, nr, atomtypes->radius);
766 nblock_bc(cr, nr, atomtypes->vol);
767 nblock_bc(cr, nr, atomtypes->surftens);
768 nblock_bc(cr, nr, atomtypes->gb_radius);
769 nblock_bc(cr, nr, atomtypes->S_hct);
772 /*! \brief Broadcasts ir and mtop from the master to all nodes in
773 * cr->mpi_comm_mygroup. */
774 static
775 void bcast_ir_mtop(const t_commrec *cr, t_inputrec *inputrec, gmx_mtop_t *mtop)
777 int i;
778 if (debug)
780 fprintf(debug, "in bc_data\n");
782 bc_inputrec(cr, inputrec);
783 if (debug)
785 fprintf(debug, "after bc_inputrec\n");
787 bc_symtab(cr, &mtop->symtab);
788 if (debug)
790 fprintf(debug, "after bc_symtab\n");
792 bc_string(cr, &mtop->symtab, &mtop->name);
793 if (debug)
795 fprintf(debug, "after bc_name\n");
798 bc_ffparams(cr, &mtop->ffparams);
800 block_bc(cr, mtop->nmoltype);
801 snew_bc(cr, mtop->moltype, mtop->nmoltype);
802 for (i = 0; i < mtop->nmoltype; i++)
804 bc_moltype(cr, &mtop->symtab, &mtop->moltype[i]);
807 block_bc(cr, mtop->bIntermolecularInteractions);
808 if (mtop->bIntermolecularInteractions)
810 snew_bc(cr, mtop->intermolecular_ilist, F_NRE);
811 bc_ilists(cr, mtop->intermolecular_ilist);
814 block_bc(cr, mtop->nmolblock);
815 snew_bc(cr, mtop->molblock, mtop->nmolblock);
816 for (i = 0; i < mtop->nmolblock; i++)
818 bc_molblock(cr, &mtop->molblock[i]);
821 block_bc(cr, mtop->natoms);
823 bc_atomtypes(cr, &mtop->atomtypes);
825 bc_block(cr, &mtop->mols);
826 bc_groups(cr, &mtop->symtab, mtop->natoms, &mtop->groups);
829 void init_parallel(t_commrec *cr, t_inputrec *inputrec,
830 gmx_mtop_t *mtop)
832 bcast_ir_mtop(cr, inputrec, mtop);