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