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 * Gromacs Runs On Most of All Computer Systems
56 #include "thread_mpi/threads.h"
63 #define MD_POLARISE (1<<2)
64 #define MD_IONIZE (1<<3)
65 #define MD_RERUN (1<<4)
66 #define MD_RERUN_VSITE (1<<5)
67 #define MD_FFSCAN (1<<6)
68 #define MD_SEPPOT (1<<7)
69 #define MD_PARTDEC (1<<9)
70 #define MD_DDBONDCHECK (1<<10)
71 #define MD_DDBONDCOMM (1<<11)
72 #define MD_CONFOUT (1<<12)
73 #define MD_REPRODUCIBLE (1<<13)
74 #define MD_READ_RNG (1<<14)
75 #define MD_APPENDFILES (1<<15)
76 #define MD_APPENDFILESSET (1<<21)
77 #define MD_KEEPANDNUMCPT (1<<16)
78 #define MD_READ_EKIN (1<<17)
79 #define MD_STARTFROMCPT (1<<18)
80 #define MD_RESETCOUNTERSHALFWAY (1<<19)
82 /* Define a number of flags to better control the information
83 * passed to compute_globals in md.c and global_stat.
86 /* We are rerunning the simulation */
87 #define CGLO_RERUNMD (1<<1)
88 /* we are computing the kinetic energy from average velocities */
89 #define CGLO_EKINAVEVEL (1<<2)
90 /* we are removing the center of mass momenta */
91 #define CGLO_STOPCM (1<<3)
92 /* bGStat is defined in do_md */
93 #define CGLO_GSTAT (1<<4)
94 /* Sum the energy terms in global computation */
95 #define CGLO_ENERGY (1<<6)
96 /* Sum the kinetic energy terms in global computation */
97 #define CGLO_TEMPERATURE (1<<7)
98 /* Sum the kinetic energy terms in global computation */
99 #define CGLO_PRESSURE (1<<8)
100 /* Sum the constraint term in global computation */
101 #define CGLO_CONSTRAINT (1<<9)
102 /* we are using an integrator that requires iteration over some steps - currently not used*/
103 #define CGLO_ITERATE (1<<10)
104 /* it is the first time we are iterating (or, only once through is required */
105 #define CGLO_FIRSTITERATE (1<<11)
106 /* Reading ekin from the trajectory */
107 #define CGLO_READEKIN (1<<12)
108 /* we need to reset the ekin rescaling factor here */
109 #define CGLO_SCALEEKIN (1<<13)
112 ddnoSEL, ddnoINTERLEAVE, ddnoPP_PME, ddnoCARTESIAN, ddnoNR
124 double time_per_step;
126 gmx_large_int_t nsteps_done;
135 gmx_bool bKeepAndNumCPT;
142 /* Variables for temporary use with the deform option,
143 * used in runner.c and md.c.
144 * (These variables should be stored in the tpx file.)
146 extern gmx_large_int_t deform_init_init_step_tpx;
147 extern matrix deform_init_box_tpx;
149 extern tMPI_Thread_mutex_t deform_init_box_mutex;
151 /* The minimum number of atoms per thread. With fewer atoms than this,
152 * the number of threads will get lowered.
154 #define MIN_ATOMS_PER_THREAD 90
158 typedef double gmx_integrator_t(FILE *log,t_commrec *cr,
159 int nfile,const t_filenm fnm[],
160 const output_env_t oenv, gmx_bool bVerbose,
161 gmx_bool bCompact, int nstglobalcomm,
162 gmx_vsite_t *vsite,gmx_constr_t constr,
164 t_inputrec *inputrec,
165 gmx_mtop_t *mtop,t_fcdata *fcd,
168 t_nrnb *nrnb,gmx_wallcycle_t wcycle,
171 int repl_ex_nst,int repl_ex_seed,
172 real cpt_period,real max_hours,
173 const char *deviceOptions,
175 gmx_runtime_t *runtime);
177 typedef struct gmx_global_stat *gmx_global_stat_t;
179 /* ROUTINES from md.c */
181 gmx_integrator_t do_md;
183 gmx_integrator_t do_md_openmm;
187 /* ROUTINES from minimize.c */
189 gmx_integrator_t do_steep;
190 /* Do steepest descents EM */
192 gmx_integrator_t do_cg;
193 /* Do conjugate gradient EM */
195 gmx_integrator_t do_lbfgs;
196 /* Do conjugate gradient L-BFGS */
198 gmx_integrator_t do_nm;
199 /* Do normal mode analysis */
201 /* ROUTINES from tpi.c */
203 gmx_integrator_t do_tpi;
204 /* Do test particle insertion */
207 /* ROUTINES from md_support.c */
209 /* return the number of steps between global communcations */
210 int check_nstglobalcomm(FILE *fplog,t_commrec *cr,
211 int nstglobalcomm,t_inputrec *ir);
213 /* check whether an 'nst'-style parameter p is a multiple of nst, and
214 set it to be one if not, with a warning. */
215 void check_nst_param(FILE *fplog,t_commrec *cr,
216 const char *desc_nst,int nst,
217 const char *desc_p,int *p);
219 /* check which of the multisim simulations has the shortest number of
220 steps and return that number of nsteps */
221 gmx_large_int_t get_multisim_nsteps(const t_commrec *cr,
222 gmx_large_int_t nsteps);
224 void rerun_parallel_comm(t_commrec *cr,t_trxframe *fr,
225 gmx_bool *bNotLastFrame);
227 /* get the conserved energy associated with the ensemble type*/
228 real compute_conserved_from_auxiliary(t_inputrec *ir, t_state *state,
231 /* reset all cycle and time counters. */
232 void reset_all_counters(FILE *fplog,t_commrec *cr,
233 gmx_large_int_t step,
234 gmx_large_int_t *step_rel,t_inputrec *ir,
235 gmx_wallcycle_t wcycle,t_nrnb *nrnb,
236 gmx_runtime_t *runtime);
240 /* ROUTINES from sim_util.c */
241 void do_pbc_first(FILE *log,matrix box,t_forcerec *fr,
242 t_graph *graph,rvec x[]);
244 void do_pbc_first_mtop(FILE *fplog,int ePBC,matrix box,
245 gmx_mtop_t *mtop,rvec x[]);
247 void do_pbc_mtop(FILE *fplog,int ePBC,matrix box,
248 gmx_mtop_t *mtop,rvec x[]);
252 /* ROUTINES from stat.c */
253 gmx_global_stat_t global_stat_init(t_inputrec *ir);
255 void global_stat_destroy(gmx_global_stat_t gs);
257 void global_stat(FILE *log,gmx_global_stat_t gs,
258 t_commrec *cr,gmx_enerdata_t *enerd,
259 tensor fvir,tensor svir,rvec mu_tot,
260 t_inputrec *inputrec,
261 gmx_ekindata_t *ekind,
262 gmx_constr_t constr,t_vcm *vcm,
264 gmx_mtop_t *top_global, t_state *state_local,
265 gmx_bool bSumEkinhOld, int flags);
266 /* Communicate statistics over cr->mpi_comm_mysim */
268 gmx_mdoutf_t *init_mdoutf(int nfile,const t_filenm fnm[],
270 const t_commrec *cr,const t_inputrec *ir,
271 const output_env_t oenv);
272 /* Returns a pointer to a data structure with all output file pointers
273 * and names required by mdrun.
276 void done_mdoutf(gmx_mdoutf_t *of);
277 /* Close all open output files and free the of pointer */
279 #define MDOF_X (1<<0)
280 #define MDOF_V (1<<1)
281 #define MDOF_F (1<<2)
282 #define MDOF_XTC (1<<3)
283 #define MDOF_CPT (1<<4)
285 void write_traj(FILE *fplog,t_commrec *cr,
288 gmx_mtop_t *top_global,
289 gmx_large_int_t step,double t,
290 t_state *state_local,t_state *state_global,
291 rvec *f_local,rvec *f_global,
292 int *n_xtc,rvec **x_xtc);
293 /* Routine that writes frames to trn, xtc and/or checkpoint.
294 * What is written is determined by the mdof_flags defined above.
295 * Data is collected to the master node only when necessary.
298 int do_per_step(gmx_large_int_t step,gmx_large_int_t nstep);
299 /* Return TRUE if io should be done */
301 int do_any_io(int step, t_inputrec *ir);
303 /* ROUTINES from sim_util.c */
305 double gmx_gettime();
307 void print_time(FILE *out, gmx_runtime_t *runtime,
308 gmx_large_int_t step,t_inputrec *ir, t_commrec *cr);
310 void runtime_start(gmx_runtime_t *runtime);
312 void runtime_end(gmx_runtime_t *runtime);
314 void runtime_upd_proc(gmx_runtime_t *runtime);
315 /* The processor time should be updated every once in a while,
316 * since on 32-bit manchines it loops after 72 minutes.
319 void print_date_and_time(FILE *log,int pid,const char *title,
320 const gmx_runtime_t *runtime);
322 void nstop_cm(FILE *log,t_commrec *cr,
323 int start,int nr_atoms,real mass[],rvec x[],rvec v[]);
325 void finish_run(FILE *log,t_commrec *cr,const char *confout,
326 t_inputrec *inputrec,
327 t_nrnb nrnb[],gmx_wallcycle_t wcycle,
328 gmx_runtime_t *runtime,
329 gmx_bool bWriteStat);
331 void calc_enervirdiff(FILE *fplog,int eDispCorr,t_forcerec *fr);
333 void calc_dispcorr(FILE *fplog,t_inputrec *ir,t_forcerec *fr,
334 gmx_large_int_t step, int natoms,
335 matrix box,real lambda,tensor pres,tensor virial,
336 real *prescorr, real *enercorr, real *dvdlcorr);
349 void check_nnodes_top(char *fn,t_topology *top);
350 /* Reset the tpr file to work with one node if necessary */
353 /* check the version */
354 void check_ir_old_tpx_versions(t_commrec *cr,FILE *fplog,
355 t_inputrec *ir,gmx_mtop_t *mtop);
357 /* Allocate and initialize node-local state entries. */
358 void set_state_entries(t_state *state,const t_inputrec *ir,int nnodes);
360 /* Broadcast the data for a simulation, and allocate node-specific settings
361 such as rng generators. */
362 void init_parallel(FILE *log, t_commrec *cr, t_inputrec *inputrec,
366 void do_constrain_first(FILE *log,gmx_constr_t constr,
367 t_inputrec *inputrec,t_mdatoms *md,
368 t_state *state,rvec *f,
369 t_graph *graph,t_commrec *cr,t_nrnb *nrnb,
370 t_forcerec *fr, gmx_localtop_t *top, tensor shake_vir);
372 void dynamic_load_balancing(gmx_bool bVerbose,t_commrec *cr,real capacity[],
373 int dimension,t_mdatoms *md,t_topology *top,
374 rvec x[],rvec v[],matrix box);
375 /* Perform load balancing, i.e. split the particles over processors
376 * based on their coordinates in the "dimension" direction.
379 int multisim_min(const gmx_multisim_t *ms,int nmin,int n);
380 /* Set an appropriate value for n across the whole multi-simulation */
382 int multisim_nstsimsync(const t_commrec *cr,
383 const t_inputrec *ir,int repl_ex_nst);
384 /* Determine the interval for inter-simulation communication */
386 void init_global_signals(globsig_t *gs,const t_commrec *cr,
387 const t_inputrec *ir,int repl_ex_nst);
388 /* Constructor for globsig_t */
390 void copy_coupling_state(t_state *statea,t_state *stateb,
391 gmx_ekindata_t *ekinda,gmx_ekindata_t *ekindb, t_grpopts* opts);
392 /* Copy stuff from state A to state B */
394 void compute_globals(FILE *fplog, gmx_global_stat_t gstat, t_commrec *cr, t_inputrec *ir,
395 t_forcerec *fr, gmx_ekindata_t *ekind,
396 t_state *state, t_state *state_global, t_mdatoms *mdatoms,
397 t_nrnb *nrnb, t_vcm *vcm, gmx_wallcycle_t wcycle,
398 gmx_enerdata_t *enerd,tensor force_vir, tensor shake_vir, tensor total_vir,
399 tensor pres, rvec mu_tot, gmx_constr_t constr,
400 globsig_t *gs,gmx_bool bInterSimGS,
401 matrix box, gmx_mtop_t *top_global, real *pcurr,
402 int natoms, gmx_bool *bSumEkinhOld, int flags);
403 /* Compute global variables during integration */
405 int mdrunner(int nthreads_requested, FILE *fplog,t_commrec *cr,int nfile,
406 const t_filenm fnm[], const output_env_t oenv, gmx_bool bVerbose,
407 gmx_bool bCompact, int nstglobalcomm, ivec ddxyz,int dd_node_order,
408 real rdd, real rconstr, const char *dddlb_opt,real dlb_scale,
409 const char *ddcsx,const char *ddcsy,const char *ddcsz,
410 int nstepout, int resetstep, int nmultisim, int repl_ex_nst,
411 int repl_ex_seed, real pforce,real cpt_period,real max_hours,
412 const char *deviceOptions, unsigned long Flags);
413 /* Driver routine, that calls the different methods */
415 void md_print_warning(const t_commrec *cr,FILE *fplog,const char *buf);
416 /* Print a warning message to stderr on the master node
417 * and to fplog if fplog!=NULL.
420 void init_md(FILE *fplog,
421 t_commrec *cr,t_inputrec *ir, const output_env_t oenv,
422 double *t,double *t0,
423 real *lambda,double *lam0,
424 t_nrnb *nrnb,gmx_mtop_t *mtop,
426 int nfile,const t_filenm fnm[],
427 gmx_mdoutf_t **outf,t_mdebin **mdebin,
428 tensor force_vir,tensor shake_vir,
430 gmx_bool *bSimAnn,t_vcm **vcm,
431 t_state *state, unsigned long Flags);
432 /* Routine in sim_util.c */
438 #endif /* _mdrun_h */