dh->derivative = derivative;
dh->nlambda = nlambda;
- snew(dh->lambda, nlambda);
+ dh->lambda.resize(nlambda);
for (i = 0; i < nlambda; i++)
{
assert(lambda);
}
- snew(dh->subblock_meta_d, dh->nlambda + 1);
+ dh->subblock_meta_d.resize(dh->nlambda + 1);
dh->ndhmax = ndhmax + 2;
- for (i = 0; i < 2; i++)
- {
- dh->bin[i] = nullptr;
- }
- snew(dh->dh, dh->ndhmax);
- snew(dh->dhf, dh->ndhmax);
+ dh->dh.resize(dh->ndhmax);
+ dh->dhf.resize(dh->ndhmax);
if (nbins <= 0 || dx < GMX_REAL_EPS * 10)
{
dh->nbins = nbins;
for (i = 0; i < dh->nhist; i++)
{
- snew(dh->bin[i], dh->nbins);
+ dh->bin[i].resize(dh->nbins);
}
}
mde_delta_h_reset(dh);
}
-static void done_mde_delta_h(t_mde_delta_h* dh)
-{
- sfree(dh->lambda);
- sfree(dh->subblock_meta_d);
- sfree(dh->dh);
- sfree(dh->dhf);
- for (int i = 0; i < dh->nhist; i++)
- {
- sfree(dh->bin[i]);
- }
-}
-
/* Add a value to the delta_h list */
static void mde_delta_h_add_dh(t_mde_delta_h* dh, double delta_h)
{
the main delta_h_coll) */
blk->sub[0].nr = 2;
blk->sub[0].type = XdrDataType::Int;
- blk->sub[0].ival = dh->subblock_meta_i;
+ blk->sub[0].ival = dh->subblock_meta_i.data();
/* subblock 2 */
for (i = 0; i < dh->nlambda; i++)
}
blk->sub[1].nr = dh->nlambda;
blk->sub[1].type = XdrDataType::Double;
- blk->sub[1].dval = dh->subblock_meta_d;
+ blk->sub[1].dval = dh->subblock_meta_d.data();
/* subblock 3 */
/* check if there's actual data to be written. */
{
dh->dhf[i] = static_cast<float>(dh->dh[i]);
}
- blk->sub[2].fval = dh->dhf;
+ blk->sub[2].fval = dh->dhf.data();
dh->written = TRUE;
}
else
dh->subblock_meta_d[1] = dh->dx;
blk->sub[0].nr = 2 + ((dh->nlambda > 1) ? dh->nlambda : 0);
blk->sub[0].type = XdrDataType::Double;
- blk->sub[0].dval = dh->subblock_meta_d;
+ blk->sub[0].dval = dh->subblock_meta_d.data();
/* subblock 2: the starting point(s) as a long integer */
dh->subblock_meta_l[0] = nhist_written;
blk->sub[1].nr = nhist_written + 3;
blk->sub[1].type = XdrDataType::Int64;
- blk->sub[1].lval = dh->subblock_meta_l;
+ blk->sub[1].lval = dh->subblock_meta_l.data();
/* subblock 3 + 4 : the histogram data */
for (i = 0; i < nhist_written; i++)
blk->sub[i + 2].nr = dh->maxbin[i] + 1; /* it's +1 because size=index+1
in C */
blk->sub[i + 2].type = XdrDataType::Int;
- blk->sub[i + 2].ival = dh->bin[i];
+ blk->sub[i + 2].ival = dh->bin[i].data();
}
}
}
/* initialize the collection*/
-void mde_delta_h_coll_init(t_mde_delta_h_coll* dhc, const t_inputrec& inputrec)
+t_mde_delta_h_coll::t_mde_delta_h_coll(const t_inputrec& inputrec)
{
int i, j, n;
double* lambda_vec;
int ndhmax = inputrec.nstenergy / inputrec.nstcalcenergy;
t_lambda* fep = inputrec.fepvals.get();
- dhc->temperature = inputrec.opts.ref_t[0]; /* only store system temperature */
- dhc->start_time = 0.;
- dhc->delta_time = inputrec.delta_t * inputrec.fepvals->nstdhdl;
- dhc->start_time_set = FALSE;
+ temperature = inputrec.opts.ref_t[0]; /* only store system temperature */
+ start_time = 0.;
+ delta_time = inputrec.delta_t * inputrec.fepvals->nstdhdl;
+ start_time_set = FALSE;
/* this is the compatibility lambda value. If it is >=0, it is valid,
and there is either an old-style lambda or a slow growth simulation. */
- dhc->start_lambda = inputrec.fepvals->init_lambda;
+ start_lambda = inputrec.fepvals->init_lambda;
/* for continuous change of lambda values */
- dhc->delta_lambda = inputrec.fepvals->delta_lambda * inputrec.fepvals->nstdhdl;
+ delta_lambda = inputrec.fepvals->delta_lambda * inputrec.fepvals->nstdhdl;
- if (dhc->start_lambda < 0)
+ if (start_lambda < 0)
{
/* create the native lambda vectors */
- dhc->lambda_index = fep->init_fep_state;
- dhc->n_lambda_vec = 0;
+ lambda_index = fep->init_fep_state;
+ n_lambda_vec = 0;
for (auto i : keysOf(fep->separate_dvdl))
{
if (fep->separate_dvdl[i])
{
- dhc->n_lambda_vec++;
+ n_lambda_vec++;
}
}
- snew(dhc->native_lambda_vec, dhc->n_lambda_vec);
- snew(dhc->native_lambda_components, dhc->n_lambda_vec);
+ native_lambda_vec.resize(n_lambda_vec);
+ native_lambda_components.resize(n_lambda_vec);
j = 0;
for (auto i : keysOf(fep->separate_dvdl))
{
if (fep->separate_dvdl[i])
{
- dhc->native_lambda_components[j] = static_cast<int>(i);
+ native_lambda_components[j] = static_cast<int>(i);
if (fep->init_fep_state >= 0 && fep->init_fep_state < fep->n_lambda)
{
- dhc->native_lambda_vec[j] = fep->all_lambda[i][fep->init_fep_state];
+ native_lambda_vec[j] = fep->all_lambda[i][fep->init_fep_state];
}
else
{
- dhc->native_lambda_vec[j] = -1;
+ native_lambda_vec[j] = -1;
}
j++;
}
else
{
/* don't allocate the meta-data subblocks for lambda vectors */
- dhc->native_lambda_vec = nullptr;
- dhc->n_lambda_vec = 0;
- dhc->native_lambda_components = nullptr;
- dhc->lambda_index = -1;
+ n_lambda_vec = 0;
+ lambda_index = -1;
}
/* allocate metadata subblocks */
- snew(dhc->subblock_d, c_subblockDNumPreEntries + dhc->n_lambda_vec);
- snew(dhc->subblock_i, c_subblockINumPreEntries + dhc->n_lambda_vec);
+ subblock_d.resize(c_subblockDNumPreEntries + n_lambda_vec);
+ subblock_i.resize(c_subblockINumPreEntries + n_lambda_vec);
/* now decide which data to write out */
- dhc->nlambda = 0;
- dhc->ndhdl = 0;
- dhc->dh_expanded_index = -1;
- dhc->dh_energy_index = -1;
- dhc->dh_pv_index = -1;
+ nlambda = 0;
+ ndhdl = 0;
+ dh_expanded_index = -1;
+ dh_energy_index = -1;
+ dh_pv_index = -1;
/* total number of raw data point collections in the sample */
- dhc->ndh = 0;
+ ndh = 0;
{
gmx_bool bExpanded = FALSE;
{
if (inputrec.fepvals->separate_dvdl[i])
{
- dhc->ndh += 1;
- dhc->ndhdl += 1;
+ ndh += 1;
+ ndhdl += 1;
}
}
}
/* add the lambdas */
- dhc->nlambda = inputrec.fepvals->lambda_stop_n - inputrec.fepvals->lambda_start_n;
- dhc->ndh += dhc->nlambda;
+ nlambda = inputrec.fepvals->lambda_stop_n - inputrec.fepvals->lambda_start_n;
+ ndh += nlambda;
/* another compatibility check */
- if (dhc->start_lambda < 0)
+ if (start_lambda < 0)
{
/* include one more for the specification of the state, by lambda or
fep_state*/
if (inputrec.expandedvals->elmcmove > LambdaMoveCalculation::No)
{
- dhc->ndh += 1;
+ ndh += 1;
bExpanded = TRUE;
}
/* whether to print energies */
if (inputrec.fepvals->edHdLPrintEnergy != FreeEnergyPrintEnergy::No)
{
- dhc->ndh += 1;
+ ndh += 1;
bEnergy = TRUE;
}
if (inputrec.epc > PressureCoupling::No)
{
- dhc->ndh += 1; /* include pressure-volume work */
+ ndh += 1; /* include pressure-volume work */
bPV = TRUE;
}
}
/* allocate them */
- snew(dhc->dh, dhc->ndh);
+ dh.resize(ndh);
/* now initialize them */
/* the order, for now, must match that of the dhdl.xvg file because of
n = 0;
if (bExpanded)
{
- dhc->dh_expanded_index = n;
- mde_delta_h_init(dhc->dh + n,
+ dh_expanded_index = n;
+ mde_delta_h_init(&dh[n],
inputrec.fepvals->dh_hist_size,
inputrec.fepvals->dh_hist_spacing,
ndhmax,
}
if (bEnergy)
{
- dhc->dh_energy_index = n;
- mde_delta_h_init(dhc->dh + n,
+ dh_energy_index = n;
+ mde_delta_h_init(&dh[n],
inputrec.fepvals->dh_hist_size,
inputrec.fepvals->dh_hist_spacing,
ndhmax,
n_lambda_components = 0;
if (fep->dhdl_derivatives == DhDlDerivativeCalculation::Yes)
{
- dhc->dh_dhdl_index = n;
+ dh_dhdl_index = n;
for (auto i : keysOf(fep->separate_dvdl))
{
if (inputrec.fepvals->separate_dvdl[i])
{
/* we give it init_lambda for compatibility */
- mde_delta_h_init(dhc->dh + n,
+ mde_delta_h_init(&dh[n],
inputrec.fepvals->dh_hist_size,
inputrec.fepvals->dh_hist_spacing,
ndhmax,
}
}
/* add the lambdas */
- dhc->dh_du_index = n;
+ dh_du_index = n;
snew(lambda_vec, n_lambda_components);
for (i = inputrec.fepvals->lambda_start_n; i < inputrec.fepvals->lambda_stop_n; i++)
{
}
}
- mde_delta_h_init(dhc->dh + n,
+ mde_delta_h_init(&dh[n],
inputrec.fepvals->dh_hist_size,
inputrec.fepvals->dh_hist_spacing,
ndhmax,
sfree(lambda_vec);
if (bPV)
{
- dhc->dh_pv_index = n;
- mde_delta_h_init(dhc->dh + n,
+ dh_pv_index = n;
+ mde_delta_h_init(&dh[n],
inputrec.fepvals->dh_hist_size,
inputrec.fepvals->dh_hist_spacing,
ndhmax,
}
}
-void done_mde_delta_h_coll(t_mde_delta_h_coll* dhc)
-{
- if (dhc == nullptr)
- {
- return;
- }
- sfree(dhc->native_lambda_vec);
- sfree(dhc->native_lambda_components);
- sfree(dhc->subblock_d);
- sfree(dhc->subblock_i);
- for (int i = 0; i < dhc->ndh; ++i)
- {
- done_mde_delta_h(&dhc->dh[i]);
- }
- sfree(dhc->dh);
- sfree(dhc);
-}
-
/* add a bunch of samples - note fep_state is double to allow for better data storage */
void mde_delta_h_coll_add_dh(t_mde_delta_h_coll* dhc,
double fep_state,
/* only allocate lambda vector component blocks if they must be written out
for backward compatibility */
- if (dhc->native_lambda_components != nullptr)
+ if (!dhc->native_lambda_components.empty())
{
add_subblocks_enxblock(blk, 2);
}
dhc->subblock_d[3] = dhc->start_lambda; /* old-style lambda at starttime */
dhc->subblock_d[4] = dhc->delta_lambda; /* lambda diff. between samples */
/* set the lambda vector components if they exist */
- if (dhc->native_lambda_components != nullptr)
+ if (!dhc->native_lambda_components.empty())
{
for (i = 0; i < dhc->n_lambda_vec; i++)
{
blk->id = enxDHCOLL;
blk->sub[0].nr = c_subblockDNumPreEntries + dhc->n_lambda_vec;
blk->sub[0].type = XdrDataType::Double;
- blk->sub[0].dval = dhc->subblock_d;
+ blk->sub[0].dval = dhc->subblock_d.data();
- if (dhc->native_lambda_components != nullptr)
+ if (!dhc->native_lambda_components.empty())
{
dhc->subblock_i[0] = dhc->lambda_index;
/* set the lambda vector component IDs if they exist */
}
blk->sub[1].nr = c_subblockINumPreEntries + dhc->n_lambda_vec;
blk->sub[1].type = XdrDataType::Int;
- blk->sub[1].ival = dhc->subblock_i;
+ blk->sub[1].ival = dhc->subblock_i.data();
}
for (i = 0; i < dhc->ndh; i++)
add_blocks_enxframe(fr, nblock);
blk = fr->block + (nblock - 1);
- mde_delta_h_handle_block(dhc->dh + i, blk);
+ mde_delta_h_handle_block(&dhc->dh[i], blk);
}
}
if (dhc->dh[i].written)
{
/* we can now throw away the data */
- mde_delta_h_reset(dhc->dh + i);
+ mde_delta_h_reset(&dhc->dh[i]);
}
}
dhc->start_time_set = FALSE;
for (int i = 0; i < dhc->ndh; i++)
{
std::vector<real>& dh = deltaH->dh[i];
- dh.resize(dhc->dh[i].ndh);
- std::copy(dhc->dh[i].dh, dhc->dh[i].dh + dhc->dh[i].ndh, dh.begin());
+ for (unsigned int j = 0; j < dhc->dh[i].ndh; j++)
+ {
+ dh.emplace_back(dhc->dh[i].dh[j]);
+ }
}
deltaH->start_time = dhc->start_time;
deltaH->start_lambda = dhc->start_lambda;
* the research papers on the package. Check out http://www.gromacs.org.
*/
-#ifndef _mdebin_bar_h
-#define _mdebin_bar_h
+#ifndef GMX_MDLIB_MDEBIN_BAR_H
+#define GMX_MDLIB_MDEBIN_BAR_H
-#include "gromacs/mdlib/energyoutput.h"
-#include "gromacs/utility/arrayref.h"
+#include "gromacs/utility/real.h"
/* The functions & data structures here describe writing
energy differences (or their histogram )for use with g_bar */
class delta_h_history_t;
struct t_enxframe;
+class energyhistory_t;
+struct t_inputrec;
+
+namespace gmx
+{
+template<typename>
+class ArrayRef;
+}
/* Data for one foreign lambda, or derivative. */
typedef struct
{
- real* dh; /* the raw energy data. */
- float* dhf; /* raw difference data -- in floats, for storage. */
- unsigned int ndh; /* number of data points */
- unsigned int ndhmax; /* the maximum number of points */
-
- int nhist; /* the number of histograms. There can either be
- 0 (for no histograms)
- 1 (for 'foreign lambda' histograms)
- 2 (for derivative histograms: there's
- a 'forward' and 'backward' histogram
- containing the minimum and maximum
- values, respectively). */
- int* bin[2]; /* the histogram(s) */
- double dx; /* the histogram spacing in kJ/mol. This is the
- same for the two histograms? */
- unsigned int nbins; /* the number of bins in the histograms*/
- int64_t x0[2]; /* the starting point in units of spacing
- of the histogram */
- unsigned int maxbin[2]; /* highest bin number with data */
-
- int type; /* the block type according to dhbtDH, etc. */
- int derivative; /* The derivative direction (as an index in the lambda
- vector) if this delta_h contains derivatives */
- double* lambda; /* lambda vector (or NULL if not applicable) */
- int nlambda; /* length of the lambda vector */
- bool written; /* whether this data has already been written out */
-
- int64_t subblock_meta_l[5]; /* metadata for an mdebin subblock for
+ std::vector<real> dh; /* the raw energy data. */
+ std::vector<float> dhf; /* raw difference data -- in floats, for storage. */
+ unsigned int ndh; /* number of data points */
+ unsigned int ndhmax; /* the maximum number of points */
+
+ int nhist; /* the number of histograms. There can either be
+ 0 (for no histograms)
+ 1 (for 'foreign lambda' histograms)
+ 2 (for derivative histograms: there's
+ a 'forward' and 'backward' histogram
+ containing the minimum and maximum
+ values, respectively). */
+ std::array<std::vector<int>, 2> bin; /* the histogram(s) */
+ double dx; /* the histogram spacing in kJ/mol. This is the
+ same for the two histograms? */
+ unsigned int nbins; /* the number of bins in the histograms*/
+ std::array<int64_t, 2> x0; /* the starting point in units of spacing
+ of the histogram */
+ std::array<unsigned int, 2> maxbin; /* highest bin number with data */
+
+ int type; /* the block type according to dhbtDH, etc. */
+ int derivative; /* The derivative direction (as an index in the lambda
+ vector) if this delta_h contains derivatives */
+ std::vector<double> lambda; /* lambda vector (or NULL if not applicable) */
+ int nlambda; /* length of the lambda vector */
+ bool written; /* whether this data has already been written out */
+
+ std::array<int64_t, 5> subblock_meta_l; /* metadata for an mdebin subblock for
I/O: for histogram counts, etc.*/
- double* subblock_meta_d; /* metadata subblock for I/O, used for
+ std::vector<double> subblock_meta_d; /* metadata subblock for I/O, used for
communicating doubles (i.e. the lambda
vector) */
- int subblock_meta_i[4]; /* metadata subblock for I/O, used for
+ std::array<int, 4> subblock_meta_i; /* metadata subblock for I/O, used for
communicating ints (i.e. derivative indices,
etc.) */
} t_mde_delta_h;
/* the type definition is in mdebin_bar.h */
struct t_mde_delta_h_coll
{
- t_mde_delta_h* dh; /* the delta h data */
- int ndh; /* the number of delta_h structures */
+ t_mde_delta_h_coll(const t_inputrec& inputrec);
+
+ std::vector<t_mde_delta_h> dh; /* the delta h data */
+ int ndh; /* the number of delta_h structures */
int nlambda; /* number of bar dU delta_h structures */
int dh_du_index; /* the delta h data (index into dh) */
double delta_lambda; /* delta lambda, for continuous motion of state */
double temperature; /* the temperature of the samples*/
- double* native_lambda_vec; /* The lambda vector describing the current
+ std::vector<double> native_lambda_vec; /* The lambda vector describing the current
lambda state if it is set (NULL otherwise) */
- int n_lambda_vec; /* the size of the native lambda vector */
- int* native_lambda_components; /* the native lambda (and by extension,
+ int n_lambda_vec; /* the size of the native lambda vector */
+ std::vector<int> native_lambda_components; /* the native lambda (and by extension,
foreign lambda) components in terms
of efptFEP, efptMASS, etc. */
- int lambda_index; /* the lambda_fep_state */
+ int lambda_index; /* the lambda_fep_state */
- double* subblock_d; /* for writing a metadata mdebin subblock for I/O */
- int* subblock_i; /* for writing a metadata mdebin subblock for I/O */
+ std::vector<double> subblock_d; /* for writing a metadata mdebin subblock for I/O */
+ std::vector<int> subblock_i; /* for writing a metadata mdebin subblock for I/O */
};
-
-/* initialize a collection of delta h histograms/sets
- dhc = the collection
- inputrec = the input record */
-
-void mde_delta_h_coll_init(t_mde_delta_h_coll* dhc, const t_inputrec& inputrec);
-
-void done_mde_delta_h_coll(t_mde_delta_h_coll* dhc);
-
/* add a bunch of samples to the delta_h collection
dhc = the collection
dhdl = the hamiltonian derivatives
/* restore the variables from an energyhistory */
void mde_delta_h_coll_restore_energyhistory(t_mde_delta_h_coll* dhc, const delta_h_history_t* deltaH);
-#endif /* _mdebin_bar_h */
+#endif /* GMX_MDLIB_MDEBIN_BAR_H */
biod.pnpi.spb.ru / alexxy / gromacs.git / commitdiff