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,2015,2016,2017 by the GROMACS development team.
7 * Copyright (c) 2018,2019,2020,2021, by the GROMACS development team, led by
8 * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
9 * and including many others, as listed in the AUTHORS file in the
10 * top-level source directory and at http://www.gromacs.org.
12 * GROMACS is free software; you can redistribute it and/or
13 * modify it under the terms of the GNU Lesser General Public License
14 * as published by the Free Software Foundation; either version 2.1
15 * of the License, or (at your option) any later version.
17 * GROMACS is distributed in the hope that it will be useful,
18 * but WITHOUT ANY WARRANTY; without even the implied warranty of
19 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
20 * Lesser General Public License for more details.
22 * You should have received a copy of the GNU Lesser General Public
23 * License along with GROMACS; if not, see
24 * http://www.gnu.org/licenses, or write to the Free Software Foundation,
25 * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
27 * If you want to redistribute modifications to GROMACS, please
28 * consider that scientific software is very special. Version
29 * control is crucial - bugs must be traceable. We will be happy to
30 * consider code for inclusion in the official distribution, but
31 * derived work must not be called official GROMACS. Details are found
32 * in the README & COPYING files - if they are missing, get the
33 * official version at http://www.gromacs.org.
35 * To help us fund GROMACS development, we humbly ask that you cite
36 * the research papers on the package. Check out http://www.gromacs.org.
39 #ifndef GMX_MDLIB_MDEBIN_BAR_H
40 #define GMX_MDLIB_MDEBIN_BAR_H
42 #include "gromacs/utility/real.h"
44 /* The functions & data structures here describe writing
45 energy differences (or their histogram )for use with g_bar */
47 class delta_h_history_t;
49 class energyhistory_t;
58 /* Data for one foreign lambda, or derivative. */
61 std::vector<real> dh; /* the raw energy data. */
62 std::vector<float> dhf; /* raw difference data -- in floats, for storage. */
63 unsigned int ndh; /* number of data points */
64 unsigned int ndhmax; /* the maximum number of points */
66 int nhist; /* the number of histograms. There can either be
68 1 (for 'foreign lambda' histograms)
69 2 (for derivative histograms: there's
70 a 'forward' and 'backward' histogram
71 containing the minimum and maximum
72 values, respectively). */
73 std::array<std::vector<int>, 2> bin; /* the histogram(s) */
74 double dx; /* the histogram spacing in kJ/mol. This is the
75 same for the two histograms? */
76 unsigned int nbins; /* the number of bins in the histograms*/
77 std::array<int64_t, 2> x0; /* the starting point in units of spacing
79 std::array<unsigned int, 2> maxbin; /* highest bin number with data */
81 int type; /* the block type according to dhbtDH, etc. */
82 int derivative; /* The derivative direction (as an index in the lambda
83 vector) if this delta_h contains derivatives */
84 std::vector<double> lambda; /* lambda vector (or NULL if not applicable) */
85 int nlambda; /* length of the lambda vector */
86 bool written; /* whether this data has already been written out */
88 std::array<int64_t, 5> subblock_meta_l; /* metadata for an mdebin subblock for
89 I/O: for histogram counts, etc.*/
90 std::vector<double> subblock_meta_d; /* metadata subblock for I/O, used for
91 communicating doubles (i.e. the lambda
93 std::array<int, 4> subblock_meta_i; /* metadata subblock for I/O, used for
94 communicating ints (i.e. derivative indices,
98 /* the type definition is in mdebin_bar.h */
99 struct t_mde_delta_h_coll
101 t_mde_delta_h_coll(const t_inputrec& inputrec);
103 std::vector<t_mde_delta_h> dh; /* the delta h data */
104 int ndh; /* the number of delta_h structures */
106 int nlambda; /* number of bar dU delta_h structures */
107 int dh_du_index; /* the delta h data (index into dh) */
109 int ndhdl; /* number of bar dU delta_h structures */
110 int dh_dhdl_index; /* the dhdl data (index into dh) */
112 int dh_energy_index; /* energy output block (index into dh) */
113 int dh_pv_index; /* pV output block (index into dh) */
114 int dh_expanded_index; /* expanded ensemble output block (index into dh) */
116 double start_time; /* start time of the current dh collection */
117 double delta_time; /* time difference between samples */
118 bool start_time_set; /* whether the start time has been set */
119 double start_lambda; /* starting lambda for continuous motion of state*/
120 double delta_lambda; /* delta lambda, for continuous motion of state */
121 double temperature; /* the temperature of the samples*/
123 std::vector<double> native_lambda_vec; /* The lambda vector describing the current
124 lambda state if it is set (NULL otherwise) */
125 int n_lambda_vec; /* the size of the native lambda vector */
126 std::vector<int> native_lambda_components; /* the native lambda (and by extension,
127 foreign lambda) components in terms
128 of efptFEP, efptMASS, etc. */
129 int lambda_index; /* the lambda_fep_state */
131 std::vector<double> subblock_d; /* for writing a metadata mdebin subblock for I/O */
132 std::vector<int> subblock_i; /* for writing a metadata mdebin subblock for I/O */
135 /* add a bunch of samples to the delta_h collection
137 dhdl = the hamiltonian derivatives
138 U = the array with energies: from enerd->enerpart_lambda.
139 time = the current simulation time.
140 current_lambda = current lambda values : primarily useful for continuous processes
141 fep_state = current fep_state
144 /* add a bunch of samples - note fep_state is double to allow for better data storage */
145 void mde_delta_h_coll_add_dh(t_mde_delta_h_coll* dhc,
149 gmx::ArrayRef<double> dhdl,
153 /* write the data associated with the du blocks collection as a collection
157 nblock = the current number of blocks */
158 void mde_delta_h_coll_handle_block(t_mde_delta_h_coll* dhc, t_enxframe* fr, int nblock);
161 /* reset the collection of delta_h buffers for a new round of
163 void mde_delta_h_coll_reset(t_mde_delta_h_coll* dhc);
166 /* set the energyhistory variables to save state */
167 void mde_delta_h_coll_update_energyhistory(const t_mde_delta_h_coll* dhc, energyhistory_t* enerhist);
169 /* restore the variables from an energyhistory */
170 void mde_delta_h_coll_restore_energyhistory(t_mde_delta_h_coll* dhc, const delta_h_history_t* deltaH);
172 #endif /* GMX_MDLIB_MDEBIN_BAR_H */