1 /* 2 * This file is part of the GROMACS molecular simulation package. 3 * 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, 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. 11 * 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. 16 * 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. 21 * 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. 26 * 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. 34 * 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. 37 */ 38 39 #ifndef _mdebin_bar_h 40 #define _mdebin_bar_h 41 42 #include "gromacs/mdlib/energyoutput.h" 43 44 /* The functions & data structures here describe writing 45 energy differences (or their histogram )for use with g_bar */ 46 47 class delta_h_history_t; 48 struct t_enxframe; 49 50 /* Data for one foreign lambda, or derivative. */ 51 typedef struct 52 { 53 real* dh; /* the raw energy data. */ 54 float* dhf; /* raw difference data -- in floats, for storage. */ 55 unsigned int ndh; /* number of data points */ 56 unsigned int ndhmax; /* the maximum number of points */ 57 58 int nhist; /* the number of histograms. There can either be 59 0 (for no histograms) 60 1 (for 'foreign lambda' histograms) 61 2 (for derivative histograms: there's 62 a 'forward' and 'backward' histogram 63 containing the minimum and maximum 64 values, respectively). */ 65 int* bin[2]; /* the histogram(s) */ 66 double dx; /* the histogram spacing in kJ/mol. This is the 67 same for the two histograms? */ 68 unsigned int nbins; /* the number of bins in the histograms*/ 69 int64_t x0[2]; /* the starting point in units of spacing 70 of the histogram */ 71 unsigned int maxbin[2]; /* highest bin number with data */ 72 73 int type; /* the block type according to dhbtDH, etc. */ 74 int derivative; /* The derivative direction (as an index in the lambda 75 vector) if this delta_h contains derivatives */ 76 double* lambda; /* lambda vector (or NULL if not applicable) */ 77 int nlambda; /* length of the lambda vector */ 78 gmx_bool written; /* whether this data has already been written out */ 79 80 int64_t subblock_meta_l[5]; /* metadata for an mdebin subblock for 81 I/O: for histogram counts, etc.*/ 82 double* subblock_meta_d; /* metadata subblock for I/O, used for 83 communicating doubles (i.e. the lambda 84 vector) */ 85 int subblock_meta_i[4]; /* metadata subblock for I/O, used for 86 communicating ints (i.e. derivative indices, 87 etc.) */ 88 } t_mde_delta_h; 89 90 /* the type definition is in mdebin_bar.h */ 91 struct t_mde_delta_h_coll 92 { 93 t_mde_delta_h* dh; /* the delta h data */ 94 int ndh; /* the number of delta_h structures */ 95 96 int nlambda; /* number of bar dU delta_h structures */ 97 t_mde_delta_h* dh_du; /* the delta h data (pointer into dh) */ 98 99 int ndhdl; /* number of bar dU delta_h structures */ 100 t_mde_delta_h* dh_dhdl; /* the dhdl data (pointer into dh) */ 101 102 t_mde_delta_h* dh_energy; /* energy output block (pointer into dh) */ 103 t_mde_delta_h* dh_pv; /* pV output block (pointer into dh) */ 104 t_mde_delta_h* dh_expanded; /* expanded ensemble output block (pointer 105 into dh) */ 106 107 double start_time; /* start time of the current dh collection */ 108 double delta_time; /* time difference between samples */ 109 gmx_bool start_time_set; /* whether the start time has been set */ 110 double start_lambda; /* starting lambda for continuous motion of state*/ 111 double delta_lambda; /* delta lambda, for continuous motion of state */ 112 double temperature; /* the temperature of the samples*/ 113 114 double* native_lambda_vec; /* The lambda vector describing the current 115 lambda state if it is set (NULL otherwise) */ 116 int n_lambda_vec; /* the size of the native lambda vector */ 117 int* native_lambda_components; /* the native lambda (and by extension, 118 foreign lambda) components in terms 119 of efptFEP, efptMASS, etc. */ 120 int lambda_index; /* the lambda_fep_state */ 121 122 double* subblock_d; /* for writing a metadata mdebin subblock for I/O */ 123 int* subblock_i; /* for writing a metadata mdebin subblock for I/O */ 124 125 double* lambda_vec_subblock; /* native lambda vector data subblock for 126 I/O */ 127 int* lambda_index_subblock; /* lambda vector index data subblock for I/O */ 128 }; 129 130 131 /* initialize a collection of delta h histograms/sets 132 dhc = the collection 133 ir = the input record */ 134 135 void mde_delta_h_coll_init(t_mde_delta_h_coll* dhc, const t_inputrec* ir); 136 137 void done_mde_delta_h_coll(t_mde_delta_h_coll* dhc); 138 139 /* add a bunch of samples to the delta_h collection 140 dhc = the collection 141 dhdl = the hamiltonian derivatives 142 U = the array with energies: from enerd->enerpart_lambda. 143 time = the current simulation time. 144 current_lambda = current lambda values : primarily useful for continuous processes 145 fep_state = current fep_state 146 */ 147 148 /* add a bunch of samples - note fep_state is double to allow for better data storage */ 149 void mde_delta_h_coll_add_dh(t_mde_delta_h_coll* dhc, 150 double fep_state, 151 double energy, 152 double pV, 153 double* dhdl, 154 double* foreign_dU, 155 double time); 156 157 /* write the data associated with the du blocks collection as a collection 158 of mdebin blocks. 159 dhc = the collection 160 fr = the enxio frame 161 nblock = the current number of blocks */ 162 void mde_delta_h_coll_handle_block(t_mde_delta_h_coll* dhc, t_enxframe* fr, int nblock); 163 164 165 /* reset the collection of delta_h buffers for a new round of 166 data gathering */ 167 void mde_delta_h_coll_reset(t_mde_delta_h_coll* dhc); 168 169 170 /* set the energyhistory variables to save state */ 171 void mde_delta_h_coll_update_energyhistory(const t_mde_delta_h_coll* dhc, energyhistory_t* enerhist); 172 173 /* restore the variables from an energyhistory */ 174 void mde_delta_h_coll_restore_energyhistory(t_mde_delta_h_coll* dhc, const delta_h_history_t* deltaH); 175 176 #endif /* _mdebin_bar_h */ 177