4 * This source code is part of
8 * GROningen MAchine for Chemical Simulations
11 * Copyright (c) 1991-2001, University of Groningen, The Netherlands
12 * This program is free software; you can redistribute it and/or
13 * modify it under the terms of the GNU General Public License
14 * as published by the Free Software Foundation; either version 2
15 * of the License, or (at your option) any later version.
17 * If you want to redistribute modifications, please consider that
18 * scientific software is very special. Version control is crucial -
19 * bugs must be traceable. We will be happy to consider code for
20 * inclusion in the official distribution, but derived work must not
21 * be called official GROMACS. Details are found in the README & COPYING
22 * files - if they are missing, get the official version at www.gromacs.org.
24 * To help us fund GROMACS development, we humbly ask that you cite
25 * the papers on the package - you can find them in the top README file.
27 * For more info, check our website at http://www.gromacs.org
30 * Gnomes, ROck Monsters And Chili Sauce
37 int n
; /* Number of terms */
38 real
*a
; /* Coeffients (V / nm ) */
39 real
*phi
; /* Phase angles */
43 int ngtc
; /* # T-Coupl groups */
44 int ngacc
; /* # Accelerate groups */
45 int ngfrz
; /* # Freeze groups */
46 int ngener
; /* # Ener groups */
47 real
*nrdf
; /* Nr of degrees of freedom in a group */
48 real
*ref_t
; /* Coupling temperature per group */
49 real
*tau_t
; /* Tau coupling time */
50 rvec
*acc
; /* Acceleration per group */
51 ivec
*nFreeze
; /* Freeze the group in each direction ? */
52 bool *eg_excl
; /* Exclusions of energy group pairs */
56 int eI
; /* Integration method */
57 int nsteps
; /* number of steps to be taken */
58 int ns_type
; /* which ns method should we use? */
59 int nstlist
; /* number of steps before pairlist is generated */
60 int ndelta
; /* number of cells per rlong */
61 bool bDomDecomp
; /* Should we do domain decomposition? */
62 int decomp_dir
; /* Direction of decomposition (may not be opt.) */
63 int nstcomm
; /* number of steps after which center of mass */
64 /* motion is removed */
65 int nstlog
; /* number of steps after which print to logfile */
66 int nstxout
; /* number of steps after which X is output */
67 int nstvout
; /* id. for V */
68 int nstfout
; /* id. for F */
69 int nstenergy
; /* number of steps after which energies printed */
70 int nstxtcout
; /* id. for compressed trj (.xtc) */
71 real init_t
; /* initial time (ps) */
72 real delta_t
; /* time step (ps) */
73 real xtcprec
; /* precision of xtc file */
74 int nkx
,nky
,nkz
; /* number of k vectors in each spatial dimension*/
75 /* for fourier methods for long range electrost.*/
76 int pme_order
; /* interpolation order for PME */
77 real ewald_rtol
; /* Real space tolerance for Ewald, determines */
78 /* the real/reciprocal space relative weight */
79 int ewald_geometry
; /* normal/3d ewald, or pseudo-2d LR corrections */
80 real epsilon_surface
; /* Epsilon for PME dipole correction */
81 bool bOptFFT
; /* optimize the fft plan at start */
82 int ePBC
; /* Type of periodic boundary conditions */
83 bool bUncStart
; /* Do not constrain the start configuration */
84 int etc
; /* temperature coupling */
85 int epc
; /* pressure coupling */
86 int epct
; /* pressure coupling type */
87 real tau_p
; /* pressure coupling time (ps) */
88 tensor ref_p
; /* reference pressure (kJ/(mol nm^3)) */
89 tensor compress
; /* compressability ((mol nm^3)/kJ) */
90 bool bSimAnn
; /* simulated annealing (SA) */
91 real zero_temp_time
; /* time at which temp becomes zero in sim. ann. */
92 real rlist
; /* short range pairlist cut-off (nm) */
93 int coulombtype
; /* Type of electrostatics treatment */
94 real rcoulomb_switch
; /* Coulomb switch range start (nm) */
95 real rcoulomb
; /* Coulomb cutoff (nm) */
96 int vdwtype
; /* Type of Van der Waals treatment */
97 real rvdw_switch
; /* Van der Waals switch range start (nm) */
98 real rvdw
; /* Van der Waals cutoff (nm) */
99 real epsilon_r
; /* relative dielectric constant */
100 int eDispCorr
; /* Perform Long range dispersion corrections */
101 real shake_tol
; /* tolerance for shake */
102 real fudgeQQ
; /* Id. for 1-4 coulomb interactions */
103 int efep
; /* free energy interpolation no/yes */
104 real init_lambda
; /* initial value for perturbation variable */
105 real delta_lambda
; /* change of lambda per time step (1/dt) */
106 real sc_alpha
; /* free energy soft-core parameter */
107 real sc_sigma
; /* free energy soft-core sigma when c6 or c12=0 */
108 real dr_fc
; /* force constant for ta_disre */
109 int eDisreWeighting
; /* type of weighting of pairs in one restraints */
110 bool bDisreMixed
; /* Use comb of time averaged and instan. viol's */
111 int nstdisreout
; /* frequency of writing pair distances to enx */
112 real dr_tau
; /* time constant for memory function in disres */
113 real orires_fc
; /* force constant for orientational restraints */
114 real orires_tau
; /* time constant for memory function in orires */
115 int nstorireout
; /* frequency of writing tr(SD) to enx */
116 real em_stepsize
; /* The stepsize for updating */
117 real em_tol
; /* The tolerance */
118 int niter
; /* Number of iterations for convergence of */
119 /* steepest descent in relax_shells */
120 real fc_stepsize
; /* Stepsize for directional minimization */
121 /* in relax_shells */
122 int nstcgsteep
; /* number of steps after which a steepest */
123 /* descents step is done while doing cg */
124 int eConstrAlg
; /* Type of constraint algorithm */
125 int nProjOrder
; /* Order of the LINCS Projection Algorithm */
126 real LincsWarnAngle
; /* If bond rotates more than %g degrees, warn */
127 bool bShakeSOR
; /* Use successive overrelaxation for shake */
128 real bd_temp
; /* Temperature for Brownian Dynamics (BD) */
129 real bd_fric
; /* Friction coefficient for BD (amu / ps) */
130 int ld_seed
; /* Random seed for SD and BD */
131 real cos_accel
; /* Acceleration for viscosity calculation */
132 int userint1
; /* User determined parameters */
140 t_grpopts opts
; /* Group options */
141 t_cosines ex
[DIM
]; /* Electric field stuff (spatial part) */
142 t_cosines et
[DIM
]; /* Electric field stuff (time part) */