1 /* -*- mode: c; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4; c-file-style: "stroustrup"; -*-
3 * This source code is part of
7 * GROningen MAchine for Chemical Simulations
10 * Written by David van der Spoel, Erik Lindahl, Berk Hess, and others.
11 * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
12 * Copyright (c) 2001-2004, The GROMACS development team,
13 * check out http://www.gromacs.org for more information.
15 * This program is free software; you can redistribute it and/or
16 * modify it under the terms of the GNU General Public License
17 * as published by the Free Software Foundation; either version 2
18 * of the License, or (at your option) any later version.
20 * If you want to redistribute modifications, please consider that
21 * scientific software is very special. Version control is crucial -
22 * bugs must be traceable. We will be happy to consider code for
23 * inclusion in the official distribution, but derived work must not
24 * be called official GROMACS. Details are found in the README & COPYING
25 * files - if they are missing, get the official version at www.gromacs.org.
27 * To help us fund GROMACS development, we humbly ask that you cite
28 * the papers on the package - you can find them in the top README file.
30 * For more info, check our website at http://www.gromacs.org
33 * GROningen Mixture of Alchemy and Childrens' Stories
35 /* This file is completely threadsafe - keep it that way! */
47 #include "gmx_fatal.h"
52 #define block_bc(cr, d) gmx_bcast( sizeof(d), &(d),(cr))
53 /* Probably the test for (nr) > 0 in the next macro is only needed
54 * on BlueGene(/L), where IBM's MPI_Bcast will segfault after
55 * dereferencing a null pointer, even when no data is to be transferred. */
56 #define nblock_bc(cr,nr,d) { if ((nr) > 0) gmx_bcast((nr)*sizeof((d)[0]), (d),(cr)); }
57 #define snew_bc(cr,d,nr) { if (!MASTER(cr)) snew((d),(nr)); }
58 /* Dirty macro with bAlloc not as an argument */
59 #define nblock_abc(cr,nr,d) { if (bAlloc) snew((d),(nr)); nblock_bc(cr,(nr),(d)); }
61 static void bc_string(const t_commrec
*cr
,t_symtab
*symtab
,char ***s
)
66 handle
= lookup_symtab(symtab
,*s
);
70 *s
= get_symtab_handle(symtab
,handle
);
74 static void bc_strings(const t_commrec
*cr
,t_symtab
*symtab
,int nr
,char ****nm
)
84 handle
[i
] = lookup_symtab(symtab
,NM
[i
]);
86 nblock_bc(cr
,nr
,handle
);
91 for (i
=0; (i
<nr
); i
++)
92 (*nm
)[i
] = get_symtab_handle(symtab
,handle
[i
]);
97 static void bc_strings_resinfo(const t_commrec
*cr
,t_symtab
*symtab
,
98 int nr
,t_resinfo
*resinfo
)
105 for(i
=0; (i
<nr
); i
++)
106 handle
[i
] = lookup_symtab(symtab
,resinfo
[i
].name
);
108 nblock_bc(cr
,nr
,handle
);
111 for (i
=0; (i
<nr
); i
++)
112 resinfo
[i
].name
= get_symtab_handle(symtab
,handle
[i
]);
117 static void bc_symtab(const t_commrec
*cr
,t_symtab
*symtab
)
122 block_bc(cr
,symtab
->nr
);
124 snew_bc(cr
,symtab
->symbuf
,1);
125 symbuf
= symtab
->symbuf
;
126 symbuf
->bufsize
= nr
;
127 snew_bc(cr
,symbuf
->buf
,nr
);
128 for (i
=0; i
<nr
; i
++) {
130 len
= strlen(symbuf
->buf
[i
]) + 1;
132 snew_bc(cr
,symbuf
->buf
[i
],len
);
133 nblock_bc(cr
,len
,symbuf
->buf
[i
]);
137 static void bc_block(const t_commrec
*cr
,t_block
*block
)
139 block_bc(cr
,block
->nr
);
140 snew_bc(cr
,block
->index
,block
->nr
+1);
141 nblock_bc(cr
,block
->nr
+1,block
->index
);
144 static void bc_blocka(const t_commrec
*cr
,t_blocka
*block
)
146 block_bc(cr
,block
->nr
);
147 snew_bc(cr
,block
->index
,block
->nr
+1);
148 nblock_bc(cr
,block
->nr
+1,block
->index
);
149 block_bc(cr
,block
->nra
);
151 snew_bc(cr
,block
->a
,block
->nra
);
152 nblock_bc(cr
,block
->nra
,block
->a
);
156 static void bc_grps(const t_commrec
*cr
,t_grps grps
[])
160 for(i
=0; (i
<egcNR
); i
++) {
161 block_bc(cr
,grps
[i
].nr
);
162 snew_bc(cr
,grps
[i
].nm_ind
,grps
[i
].nr
);
163 nblock_bc(cr
,grps
[i
].nr
,grps
[i
].nm_ind
);
167 static void bc_atoms(const t_commrec
*cr
,t_symtab
*symtab
,t_atoms
*atoms
)
171 block_bc(cr
,atoms
->nr
);
172 snew_bc(cr
,atoms
->atom
,atoms
->nr
);
173 nblock_bc(cr
,atoms
->nr
,atoms
->atom
);
174 bc_strings(cr
,symtab
,atoms
->nr
,&atoms
->atomname
);
175 block_bc(cr
,atoms
->nres
);
176 snew_bc(cr
,atoms
->resinfo
,atoms
->nres
);
177 nblock_bc(cr
,atoms
->nres
,atoms
->resinfo
);
178 bc_strings_resinfo(cr
,symtab
,atoms
->nres
,atoms
->resinfo
);
179 /* QMMM requires atomtypes to be known on all nodes as well */
180 bc_strings(cr
,symtab
,atoms
->nr
,&atoms
->atomtype
);
181 bc_strings(cr
,symtab
,atoms
->nr
,&atoms
->atomtypeB
);
184 static void bc_groups(const t_commrec
*cr
,t_symtab
*symtab
,
185 int natoms
,gmx_groups_t
*groups
)
190 bc_grps(cr
,groups
->grps
);
191 block_bc(cr
,groups
->ngrpname
);
192 bc_strings(cr
,symtab
,groups
->ngrpname
,&groups
->grpname
);
193 for(g
=0; g
<egcNR
; g
++) {
195 if (groups
->grpnr
[g
]) {
203 groups
->grpnr
[g
] = NULL
;
205 snew_bc(cr
,groups
->grpnr
[g
],n
);
206 nblock_bc(cr
,n
,groups
->grpnr
[g
]);
209 if (debug
) fprintf(debug
,"after bc_groups\n");
212 void bcast_state_setup(const t_commrec
*cr
,t_state
*state
)
214 block_bc(cr
,state
->natoms
);
215 block_bc(cr
,state
->ngtc
);
216 block_bc(cr
,state
->nnhpres
);
217 block_bc(cr
,state
->nhchainlength
);
218 block_bc(cr
,state
->nrng
);
219 block_bc(cr
,state
->nrngi
);
220 block_bc(cr
,state
->flags
);
223 void bcast_state(const t_commrec
*cr
,t_state
*state
,gmx_bool bAlloc
)
227 bcast_state_setup(cr
,state
);
229 nnht
= (state
->ngtc
)*(state
->nhchainlength
);
230 nnhtp
= (state
->nnhpres
)*(state
->nhchainlength
);
236 state
->nalloc
= state
->natoms
;
238 for(i
=0; i
<estNR
; i
++) {
239 if (state
->flags
& (1<<i
)) {
241 case estLAMBDA
: block_bc(cr
,state
->lambda
); break;
242 case estBOX
: block_bc(cr
,state
->box
); break;
243 case estBOX_REL
: block_bc(cr
,state
->box_rel
); break;
244 case estBOXV
: block_bc(cr
,state
->boxv
); break;
245 case estPRES_PREV
: block_bc(cr
,state
->pres_prev
); break;
246 case estSVIR_PREV
: block_bc(cr
,state
->svir_prev
); break;
247 case estFVIR_PREV
: block_bc(cr
,state
->fvir_prev
); break;
248 case estNH_XI
: nblock_abc(cr
,nnht
,state
->nosehoover_xi
); break;
249 case estNH_VXI
: nblock_abc(cr
,nnht
,state
->nosehoover_vxi
); break;
250 case estNHPRES_XI
: nblock_abc(cr
,nnhtp
,state
->nhpres_xi
); break;
251 case estNHPRES_VXI
: nblock_abc(cr
,nnhtp
,state
->nhpres_vxi
); break;
252 case estTC_INT
: nblock_abc(cr
,state
->ngtc
,state
->therm_integral
); break;
253 case estVETA
: block_bc(cr
,state
->veta
); break;
254 case estVOL0
: block_bc(cr
,state
->vol0
); break;
255 case estX
: nblock_abc(cr
,state
->natoms
,state
->x
); break;
256 case estV
: nblock_abc(cr
,state
->natoms
,state
->v
); break;
257 case estSDX
: nblock_abc(cr
,state
->natoms
,state
->sd_X
); break;
258 case estCGP
: nblock_abc(cr
,state
->natoms
,state
->cg_p
); break;
259 case estLD_RNG
: if(state
->nrngi
== 1) nblock_abc(cr
,state
->nrng
,state
->ld_rng
); break;
260 case estLD_RNGI
: if(state
->nrngi
== 1) nblock_abc(cr
,state
->nrngi
,state
->ld_rngi
); break;
261 case estDISRE_INITF
: block_bc(cr
,state
->hist
.disre_initf
); break;
262 case estDISRE_RM3TAV
:
263 block_bc(cr
,state
->hist
.ndisrepairs
);
264 nblock_abc(cr
,state
->hist
.ndisrepairs
,state
->hist
.disre_rm3tav
);
266 case estORIRE_INITF
: block_bc(cr
,state
->hist
.orire_initf
); break;
268 block_bc(cr
,state
->hist
.norire_Dtav
);
269 nblock_abc(cr
,state
->hist
.norire_Dtav
,state
->hist
.orire_Dtav
);
273 "Communication is not implemented for %s in bcast_state",
280 static void bc_ilists(const t_commrec
*cr
,t_ilist
*ilist
)
284 /* Here we only communicate the non-zero length ilists */
286 for(ftype
=0; ftype
<F_NRE
; ftype
++) {
287 if (ilist
[ftype
].nr
> 0) {
289 block_bc(cr
,ilist
[ftype
].nr
);
290 nblock_bc(cr
,ilist
[ftype
].nr
,ilist
[ftype
].iatoms
);
296 for(ftype
=0; ftype
<F_NRE
; ftype
++) {
302 block_bc(cr
,ilist
[ftype
].nr
);
303 snew_bc(cr
,ilist
[ftype
].iatoms
,ilist
[ftype
].nr
);
304 nblock_bc(cr
,ilist
[ftype
].nr
,ilist
[ftype
].iatoms
);
306 } while (ftype
>= 0);
309 if (debug
) fprintf(debug
,"after bc_ilists\n");
312 static void bc_idef(const t_commrec
*cr
,t_idef
*idef
)
314 block_bc(cr
,idef
->ntypes
);
315 block_bc(cr
,idef
->atnr
);
316 snew_bc(cr
,idef
->functype
,idef
->ntypes
);
317 snew_bc(cr
,idef
->iparams
,idef
->ntypes
);
318 nblock_bc(cr
,idef
->ntypes
,idef
->functype
);
319 nblock_bc(cr
,idef
->ntypes
,idef
->iparams
);
320 block_bc(cr
,idef
->fudgeQQ
);
321 bc_ilists(cr
,idef
->il
);
322 block_bc(cr
,idef
->ilsort
);
325 static void bc_cmap(const t_commrec
*cr
, gmx_cmap_t
*cmap_grid
)
329 block_bc(cr
,cmap_grid
->ngrid
);
330 block_bc(cr
,cmap_grid
->grid_spacing
);
332 ngrid
= cmap_grid
->ngrid
;
333 nelem
= cmap_grid
->grid_spacing
* cmap_grid
->grid_spacing
;
337 snew_bc(cr
,cmap_grid
->cmapdata
,ngrid
);
341 snew_bc(cr
,cmap_grid
->cmapdata
[i
].cmap
,4*nelem
);
342 nblock_bc(cr
,4*nelem
,cmap_grid
->cmapdata
[i
].cmap
);
347 static void bc_ffparams(const t_commrec
*cr
,gmx_ffparams_t
*ffp
)
351 block_bc(cr
,ffp
->ntypes
);
352 block_bc(cr
,ffp
->atnr
);
353 snew_bc(cr
,ffp
->functype
,ffp
->ntypes
);
354 snew_bc(cr
,ffp
->iparams
,ffp
->ntypes
);
355 nblock_bc(cr
,ffp
->ntypes
,ffp
->functype
);
356 nblock_bc(cr
,ffp
->ntypes
,ffp
->iparams
);
357 block_bc(cr
,ffp
->reppow
);
358 block_bc(cr
,ffp
->fudgeQQ
);
359 bc_cmap(cr
,&ffp
->cmap_grid
);
362 static void bc_grpopts(const t_commrec
*cr
,t_grpopts
*g
)
366 block_bc(cr
,g
->ngtc
);
367 block_bc(cr
,g
->ngacc
);
368 block_bc(cr
,g
->ngfrz
);
369 block_bc(cr
,g
->ngener
);
370 snew_bc(cr
,g
->nrdf
,g
->ngtc
);
371 snew_bc(cr
,g
->tau_t
,g
->ngtc
);
372 snew_bc(cr
,g
->ref_t
,g
->ngtc
);
373 snew_bc(cr
,g
->acc
,g
->ngacc
);
374 snew_bc(cr
,g
->nFreeze
,g
->ngfrz
);
375 snew_bc(cr
,g
->egp_flags
,g
->ngener
*g
->ngener
);
377 nblock_bc(cr
,g
->ngtc
,g
->nrdf
);
378 nblock_bc(cr
,g
->ngtc
,g
->tau_t
);
379 nblock_bc(cr
,g
->ngtc
,g
->ref_t
);
380 nblock_bc(cr
,g
->ngacc
,g
->acc
);
381 nblock_bc(cr
,g
->ngfrz
,g
->nFreeze
);
382 nblock_bc(cr
,g
->ngener
*g
->ngener
,g
->egp_flags
);
383 snew_bc(cr
,g
->annealing
,g
->ngtc
);
384 snew_bc(cr
,g
->anneal_npoints
,g
->ngtc
);
385 snew_bc(cr
,g
->anneal_time
,g
->ngtc
);
386 snew_bc(cr
,g
->anneal_temp
,g
->ngtc
);
387 nblock_bc(cr
,g
->ngtc
,g
->annealing
);
388 nblock_bc(cr
,g
->ngtc
,g
->anneal_npoints
);
389 for(i
=0;(i
<g
->ngtc
); i
++) {
390 n
= g
->anneal_npoints
[i
];
392 snew_bc(cr
,g
->anneal_time
[i
],n
);
393 snew_bc(cr
,g
->anneal_temp
[i
],n
);
394 nblock_bc(cr
,n
,g
->anneal_time
[i
]);
395 nblock_bc(cr
,n
,g
->anneal_temp
[i
]);
399 /* QMMM stuff, see inputrec */
400 block_bc(cr
,g
->ngQM
);
401 snew_bc(cr
,g
->QMmethod
,g
->ngQM
);
402 snew_bc(cr
,g
->QMbasis
,g
->ngQM
);
403 snew_bc(cr
,g
->QMcharge
,g
->ngQM
);
404 snew_bc(cr
,g
->QMmult
,g
->ngQM
);
405 snew_bc(cr
,g
->bSH
,g
->ngQM
);
406 snew_bc(cr
,g
->CASorbitals
,g
->ngQM
);
407 snew_bc(cr
,g
->CASelectrons
,g
->ngQM
);
408 snew_bc(cr
,g
->SAon
,g
->ngQM
);
409 snew_bc(cr
,g
->SAoff
,g
->ngQM
);
410 snew_bc(cr
,g
->SAsteps
,g
->ngQM
);
414 nblock_bc(cr
,g
->ngQM
,g
->QMmethod
);
415 nblock_bc(cr
,g
->ngQM
,g
->QMbasis
);
416 nblock_bc(cr
,g
->ngQM
,g
->QMcharge
);
417 nblock_bc(cr
,g
->ngQM
,g
->QMmult
);
418 nblock_bc(cr
,g
->ngQM
,g
->bSH
);
419 nblock_bc(cr
,g
->ngQM
,g
->CASorbitals
);
420 nblock_bc(cr
,g
->ngQM
,g
->CASelectrons
);
421 nblock_bc(cr
,g
->ngQM
,g
->SAon
);
422 nblock_bc(cr
,g
->ngQM
,g
->SAoff
);
423 nblock_bc(cr
,g
->ngQM
,g
->SAsteps
);
424 /* end of QMMM stuff */
428 static void bc_cosines(const t_commrec
*cr
,t_cosines
*cs
)
431 snew_bc(cr
,cs
->a
,cs
->n
);
432 snew_bc(cr
,cs
->phi
,cs
->n
);
434 nblock_bc(cr
,cs
->n
,cs
->a
);
435 nblock_bc(cr
,cs
->n
,cs
->phi
);
439 static void bc_pullgrp(const t_commrec
*cr
,t_pullgrp
*pgrp
)
443 snew_bc(cr
,pgrp
->ind
,pgrp
->nat
);
444 nblock_bc(cr
,pgrp
->nat
,pgrp
->ind
);
446 if (pgrp
->nweight
> 0) {
447 snew_bc(cr
,pgrp
->weight
,pgrp
->nweight
);
448 nblock_bc(cr
,pgrp
->nweight
,pgrp
->weight
);
452 static void bc_pull(const t_commrec
*cr
,t_pull
*pull
)
457 snew_bc(cr
,pull
->grp
,pull
->ngrp
+1);
458 for(g
=0; g
<pull
->ngrp
+1; g
++)
460 bc_pullgrp(cr
,&pull
->grp
[g
]);
464 static void bc_rotgrp(const t_commrec
*cr
,t_rotgrp
*rotg
)
468 snew_bc(cr
,rotg
->ind
,rotg
->nat
);
469 nblock_bc(cr
,rotg
->nat
,rotg
->ind
);
470 snew_bc(cr
,rotg
->x_ref
,rotg
->nat
);
471 nblock_bc(cr
,rotg
->nat
,rotg
->x_ref
);
475 static void bc_rot(const t_commrec
*cr
,t_rot
*rot
)
480 snew_bc(cr
,rot
->grp
,rot
->ngrp
);
481 for(g
=0; g
<rot
->ngrp
; g
++)
482 bc_rotgrp(cr
,&rot
->grp
[g
]);
485 static void bc_inputrec(const t_commrec
*cr
,t_inputrec
*inputrec
)
487 gmx_bool bAlloc
=TRUE
;
490 block_bc(cr
,*inputrec
);
491 snew_bc(cr
,inputrec
->flambda
,inputrec
->n_flambda
);
492 nblock_bc(cr
,inputrec
->n_flambda
,inputrec
->flambda
);
493 bc_grpopts(cr
,&(inputrec
->opts
));
494 if (inputrec
->ePull
!= epullNO
) {
495 snew_bc(cr
,inputrec
->pull
,1);
496 bc_pull(cr
,inputrec
->pull
);
498 if (inputrec
->bRot
) {
499 snew_bc(cr
,inputrec
->rot
,1);
500 bc_rot(cr
,inputrec
->rot
);
502 for(i
=0; (i
<DIM
); i
++) {
503 bc_cosines(cr
,&(inputrec
->ex
[i
]));
504 bc_cosines(cr
,&(inputrec
->et
[i
]));
506 if (inputrec
->n_adress_tf_grps
> 0){
507 snew_bc(cr
,inputrec
->adress_tf_table_index
,inputrec
->n_adress_tf_grps
);
508 nblock_bc(cr
,inputrec
->n_adress_tf_grps
,inputrec
->adress_tf_table_index
);
510 if (inputrec
->n_energy_grps
> 0){
511 snew_bc(cr
,inputrec
->adress_group_explicit
,inputrec
->n_energy_grps
);
512 nblock_bc(cr
,inputrec
->n_energy_grps
,inputrec
->adress_group_explicit
);
516 static void bc_moltype(const t_commrec
*cr
,t_symtab
*symtab
,
517 gmx_moltype_t
*moltype
)
519 bc_string(cr
,symtab
,&moltype
->name
);
520 bc_atoms(cr
,symtab
,&moltype
->atoms
);
521 if (debug
) fprintf(debug
,"after bc_atoms\n");
523 bc_ilists(cr
,moltype
->ilist
);
524 bc_block(cr
,&moltype
->cgs
);
525 bc_blocka(cr
,&moltype
->excls
);
528 static void bc_molblock(const t_commrec
*cr
,gmx_molblock_t
*molb
)
530 gmx_bool bAlloc
=TRUE
;
532 block_bc(cr
,molb
->type
);
533 block_bc(cr
,molb
->nmol
);
534 block_bc(cr
,molb
->natoms_mol
);
535 block_bc(cr
,molb
->nposres_xA
);
536 if (molb
->nposres_xA
> 0) {
537 snew_bc(cr
,molb
->posres_xA
,molb
->nposres_xA
);
538 nblock_bc(cr
,molb
->nposres_xA
*DIM
,molb
->posres_xA
[0]);
540 block_bc(cr
,molb
->nposres_xB
);
541 if (molb
->nposres_xB
> 0) {
542 snew_bc(cr
,molb
->posres_xB
,molb
->nposres_xB
);
543 nblock_bc(cr
,molb
->nposres_xB
*DIM
,molb
->posres_xB
[0]);
545 if (debug
) fprintf(debug
,"after bc_molblock\n");
548 static void bc_atomtypes(const t_commrec
*cr
, t_atomtypes
*atomtypes
)
552 block_bc(cr
,atomtypes
->nr
);
556 snew_bc(cr
,atomtypes
->radius
,nr
);
557 snew_bc(cr
,atomtypes
->vol
,nr
);
558 snew_bc(cr
,atomtypes
->surftens
,nr
);
559 snew_bc(cr
,atomtypes
->gb_radius
,nr
);
560 snew_bc(cr
,atomtypes
->S_hct
,nr
);
562 nblock_bc(cr
,nr
,atomtypes
->radius
);
563 nblock_bc(cr
,nr
,atomtypes
->vol
);
564 nblock_bc(cr
,nr
,atomtypes
->surftens
);
565 nblock_bc(cr
,nr
,atomtypes
->gb_radius
);
566 nblock_bc(cr
,nr
,atomtypes
->S_hct
);
570 void bcast_ir_mtop(const t_commrec
*cr
,t_inputrec
*inputrec
,gmx_mtop_t
*mtop
)
573 if (debug
) fprintf(debug
,"in bc_data\n");
574 bc_inputrec(cr
,inputrec
);
575 if (debug
) fprintf(debug
,"after bc_inputrec\n");
576 bc_symtab(cr
,&mtop
->symtab
);
577 if (debug
) fprintf(debug
,"after bc_symtab\n");
578 bc_string(cr
,&mtop
->symtab
,&mtop
->name
);
579 if (debug
) fprintf(debug
,"after bc_name\n");
581 bc_ffparams(cr
,&mtop
->ffparams
);
583 block_bc(cr
,mtop
->nmoltype
);
584 snew_bc(cr
,mtop
->moltype
,mtop
->nmoltype
);
585 for(i
=0; i
<mtop
->nmoltype
; i
++) {
586 bc_moltype(cr
,&mtop
->symtab
,&mtop
->moltype
[i
]);
589 block_bc(cr
,mtop
->nmolblock
);
590 snew_bc(cr
,mtop
->molblock
,mtop
->nmolblock
);
591 for(i
=0; i
<mtop
->nmolblock
; i
++) {
592 bc_molblock(cr
,&mtop
->molblock
[i
]);
595 block_bc(cr
,mtop
->natoms
);
597 bc_atomtypes(cr
,&mtop
->atomtypes
);
599 bc_block(cr
,&mtop
->mols
);
600 bc_groups(cr
,&mtop
->symtab
,mtop
->natoms
,&mtop
->groups
);