1 /* -*- mode: c; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4; c-file-style: "stroustrup"; -*-
4 * This source code is part of
8 * GROningen MAchine for Chemical Simulations
11 * Written by David van der Spoel, Erik Lindahl, Berk Hess, and others.
12 * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
13 * Copyright (c) 2001-2004, The GROMACS development team,
14 * check out http://www.gromacs.org for more information.
16 * This program is free software; you can redistribute it and/or
17 * modify it under the terms of the GNU General Public License
18 * as published by the Free Software Foundation; either version 2
19 * of the License, or (at your option) any later version.
21 * If you want to redistribute modifications, please consider that
22 * scientific software is very special. Version control is crucial -
23 * bugs must be traceable. We will be happy to consider code for
24 * inclusion in the official distribution, but derived work must not
25 * be called official GROMACS. Details are found in the README & COPYING
26 * files - if they are missing, get the official version at www.gromacs.org.
28 * To help us fund GROMACS development, we humbly ask that you cite
29 * the papers on the package - you can find them in the top README file.
31 * For more info, check our website at http://www.gromacs.org
34 * Gallium Rubidium Oxygen Manganese Argon Carbon Silicon
43 #if ((defined WIN32 || defined _WIN32 || defined WIN64 || defined _WIN64) && !defined __CYGWIN__ && !defined __CYGWIN32__)
64 #include "mpelogging.h"
70 #include "checkpoint.h"
71 #include "mtop_util.h"
72 #include "sighandler.h"
75 #include "pull_rotation.h"
76 #include "md_openmm.h"
90 #include "md_openmm.h"
95 gmx_integrator_t
*func
;
98 /* The array should match the eI array in include/types/enums.h */
99 #ifdef GMX_OPENMM /* FIXME do_md_openmm needs fixing */
100 const gmx_intp_t integrator
[eiNR
] = { {do_md_openmm
}, {do_md_openmm
}, {do_md_openmm
}, {do_md_openmm
}, {do_md_openmm
}, {do_md_openmm
}, {do_md_openmm
}, {do_md_openmm
}, {do_md_openmm
}, {do_md_openmm
}, {do_md_openmm
},{do_md_openmm
}};
102 const gmx_intp_t integrator
[eiNR
] = { {do_md
}, {do_steep
}, {do_cg
}, {do_md
}, {do_md
}, {do_nm
}, {do_lbfgs
}, {do_tpi
}, {do_tpi
}, {do_md
}, {do_md
},{do_md
}};
105 gmx_large_int_t deform_init_init_step_tpx
;
106 matrix deform_init_box_tpx
;
108 tMPI_Thread_mutex_t deform_init_box_mutex
=TMPI_THREAD_MUTEX_INITIALIZER
;
113 struct mdrunner_arglist
127 const char *dddlb_opt
;
140 const char *deviceOptions
;
142 int ret
; /* return value */
146 /* The function used for spawning threads. Extracts the mdrunner()
147 arguments from its one argument and calls mdrunner(), after making
149 static void mdrunner_start_fn(void *arg
)
151 struct mdrunner_arglist
*mda
=(struct mdrunner_arglist
*)arg
;
152 struct mdrunner_arglist mc
=*mda
; /* copy the arg list to make sure
153 that it's thread-local. This doesn't
154 copy pointed-to items, of course,
155 but those are all const. */
156 t_commrec
*cr
; /* we need a local version of this */
160 fnm
= dup_tfn(mc
.nfile
, mc
.fnm
);
162 cr
= init_par_threads(mc
.cr
);
169 mda
->ret
=mdrunner(cr
->nnodes
, fplog
, cr
, mc
.nfile
, fnm
, mc
.oenv
,
170 mc
.bVerbose
, mc
.bCompact
, mc
.nstglobalcomm
,
171 mc
.ddxyz
, mc
.dd_node_order
, mc
.rdd
,
172 mc
.rconstr
, mc
.dddlb_opt
, mc
.dlb_scale
,
173 mc
.ddcsx
, mc
.ddcsy
, mc
.ddcsz
, mc
.nstepout
, mc
.resetstep
,
174 mc
.nmultisim
, mc
.repl_ex_nst
, mc
.repl_ex_seed
, mc
.pforce
,
175 mc
.cpt_period
, mc
.max_hours
, mc
.deviceOptions
, mc
.Flags
);
178 /* called by mdrunner() to start a specific number of threads (including
179 the main thread) for thread-parallel runs. This in turn calls mdrunner()
181 All options besides nthreads are the same as for mdrunner(). */
182 static t_commrec
*mdrunner_start_threads(int nthreads
,
183 FILE *fplog
,t_commrec
*cr
,int nfile
,
184 const t_filenm fnm
[], const output_env_t oenv
, gmx_bool bVerbose
,
185 gmx_bool bCompact
, int nstglobalcomm
,
186 ivec ddxyz
,int dd_node_order
,real rdd
,real rconstr
,
187 const char *dddlb_opt
,real dlb_scale
,
188 const char *ddcsx
,const char *ddcsy
,const char *ddcsz
,
189 int nstepout
,int resetstep
,int nmultisim
,int repl_ex_nst
,
190 int repl_ex_seed
, real pforce
,real cpt_period
, real max_hours
,
191 const char *deviceOptions
, unsigned long Flags
)
194 struct mdrunner_arglist
*mda
;
195 t_commrec
*crn
; /* the new commrec */
198 /* first check whether we even need to start tMPI */
202 /* a few small, one-time, almost unavoidable memory leaks: */
204 fnmn
=dup_tfn(nfile
, fnm
);
206 /* fill the data structure to pass as void pointer to thread start fn */
212 mda
->bVerbose
=bVerbose
;
213 mda
->bCompact
=bCompact
;
214 mda
->nstglobalcomm
=nstglobalcomm
;
215 mda
->ddxyz
[XX
]=ddxyz
[XX
];
216 mda
->ddxyz
[YY
]=ddxyz
[YY
];
217 mda
->ddxyz
[ZZ
]=ddxyz
[ZZ
];
218 mda
->dd_node_order
=dd_node_order
;
220 mda
->rconstr
=rconstr
;
221 mda
->dddlb_opt
=dddlb_opt
;
222 mda
->dlb_scale
=dlb_scale
;
226 mda
->nstepout
=nstepout
;
227 mda
->resetstep
=resetstep
;
228 mda
->nmultisim
=nmultisim
;
229 mda
->repl_ex_nst
=repl_ex_nst
;
230 mda
->repl_ex_seed
=repl_ex_seed
;
232 mda
->cpt_period
=cpt_period
;
233 mda
->max_hours
=max_hours
;
234 mda
->deviceOptions
=deviceOptions
;
237 fprintf(stderr
, "Starting %d threads\n",nthreads
);
239 /* now spawn new threads that start mdrunner_start_fn(), while
240 the main thread returns */
241 ret
=tMPI_Init_fn(TRUE
, nthreads
, mdrunner_start_fn
, (void*)(mda
) );
242 if (ret
!=TMPI_SUCCESS
)
245 /* make a new comm_rec to reflect the new situation */
246 crn
=init_par_threads(cr
);
251 /* get the number of threads based on how many there were requested,
252 which algorithms we're using, and how many particles there are. */
253 static int get_nthreads(int nthreads_requested
, t_inputrec
*inputrec
,
256 int nthreads
,nthreads_new
;
257 int min_atoms_per_thread
;
260 nthreads
= nthreads_requested
;
262 /* determine # of hardware threads. */
263 if (nthreads_requested
< 1)
265 if ((env
= getenv("GMX_MAX_THREADS")) != NULL
)
268 sscanf(env
,"%d",&nthreads
);
271 gmx_fatal(FARGS
,"GMX_MAX_THREADS (%d) should be larger than 0",
277 nthreads
= tMPI_Get_recommended_nthreads();
281 if (inputrec
->eI
== eiNM
|| EI_TPI(inputrec
->eI
))
283 /* Steps are divided over the nodes iso splitting the atoms */
284 min_atoms_per_thread
= 0;
288 min_atoms_per_thread
= MIN_ATOMS_PER_THREAD
;
291 /* Check if an algorithm does not support parallel simulation. */
293 ( inputrec
->eI
== eiLBFGS
||
294 inputrec
->coulombtype
== eelEWALD
) )
296 fprintf(stderr
,"\nThe integration or electrostatics algorithm doesn't support parallel runs. Not starting any threads.\n");
299 else if (nthreads_requested
< 1 &&
300 mtop
->natoms
/nthreads
< min_atoms_per_thread
)
302 /* the thread number was chosen automatically, but there are too many
303 threads (too few atoms per thread) */
304 nthreads_new
= max(1,mtop
->natoms
/min_atoms_per_thread
);
306 if (nthreads_new
> 8 || (nthreads
== 8 && nthreads_new
> 4))
308 /* Use only multiples of 4 above 8 threads
309 * or with an 8-core processor
310 * (to avoid 6 threads on 8 core processors with 4 real cores).
312 nthreads_new
= (nthreads_new
/4)*4;
314 else if (nthreads_new
> 4)
316 /* Avoid 5 or 7 threads */
317 nthreads_new
= (nthreads_new
/2)*2;
320 nthreads
= nthreads_new
;
322 fprintf(stderr
,"\n");
323 fprintf(stderr
,"NOTE: Parallelization is limited by the small number of atoms,\n");
324 fprintf(stderr
," only starting %d threads.\n",nthreads
);
325 fprintf(stderr
," You can use the -nt option to optimize the number of threads.\n\n");
332 int mdrunner(int nthreads_requested
, FILE *fplog
,t_commrec
*cr
,int nfile
,
333 const t_filenm fnm
[], const output_env_t oenv
, gmx_bool bVerbose
,
334 gmx_bool bCompact
, int nstglobalcomm
,
335 ivec ddxyz
,int dd_node_order
,real rdd
,real rconstr
,
336 const char *dddlb_opt
,real dlb_scale
,
337 const char *ddcsx
,const char *ddcsy
,const char *ddcsz
,
338 int nstepout
,int resetstep
,int nmultisim
,int repl_ex_nst
,
339 int repl_ex_seed
, real pforce
,real cpt_period
,real max_hours
,
340 const char *deviceOptions
, unsigned long Flags
)
342 double nodetime
=0,realtime
;
343 t_inputrec
*inputrec
;
346 gmx_ddbox_t ddbox
={0};
347 int npme_major
,npme_minor
;
350 gmx_mtop_t
*mtop
=NULL
;
351 t_mdatoms
*mdatoms
=NULL
;
355 gmx_pme_t
*pmedata
=NULL
;
356 gmx_vsite_t
*vsite
=NULL
;
358 int i
,m
,nChargePerturbed
=-1,status
,nalloc
;
360 gmx_wallcycle_t wcycle
;
361 gmx_bool bReadRNG
,bReadEkin
;
363 gmx_runtime_t runtime
;
365 gmx_large_int_t reset_counters
;
367 t_commrec
*cr_old
=cr
;
370 /* CAUTION: threads may be started later on in this function, so
371 cr doesn't reflect the final parallel state right now */
375 if (bVerbose
&& SIMMASTER(cr
))
377 fprintf(stderr
,"Getting Loaded...\n");
380 if (Flags
& MD_APPENDFILES
)
388 /* Read (nearly) all data required for the simulation */
389 read_tpx_state(ftp2fn(efTPX
,nfile
,fnm
),inputrec
,state
,NULL
,mtop
);
391 /* NOW the threads will be started: */
393 nthreads
= get_nthreads(nthreads_requested
, inputrec
, mtop
);
397 /* now start the threads. */
398 cr
=mdrunner_start_threads(nthreads
, fplog
, cr_old
, nfile
, fnm
,
399 oenv
, bVerbose
, bCompact
, nstglobalcomm
,
400 ddxyz
, dd_node_order
, rdd
, rconstr
,
401 dddlb_opt
, dlb_scale
, ddcsx
, ddcsy
, ddcsz
,
402 nstepout
, resetstep
, nmultisim
,
403 repl_ex_nst
, repl_ex_seed
, pforce
,
404 cpt_period
, max_hours
, deviceOptions
,
406 /* the main thread continues here with a new cr. We don't deallocate
407 the old cr because other threads may still be reading it. */
410 gmx_comm("Failed to spawn threads");
415 /* END OF CAUTION: cr is now reliable */
419 /* now broadcast everything to the non-master nodes/threads: */
420 init_parallel(fplog
, cr
, inputrec
, mtop
);
424 pr_inputrec(fplog
,0,"Input Parameters",inputrec
,FALSE
);
427 /* now make sure the state is initialized and propagated */
428 set_state_entries(state
,inputrec
,cr
->nnodes
);
430 /* A parallel command line option consistency check that we can
431 only do after any threads have started. */
433 (ddxyz
[XX
] > 1 || ddxyz
[YY
] > 1 || ddxyz
[ZZ
] > 1 || cr
->npmenodes
> 0))
436 "The -dd or -npme option request a parallel simulation, "
438 "but mdrun was compiled without threads or MPI enabled"
441 "but the number of threads (option -nt) is 1"
443 "but mdrun was not started through mpirun/mpiexec or only one process was requested through mpirun/mpiexec"
449 if (can_use_allvsall(inputrec
,mtop
,TRUE
,cr
,fplog
))
451 /* All-vs-all loops do not work with domain decomposition */
455 if (!EEL_PME(inputrec
->coulombtype
) || (Flags
& MD_PARTDEC
))
461 fcRegisterSteps(inputrec
->nsteps
,inputrec
->init_step
);
464 /* NMR restraints must be initialized before load_checkpoint,
465 * since with time averaging the history is added to t_state.
466 * For proper consistency check we therefore need to extend
468 * So the PME-only nodes (if present) will also initialize
469 * the distance restraints.
473 /* This needs to be called before read_checkpoint to extend the state */
474 init_disres(fplog
,mtop
,inputrec
,cr
,Flags
& MD_PARTDEC
,fcd
,state
);
476 if (gmx_mtop_ftype_count(mtop
,F_ORIRES
) > 0)
478 if (PAR(cr
) && !(Flags
& MD_PARTDEC
))
480 gmx_fatal(FARGS
,"Orientation restraints do not work (yet) with domain decomposition, use particle decomposition (mdrun option -pd)");
482 /* Orientation restraints */
485 init_orires(fplog
,mtop
,state
->x
,inputrec
,cr
->ms
,&(fcd
->orires
),
490 if (DEFORM(*inputrec
))
492 /* Store the deform reference box before reading the checkpoint */
495 copy_mat(state
->box
,box
);
499 gmx_bcast(sizeof(box
),box
,cr
);
501 /* Because we do not have the update struct available yet
502 * in which the reference values should be stored,
503 * we store them temporarily in static variables.
504 * This should be thread safe, since they are only written once
505 * and with identical values.
508 tMPI_Thread_mutex_lock(&deform_init_box_mutex
);
510 deform_init_init_step_tpx
= inputrec
->init_step
;
511 copy_mat(box
,deform_init_box_tpx
);
513 tMPI_Thread_mutex_unlock(&deform_init_box_mutex
);
517 if (opt2bSet("-cpi",nfile
,fnm
))
519 /* Check if checkpoint file exists before doing continuation.
520 * This way we can use identical input options for the first and subsequent runs...
522 if( gmx_fexist_master(opt2fn_master("-cpi",nfile
,fnm
,cr
),cr
) )
524 load_checkpoint(opt2fn_master("-cpi",nfile
,fnm
,cr
),&fplog
,
525 cr
,Flags
& MD_PARTDEC
,ddxyz
,
526 inputrec
,state
,&bReadRNG
,&bReadEkin
,
527 (Flags
& MD_APPENDFILES
));
531 Flags
|= MD_READ_RNG
;
535 Flags
|= MD_READ_EKIN
;
540 if (((MASTER(cr
) || (Flags
& MD_SEPPOT
)) && (Flags
& MD_APPENDFILES
))
542 /* With thread MPI only the master node/thread exists in mdrun.c,
543 * therefore non-master nodes need to open the "seppot" log file here.
545 || (!MASTER(cr
) && (Flags
& MD_SEPPOT
))
549 gmx_log_open(ftp2fn(efLOG
,nfile
,fnm
),cr
,!(Flags
& MD_SEPPOT
),
555 copy_mat(state
->box
,box
);
560 gmx_bcast(sizeof(box
),box
,cr
);
563 /* Essential dynamics */
564 if (opt2bSet("-ei",nfile
,fnm
))
566 /* Open input and output files, allocate space for ED data structure */
567 ed
= ed_open(nfile
,fnm
,Flags
,cr
);
570 if (bVerbose
&& SIMMASTER(cr
))
572 fprintf(stderr
,"Loaded with Money\n\n");
575 if (PAR(cr
) && !((Flags
& MD_PARTDEC
) ||
576 EI_TPI(inputrec
->eI
) ||
577 inputrec
->eI
== eiNM
))
579 cr
->dd
= init_domain_decomposition(fplog
,cr
,Flags
,ddxyz
,rdd
,rconstr
,
584 &ddbox
,&npme_major
,&npme_minor
);
586 make_dd_communicators(fplog
,cr
,dd_node_order
);
588 /* Set overallocation to avoid frequent reallocation of arrays */
589 set_over_alloc_dd(TRUE
);
593 /* PME, if used, is done on all nodes with 1D decomposition */
595 cr
->duty
= (DUTY_PP
| DUTY_PME
);
598 if (!EI_TPI(inputrec
->eI
))
600 npme_major
= cr
->nnodes
;
603 if (inputrec
->ePBC
== epbcSCREW
)
606 "pbc=%s is only implemented with domain decomposition",
607 epbc_names
[inputrec
->ePBC
]);
613 /* After possible communicator splitting in make_dd_communicators.
614 * we can set up the intra/inter node communication.
616 gmx_setup_nodecomm(fplog
,cr
);
619 wcycle
= wallcycle_init(fplog
,resetstep
,cr
);
622 /* Master synchronizes its value of reset_counters with all nodes
623 * including PME only nodes */
624 reset_counters
= wcycle_get_reset_counters(wcycle
);
625 gmx_bcast_sim(sizeof(reset_counters
),&reset_counters
,cr
);
626 wcycle_set_reset_counters(wcycle
, reset_counters
);
631 if (cr
->duty
& DUTY_PP
)
633 /* For domain decomposition we allocate dynamically
634 * in dd_partition_system.
636 if (DOMAINDECOMP(cr
))
638 bcast_state_setup(cr
,state
);
648 bcast_state(cr
,state
,TRUE
);
652 /* Dihedral Restraints */
653 if (gmx_mtop_ftype_count(mtop
,F_DIHRES
) > 0)
655 init_dihres(fplog
,mtop
,inputrec
,fcd
);
658 /* Initiate forcerecord */
660 init_forcerec(fplog
,oenv
,fr
,fcd
,inputrec
,mtop
,cr
,box
,FALSE
,
661 opt2fn("-table",nfile
,fnm
),
662 opt2fn("-tablep",nfile
,fnm
),
663 opt2fn("-tableb",nfile
,fnm
),FALSE
,pforce
);
665 /* version for PCA_NOT_READ_NODE (see md.c) */
666 /*init_forcerec(fplog,fr,fcd,inputrec,mtop,cr,box,FALSE,
667 "nofile","nofile","nofile",FALSE,pforce);
669 fr
->bSepDVDL
= ((Flags
& MD_SEPPOT
) == MD_SEPPOT
);
671 /* Initialize QM-MM */
674 init_QMMMrec(cr
,box
,mtop
,inputrec
,fr
);
677 /* Initialize the mdatoms structure.
678 * mdatoms is not filled with atom data,
679 * as this can not be done now with domain decomposition.
681 mdatoms
= init_mdatoms(fplog
,mtop
,inputrec
->efep
!=efepNO
);
683 /* Initialize the virtual site communication */
684 vsite
= init_vsite(mtop
,cr
);
686 calc_shifts(box
,fr
->shift_vec
);
688 /* With periodic molecules the charge groups should be whole at start up
689 * and the virtual sites should not be far from their proper positions.
691 if (!inputrec
->bContinuation
&& MASTER(cr
) &&
692 !(inputrec
->ePBC
!= epbcNONE
&& inputrec
->bPeriodicMols
))
694 /* Make molecules whole at start of run */
695 if (fr
->ePBC
!= epbcNONE
)
697 do_pbc_first_mtop(fplog
,inputrec
->ePBC
,box
,mtop
,state
->x
);
701 /* Correct initial vsite positions are required
702 * for the initial distribution in the domain decomposition
703 * and for the initial shell prediction.
705 construct_vsites_mtop(fplog
,vsite
,mtop
,state
->x
);
709 /* Initiate PPPM if necessary */
710 if (fr
->eeltype
== eelPPPM
)
712 if (mdatoms
->nChargePerturbed
)
714 gmx_fatal(FARGS
,"Free energy with %s is not implemented",
715 eel_names
[fr
->eeltype
]);
717 status
= gmx_pppm_init(fplog
,cr
,oenv
,FALSE
,TRUE
,box
,
718 getenv("GMXGHAT"),inputrec
, (Flags
& MD_REPRODUCIBLE
));
721 gmx_fatal(FARGS
,"Error %d initializing PPPM",status
);
725 if (EEL_PME(fr
->eeltype
))
727 ewaldcoeff
= fr
->ewaldcoeff
;
728 pmedata
= &fr
->pmedata
;
737 /* This is a PME only node */
739 /* We don't need the state */
742 ewaldcoeff
= calc_ewaldcoeff(inputrec
->rcoulomb
, inputrec
->ewald_rtol
);
746 /* Initiate PME if necessary,
747 * either on all nodes or on dedicated PME nodes only. */
748 if (EEL_PME(inputrec
->coulombtype
))
752 nChargePerturbed
= mdatoms
->nChargePerturbed
;
754 if (cr
->npmenodes
> 0)
756 /* The PME only nodes need to know nChargePerturbed */
757 gmx_bcast_sim(sizeof(nChargePerturbed
),&nChargePerturbed
,cr
);
759 if (cr
->duty
& DUTY_PME
)
761 status
= gmx_pme_init(pmedata
,cr
,npme_major
,npme_minor
,inputrec
,
762 mtop
? mtop
->natoms
: 0,nChargePerturbed
,
763 (Flags
& MD_REPRODUCIBLE
));
766 gmx_fatal(FARGS
,"Error %d initializing PME",status
);
772 if (integrator
[inputrec
->eI
].func
== do_md
775 integrator
[inputrec
->eI
].func
== do_md_openmm
779 /* Turn on signal handling on all nodes */
781 * (A user signal from the PME nodes (if any)
782 * is communicated to the PP nodes.
784 signal_handler_install();
787 if (cr
->duty
& DUTY_PP
)
789 if (inputrec
->ePull
!= epullNO
)
791 /* Initialize pull code */
792 init_pull(fplog
,inputrec
,nfile
,fnm
,mtop
,cr
,oenv
,
793 EI_DYNAMICS(inputrec
->eI
) && MASTER(cr
),Flags
);
798 /* Initialize enforced rotation code */
799 init_rot(fplog
,inputrec
,nfile
,fnm
,cr
,state
->x
,state
->box
,mtop
,oenv
,
803 constr
= init_constraints(fplog
,mtop
,inputrec
,ed
,state
,cr
);
805 if (DOMAINDECOMP(cr
))
807 dd_init_bondeds(fplog
,cr
->dd
,mtop
,vsite
,constr
,inputrec
,
808 Flags
& MD_DDBONDCHECK
,fr
->cginfo_mb
);
810 set_dd_parameters(fplog
,cr
->dd
,dlb_scale
,inputrec
,fr
,&ddbox
);
812 setup_dd_grid(fplog
,cr
->dd
);
815 /* Now do whatever the user wants us to do (how flexible...) */
816 integrator
[inputrec
->eI
].func(fplog
,cr
,nfile
,fnm
,
817 oenv
,bVerbose
,bCompact
,
820 nstepout
,inputrec
,mtop
,
822 mdatoms
,nrnb
,wcycle
,ed
,fr
,
823 repl_ex_nst
,repl_ex_seed
,
824 cpt_period
,max_hours
,
829 if (inputrec
->ePull
!= epullNO
)
831 finish_pull(fplog
,inputrec
->pull
);
836 finish_rot(fplog
,inputrec
->rot
);
843 gmx_pmeonly(*pmedata
,cr
,nrnb
,wcycle
,ewaldcoeff
,FALSE
,inputrec
);
846 if (EI_DYNAMICS(inputrec
->eI
) || EI_TPI(inputrec
->eI
))
848 /* Some timing stats */
851 if (runtime
.proc
== 0)
853 runtime
.proc
= runtime
.real
;
862 wallcycle_stop(wcycle
,ewcRUN
);
864 /* Finish up, write some stuff
865 * if rerunMD, don't write last frame again
867 finish_run(fplog
,cr
,ftp2fn(efSTO
,nfile
,fnm
),
868 inputrec
,nrnb
,wcycle
,&runtime
,
869 EI_DYNAMICS(inputrec
->eI
) && !MULTISIM(cr
));
871 /* Does what it says */
872 print_date_and_time(fplog
,cr
->nodeid
,"Finished mdrun",&runtime
);
874 /* Close logfile already here if we were appending to it */
875 if (MASTER(cr
) && (Flags
& MD_APPENDFILES
))
877 gmx_log_close(fplog
);
880 rc
=(int)gmx_get_stop_condition();
883 /* we need to join all threads. The sub-threads join when they
884 exit this function, but the master thread needs to be told to
886 if (PAR(cr
) && MASTER(cr
))
895 void md_print_warning(const t_commrec
*cr
,FILE *fplog
,const char *buf
)
899 fprintf(stderr
,"\n%s\n",buf
);
903 fprintf(fplog
,"\n%s\n",buf
);