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 * GRoups of Organic Molecules in ACtion for Science
46 * The t_state struct should contain all the (possibly) non-static
47 * information required to define the state of the system.
48 * Currently the random seeds for SD and BD are missing.
51 /* for now, define the length of the NH chains here */
52 #define NHCHAINLENGTH 10
53 #define MAXLAMBDAS 1024
55 /* These enums are used in flags as (1<<est...).
56 * The order of these enums should not be changed,
57 * since that affects the checkpoint (.cpt) file format.
60 estBOX, estBOX_REL, estBOXV, estPRES_PREV, estNH_XI, estTC_INT,
61 estX, estV, estSDX, estCGP, estLD_RNG, estLD_RNGI,
62 estDISRE_INITF, estDISRE_RM3TAV,
63 estORIRE_INITF, estORIRE_DTAV,
64 estSVIR_PREV, estNH_VXI, estVETA, estVOL0, estNHPRES_XI, estNHPRES_VXI, estFVIR_PREV,
65 estFEPSTATE, estMC_RNG, estMC_RNGI,
68 #define EST_DISTR(e) (!(((e) >= estLAMBDA && (e) <= estTC_INT) || ((e) >= estSVIR_PREV && (e) <= estMC_RNGI)))
70 /* The names of the state entries, defined in src/gmxlib/checkpoint.c */
71 extern const char *est_names[estNR];
75 real disre_initf; /* The scaling factor for initializing the time av. */
76 int ndisrepairs; /* The number of distance restraints */
77 real *disre_rm3tav; /* The r^-3 time averaged pair distances */
78 real orire_initf; /* The scaling factor for initializing the time av. */
79 int norire_Dtav; /* The number of matrix element in dtav (npair*5) */
80 real *orire_Dtav; /* The time averaged orientation tensors */
83 /* Struct used for checkpointing only.
84 * This struct would not be required with unlimited precision.
85 * But because of limited precision, the COM motion removal implementation
86 * can cause the kinetic energy in the MD loop to differ by a few bits from
87 * the kinetic energy one would determine from state.v.
97 double *ekinscalef_nhc;
98 double *ekinscaleh_nhc;
104 /* energy history for delta_h histograms */
107 int nndh; /* the number of energy difference lists */
108 int *ndh; /* the number in each energy difference list */
109 real **dh; /* the energy difference lists */
111 double start_time; /* the start time of these energy diff blocks */
112 double start_lambda; /* lambda at start time */
114 gmx_bool start_lambda_set; /* whether the lambda value is set. Here
115 For backward-compatibility. */
120 int nlambda; /* total number of lambda states - for history*/
122 gmx_bool bEquil; /* reached equilibration */
123 int *n_at_lam; /* number of points observed at each lambda */
124 real *wl_histo; /* histogram for WL flatness determination */
125 real wl_delta; /* current wang-landau delta */
127 real *sum_weights; /* weights of the states */
128 real *sum_dg; /* free energies of the states -- not actually used for weighting, but informational */
129 real *sum_minvar; /* corrections to weights for minimum variance */
130 real *sum_variance; /* variances of the states */
132 real **accum_p; /* accumulated bennett weights for n+1 */
133 real **accum_m; /* accumulated bennett weights for n-1 */
134 real **accum_p2; /* accumulated squared bennett weights for n+1 */
135 real **accum_m2; /* accumulated squared bennett weights for n-1 */
137 real **Tij; /* transition matrix */
138 real **Tij_empirical; /* Empirical transition matrix */
143 gmx_large_int_t nsteps; /* The number of steps in the history */
144 gmx_large_int_t nsum; /* The nr. of steps in the ener_ave and ener_sum */
145 double * ener_ave; /* Energy term history sum to get fluctuations */
146 double * ener_sum; /* Energy term history sum to get fluctuations */
147 int nener; /* Number of energy terms in two previous arrays */
148 gmx_large_int_t nsteps_sim; /* The number of steps in ener_sum_sim */
149 gmx_large_int_t nsum_sim; /* The number of frames in ener_sum_sim */
150 double * ener_sum_sim; /* Energy term history sum of the whole sim */
152 delta_h_history_t *dht; /* The BAR energy differences */
161 int nhchainlength; /* number of nose-hoover chains */
164 int flags; /* Flags telling which entries are present */
165 int fep_state; /* indicates which of the alchemical states we are in */
166 real *lambda; /* lambda vector */
167 matrix box; /* box vector coordinates */
168 matrix box_rel; /* Relitaive box vectors to preserve shape */
169 matrix boxv; /* box velocitites for Parrinello-Rahman pcoupl */
170 matrix pres_prev; /* Pressure of the previous step for pcoupl */
171 matrix svir_prev; /* Shake virial for previous step for pcoupl */
172 matrix fvir_prev; /* Force virial of the previous step for pcoupl */
173 double *nosehoover_xi; /* for Nose-Hoover tcoupl (ngtc) */
174 double *nosehoover_vxi; /* for N-H tcoupl (ngtc) */
175 double *nhpres_xi; /* for Nose-Hoover pcoupl for barostat */
176 double *nhpres_vxi; /* for Nose-Hoover pcoupl for barostat */
177 double *therm_integral; /* for N-H/V-rescale tcoupl (ngtc) */
178 real veta; /* trotter based isotropic P-coupling */
179 real vol0; /* initial volume,required for computing NPT conserverd quantity */
180 int nalloc; /* Allocation size for x, v and sd_x when !=NULL*/
181 rvec *x; /* the coordinates (natoms) */
182 rvec *v; /* the velocities (natoms) */
183 rvec *sd_X; /* random part of the x update for stoch. dyn. */
184 rvec *cg_p; /* p vector for conjugate gradient minimization */
186 unsigned int *ld_rng; /* RNG random state */
187 int *ld_rngi; /* RNG index */
189 int nmcrng; /* number of RNG states */
190 unsigned int *mc_rng; /* lambda MC RNG random state */
191 int *mc_rngi; /* lambda MC RNG index */
193 history_t hist; /* Time history for restraints */
195 ekinstate_t ekinstate; /* The state of the kinetic energy data */
197 energyhistory_t enerhist; /* Energy history for statistics */
198 df_history_t dfhist; /*Free energy history for free energy analysis */
200 int ddp_count; /* The DD partitioning count for this state */
201 int ddp_count_cg_gl; /* The DD part. count for index_gl */
202 int ncg_gl; /* The number of local charge groups */
203 int *cg_gl; /* The global cg number of the local cgs */
204 int cg_gl_nalloc; /* Allocation size of cg_gl; */
209 double *Qinv; /* inverse mass of thermostat -- computed from inputs, but a good place to store */
210 double *QPinv; /* inverse mass of thermostat for barostat -- computed from inputs, but a good place to store */
211 double Winv; /* Pressure mass inverse -- computed, not input, but a good place to store. Need to make a matrix later */
212 tensor Winvm; /* inverse pressure mass tensor, computed */
231 #endif /* _state_h_ */