2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
5 * Copyright (c) 2001-2004, The GROMACS development team.
6 * Copyright (c) 2013,2014, by the GROMACS development team, led by
7 * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
8 * and including many others, as listed in the AUTHORS file in the
9 * top-level source directory and at http://www.gromacs.org.
11 * GROMACS is free software; you can redistribute it and/or
12 * modify it under the terms of the GNU Lesser General Public License
13 * as published by the Free Software Foundation; either version 2.1
14 * of the License, or (at your option) any later version.
16 * GROMACS is distributed in the hope that it will be useful,
17 * but WITHOUT ANY WARRANTY; without even the implied warranty of
18 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
19 * Lesser General Public License for more details.
21 * You should have received a copy of the GNU Lesser General Public
22 * License along with GROMACS; if not, see
23 * http://www.gnu.org/licenses, or write to the Free Software Foundation,
24 * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
26 * If you want to redistribute modifications to GROMACS, please
27 * consider that scientific software is very special. Version
28 * control is crucial - bugs must be traceable. We will be happy to
29 * consider code for inclusion in the official distribution, but
30 * derived work must not be called official GROMACS. Details are found
31 * in the README & COPYING files - if they are missing, get the
32 * official version at http://www.gromacs.org.
34 * To help us fund GROMACS development, we humbly ask that you cite
35 * the research papers on the package. Check out http://www.gromacs.org.
41 #include "gromacs/legacyheaders/types/simple.h"
42 #include "gromacs/swap/enums.h"
49 * The t_state struct should contain all the (possibly) non-static
50 * information required to define the state of the system.
51 * Currently the random seeds for SD and BD are missing.
54 /* These enums are used in flags as (1<<est...).
55 * The order of these enums should not be changed,
56 * 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,
69 #define EST_DISTR(e) (!(((e) >= estLAMBDA && (e) <= estTC_INT) || ((e) >= estSVIR_PREV && (e) <= estMC_RNGI)))
71 /* The names of the state entries, defined in src/gmxlib/checkpoint.c */
72 extern const char *est_names[estNR];
76 real disre_initf; /* The scaling factor for initializing the time av. */
77 int ndisrepairs; /* The number of distance restraints */
78 real *disre_rm3tav; /* The r^-3 time averaged pair distances */
79 real orire_initf; /* The scaling factor for initializing the time av. */
80 int norire_Dtav; /* The number of matrix element in dtav (npair*5) */
81 real *orire_Dtav; /* The time averaged orientation tensors */
84 /* Struct used for checkpointing only.
85 * This struct would not be required with unlimited precision.
86 * But because of limited precision, the COM motion removal implementation
87 * can cause the kinetic energy in the MD loop to differ by a few bits from
88 * the kinetic energy one would determine from state.v.
98 double *ekinscalef_nhc;
99 double *ekinscaleh_nhc;
105 /* energy history for delta_h histograms */
108 int nndh; /* the number of energy difference lists */
109 int *ndh; /* the number in each energy difference list */
110 real **dh; /* the energy difference lists */
112 double start_time; /* the start time of these energy diff blocks */
113 double start_lambda; /* lambda at start time */
115 gmx_bool start_lambda_set; /* whether the lambda value is set. Here
116 For backward-compatibility. */
121 int nlambda; /* total number of lambda states - for history*/
123 gmx_bool bEquil; /* Have we reached equilibration */
124 int *n_at_lam; /* number of points observed at each lambda */
125 real *wl_histo; /* histogram for WL flatness determination */
126 real wl_delta; /* current wang-landau delta */
128 real *sum_weights; /* weights of the states */
129 real *sum_dg; /* free energies of the states -- not actually used for weighting, but informational */
130 real *sum_minvar; /* corrections to weights for minimum variance */
131 real *sum_variance; /* variances of the states */
133 real **accum_p; /* accumulated bennett weights for n+1 */
134 real **accum_m; /* accumulated bennett weights for n-1 */
135 real **accum_p2; /* accumulated squared bennett weights for n+1 */
136 real **accum_m2; /* accumulated squared bennett weights for n-1 */
138 real **Tij; /* transition matrix */
139 real **Tij_empirical; /* Empirical transition matrix */
145 gmx_int64_t nsteps; /* The number of steps in the history */
146 gmx_int64_t nsum; /* The nr. of steps in the ener_ave and ener_sum */
147 double * ener_ave; /* Energy term history sum to get fluctuations */
148 double * ener_sum; /* Energy term history sum to get fluctuations */
149 int nener; /* Number of energy terms in two previous arrays */
150 gmx_int64_t nsteps_sim; /* The number of steps in ener_sum_sim */
151 gmx_int64_t nsum_sim; /* The number of frames in ener_sum_sim */
152 double * ener_sum_sim; /* Energy term history sum of the whole sim */
154 delta_h_history_t *dht; /* The BAR energy differences */
160 /* If one uses essential dynamics or flooding on a group of atoms from
161 * more than one molecule, we cannot make this group whole with
162 * do_pbc_first_mtop(). We assume that the ED group has the correct PBC
163 * representation at the beginning of the simulation and keep track
164 * of the shifts to always get it into that representation.
165 * For proper restarts from a checkpoint we store the positions of the
166 * reference group at the time of checkpoint writing */
167 gmx_bool bFromCpt; /* Did we start from a checkpoint file? */
168 int nED; /* No. of ED/Flooding data sets, if <1 no ED */
169 int *nref; /* No. of atoms in i'th reference structure */
170 int *nav; /* Same for average structure */
171 rvec **old_sref; /* Positions of the reference atoms
172 at the last time step (with correct PBC
174 rvec **old_sref_p; /* Pointer to these positions */
175 rvec **old_sav; /* Same for the average positions */
183 int eSwapCoords; /* Swapping along x, y, or z-direction? */
184 int nat_req[eCompNR][eIonNR]; /* Requested ion numbers per type an comp. */
185 int *nat_req_p[eCompNR][eIonNR]; /* Pointer to this data (for .cpt writing) */
186 int nAverage; /* Use average over this many swap attempt
187 steps when determining the ion counts */
188 int inflow_netto[eCompNR][eIonNR]; /* Flux determined from the # of swaps */
189 int *inflow_netto_p[eCompNR][eIonNR]; /* Pointer to this data */
190 int *nat_past[eCompNR][eIonNR]; /* Array with nAverage entries for history */
191 int *nat_past_p[eCompNR][eIonNR]; /* Pointer points to the first entry only */
193 /* Channel flux detection, this is counting only and has no influence on whether swaps
194 * are performed or not: */
195 int fluxfromAtoB[eCompNR][eIonNR]; /* Flux determined from the split cylinders */
196 int *fluxfromAtoB_p[eCompNR][eIonNR]; /* Pointer to this data */
197 int *fluxleak; /* Flux not going through any channel */
198 int nions; /* Size of the following arrays */
199 unsigned char *comp_from; /* Ion came from which compartment? */
200 unsigned char *channel_label; /* Through which channel did this ion pass? */
202 /* To also make multimeric channel proteins whole, we save the last whole configuration of
203 * the channels in the checkpoint file. If we have no checkpoint file, we assume that the
204 * starting configuration hast the correct PBC representation after making the individual
206 gmx_bool bFromCpt; /* Did we started from a checkpoint file? */
207 int nat[eChanNR]; /* Size of xc_old_whole, i.e. the number of
208 atoms in each channel */
209 rvec *xc_old_whole[eChanNR]; /* Last known whole positions of the two
210 channels (important for multimeric ch.!) */
211 rvec **xc_old_whole_p[eChanNR]; /* Pointer to these positions */
221 int nhchainlength; /* number of nose-hoover chains */
222 int flags; /* Flags telling which entries are present */
223 int fep_state; /* indicates which of the alchemical states we are in */
224 real *lambda; /* lambda vector */
225 matrix box; /* box vector coordinates */
226 matrix box_rel; /* Relitaive box vectors to preserve shape */
227 matrix boxv; /* box velocitites for Parrinello-Rahman pcoupl */
228 matrix pres_prev; /* Pressure of the previous step for pcoupl */
229 matrix svir_prev; /* Shake virial for previous step for pcoupl */
230 matrix fvir_prev; /* Force virial of the previous step for pcoupl */
231 double *nosehoover_xi; /* for Nose-Hoover tcoupl (ngtc) */
232 double *nosehoover_vxi; /* for N-H tcoupl (ngtc) */
233 double *nhpres_xi; /* for Nose-Hoover pcoupl for barostat */
234 double *nhpres_vxi; /* for Nose-Hoover pcoupl for barostat */
235 double *therm_integral; /* for N-H/V-rescale tcoupl (ngtc) */
236 real veta; /* trotter based isotropic P-coupling */
237 real vol0; /* initial volume,required for computing NPT conserverd quantity */
238 int nalloc; /* Allocation size for x, v and sd_x when !=NULL*/
239 rvec *x; /* the coordinates (natoms) */
240 rvec *v; /* the velocities (natoms) */
241 rvec *sd_X; /* random part of the x update for stoch. dyn. */
242 rvec *cg_p; /* p vector for conjugate gradient minimization */
244 history_t hist; /* Time history for restraints */
246 ekinstate_t ekinstate; /* The state of the kinetic energy data */
248 energyhistory_t enerhist; /* Energy history for statistics */
249 swapstate_t swapstate; /* Position swapping */
250 df_history_t dfhist; /*Free energy history for free energy analysis */
251 edsamstate_t edsamstate; /* Essential dynamics / flooding history */
253 int ddp_count; /* The DD partitioning count for this state */
254 int ddp_count_cg_gl; /* The DD part. count for index_gl */
255 int ncg_gl; /* The number of local charge groups */
256 int *cg_gl; /* The global cg number of the local cgs */
257 int cg_gl_nalloc; /* Allocation size of cg_gl; */
262 double *Qinv; /* inverse mass of thermostat -- computed from inputs, but a good place to store */
263 double *QPinv; /* inverse mass of thermostat for barostat -- computed from inputs, but a good place to store */
264 double Winv; /* Pressure mass inverse -- computed, not input, but a good place to store. Need to make a matrix later */
265 tensor Winvm; /* inverse pressure mass tensor, computed */
284 #endif /* _state_h_ */