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! */
40 #include "broadcaststructs.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/inmemoryserializer.h"
57 #include "gromacs/utility/keyvaluetree.h"
58 #include "gromacs/utility/keyvaluetreeserializer.h"
59 #include "gromacs/utility/smalloc.h"
61 static void bc_cstring(const t_commrec
*cr
, char **s
)
65 if (MASTER(cr
) && *s
!= nullptr)
67 /* Size of the char buffer is string length + 1 for '\0' */
68 size
= strlen(*s
) + 1;
77 nblock_bc(cr
, size
, *s
);
79 else if (!MASTER(cr
) && *s
!= nullptr)
86 static void bc_string(const t_commrec
*cr
, t_symtab
*symtab
, char ***s
)
92 handle
= lookup_symtab(symtab
, *s
);
97 *s
= get_symtab_handle(symtab
, handle
);
101 static void bc_strings(const t_commrec
*cr
, t_symtab
*symtab
, int nr
, char ****nm
)
109 for (i
= 0; (i
< nr
); i
++)
111 handle
[i
] = lookup_symtab(symtab
, (*nm
)[i
]);
114 nblock_bc(cr
, nr
, handle
);
118 snew_bc(cr
, *nm
, nr
);
119 for (i
= 0; (i
< nr
); i
++)
121 (*nm
)[i
] = get_symtab_handle(symtab
, handle
[i
]);
127 static void bc_strings_resinfo(const t_commrec
*cr
, t_symtab
*symtab
,
128 int nr
, t_resinfo
*resinfo
)
136 for (i
= 0; (i
< nr
); i
++)
138 handle
[i
] = lookup_symtab(symtab
, resinfo
[i
].name
);
141 nblock_bc(cr
, nr
, handle
);
145 for (i
= 0; (i
< nr
); i
++)
147 resinfo
[i
].name
= get_symtab_handle(symtab
, handle
[i
]);
153 static void bc_symtab(const t_commrec
*cr
, t_symtab
*symtab
)
158 block_bc(cr
, symtab
->nr
);
160 snew_bc(cr
, symtab
->symbuf
, 1);
161 symbuf
= symtab
->symbuf
;
162 symbuf
->bufsize
= nr
;
163 snew_bc(cr
, symbuf
->buf
, nr
);
164 for (i
= 0; i
< nr
; i
++)
168 len
= strlen(symbuf
->buf
[i
]) + 1;
171 snew_bc(cr
, symbuf
->buf
[i
], len
);
172 nblock_bc(cr
, len
, symbuf
->buf
[i
]);
176 static void bc_block(const t_commrec
*cr
, t_block
*block
)
178 block_bc(cr
, block
->nr
);
179 snew_bc(cr
, block
->index
, block
->nr
+1);
180 nblock_bc(cr
, block
->nr
+1, block
->index
);
183 static void bc_blocka(const t_commrec
*cr
, t_blocka
*block
)
185 block_bc(cr
, block
->nr
);
186 snew_bc(cr
, block
->index
, block
->nr
+1);
187 nblock_bc(cr
, block
->nr
+1, block
->index
);
188 block_bc(cr
, block
->nra
);
191 snew_bc(cr
, block
->a
, block
->nra
);
192 nblock_bc(cr
, block
->nra
, block
->a
);
196 static void bc_grps(const t_commrec
*cr
, t_grps grps
[])
200 for (i
= 0; (i
< egcNR
); i
++)
202 block_bc(cr
, grps
[i
].nr
);
203 snew_bc(cr
, grps
[i
].nm_ind
, grps
[i
].nr
);
204 nblock_bc(cr
, grps
[i
].nr
, grps
[i
].nm_ind
);
208 static void bc_atoms(const t_commrec
*cr
, t_symtab
*symtab
, t_atoms
*atoms
)
210 block_bc(cr
, atoms
->nr
);
211 snew_bc(cr
, atoms
->atom
, atoms
->nr
);
212 nblock_bc(cr
, atoms
->nr
, atoms
->atom
);
213 bc_strings(cr
, symtab
, atoms
->nr
, &atoms
->atomname
);
214 block_bc(cr
, atoms
->nres
);
215 snew_bc(cr
, atoms
->resinfo
, atoms
->nres
);
216 nblock_bc(cr
, atoms
->nres
, atoms
->resinfo
);
217 bc_strings_resinfo(cr
, symtab
, atoms
->nres
, atoms
->resinfo
);
218 /* QMMM requires atomtypes to be known on all nodes as well */
219 bc_strings(cr
, symtab
, atoms
->nr
, &atoms
->atomtype
);
220 bc_strings(cr
, symtab
, atoms
->nr
, &atoms
->atomtypeB
);
223 static void bc_groups(const t_commrec
*cr
, t_symtab
*symtab
,
224 int natoms
, gmx_groups_t
*groups
)
228 bc_grps(cr
, groups
->grps
);
229 block_bc(cr
, groups
->ngrpname
);
230 bc_strings(cr
, symtab
, groups
->ngrpname
, &groups
->grpname
);
231 for (g
= 0; g
< egcNR
; g
++)
235 if (groups
->grpnr
[g
])
247 groups
->grpnr
[g
] = nullptr;
251 snew_bc(cr
, groups
->grpnr
[g
], n
);
252 nblock_bc(cr
, n
, groups
->grpnr
[g
]);
257 fprintf(debug
, "after bc_groups\n");
261 static void bcastPaddedRVecVector(const t_commrec
*cr
, PaddedRVecVector
*v
, unsigned int n
)
263 /* We need to allocate one element extra, since we might use
264 * (unaligned) 4-wide SIMD loads to access rvec entries.
267 nblock_bc(cr
, n
, as_rvec_array(v
->data()));
270 void bcast_state(const t_commrec
*cr
, t_state
*state
)
274 if (!PAR(cr
) || (cr
->nnodes
- cr
->npmenodes
<= 1))
279 /* Broadcasts the state sizes and flags from the master to all nodes
280 * in cr->mpi_comm_mygroup. The arrays are not broadcasted. */
281 block_bc(cr
, state
->natoms
);
282 block_bc(cr
, state
->ngtc
);
283 block_bc(cr
, state
->nnhpres
);
284 block_bc(cr
, state
->nhchainlength
);
285 block_bc(cr
, state
->flags
);
289 /* We allocate dynamically in dd_partition_system. */
292 /* The code below is reachable only by TPI and NM, so it is not
293 tested by anything. */
295 nnht
= (state
->ngtc
)*(state
->nhchainlength
);
296 nnhtp
= (state
->nnhpres
)*(state
->nhchainlength
);
298 for (i
= 0; i
< estNR
; i
++)
300 if (state
->flags
& (1<<i
))
304 case estLAMBDA
: nblock_bc(cr
, efptNR
, state
->lambda
.data()); break;
305 case estFEPSTATE
: block_bc(cr
, state
->fep_state
); break;
306 case estBOX
: block_bc(cr
, state
->box
); break;
307 case estBOX_REL
: block_bc(cr
, state
->box_rel
); break;
308 case estBOXV
: block_bc(cr
, state
->boxv
); break;
309 case estPRES_PREV
: block_bc(cr
, state
->pres_prev
); break;
310 case estSVIR_PREV
: block_bc(cr
, state
->svir_prev
); break;
311 case estFVIR_PREV
: block_bc(cr
, state
->fvir_prev
); break;
312 case estNH_XI
: nblock_abc(cr
, nnht
, &state
->nosehoover_xi
); break;
313 case estNH_VXI
: nblock_abc(cr
, nnht
, &state
->nosehoover_vxi
); break;
314 case estNHPRES_XI
: nblock_abc(cr
, nnhtp
, &state
->nhpres_xi
); break;
315 case estNHPRES_VXI
: nblock_abc(cr
, nnhtp
, &state
->nhpres_vxi
); break;
316 case estTC_INT
: nblock_abc(cr
, state
->ngtc
, &state
->therm_integral
); break;
317 case estVETA
: block_bc(cr
, state
->veta
); break;
318 case estVOL0
: block_bc(cr
, state
->vol0
); break;
319 case estX
: bcastPaddedRVecVector(cr
, &state
->x
, state
->natoms
); break;
320 case estV
: bcastPaddedRVecVector(cr
, &state
->v
, state
->natoms
); break;
321 case estCGP
: bcastPaddedRVecVector(cr
, &state
->cg_p
, state
->natoms
); break;
322 case estDISRE_INITF
: block_bc(cr
, state
->hist
.disre_initf
); break;
323 case estDISRE_RM3TAV
:
324 block_bc(cr
, state
->hist
.ndisrepairs
);
325 nblock_abc(cr
, state
->hist
.ndisrepairs
, &state
->hist
.disre_rm3tav
);
327 case estORIRE_INITF
: block_bc(cr
, state
->hist
.orire_initf
); break;
329 block_bc(cr
, state
->hist
.norire_Dtav
);
330 nblock_abc(cr
, state
->hist
.norire_Dtav
, &state
->hist
.orire_Dtav
);
334 "Communication is not implemented for %s in bcast_state",
341 static void bc_ilists(const t_commrec
*cr
, t_ilist
*ilist
)
345 /* Here we only communicate the non-zero length ilists */
348 for (ftype
= 0; ftype
< F_NRE
; ftype
++)
350 if (ilist
[ftype
].nr
> 0)
353 block_bc(cr
, ilist
[ftype
].nr
);
354 nblock_bc(cr
, ilist
[ftype
].nr
, ilist
[ftype
].iatoms
);
362 for (ftype
= 0; ftype
< F_NRE
; ftype
++)
371 block_bc(cr
, ilist
[ftype
].nr
);
372 snew_bc(cr
, ilist
[ftype
].iatoms
, ilist
[ftype
].nr
);
373 nblock_bc(cr
, ilist
[ftype
].nr
, ilist
[ftype
].iatoms
);
381 fprintf(debug
, "after bc_ilists\n");
385 static void bc_cmap(const t_commrec
*cr
, gmx_cmap_t
*cmap_grid
)
389 block_bc(cr
, cmap_grid
->ngrid
);
390 block_bc(cr
, cmap_grid
->grid_spacing
);
392 ngrid
= cmap_grid
->ngrid
;
393 nelem
= cmap_grid
->grid_spacing
* cmap_grid
->grid_spacing
;
397 snew_bc(cr
, cmap_grid
->cmapdata
, ngrid
);
399 for (i
= 0; i
< ngrid
; i
++)
401 snew_bc(cr
, cmap_grid
->cmapdata
[i
].cmap
, 4*nelem
);
402 nblock_bc(cr
, 4*nelem
, cmap_grid
->cmapdata
[i
].cmap
);
407 static void bc_ffparams(const t_commrec
*cr
, gmx_ffparams_t
*ffp
)
409 block_bc(cr
, ffp
->ntypes
);
410 block_bc(cr
, ffp
->atnr
);
411 snew_bc(cr
, ffp
->functype
, ffp
->ntypes
);
412 snew_bc(cr
, ffp
->iparams
, ffp
->ntypes
);
413 nblock_bc(cr
, ffp
->ntypes
, ffp
->functype
);
414 nblock_bc(cr
, ffp
->ntypes
, ffp
->iparams
);
415 block_bc(cr
, ffp
->reppow
);
416 block_bc(cr
, ffp
->fudgeQQ
);
417 bc_cmap(cr
, &ffp
->cmap_grid
);
420 static void bc_grpopts(const t_commrec
*cr
, t_grpopts
*g
)
424 block_bc(cr
, g
->ngtc
);
425 block_bc(cr
, g
->ngacc
);
426 block_bc(cr
, g
->ngfrz
);
427 block_bc(cr
, g
->ngener
);
428 snew_bc(cr
, g
->nrdf
, g
->ngtc
);
429 snew_bc(cr
, g
->tau_t
, g
->ngtc
);
430 snew_bc(cr
, g
->ref_t
, g
->ngtc
);
431 snew_bc(cr
, g
->acc
, g
->ngacc
);
432 snew_bc(cr
, g
->nFreeze
, g
->ngfrz
);
433 snew_bc(cr
, g
->egp_flags
, g
->ngener
*g
->ngener
);
435 nblock_bc(cr
, g
->ngtc
, g
->nrdf
);
436 nblock_bc(cr
, g
->ngtc
, g
->tau_t
);
437 nblock_bc(cr
, g
->ngtc
, g
->ref_t
);
438 nblock_bc(cr
, g
->ngacc
, g
->acc
);
439 nblock_bc(cr
, g
->ngfrz
, g
->nFreeze
);
440 nblock_bc(cr
, g
->ngener
*g
->ngener
, g
->egp_flags
);
441 snew_bc(cr
, g
->annealing
, g
->ngtc
);
442 snew_bc(cr
, g
->anneal_npoints
, g
->ngtc
);
443 snew_bc(cr
, g
->anneal_time
, g
->ngtc
);
444 snew_bc(cr
, g
->anneal_temp
, g
->ngtc
);
445 nblock_bc(cr
, g
->ngtc
, g
->annealing
);
446 nblock_bc(cr
, g
->ngtc
, g
->anneal_npoints
);
447 for (i
= 0; (i
< g
->ngtc
); i
++)
449 n
= g
->anneal_npoints
[i
];
452 snew_bc(cr
, g
->anneal_time
[i
], n
);
453 snew_bc(cr
, g
->anneal_temp
[i
], n
);
454 nblock_bc(cr
, n
, g
->anneal_time
[i
]);
455 nblock_bc(cr
, n
, g
->anneal_temp
[i
]);
459 /* QMMM stuff, see inputrec */
460 block_bc(cr
, g
->ngQM
);
461 snew_bc(cr
, g
->QMmethod
, g
->ngQM
);
462 snew_bc(cr
, g
->QMbasis
, g
->ngQM
);
463 snew_bc(cr
, g
->QMcharge
, g
->ngQM
);
464 snew_bc(cr
, g
->QMmult
, g
->ngQM
);
465 snew_bc(cr
, g
->bSH
, g
->ngQM
);
466 snew_bc(cr
, g
->CASorbitals
, g
->ngQM
);
467 snew_bc(cr
, g
->CASelectrons
, g
->ngQM
);
468 snew_bc(cr
, g
->SAon
, g
->ngQM
);
469 snew_bc(cr
, g
->SAoff
, g
->ngQM
);
470 snew_bc(cr
, g
->SAsteps
, g
->ngQM
);
474 nblock_bc(cr
, g
->ngQM
, g
->QMmethod
);
475 nblock_bc(cr
, g
->ngQM
, g
->QMbasis
);
476 nblock_bc(cr
, g
->ngQM
, g
->QMcharge
);
477 nblock_bc(cr
, g
->ngQM
, g
->QMmult
);
478 nblock_bc(cr
, g
->ngQM
, g
->bSH
);
479 nblock_bc(cr
, g
->ngQM
, g
->CASorbitals
);
480 nblock_bc(cr
, g
->ngQM
, g
->CASelectrons
);
481 nblock_bc(cr
, g
->ngQM
, g
->SAon
);
482 nblock_bc(cr
, g
->ngQM
, g
->SAoff
);
483 nblock_bc(cr
, g
->ngQM
, g
->SAsteps
);
484 /* end of QMMM stuff */
488 static void bc_pull_group(const t_commrec
*cr
, t_pull_group
*pgrp
)
493 snew_bc(cr
, pgrp
->ind
, pgrp
->nat
);
494 nblock_bc(cr
, pgrp
->nat
, pgrp
->ind
);
496 if (pgrp
->nweight
> 0)
498 snew_bc(cr
, pgrp
->weight
, pgrp
->nweight
);
499 nblock_bc(cr
, pgrp
->nweight
, pgrp
->weight
);
503 static void bc_pull(const t_commrec
*cr
, pull_params_t
*pull
)
508 snew_bc(cr
, pull
->group
, pull
->ngroup
);
509 for (g
= 0; g
< pull
->ngroup
; g
++)
511 bc_pull_group(cr
, &pull
->group
[g
]);
513 snew_bc(cr
, pull
->coord
, pull
->ncoord
);
514 nblock_bc(cr
, pull
->ncoord
, pull
->coord
);
515 for (int c
= 0; c
< pull
->ncoord
; c
++)
519 pull
->coord
[c
].externalPotentialProvider
= nullptr;
521 if (pull
->coord
[c
].eType
== epullEXTERNAL
)
523 bc_cstring(cr
, &pull
->coord
[c
].externalPotentialProvider
);
528 static void bc_rotgrp(const t_commrec
*cr
, t_rotgrp
*rotg
)
533 snew_bc(cr
, rotg
->ind
, rotg
->nat
);
534 nblock_bc(cr
, rotg
->nat
, rotg
->ind
);
535 snew_bc(cr
, rotg
->x_ref
, rotg
->nat
);
536 nblock_bc(cr
, rotg
->nat
, rotg
->x_ref
);
540 static void bc_rot(const t_commrec
*cr
, t_rot
*rot
)
545 snew_bc(cr
, rot
->grp
, rot
->ngrp
);
546 for (g
= 0; g
< rot
->ngrp
; g
++)
548 bc_rotgrp(cr
, &rot
->grp
[g
]);
552 static void bc_imd(const t_commrec
*cr
, t_IMD
*imd
)
555 snew_bc(cr
, imd
->ind
, imd
->nat
);
556 nblock_bc(cr
, imd
->nat
, imd
->ind
);
559 static void bc_fepvals(const t_commrec
*cr
, t_lambda
*fep
)
563 block_bc(cr
, fep
->nstdhdl
);
564 block_bc(cr
, fep
->init_lambda
);
565 block_bc(cr
, fep
->init_fep_state
);
566 block_bc(cr
, fep
->delta_lambda
);
567 block_bc(cr
, fep
->edHdLPrintEnergy
);
568 block_bc(cr
, fep
->n_lambda
);
569 if (fep
->n_lambda
> 0)
571 snew_bc(cr
, fep
->all_lambda
, efptNR
);
572 nblock_bc(cr
, efptNR
, fep
->all_lambda
);
573 for (i
= 0; i
< efptNR
; i
++)
575 snew_bc(cr
, fep
->all_lambda
[i
], fep
->n_lambda
);
576 nblock_bc(cr
, fep
->n_lambda
, fep
->all_lambda
[i
]);
579 block_bc(cr
, fep
->sc_alpha
);
580 block_bc(cr
, fep
->sc_power
);
581 block_bc(cr
, fep
->sc_r_power
);
582 block_bc(cr
, fep
->sc_sigma
);
583 block_bc(cr
, fep
->sc_sigma_min
);
584 block_bc(cr
, fep
->bScCoul
);
585 nblock_bc(cr
, efptNR
, &(fep
->separate_dvdl
[0]));
586 block_bc(cr
, fep
->dhdl_derivatives
);
587 block_bc(cr
, fep
->dh_hist_size
);
588 block_bc(cr
, fep
->dh_hist_spacing
);
591 fprintf(debug
, "after bc_fepvals\n");
595 static void bc_expandedvals(const t_commrec
*cr
, t_expanded
*expand
, int n_lambda
)
597 block_bc(cr
, expand
->nstexpanded
);
598 block_bc(cr
, expand
->elamstats
);
599 block_bc(cr
, expand
->elmcmove
);
600 block_bc(cr
, expand
->elmceq
);
601 block_bc(cr
, expand
->equil_n_at_lam
);
602 block_bc(cr
, expand
->equil_wl_delta
);
603 block_bc(cr
, expand
->equil_ratio
);
604 block_bc(cr
, expand
->equil_steps
);
605 block_bc(cr
, expand
->equil_samples
);
606 block_bc(cr
, expand
->lmc_seed
);
607 block_bc(cr
, expand
->minvar
);
608 block_bc(cr
, expand
->minvar_const
);
609 block_bc(cr
, expand
->c_range
);
610 block_bc(cr
, expand
->bSymmetrizedTMatrix
);
611 block_bc(cr
, expand
->nstTij
);
612 block_bc(cr
, expand
->lmc_repeats
);
613 block_bc(cr
, expand
->lmc_forced_nstart
);
614 block_bc(cr
, expand
->gibbsdeltalam
);
615 block_bc(cr
, expand
->wl_scale
);
616 block_bc(cr
, expand
->wl_ratio
);
617 block_bc(cr
, expand
->init_wl_delta
);
618 block_bc(cr
, expand
->bInit_weights
);
619 snew_bc(cr
, expand
->init_lambda_weights
, n_lambda
);
620 nblock_bc(cr
, n_lambda
, expand
->init_lambda_weights
);
621 block_bc(cr
, expand
->mc_temp
);
624 fprintf(debug
, "after bc_expandedvals\n");
628 static void bc_simtempvals(const t_commrec
*cr
, t_simtemp
*simtemp
, int n_lambda
)
630 block_bc(cr
, simtemp
->simtemp_low
);
631 block_bc(cr
, simtemp
->simtemp_high
);
632 block_bc(cr
, simtemp
->eSimTempScale
);
633 snew_bc(cr
, simtemp
->temperatures
, n_lambda
);
634 nblock_bc(cr
, n_lambda
, simtemp
->temperatures
);
637 fprintf(debug
, "after bc_simtempvals\n");
642 static void bc_swapions(const t_commrec
*cr
, t_swapcoords
*swap
)
646 /* Broadcast atom indices for split groups, solvent group, and for all user-defined swap groups */
647 snew_bc(cr
, swap
->grp
, swap
->ngrp
);
648 for (int i
= 0; i
< swap
->ngrp
; i
++)
650 t_swapGroup
*g
= &swap
->grp
[i
];
653 snew_bc(cr
, g
->ind
, g
->nat
);
654 nblock_bc(cr
, g
->nat
, g
->ind
);
659 len
= strlen(g
->molname
);
662 snew_bc(cr
, g
->molname
, len
);
663 nblock_bc(cr
, len
, g
->molname
);
668 static void bc_inputrec(const t_commrec
*cr
, t_inputrec
*inputrec
)
670 // Note that this overwrites pointers in inputrec, so all pointer fields
671 // Must be initialized separately below.
672 block_bc(cr
, *inputrec
);
675 gmx::InMemorySerializer serializer
;
676 gmx::serializeKeyValueTree(*inputrec
->params
, &serializer
);
677 std::vector
<char> buffer
= serializer
.finishAndGetBuffer();
678 size_t size
= buffer
.size();
680 nblock_bc(cr
, size
, buffer
.data());
684 // block_bc() above overwrites the old pointer, so set it to a
685 // reasonable value in case code below throws.
686 // cppcheck-suppress redundantAssignment
687 inputrec
->params
= nullptr;
688 std::vector
<char> buffer
;
691 nblock_abc(cr
, size
, &buffer
);
692 gmx::InMemoryDeserializer
serializer(buffer
);
693 // cppcheck-suppress redundantAssignment
694 inputrec
->params
= new gmx::KeyValueTreeObject(
695 gmx::deserializeKeyValueTree(&serializer
));
698 bc_grpopts(cr
, &(inputrec
->opts
));
700 /* even if efep is efepNO, we need to initialize to make sure that
701 * n_lambda is set to zero */
703 snew_bc(cr
, inputrec
->fepvals
, 1);
704 if (inputrec
->efep
!= efepNO
|| inputrec
->bSimTemp
)
706 bc_fepvals(cr
, inputrec
->fepvals
);
708 /* need to initialize this as well because of data checked for in the logic */
709 snew_bc(cr
, inputrec
->expandedvals
, 1);
710 if (inputrec
->bExpanded
)
712 bc_expandedvals(cr
, inputrec
->expandedvals
, inputrec
->fepvals
->n_lambda
);
714 snew_bc(cr
, inputrec
->simtempvals
, 1);
715 if (inputrec
->bSimTemp
)
717 bc_simtempvals(cr
, inputrec
->simtempvals
, inputrec
->fepvals
->n_lambda
);
721 snew_bc(cr
, inputrec
->pull
, 1);
722 bc_pull(cr
, inputrec
->pull
);
726 snew_bc(cr
, inputrec
->rot
, 1);
727 bc_rot(cr
, inputrec
->rot
);
731 snew_bc(cr
, inputrec
->imd
, 1);
732 bc_imd(cr
, inputrec
->imd
);
734 if (inputrec
->eSwapCoords
!= eswapNO
)
736 snew_bc(cr
, inputrec
->swap
, 1);
737 bc_swapions(cr
, inputrec
->swap
);
741 static void bc_moltype(const t_commrec
*cr
, t_symtab
*symtab
,
742 gmx_moltype_t
*moltype
)
744 bc_string(cr
, symtab
, &moltype
->name
);
745 bc_atoms(cr
, symtab
, &moltype
->atoms
);
748 fprintf(debug
, "after bc_atoms\n");
751 bc_ilists(cr
, moltype
->ilist
);
752 bc_block(cr
, &moltype
->cgs
);
753 bc_blocka(cr
, &moltype
->excls
);
756 static void bc_molblock(const t_commrec
*cr
, gmx_molblock_t
*molb
)
759 if (molb
->nposres_xA
> 0)
761 snew_bc(cr
, molb
->posres_xA
, molb
->nposres_xA
);
762 nblock_bc(cr
, molb
->nposres_xA
*DIM
, molb
->posres_xA
[0]);
764 if (molb
->nposres_xB
> 0)
766 snew_bc(cr
, molb
->posres_xB
, molb
->nposres_xB
);
767 nblock_bc(cr
, molb
->nposres_xB
*DIM
, molb
->posres_xB
[0]);
771 fprintf(debug
, "after bc_molblock\n");
775 static void bc_atomtypes(const t_commrec
*cr
, t_atomtypes
*atomtypes
)
779 block_bc(cr
, atomtypes
->nr
);
783 snew_bc(cr
, atomtypes
->radius
, nr
);
784 snew_bc(cr
, atomtypes
->vol
, nr
);
785 snew_bc(cr
, atomtypes
->surftens
, nr
);
786 snew_bc(cr
, atomtypes
->gb_radius
, nr
);
787 snew_bc(cr
, atomtypes
->S_hct
, nr
);
789 nblock_bc(cr
, nr
, atomtypes
->radius
);
790 nblock_bc(cr
, nr
, atomtypes
->vol
);
791 nblock_bc(cr
, nr
, atomtypes
->surftens
);
792 nblock_bc(cr
, nr
, atomtypes
->gb_radius
);
793 nblock_bc(cr
, nr
, atomtypes
->S_hct
);
796 /*! \brief Broadcasts ir and mtop from the master to all nodes in
797 * cr->mpi_comm_mygroup. */
799 void bcast_ir_mtop(const t_commrec
*cr
, t_inputrec
*inputrec
, gmx_mtop_t
*mtop
)
804 fprintf(debug
, "in bc_data\n");
806 bc_inputrec(cr
, inputrec
);
809 fprintf(debug
, "after bc_inputrec\n");
811 bc_symtab(cr
, &mtop
->symtab
);
814 fprintf(debug
, "after bc_symtab\n");
816 bc_string(cr
, &mtop
->symtab
, &mtop
->name
);
819 fprintf(debug
, "after bc_name\n");
822 bc_ffparams(cr
, &mtop
->ffparams
);
824 block_bc(cr
, mtop
->nmoltype
);
825 snew_bc(cr
, mtop
->moltype
, mtop
->nmoltype
);
826 for (i
= 0; i
< mtop
->nmoltype
; i
++)
828 bc_moltype(cr
, &mtop
->symtab
, &mtop
->moltype
[i
]);
831 block_bc(cr
, mtop
->bIntermolecularInteractions
);
832 if (mtop
->bIntermolecularInteractions
)
834 snew_bc(cr
, mtop
->intermolecular_ilist
, F_NRE
);
835 bc_ilists(cr
, mtop
->intermolecular_ilist
);
838 block_bc(cr
, mtop
->nmolblock
);
839 snew_bc(cr
, mtop
->molblock
, mtop
->nmolblock
);
840 for (i
= 0; i
< mtop
->nmolblock
; i
++)
842 bc_molblock(cr
, &mtop
->molblock
[i
]);
845 block_bc(cr
, mtop
->natoms
);
847 bc_atomtypes(cr
, &mtop
->atomtypes
);
849 bc_block(cr
, &mtop
->mols
);
850 bc_groups(cr
, &mtop
->symtab
, mtop
->natoms
, &mtop
->groups
);
853 void init_parallel(t_commrec
*cr
, t_inputrec
*inputrec
,
856 bcast_ir_mtop(cr
, inputrec
, mtop
);