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, 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! */
42 #include "gromacs/gmxlib/main.h"
43 #include "gromacs/gmxlib/network.h"
44 #include "gromacs/legacyheaders/types/commrec.h"
45 #include "gromacs/legacyheaders/types/enums.h"
46 #include "gromacs/math/vec.h"
47 #include "gromacs/mdlib/mdrun.h"
48 #include "gromacs/mdlib/tgroup.h"
49 #include "gromacs/mdtypes/inputrec.h"
50 #include "gromacs/mdtypes/state.h"
51 #include "gromacs/topology/symtab.h"
52 #include "gromacs/topology/topology.h"
53 #include "gromacs/utility/fatalerror.h"
54 #include "gromacs/utility/smalloc.h"
56 #define block_bc(cr, d) gmx_bcast( sizeof(d), &(d), (cr))
57 /* Probably the test for (nr) > 0 in the next macro is only needed
58 * on BlueGene(/L), where IBM's MPI_Bcast will segfault after
59 * dereferencing a null pointer, even when no data is to be transferred. */
60 #define nblock_bc(cr, nr, d) { if ((nr) > 0) {gmx_bcast((nr)*sizeof((d)[0]), (d), (cr)); }}
61 #define snew_bc(cr, d, nr) { if (!MASTER(cr)) {snew((d), (nr)); }}
62 /* Dirty macro with bAlloc not as an argument */
63 #define nblock_abc(cr, nr, d) { if (bAlloc) {snew((d), (nr)); } nblock_bc(cr, (nr), (d)); }
65 static void bc_string(const t_commrec
*cr
, t_symtab
*symtab
, char ***s
)
71 handle
= lookup_symtab(symtab
, *s
);
76 *s
= get_symtab_handle(symtab
, handle
);
80 static void bc_strings(const t_commrec
*cr
, t_symtab
*symtab
, int nr
, char ****nm
)
88 for (i
= 0; (i
< nr
); i
++)
90 handle
[i
] = lookup_symtab(symtab
, (*nm
)[i
]);
93 nblock_bc(cr
, nr
, handle
);
98 for (i
= 0; (i
< nr
); i
++)
100 (*nm
)[i
] = get_symtab_handle(symtab
, handle
[i
]);
106 static void bc_strings_resinfo(const t_commrec
*cr
, t_symtab
*symtab
,
107 int nr
, t_resinfo
*resinfo
)
115 for (i
= 0; (i
< nr
); i
++)
117 handle
[i
] = lookup_symtab(symtab
, resinfo
[i
].name
);
120 nblock_bc(cr
, nr
, handle
);
124 for (i
= 0; (i
< nr
); i
++)
126 resinfo
[i
].name
= get_symtab_handle(symtab
, handle
[i
]);
132 static void bc_symtab(const t_commrec
*cr
, t_symtab
*symtab
)
137 block_bc(cr
, symtab
->nr
);
139 snew_bc(cr
, symtab
->symbuf
, 1);
140 symbuf
= symtab
->symbuf
;
141 symbuf
->bufsize
= nr
;
142 snew_bc(cr
, symbuf
->buf
, nr
);
143 for (i
= 0; i
< nr
; i
++)
147 len
= strlen(symbuf
->buf
[i
]) + 1;
150 snew_bc(cr
, symbuf
->buf
[i
], len
);
151 nblock_bc(cr
, len
, symbuf
->buf
[i
]);
155 static void bc_block(const t_commrec
*cr
, t_block
*block
)
157 block_bc(cr
, block
->nr
);
158 snew_bc(cr
, block
->index
, block
->nr
+1);
159 nblock_bc(cr
, block
->nr
+1, block
->index
);
162 static void bc_blocka(const t_commrec
*cr
, t_blocka
*block
)
164 block_bc(cr
, block
->nr
);
165 snew_bc(cr
, block
->index
, block
->nr
+1);
166 nblock_bc(cr
, block
->nr
+1, block
->index
);
167 block_bc(cr
, block
->nra
);
170 snew_bc(cr
, block
->a
, block
->nra
);
171 nblock_bc(cr
, block
->nra
, block
->a
);
175 static void bc_grps(const t_commrec
*cr
, t_grps grps
[])
179 for (i
= 0; (i
< egcNR
); i
++)
181 block_bc(cr
, grps
[i
].nr
);
182 snew_bc(cr
, grps
[i
].nm_ind
, grps
[i
].nr
);
183 nblock_bc(cr
, grps
[i
].nr
, grps
[i
].nm_ind
);
187 static void bc_atoms(const t_commrec
*cr
, t_symtab
*symtab
, t_atoms
*atoms
)
189 block_bc(cr
, atoms
->nr
);
190 snew_bc(cr
, atoms
->atom
, atoms
->nr
);
191 nblock_bc(cr
, atoms
->nr
, atoms
->atom
);
192 bc_strings(cr
, symtab
, atoms
->nr
, &atoms
->atomname
);
193 block_bc(cr
, atoms
->nres
);
194 snew_bc(cr
, atoms
->resinfo
, atoms
->nres
);
195 nblock_bc(cr
, atoms
->nres
, atoms
->resinfo
);
196 bc_strings_resinfo(cr
, symtab
, atoms
->nres
, atoms
->resinfo
);
197 /* QMMM requires atomtypes to be known on all nodes as well */
198 bc_strings(cr
, symtab
, atoms
->nr
, &atoms
->atomtype
);
199 bc_strings(cr
, symtab
, atoms
->nr
, &atoms
->atomtypeB
);
202 static void bc_groups(const t_commrec
*cr
, t_symtab
*symtab
,
203 int natoms
, gmx_groups_t
*groups
)
207 bc_grps(cr
, groups
->grps
);
208 block_bc(cr
, groups
->ngrpname
);
209 bc_strings(cr
, symtab
, groups
->ngrpname
, &groups
->grpname
);
210 for (g
= 0; g
< egcNR
; g
++)
214 if (groups
->grpnr
[g
])
226 groups
->grpnr
[g
] = NULL
;
230 snew_bc(cr
, groups
->grpnr
[g
], n
);
231 nblock_bc(cr
, n
, groups
->grpnr
[g
]);
236 fprintf(debug
, "after bc_groups\n");
240 void bcast_state(const t_commrec
*cr
, t_state
*state
)
245 if (!PAR(cr
) || (cr
->nnodes
- cr
->npmenodes
<= 1))
250 /* Broadcasts the state sizes and flags from the master to all nodes
251 * in cr->mpi_comm_mygroup. The arrays are not broadcasted. */
252 block_bc(cr
, state
->natoms
);
253 block_bc(cr
, state
->ngtc
);
254 block_bc(cr
, state
->nnhpres
);
255 block_bc(cr
, state
->nhchainlength
);
256 block_bc(cr
, state
->flags
);
257 if (state
->lambda
== NULL
)
259 snew_bc(cr
, state
->lambda
, efptNR
)
264 /* We allocate dynamically in dd_partition_system. */
267 /* The code below is reachable only by TPI and NM, so it is not
268 tested by anything. */
270 nnht
= (state
->ngtc
)*(state
->nhchainlength
);
271 nnhtp
= (state
->nnhpres
)*(state
->nhchainlength
);
273 /* We still need to allocate the arrays in state for non-master
274 * ranks, which is done (implicitly via bAlloc) in the dirty,
275 * dirty nblock_abc macro. */
276 bAlloc
= !MASTER(cr
);
279 state
->nalloc
= state
->natoms
;
281 for (i
= 0; i
< estNR
; i
++)
283 if (state
->flags
& (1<<i
))
287 case estLAMBDA
: nblock_bc(cr
, efptNR
, state
->lambda
); break;
288 case estFEPSTATE
: block_bc(cr
, state
->fep_state
); break;
289 case estBOX
: block_bc(cr
, state
->box
); break;
290 case estBOX_REL
: block_bc(cr
, state
->box_rel
); break;
291 case estBOXV
: block_bc(cr
, state
->boxv
); break;
292 case estPRES_PREV
: block_bc(cr
, state
->pres_prev
); break;
293 case estSVIR_PREV
: block_bc(cr
, state
->svir_prev
); break;
294 case estFVIR_PREV
: block_bc(cr
, state
->fvir_prev
); break;
295 case estNH_XI
: nblock_abc(cr
, nnht
, state
->nosehoover_xi
); break;
296 case estNH_VXI
: nblock_abc(cr
, nnht
, state
->nosehoover_vxi
); break;
297 case estNHPRES_XI
: nblock_abc(cr
, nnhtp
, state
->nhpres_xi
); break;
298 case estNHPRES_VXI
: nblock_abc(cr
, nnhtp
, state
->nhpres_vxi
); break;
299 case estTC_INT
: nblock_abc(cr
, state
->ngtc
, state
->therm_integral
); break;
300 case estVETA
: block_bc(cr
, state
->veta
); break;
301 case estVOL0
: block_bc(cr
, state
->vol0
); break;
302 case estX
: nblock_abc(cr
, state
->natoms
, state
->x
); break;
303 case estV
: nblock_abc(cr
, state
->natoms
, state
->v
); break;
304 case estSDX
: nblock_abc(cr
, state
->natoms
, state
->sd_X
); break;
305 case estCGP
: nblock_abc(cr
, state
->natoms
, state
->cg_p
); break;
306 case estDISRE_INITF
: block_bc(cr
, state
->hist
.disre_initf
); break;
307 case estDISRE_RM3TAV
:
308 block_bc(cr
, state
->hist
.ndisrepairs
);
309 nblock_abc(cr
, state
->hist
.ndisrepairs
, state
->hist
.disre_rm3tav
);
311 case estORIRE_INITF
: block_bc(cr
, state
->hist
.orire_initf
); break;
313 block_bc(cr
, state
->hist
.norire_Dtav
);
314 nblock_abc(cr
, state
->hist
.norire_Dtav
, state
->hist
.orire_Dtav
);
318 "Communication is not implemented for %s in bcast_state",
325 static void bc_ilists(const t_commrec
*cr
, t_ilist
*ilist
)
329 /* Here we only communicate the non-zero length ilists */
332 for (ftype
= 0; ftype
< F_NRE
; ftype
++)
334 if (ilist
[ftype
].nr
> 0)
337 block_bc(cr
, ilist
[ftype
].nr
);
338 nblock_bc(cr
, ilist
[ftype
].nr
, ilist
[ftype
].iatoms
);
346 for (ftype
= 0; ftype
< F_NRE
; ftype
++)
355 block_bc(cr
, ilist
[ftype
].nr
);
356 snew_bc(cr
, ilist
[ftype
].iatoms
, ilist
[ftype
].nr
);
357 nblock_bc(cr
, ilist
[ftype
].nr
, ilist
[ftype
].iatoms
);
365 fprintf(debug
, "after bc_ilists\n");
369 static void bc_cmap(const t_commrec
*cr
, gmx_cmap_t
*cmap_grid
)
373 block_bc(cr
, cmap_grid
->ngrid
);
374 block_bc(cr
, cmap_grid
->grid_spacing
);
376 ngrid
= cmap_grid
->ngrid
;
377 nelem
= cmap_grid
->grid_spacing
* cmap_grid
->grid_spacing
;
381 snew_bc(cr
, cmap_grid
->cmapdata
, ngrid
);
383 for (i
= 0; i
< ngrid
; i
++)
385 snew_bc(cr
, cmap_grid
->cmapdata
[i
].cmap
, 4*nelem
);
386 nblock_bc(cr
, 4*nelem
, cmap_grid
->cmapdata
[i
].cmap
);
391 static void bc_ffparams(const t_commrec
*cr
, gmx_ffparams_t
*ffp
)
393 block_bc(cr
, ffp
->ntypes
);
394 block_bc(cr
, ffp
->atnr
);
395 snew_bc(cr
, ffp
->functype
, ffp
->ntypes
);
396 snew_bc(cr
, ffp
->iparams
, ffp
->ntypes
);
397 nblock_bc(cr
, ffp
->ntypes
, ffp
->functype
);
398 nblock_bc(cr
, ffp
->ntypes
, ffp
->iparams
);
399 block_bc(cr
, ffp
->reppow
);
400 block_bc(cr
, ffp
->fudgeQQ
);
401 bc_cmap(cr
, &ffp
->cmap_grid
);
404 static void bc_grpopts(const t_commrec
*cr
, t_grpopts
*g
)
408 block_bc(cr
, g
->ngtc
);
409 block_bc(cr
, g
->ngacc
);
410 block_bc(cr
, g
->ngfrz
);
411 block_bc(cr
, g
->ngener
);
412 snew_bc(cr
, g
->nrdf
, g
->ngtc
);
413 snew_bc(cr
, g
->tau_t
, g
->ngtc
);
414 snew_bc(cr
, g
->ref_t
, g
->ngtc
);
415 snew_bc(cr
, g
->acc
, g
->ngacc
);
416 snew_bc(cr
, g
->nFreeze
, g
->ngfrz
);
417 snew_bc(cr
, g
->egp_flags
, g
->ngener
*g
->ngener
);
419 nblock_bc(cr
, g
->ngtc
, g
->nrdf
);
420 nblock_bc(cr
, g
->ngtc
, g
->tau_t
);
421 nblock_bc(cr
, g
->ngtc
, g
->ref_t
);
422 nblock_bc(cr
, g
->ngacc
, g
->acc
);
423 nblock_bc(cr
, g
->ngfrz
, g
->nFreeze
);
424 nblock_bc(cr
, g
->ngener
*g
->ngener
, g
->egp_flags
);
425 snew_bc(cr
, g
->annealing
, g
->ngtc
);
426 snew_bc(cr
, g
->anneal_npoints
, g
->ngtc
);
427 snew_bc(cr
, g
->anneal_time
, g
->ngtc
);
428 snew_bc(cr
, g
->anneal_temp
, g
->ngtc
);
429 nblock_bc(cr
, g
->ngtc
, g
->annealing
);
430 nblock_bc(cr
, g
->ngtc
, g
->anneal_npoints
);
431 for (i
= 0; (i
< g
->ngtc
); i
++)
433 n
= g
->anneal_npoints
[i
];
436 snew_bc(cr
, g
->anneal_time
[i
], n
);
437 snew_bc(cr
, g
->anneal_temp
[i
], n
);
438 nblock_bc(cr
, n
, g
->anneal_time
[i
]);
439 nblock_bc(cr
, n
, g
->anneal_temp
[i
]);
443 /* QMMM stuff, see inputrec */
444 block_bc(cr
, g
->ngQM
);
445 snew_bc(cr
, g
->QMmethod
, g
->ngQM
);
446 snew_bc(cr
, g
->QMbasis
, g
->ngQM
);
447 snew_bc(cr
, g
->QMcharge
, g
->ngQM
);
448 snew_bc(cr
, g
->QMmult
, g
->ngQM
);
449 snew_bc(cr
, g
->bSH
, g
->ngQM
);
450 snew_bc(cr
, g
->CASorbitals
, g
->ngQM
);
451 snew_bc(cr
, g
->CASelectrons
, g
->ngQM
);
452 snew_bc(cr
, g
->SAon
, g
->ngQM
);
453 snew_bc(cr
, g
->SAoff
, g
->ngQM
);
454 snew_bc(cr
, g
->SAsteps
, g
->ngQM
);
458 nblock_bc(cr
, g
->ngQM
, g
->QMmethod
);
459 nblock_bc(cr
, g
->ngQM
, g
->QMbasis
);
460 nblock_bc(cr
, g
->ngQM
, g
->QMcharge
);
461 nblock_bc(cr
, g
->ngQM
, g
->QMmult
);
462 nblock_bc(cr
, g
->ngQM
, g
->bSH
);
463 nblock_bc(cr
, g
->ngQM
, g
->CASorbitals
);
464 nblock_bc(cr
, g
->ngQM
, g
->CASelectrons
);
465 nblock_bc(cr
, g
->ngQM
, g
->SAon
);
466 nblock_bc(cr
, g
->ngQM
, g
->SAoff
);
467 nblock_bc(cr
, g
->ngQM
, g
->SAsteps
);
468 /* end of QMMM stuff */
472 static void bc_cosines(const t_commrec
*cr
, t_cosines
*cs
)
475 snew_bc(cr
, cs
->a
, cs
->n
);
476 snew_bc(cr
, cs
->phi
, cs
->n
);
479 nblock_bc(cr
, cs
->n
, cs
->a
);
480 nblock_bc(cr
, cs
->n
, cs
->phi
);
484 static void bc_pull_group(const t_commrec
*cr
, t_pull_group
*pgrp
)
489 snew_bc(cr
, pgrp
->ind
, pgrp
->nat
);
490 nblock_bc(cr
, pgrp
->nat
, pgrp
->ind
);
492 if (pgrp
->nweight
> 0)
494 snew_bc(cr
, pgrp
->weight
, pgrp
->nweight
);
495 nblock_bc(cr
, pgrp
->nweight
, pgrp
->weight
);
499 static void bc_pull(const t_commrec
*cr
, pull_params_t
*pull
)
504 snew_bc(cr
, pull
->group
, pull
->ngroup
);
505 for (g
= 0; g
< pull
->ngroup
; g
++)
507 bc_pull_group(cr
, &pull
->group
[g
]);
509 snew_bc(cr
, pull
->coord
, pull
->ncoord
);
510 nblock_bc(cr
, pull
->ncoord
, pull
->coord
);
513 static void bc_rotgrp(const t_commrec
*cr
, t_rotgrp
*rotg
)
518 snew_bc(cr
, rotg
->ind
, rotg
->nat
);
519 nblock_bc(cr
, rotg
->nat
, rotg
->ind
);
520 snew_bc(cr
, rotg
->x_ref
, rotg
->nat
);
521 nblock_bc(cr
, rotg
->nat
, rotg
->x_ref
);
525 static void bc_rot(const t_commrec
*cr
, t_rot
*rot
)
530 snew_bc(cr
, rot
->grp
, rot
->ngrp
);
531 for (g
= 0; g
< rot
->ngrp
; g
++)
533 bc_rotgrp(cr
, &rot
->grp
[g
]);
537 static void bc_imd(const t_commrec
*cr
, t_IMD
*imd
)
540 snew_bc(cr
, imd
->ind
, imd
->nat
);
541 nblock_bc(cr
, imd
->nat
, imd
->ind
);
544 static void bc_fepvals(const t_commrec
*cr
, t_lambda
*fep
)
548 block_bc(cr
, fep
->nstdhdl
);
549 block_bc(cr
, fep
->init_lambda
);
550 block_bc(cr
, fep
->init_fep_state
);
551 block_bc(cr
, fep
->delta_lambda
);
552 block_bc(cr
, fep
->edHdLPrintEnergy
);
553 block_bc(cr
, fep
->n_lambda
);
554 if (fep
->n_lambda
> 0)
556 snew_bc(cr
, fep
->all_lambda
, efptNR
);
557 nblock_bc(cr
, efptNR
, fep
->all_lambda
);
558 for (i
= 0; i
< efptNR
; i
++)
560 snew_bc(cr
, fep
->all_lambda
[i
], fep
->n_lambda
);
561 nblock_bc(cr
, fep
->n_lambda
, fep
->all_lambda
[i
]);
564 block_bc(cr
, fep
->sc_alpha
);
565 block_bc(cr
, fep
->sc_power
);
566 block_bc(cr
, fep
->sc_r_power
);
567 block_bc(cr
, fep
->sc_sigma
);
568 block_bc(cr
, fep
->sc_sigma_min
);
569 block_bc(cr
, fep
->bScCoul
);
570 nblock_bc(cr
, efptNR
, &(fep
->separate_dvdl
[0]));
571 block_bc(cr
, fep
->dhdl_derivatives
);
572 block_bc(cr
, fep
->dh_hist_size
);
573 block_bc(cr
, fep
->dh_hist_spacing
);
576 fprintf(debug
, "after bc_fepvals\n");
580 static void bc_expandedvals(const t_commrec
*cr
, t_expanded
*expand
, int n_lambda
)
582 block_bc(cr
, expand
->nstexpanded
);
583 block_bc(cr
, expand
->elamstats
);
584 block_bc(cr
, expand
->elmcmove
);
585 block_bc(cr
, expand
->elmceq
);
586 block_bc(cr
, expand
->equil_n_at_lam
);
587 block_bc(cr
, expand
->equil_wl_delta
);
588 block_bc(cr
, expand
->equil_ratio
);
589 block_bc(cr
, expand
->equil_steps
);
590 block_bc(cr
, expand
->equil_samples
);
591 block_bc(cr
, expand
->lmc_seed
);
592 block_bc(cr
, expand
->minvar
);
593 block_bc(cr
, expand
->minvar_const
);
594 block_bc(cr
, expand
->c_range
);
595 block_bc(cr
, expand
->bSymmetrizedTMatrix
);
596 block_bc(cr
, expand
->nstTij
);
597 block_bc(cr
, expand
->lmc_repeats
);
598 block_bc(cr
, expand
->lmc_forced_nstart
);
599 block_bc(cr
, expand
->gibbsdeltalam
);
600 block_bc(cr
, expand
->wl_scale
);
601 block_bc(cr
, expand
->wl_ratio
);
602 block_bc(cr
, expand
->init_wl_delta
);
603 block_bc(cr
, expand
->bInit_weights
);
604 snew_bc(cr
, expand
->init_lambda_weights
, n_lambda
);
605 nblock_bc(cr
, n_lambda
, expand
->init_lambda_weights
);
606 block_bc(cr
, expand
->mc_temp
);
609 fprintf(debug
, "after bc_expandedvals\n");
613 static void bc_simtempvals(const t_commrec
*cr
, t_simtemp
*simtemp
, int n_lambda
)
615 block_bc(cr
, simtemp
->simtemp_low
);
616 block_bc(cr
, simtemp
->simtemp_high
);
617 block_bc(cr
, simtemp
->eSimTempScale
);
618 snew_bc(cr
, simtemp
->temperatures
, n_lambda
);
619 nblock_bc(cr
, n_lambda
, simtemp
->temperatures
);
622 fprintf(debug
, "after bc_simtempvals\n");
627 static void bc_swapions(const t_commrec
*cr
, t_swapcoords
*swap
)
631 /* Broadcast atom indices for split groups, solvent group, and for all user-defined swap groups */
632 snew_bc(cr
, swap
->grp
, swap
->ngrp
);
633 for (int i
= 0; i
< swap
->ngrp
; i
++)
635 t_swapGroup
*g
= &swap
->grp
[i
];
638 snew_bc(cr
, g
->ind
, g
->nat
);
639 nblock_bc(cr
, g
->nat
, g
->ind
);
644 len
= strlen(g
->molname
);
647 snew_bc(cr
, g
->molname
, len
);
648 nblock_bc(cr
, len
, g
->molname
);
653 static void bc_inputrec(const t_commrec
*cr
, t_inputrec
*inputrec
)
657 block_bc(cr
, *inputrec
);
659 bc_grpopts(cr
, &(inputrec
->opts
));
661 /* even if efep is efepNO, we need to initialize to make sure that
662 * n_lambda is set to zero */
664 snew_bc(cr
, inputrec
->fepvals
, 1);
665 if (inputrec
->efep
!= efepNO
|| inputrec
->bSimTemp
)
667 bc_fepvals(cr
, inputrec
->fepvals
);
669 /* need to initialize this as well because of data checked for in the logic */
670 snew_bc(cr
, inputrec
->expandedvals
, 1);
671 if (inputrec
->bExpanded
)
673 bc_expandedvals(cr
, inputrec
->expandedvals
, inputrec
->fepvals
->n_lambda
);
675 snew_bc(cr
, inputrec
->simtempvals
, 1);
676 if (inputrec
->bSimTemp
)
678 bc_simtempvals(cr
, inputrec
->simtempvals
, inputrec
->fepvals
->n_lambda
);
682 snew_bc(cr
, inputrec
->pull
, 1);
683 bc_pull(cr
, inputrec
->pull
);
687 snew_bc(cr
, inputrec
->rot
, 1);
688 bc_rot(cr
, inputrec
->rot
);
692 snew_bc(cr
, inputrec
->imd
, 1);
693 bc_imd(cr
, inputrec
->imd
);
695 for (i
= 0; (i
< DIM
); i
++)
697 bc_cosines(cr
, &(inputrec
->ex
[i
]));
698 bc_cosines(cr
, &(inputrec
->et
[i
]));
700 if (inputrec
->eSwapCoords
!= eswapNO
)
702 snew_bc(cr
, inputrec
->swap
, 1);
703 bc_swapions(cr
, inputrec
->swap
);
707 static void bc_moltype(const t_commrec
*cr
, t_symtab
*symtab
,
708 gmx_moltype_t
*moltype
)
710 bc_string(cr
, symtab
, &moltype
->name
);
711 bc_atoms(cr
, symtab
, &moltype
->atoms
);
714 fprintf(debug
, "after bc_atoms\n");
717 bc_ilists(cr
, moltype
->ilist
);
718 bc_block(cr
, &moltype
->cgs
);
719 bc_blocka(cr
, &moltype
->excls
);
722 static void bc_molblock(const t_commrec
*cr
, gmx_molblock_t
*molb
)
724 block_bc(cr
, molb
->type
);
725 block_bc(cr
, molb
->nmol
);
726 block_bc(cr
, molb
->natoms_mol
);
727 block_bc(cr
, molb
->nposres_xA
);
728 if (molb
->nposres_xA
> 0)
730 snew_bc(cr
, molb
->posres_xA
, molb
->nposres_xA
);
731 nblock_bc(cr
, molb
->nposres_xA
*DIM
, molb
->posres_xA
[0]);
733 block_bc(cr
, molb
->nposres_xB
);
734 if (molb
->nposres_xB
> 0)
736 snew_bc(cr
, molb
->posres_xB
, molb
->nposres_xB
);
737 nblock_bc(cr
, molb
->nposres_xB
*DIM
, molb
->posres_xB
[0]);
741 fprintf(debug
, "after bc_molblock\n");
745 static void bc_atomtypes(const t_commrec
*cr
, t_atomtypes
*atomtypes
)
749 block_bc(cr
, atomtypes
->nr
);
753 snew_bc(cr
, atomtypes
->radius
, nr
);
754 snew_bc(cr
, atomtypes
->vol
, nr
);
755 snew_bc(cr
, atomtypes
->surftens
, nr
);
756 snew_bc(cr
, atomtypes
->gb_radius
, nr
);
757 snew_bc(cr
, atomtypes
->S_hct
, nr
);
759 nblock_bc(cr
, nr
, atomtypes
->radius
);
760 nblock_bc(cr
, nr
, atomtypes
->vol
);
761 nblock_bc(cr
, nr
, atomtypes
->surftens
);
762 nblock_bc(cr
, nr
, atomtypes
->gb_radius
);
763 nblock_bc(cr
, nr
, atomtypes
->S_hct
);
766 /*! \brief Broadcasts ir and mtop from the master to all nodes in
767 * cr->mpi_comm_mygroup. */
769 void bcast_ir_mtop(const t_commrec
*cr
, t_inputrec
*inputrec
, gmx_mtop_t
*mtop
)
774 fprintf(debug
, "in bc_data\n");
776 bc_inputrec(cr
, inputrec
);
779 fprintf(debug
, "after bc_inputrec\n");
781 bc_symtab(cr
, &mtop
->symtab
);
784 fprintf(debug
, "after bc_symtab\n");
786 bc_string(cr
, &mtop
->symtab
, &mtop
->name
);
789 fprintf(debug
, "after bc_name\n");
792 bc_ffparams(cr
, &mtop
->ffparams
);
794 block_bc(cr
, mtop
->nmoltype
);
795 snew_bc(cr
, mtop
->moltype
, mtop
->nmoltype
);
796 for (i
= 0; i
< mtop
->nmoltype
; i
++)
798 bc_moltype(cr
, &mtop
->symtab
, &mtop
->moltype
[i
]);
801 block_bc(cr
, mtop
->bIntermolecularInteractions
);
802 if (mtop
->bIntermolecularInteractions
)
804 snew_bc(cr
, mtop
->intermolecular_ilist
, F_NRE
);
805 bc_ilists(cr
, mtop
->intermolecular_ilist
);
808 block_bc(cr
, mtop
->nmolblock
);
809 snew_bc(cr
, mtop
->molblock
, mtop
->nmolblock
);
810 for (i
= 0; i
< mtop
->nmolblock
; i
++)
812 bc_molblock(cr
, &mtop
->molblock
[i
]);
815 block_bc(cr
, mtop
->natoms
);
817 bc_atomtypes(cr
, &mtop
->atomtypes
);
819 bc_block(cr
, &mtop
->mols
);
820 bc_groups(cr
, &mtop
->symtab
, mtop
->natoms
, &mtop
->groups
);
823 void init_parallel(t_commrec
*cr
, t_inputrec
*inputrec
,
826 bcast_ir_mtop(cr
, inputrec
, mtop
);