Changed raw arrays to std::vector
Moved init functionality into constructor
Made owning gmx_update_t use a unique_ptr
Changed functions using struct to take constant references
Change-Id: I3664d9cc2b7e89eba57ad36a89e9c44b7895a7ef
void andersen_tcoupl(const t_inputrec *ir, int64_t step,
const t_commrec *cr, const t_mdatoms *md,
gmx::ArrayRef<gmx::RVec> v,
void andersen_tcoupl(const t_inputrec *ir, int64_t step,
const t_commrec *cr, const t_mdatoms *md,
gmx::ArrayRef<gmx::RVec> v,
- real rate, const gmx_bool *randomize, const real *boltzfac)
+ real rate, const std::vector<bool> &randomize,
+ gmx::ArrayRef<const real> boltzfac)
{
const int *gatindex = (DOMAINDECOMP(cr) ? cr->dd->globalAtomIndices.data() : nullptr);
int i;
{
const int *gatindex = (DOMAINDECOMP(cr) ? cr->dd->globalAtomIndices.data() : nullptr);
int i;
+#include "gromacs/compat/make_unique.h"
#include "gromacs/domdec/domdec_struct.h"
#include "gromacs/fileio/confio.h"
#include "gromacs/gmxlib/network.h"
#include "gromacs/domdec/domdec_struct.h"
#include "gromacs/fileio/confio.h"
#include "gromacs/gmxlib/network.h"
real V;
} gmx_sd_sigma_t;
real V;
} gmx_sd_sigma_t;
+ std::vector<real> bd_rf;
- gmx_sd_const_t *sdc;
- gmx_sd_sigma_t *sdsig;
+ std::vector<gmx_sd_const_t> sdc;
+ std::vector<gmx_sd_sigma_t> sdsig;
/* andersen temperature control stuff */
/* andersen temperature control stuff */
- gmx_bool *randomize_group;
- real *boltzfac;
-} gmx_stochd_t;
+ std::vector<bool> randomize_group;
+ std::vector<real> boltzfac;
+
+ gmx_stochd_t(const t_inputrec *ir);
+};
+ std::unique_ptr<gmx_stochd_t> sd;
/* xprime for constraint algorithms */
/* xprime for constraint algorithms */
- PaddedVector<gmx::RVec> xp;
+ PaddedVector<gmx::RVec> xp;
/* Variables for the deform algorithm */
int64_t deformref_step;
/* Variables for the deform algorithm */
int64_t deformref_step;
}
} /* do_update_vv_pos */
}
} /* do_update_vv_pos */
-static gmx_stochd_t *init_stochd(const t_inputrec *ir)
+gmx_stochd_t::gmx_stochd_t(const t_inputrec *ir)
- gmx_stochd_t *sd;
-
- snew(sd, 1);
-
const t_grpopts *opts = &ir->opts;
const t_grpopts *opts = &ir->opts;
+ const int ngtc = opts->ngtc;
}
else if (EI_SD(ir->eI))
{
}
else if (EI_SD(ir->eI))
{
- snew(sd->sdc, ngtc);
- snew(sd->sdsig, ngtc);
-
- gmx_sd_const_t *sdc = sd->sdc;
+ sdc.resize(ngtc);
+ sdsig.resize(ngtc);
for (int gt = 0; gt < ngtc; gt++)
{
for (int gt = 0; gt < ngtc; gt++)
{
}
else if (ETC_ANDERSEN(ir->etc))
{
}
else if (ETC_ANDERSEN(ir->etc))
{
- snew(sd->randomize_group, ngtc);
- snew(sd->boltzfac, ngtc);
+ randomize_group.resize(ngtc);
+ boltzfac.resize(ngtc);
/* for now, assume that all groups, if randomized, are randomized at the same rate, i.e. tau_t is the same. */
/* since constraint groups don't necessarily match up with temperature groups! This is checked in readir.c */
/* for now, assume that all groups, if randomized, are randomized at the same rate, i.e. tau_t is the same. */
/* since constraint groups don't necessarily match up with temperature groups! This is checked in readir.c */
real reft = std::max<real>(0, opts->ref_t[gt]);
if ((opts->tau_t[gt] > 0) && (reft > 0)) /* tau_t or ref_t = 0 means that no randomization is done */
{
real reft = std::max<real>(0, opts->ref_t[gt]);
if ((opts->tau_t[gt] > 0) && (reft > 0)) /* tau_t or ref_t = 0 means that no randomization is done */
{
- sd->randomize_group[gt] = TRUE;
- sd->boltzfac[gt] = BOLTZ*opts->ref_t[gt];
+ randomize_group[gt] = true;
+ boltzfac[gt] = BOLTZ*opts->ref_t[gt];
- sd->randomize_group[gt] = FALSE;
+ randomize_group[gt] = false;
}
void update_temperature_constants(gmx_update_t *upd, const t_inputrec *ir)
}
void update_temperature_constants(gmx_update_t *upd, const t_inputrec *ir)
{
gmx_update_t *upd = new(gmx_update_t);
{
gmx_update_t *upd = new(gmx_update_t);
- if (ir->eI == eiBD || EI_SD(ir->eI) || ir->etc == etcVRESCALE || ETC_ANDERSEN(ir->etc))
- {
- upd->sd = init_stochd(ir);
- }
+ upd->sd = gmx::compat::make_unique<gmx_stochd_t>(ir);
update_temperature_constants(upd, ir);
update_temperature_constants(upd, ir);
* two with only one contribution, and one with both contributions. */
template <SDUpdate updateType>
static void
* two with only one contribution, and one with both contributions. */
template <SDUpdate updateType>
static void
-doSDUpdateGeneral(gmx_stochd_t *sd,
+doSDUpdateGeneral(const gmx_stochd_t &sd,
int start, int nrend, real dt,
const rvec accel[], const ivec nFreeze[],
const real invmass[], const unsigned short ptype[],
int start, int nrend, real dt,
const rvec accel[], const ivec nFreeze[],
const real invmass[], const unsigned short ptype[],
gmx::ThreeFry2x64<0> rng(seed, gmx::RandomDomain::UpdateCoordinates);
gmx::TabulatedNormalDistribution<real, 14> dist;
gmx::ThreeFry2x64<0> rng(seed, gmx::RandomDomain::UpdateCoordinates);
gmx::TabulatedNormalDistribution<real, 14> dist;
- gmx_sd_const_t *sdc = sd->sdc;
- gmx_sd_sigma_t *sig = sd->sdsig;
-
for (int n = start; n < nrend; n++)
{
int globalAtomIndex = gatindex ? gatindex[n] : n;
for (int n = start; n < nrend; n++)
{
int globalAtomIndex = gatindex ? gatindex[n] : n;
else if (updateType == SDUpdate::FrictionAndNoiseOnly)
{
real vn = v[n][d];
else if (updateType == SDUpdate::FrictionAndNoiseOnly)
{
real vn = v[n][d];
- v[n][d] = (vn*sdc[temperatureGroup].em +
- invsqrtMass*sig[temperatureGroup].V*dist(rng));
+ v[n][d] = (vn*sd.sdc[temperatureGroup].em +
+ invsqrtMass*sd.sdsig[temperatureGroup].V*dist(rng));
// The previous phase already updated the
// positions with a full v*dt term that must
// now be half removed.
// The previous phase already updated the
// positions with a full v*dt term that must
// now be half removed.
else
{
real vn = v[n][d] + (inverseMass*f[n][d] + accel[accelerationGroup][d])*dt;
else
{
real vn = v[n][d] + (inverseMass*f[n][d] + accel[accelerationGroup][d])*dt;
- v[n][d] = (vn*sdc[temperatureGroup].em +
- invsqrtMass*sig[temperatureGroup].V*dist(rng));
+ v[n][d] = (vn*sd.sdc[temperatureGroup].em +
+ invsqrtMass*sd.sdsig[temperatureGroup].V*dist(rng));
// Here we include half of the friction+noise
// update of v into the position update.
xprime[n][d] = x[n][d] + 0.5*(vn + v[n][d])*dt;
// Here we include half of the friction+noise
// update of v into the position update.
xprime[n][d] = x[n][d] + 0.5*(vn + v[n][d])*dt;
getThreadAtomRange(nth, th, homenr, &start_th, &end_th);
doSDUpdateGeneral<SDUpdate::FrictionAndNoiseOnly>
getThreadAtomRange(nth, th, homenr, &start_th, &end_th);
doSDUpdateGeneral<SDUpdate::FrictionAndNoiseOnly>
start_th, end_th, dt,
inputrec->opts.acc, inputrec->opts.nFreeze,
md->invmass, md->ptype,
start_th, end_th, dt,
inputrec->opts.acc, inputrec->opts.nFreeze,
md->invmass, md->ptype,
{
// With constraints, the SD update is done in 2 parts
doSDUpdateGeneral<SDUpdate::ForcesOnly>
{
// With constraints, the SD update is done in 2 parts
doSDUpdateGeneral<SDUpdate::ForcesOnly>
start_th, end_th, dt,
inputrec->opts.acc, inputrec->opts.nFreeze,
md->invmass, md->ptype,
start_th, end_th, dt,
inputrec->opts.acc, inputrec->opts.nFreeze,
md->invmass, md->ptype,
else
{
doSDUpdateGeneral<SDUpdate::Combined>
else
{
doSDUpdateGeneral<SDUpdate::Combined>
start_th, end_th, dt,
inputrec->opts.acc, inputrec->opts.nFreeze,
md->invmass, md->ptype,
start_th, end_th, dt,
inputrec->opts.acc, inputrec->opts.nFreeze,
md->invmass, md->ptype,
md->cFREEZE, md->cTC,
x_rvec, xp_rvec, v_rvec, f_rvec,
inputrec->bd_fric,
md->cFREEZE, md->cTC,
x_rvec, xp_rvec, v_rvec, f_rvec,
inputrec->bd_fric,
step, inputrec->ld_seed, DOMAINDECOMP(cr) ? cr->dd->globalAtomIndices.data() : nullptr);
break;
case (eiVV):
step, inputrec->ld_seed, DOMAINDECOMP(cr) ? cr->dd->globalAtomIndices.data() : nullptr);
break;
case (eiVV):
if ((ir->etc == etcANDERSEN) || do_per_step(step, roundToInt(1.0/rate)))
{
andersen_tcoupl(ir, step, cr, md, v, rate,
if ((ir->etc == etcANDERSEN) || do_per_step(step, roundToInt(1.0/rate)))
{
andersen_tcoupl(ir, step, cr, md, v, rate,
- upd->sd->randomize_group, upd->sd->boltzfac);
+ upd->sd->randomize_group,
+ upd->sd->boltzfac);
return TRUE;
}
return FALSE;
return TRUE;
}
return FALSE;
void andersen_tcoupl(const t_inputrec *ir, int64_t step,
const t_commrec *cr, const t_mdatoms *md,
gmx::ArrayRef<gmx::RVec> v,
void andersen_tcoupl(const t_inputrec *ir, int64_t step,
const t_commrec *cr, const t_mdatoms *md,
gmx::ArrayRef<gmx::RVec> v,
- real rate, const gmx_bool *randomize, const real *boltzfac);
+ real rate, const std::vector<bool> &randomize,
+ gmx::ArrayRef<const real> boltzfac);
void nosehoover_tcoupl(const t_grpopts *opts, const gmx_ekindata_t *ekind, real dt,
double xi[], double vxi[], const t_extmass *MassQ);
void nosehoover_tcoupl(const t_grpopts *opts, const gmx_ekindata_t *ekind, real dt,
double xi[], double vxi[], const t_extmass *MassQ);