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
, int numAtoms
)
263 v
->resize(gmx::paddedRVecVectorSize(numAtoms
));
264 nblock_bc(cr
, numAtoms
, as_rvec_array(v
->data()));
267 void broadcastStateWithoutDynamics(const t_commrec
*cr
, t_state
*state
)
269 GMX_RELEASE_ASSERT(!DOMAINDECOMP(cr
), "broadcastStateWithoutDynamics should only be used for special cases without domain decomposition");
276 /* Broadcasts the state sizes and flags from the master to all ranks
277 * in cr->mpi_comm_mygroup.
279 block_bc(cr
, state
->natoms
);
280 block_bc(cr
, state
->flags
);
282 for (int i
= 0; i
< estNR
; i
++)
284 if (state
->flags
& (1 << i
))
289 nblock_bc(cr
, efptNR
, state
->lambda
.data());
292 block_bc(cr
, state
->fep_state
);
295 block_bc(cr
, state
->box
);
298 bcastPaddedRVecVector(cr
, &state
->x
, state
->natoms
);
301 GMX_RELEASE_ASSERT(false, "The state has a dynamic entry, while no dynamic entries should be present");
308 static void bc_ilists(const t_commrec
*cr
, t_ilist
*ilist
)
312 /* Here we only communicate the non-zero length ilists */
315 for (ftype
= 0; ftype
< F_NRE
; ftype
++)
317 if (ilist
[ftype
].nr
> 0)
320 block_bc(cr
, ilist
[ftype
].nr
);
321 nblock_bc(cr
, ilist
[ftype
].nr
, ilist
[ftype
].iatoms
);
329 for (ftype
= 0; ftype
< F_NRE
; ftype
++)
338 block_bc(cr
, ilist
[ftype
].nr
);
339 snew_bc(cr
, ilist
[ftype
].iatoms
, ilist
[ftype
].nr
);
340 nblock_bc(cr
, ilist
[ftype
].nr
, ilist
[ftype
].iatoms
);
348 fprintf(debug
, "after bc_ilists\n");
352 static void bc_cmap(const t_commrec
*cr
, gmx_cmap_t
*cmap_grid
)
356 block_bc(cr
, cmap_grid
->ngrid
);
357 block_bc(cr
, cmap_grid
->grid_spacing
);
359 ngrid
= cmap_grid
->ngrid
;
360 nelem
= cmap_grid
->grid_spacing
* cmap_grid
->grid_spacing
;
364 snew_bc(cr
, cmap_grid
->cmapdata
, ngrid
);
366 for (i
= 0; i
< ngrid
; i
++)
368 snew_bc(cr
, cmap_grid
->cmapdata
[i
].cmap
, 4*nelem
);
369 nblock_bc(cr
, 4*nelem
, cmap_grid
->cmapdata
[i
].cmap
);
374 static void bc_ffparams(const t_commrec
*cr
, gmx_ffparams_t
*ffp
)
376 block_bc(cr
, ffp
->ntypes
);
377 block_bc(cr
, ffp
->atnr
);
378 snew_bc(cr
, ffp
->functype
, ffp
->ntypes
);
379 snew_bc(cr
, ffp
->iparams
, ffp
->ntypes
);
380 nblock_bc(cr
, ffp
->ntypes
, ffp
->functype
);
381 nblock_bc(cr
, ffp
->ntypes
, ffp
->iparams
);
382 block_bc(cr
, ffp
->reppow
);
383 block_bc(cr
, ffp
->fudgeQQ
);
384 bc_cmap(cr
, &ffp
->cmap_grid
);
387 static void bc_grpopts(const t_commrec
*cr
, t_grpopts
*g
)
391 block_bc(cr
, g
->ngtc
);
392 block_bc(cr
, g
->ngacc
);
393 block_bc(cr
, g
->ngfrz
);
394 block_bc(cr
, g
->ngener
);
395 snew_bc(cr
, g
->nrdf
, g
->ngtc
);
396 snew_bc(cr
, g
->tau_t
, g
->ngtc
);
397 snew_bc(cr
, g
->ref_t
, g
->ngtc
);
398 snew_bc(cr
, g
->acc
, g
->ngacc
);
399 snew_bc(cr
, g
->nFreeze
, g
->ngfrz
);
400 snew_bc(cr
, g
->egp_flags
, g
->ngener
*g
->ngener
);
402 nblock_bc(cr
, g
->ngtc
, g
->nrdf
);
403 nblock_bc(cr
, g
->ngtc
, g
->tau_t
);
404 nblock_bc(cr
, g
->ngtc
, g
->ref_t
);
405 nblock_bc(cr
, g
->ngacc
, g
->acc
);
406 nblock_bc(cr
, g
->ngfrz
, g
->nFreeze
);
407 nblock_bc(cr
, g
->ngener
*g
->ngener
, g
->egp_flags
);
408 snew_bc(cr
, g
->annealing
, g
->ngtc
);
409 snew_bc(cr
, g
->anneal_npoints
, g
->ngtc
);
410 snew_bc(cr
, g
->anneal_time
, g
->ngtc
);
411 snew_bc(cr
, g
->anneal_temp
, g
->ngtc
);
412 nblock_bc(cr
, g
->ngtc
, g
->annealing
);
413 nblock_bc(cr
, g
->ngtc
, g
->anneal_npoints
);
414 for (i
= 0; (i
< g
->ngtc
); i
++)
416 n
= g
->anneal_npoints
[i
];
419 snew_bc(cr
, g
->anneal_time
[i
], n
);
420 snew_bc(cr
, g
->anneal_temp
[i
], n
);
421 nblock_bc(cr
, n
, g
->anneal_time
[i
]);
422 nblock_bc(cr
, n
, g
->anneal_temp
[i
]);
426 /* QMMM stuff, see inputrec */
427 block_bc(cr
, g
->ngQM
);
428 snew_bc(cr
, g
->QMmethod
, g
->ngQM
);
429 snew_bc(cr
, g
->QMbasis
, g
->ngQM
);
430 snew_bc(cr
, g
->QMcharge
, g
->ngQM
);
431 snew_bc(cr
, g
->QMmult
, g
->ngQM
);
432 snew_bc(cr
, g
->bSH
, g
->ngQM
);
433 snew_bc(cr
, g
->CASorbitals
, g
->ngQM
);
434 snew_bc(cr
, g
->CASelectrons
, g
->ngQM
);
435 snew_bc(cr
, g
->SAon
, g
->ngQM
);
436 snew_bc(cr
, g
->SAoff
, g
->ngQM
);
437 snew_bc(cr
, g
->SAsteps
, g
->ngQM
);
441 nblock_bc(cr
, g
->ngQM
, g
->QMmethod
);
442 nblock_bc(cr
, g
->ngQM
, g
->QMbasis
);
443 nblock_bc(cr
, g
->ngQM
, g
->QMcharge
);
444 nblock_bc(cr
, g
->ngQM
, g
->QMmult
);
445 nblock_bc(cr
, g
->ngQM
, g
->bSH
);
446 nblock_bc(cr
, g
->ngQM
, g
->CASorbitals
);
447 nblock_bc(cr
, g
->ngQM
, g
->CASelectrons
);
448 nblock_bc(cr
, g
->ngQM
, g
->SAon
);
449 nblock_bc(cr
, g
->ngQM
, g
->SAoff
);
450 nblock_bc(cr
, g
->ngQM
, g
->SAsteps
);
451 /* end of QMMM stuff */
455 static void bc_pull_group(const t_commrec
*cr
, t_pull_group
*pgrp
)
460 snew_bc(cr
, pgrp
->ind
, pgrp
->nat
);
461 nblock_bc(cr
, pgrp
->nat
, pgrp
->ind
);
463 if (pgrp
->nweight
> 0)
465 snew_bc(cr
, pgrp
->weight
, pgrp
->nweight
);
466 nblock_bc(cr
, pgrp
->nweight
, pgrp
->weight
);
470 static void bc_pull(const t_commrec
*cr
, pull_params_t
*pull
)
475 snew_bc(cr
, pull
->group
, pull
->ngroup
);
476 for (g
= 0; g
< pull
->ngroup
; g
++)
478 bc_pull_group(cr
, &pull
->group
[g
]);
480 snew_bc(cr
, pull
->coord
, pull
->ncoord
);
481 nblock_bc(cr
, pull
->ncoord
, pull
->coord
);
482 for (int c
= 0; c
< pull
->ncoord
; c
++)
486 pull
->coord
[c
].externalPotentialProvider
= nullptr;
488 if (pull
->coord
[c
].eType
== epullEXTERNAL
)
490 bc_cstring(cr
, &pull
->coord
[c
].externalPotentialProvider
);
495 static void bc_rotgrp(const t_commrec
*cr
, t_rotgrp
*rotg
)
500 snew_bc(cr
, rotg
->ind
, rotg
->nat
);
501 nblock_bc(cr
, rotg
->nat
, rotg
->ind
);
502 snew_bc(cr
, rotg
->x_ref
, rotg
->nat
);
503 nblock_bc(cr
, rotg
->nat
, rotg
->x_ref
);
507 static void bc_rot(const t_commrec
*cr
, t_rot
*rot
)
512 snew_bc(cr
, rot
->grp
, rot
->ngrp
);
513 for (g
= 0; g
< rot
->ngrp
; g
++)
515 bc_rotgrp(cr
, &rot
->grp
[g
]);
519 static void bc_imd(const t_commrec
*cr
, t_IMD
*imd
)
522 snew_bc(cr
, imd
->ind
, imd
->nat
);
523 nblock_bc(cr
, imd
->nat
, imd
->ind
);
526 static void bc_fepvals(const t_commrec
*cr
, t_lambda
*fep
)
530 block_bc(cr
, fep
->nstdhdl
);
531 block_bc(cr
, fep
->init_lambda
);
532 block_bc(cr
, fep
->init_fep_state
);
533 block_bc(cr
, fep
->delta_lambda
);
534 block_bc(cr
, fep
->edHdLPrintEnergy
);
535 block_bc(cr
, fep
->n_lambda
);
536 if (fep
->n_lambda
> 0)
538 snew_bc(cr
, fep
->all_lambda
, efptNR
);
539 nblock_bc(cr
, efptNR
, fep
->all_lambda
);
540 for (i
= 0; i
< efptNR
; i
++)
542 snew_bc(cr
, fep
->all_lambda
[i
], fep
->n_lambda
);
543 nblock_bc(cr
, fep
->n_lambda
, fep
->all_lambda
[i
]);
546 block_bc(cr
, fep
->sc_alpha
);
547 block_bc(cr
, fep
->sc_power
);
548 block_bc(cr
, fep
->sc_r_power
);
549 block_bc(cr
, fep
->sc_sigma
);
550 block_bc(cr
, fep
->sc_sigma_min
);
551 block_bc(cr
, fep
->bScCoul
);
552 nblock_bc(cr
, efptNR
, &(fep
->separate_dvdl
[0]));
553 block_bc(cr
, fep
->dhdl_derivatives
);
554 block_bc(cr
, fep
->dh_hist_size
);
555 block_bc(cr
, fep
->dh_hist_spacing
);
558 fprintf(debug
, "after bc_fepvals\n");
562 static void bc_expandedvals(const t_commrec
*cr
, t_expanded
*expand
, int n_lambda
)
564 block_bc(cr
, expand
->nstexpanded
);
565 block_bc(cr
, expand
->elamstats
);
566 block_bc(cr
, expand
->elmcmove
);
567 block_bc(cr
, expand
->elmceq
);
568 block_bc(cr
, expand
->equil_n_at_lam
);
569 block_bc(cr
, expand
->equil_wl_delta
);
570 block_bc(cr
, expand
->equil_ratio
);
571 block_bc(cr
, expand
->equil_steps
);
572 block_bc(cr
, expand
->equil_samples
);
573 block_bc(cr
, expand
->lmc_seed
);
574 block_bc(cr
, expand
->minvar
);
575 block_bc(cr
, expand
->minvar_const
);
576 block_bc(cr
, expand
->c_range
);
577 block_bc(cr
, expand
->bSymmetrizedTMatrix
);
578 block_bc(cr
, expand
->nstTij
);
579 block_bc(cr
, expand
->lmc_repeats
);
580 block_bc(cr
, expand
->lmc_forced_nstart
);
581 block_bc(cr
, expand
->gibbsdeltalam
);
582 block_bc(cr
, expand
->wl_scale
);
583 block_bc(cr
, expand
->wl_ratio
);
584 block_bc(cr
, expand
->init_wl_delta
);
585 block_bc(cr
, expand
->bInit_weights
);
586 snew_bc(cr
, expand
->init_lambda_weights
, n_lambda
);
587 nblock_bc(cr
, n_lambda
, expand
->init_lambda_weights
);
588 block_bc(cr
, expand
->mc_temp
);
591 fprintf(debug
, "after bc_expandedvals\n");
595 static void bc_simtempvals(const t_commrec
*cr
, t_simtemp
*simtemp
, int n_lambda
)
597 block_bc(cr
, simtemp
->simtemp_low
);
598 block_bc(cr
, simtemp
->simtemp_high
);
599 block_bc(cr
, simtemp
->eSimTempScale
);
600 snew_bc(cr
, simtemp
->temperatures
, n_lambda
);
601 nblock_bc(cr
, n_lambda
, simtemp
->temperatures
);
604 fprintf(debug
, "after bc_simtempvals\n");
609 static void bc_swapions(const t_commrec
*cr
, t_swapcoords
*swap
)
613 /* Broadcast atom indices for split groups, solvent group, and for all user-defined swap groups */
614 snew_bc(cr
, swap
->grp
, swap
->ngrp
);
615 for (int i
= 0; i
< swap
->ngrp
; i
++)
617 t_swapGroup
*g
= &swap
->grp
[i
];
620 snew_bc(cr
, g
->ind
, g
->nat
);
621 nblock_bc(cr
, g
->nat
, g
->ind
);
626 len
= strlen(g
->molname
);
629 snew_bc(cr
, g
->molname
, len
);
630 nblock_bc(cr
, len
, g
->molname
);
635 static void bc_inputrec(const t_commrec
*cr
, t_inputrec
*inputrec
)
637 // Note that this overwrites pointers in inputrec, so all pointer fields
638 // Must be initialized separately below.
639 block_bc(cr
, *inputrec
);
642 gmx::InMemorySerializer serializer
;
643 gmx::serializeKeyValueTree(*inputrec
->params
, &serializer
);
644 std::vector
<char> buffer
= serializer
.finishAndGetBuffer();
645 size_t size
= buffer
.size();
647 nblock_bc(cr
, size
, buffer
.data());
651 // block_bc() above overwrites the old pointer, so set it to a
652 // reasonable value in case code below throws.
653 // cppcheck-suppress redundantAssignment
654 inputrec
->params
= nullptr;
655 std::vector
<char> buffer
;
658 nblock_abc(cr
, size
, &buffer
);
659 gmx::InMemoryDeserializer
serializer(buffer
);
660 // cppcheck-suppress redundantAssignment
661 inputrec
->params
= new gmx::KeyValueTreeObject(
662 gmx::deserializeKeyValueTree(&serializer
));
665 bc_grpopts(cr
, &(inputrec
->opts
));
667 /* even if efep is efepNO, we need to initialize to make sure that
668 * n_lambda is set to zero */
670 snew_bc(cr
, inputrec
->fepvals
, 1);
671 if (inputrec
->efep
!= efepNO
|| inputrec
->bSimTemp
)
673 bc_fepvals(cr
, inputrec
->fepvals
);
675 /* need to initialize this as well because of data checked for in the logic */
676 snew_bc(cr
, inputrec
->expandedvals
, 1);
677 if (inputrec
->bExpanded
)
679 bc_expandedvals(cr
, inputrec
->expandedvals
, inputrec
->fepvals
->n_lambda
);
681 snew_bc(cr
, inputrec
->simtempvals
, 1);
682 if (inputrec
->bSimTemp
)
684 bc_simtempvals(cr
, inputrec
->simtempvals
, inputrec
->fepvals
->n_lambda
);
688 snew_bc(cr
, inputrec
->pull
, 1);
689 bc_pull(cr
, inputrec
->pull
);
693 snew_bc(cr
, inputrec
->rot
, 1);
694 bc_rot(cr
, inputrec
->rot
);
698 snew_bc(cr
, inputrec
->imd
, 1);
699 bc_imd(cr
, inputrec
->imd
);
701 if (inputrec
->eSwapCoords
!= eswapNO
)
703 snew_bc(cr
, inputrec
->swap
, 1);
704 bc_swapions(cr
, inputrec
->swap
);
708 static void bc_moltype(const t_commrec
*cr
, t_symtab
*symtab
,
709 gmx_moltype_t
*moltype
)
711 bc_string(cr
, symtab
, &moltype
->name
);
712 bc_atoms(cr
, symtab
, &moltype
->atoms
);
715 fprintf(debug
, "after bc_atoms\n");
718 bc_ilists(cr
, moltype
->ilist
);
719 bc_block(cr
, &moltype
->cgs
);
720 bc_blocka(cr
, &moltype
->excls
);
723 static void bc_molblock(const t_commrec
*cr
, gmx_molblock_t
*molb
)
726 if (molb
->nposres_xA
> 0)
728 snew_bc(cr
, molb
->posres_xA
, molb
->nposres_xA
);
729 nblock_bc(cr
, molb
->nposres_xA
*DIM
, molb
->posres_xA
[0]);
731 if (molb
->nposres_xB
> 0)
733 snew_bc(cr
, molb
->posres_xB
, molb
->nposres_xB
);
734 nblock_bc(cr
, molb
->nposres_xB
*DIM
, molb
->posres_xB
[0]);
738 fprintf(debug
, "after bc_molblock\n");
742 static void bc_atomtypes(const t_commrec
*cr
, t_atomtypes
*atomtypes
)
746 block_bc(cr
, atomtypes
->nr
);
750 snew_bc(cr
, atomtypes
->radius
, nr
);
751 snew_bc(cr
, atomtypes
->vol
, nr
);
752 snew_bc(cr
, atomtypes
->surftens
, nr
);
753 snew_bc(cr
, atomtypes
->gb_radius
, nr
);
754 snew_bc(cr
, atomtypes
->S_hct
, nr
);
756 nblock_bc(cr
, nr
, atomtypes
->radius
);
757 nblock_bc(cr
, nr
, atomtypes
->vol
);
758 nblock_bc(cr
, nr
, atomtypes
->surftens
);
759 nblock_bc(cr
, nr
, atomtypes
->gb_radius
);
760 nblock_bc(cr
, nr
, atomtypes
->S_hct
);
763 /*! \brief Broadcasts ir and mtop from the master to all nodes in
764 * cr->mpi_comm_mygroup. */
766 void bcast_ir_mtop(const t_commrec
*cr
, t_inputrec
*inputrec
, gmx_mtop_t
*mtop
)
771 fprintf(debug
, "in bc_data\n");
773 bc_inputrec(cr
, inputrec
);
776 fprintf(debug
, "after bc_inputrec\n");
778 bc_symtab(cr
, &mtop
->symtab
);
781 fprintf(debug
, "after bc_symtab\n");
783 bc_string(cr
, &mtop
->symtab
, &mtop
->name
);
786 fprintf(debug
, "after bc_name\n");
789 bc_ffparams(cr
, &mtop
->ffparams
);
791 block_bc(cr
, mtop
->nmoltype
);
792 snew_bc(cr
, mtop
->moltype
, mtop
->nmoltype
);
793 for (i
= 0; i
< mtop
->nmoltype
; i
++)
795 bc_moltype(cr
, &mtop
->symtab
, &mtop
->moltype
[i
]);
798 block_bc(cr
, mtop
->bIntermolecularInteractions
);
799 if (mtop
->bIntermolecularInteractions
)
801 snew_bc(cr
, mtop
->intermolecular_ilist
, F_NRE
);
802 bc_ilists(cr
, mtop
->intermolecular_ilist
);
805 block_bc(cr
, mtop
->nmolblock
);
806 snew_bc(cr
, mtop
->molblock
, mtop
->nmolblock
);
807 for (i
= 0; i
< mtop
->nmolblock
; i
++)
809 bc_molblock(cr
, &mtop
->molblock
[i
]);
812 block_bc(cr
, mtop
->natoms
);
814 bc_atomtypes(cr
, &mtop
->atomtypes
);
816 bc_block(cr
, &mtop
->mols
);
817 bc_groups(cr
, &mtop
->symtab
, mtop
->natoms
, &mtop
->groups
);
820 void init_parallel(t_commrec
*cr
, t_inputrec
*inputrec
,
823 bcast_ir_mtop(cr
, inputrec
, mtop
);