* Copyright (c) 1991-2000, University of Groningen, The Netherlands.
* Copyright (c) 2001-2004, The GROMACS development team.
* Copyright (c) 2013,2014,2015,2016,2017 by the GROMACS development team.
- * Copyright (c) 2018,2019,2020, by the GROMACS development team, led by
+ * Copyright (c) 2018,2019,2020,2021, by the GROMACS development team, led by
* Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
* and including many others, as listed in the AUTHORS file in the
* top-level source directory and at http://www.gromacs.org.
enum class AccelerationType
{
none,
- group,
cosine
};
* \param[in] dt The time step.
* \param[in] dtPressureCouple Time step for pressure coupling, is 0 when no pressure
* coupling should be applied at this step.
- * \param[in] accel Acceleration per group.
* \param[in] md Atom properties.
* \param[in] ekind Kinetic energy data.
* \param[in] box The box dimensions.
bool doNoseHoover,
real dt,
real dtPressureCouple,
- const rvec* accel,
const t_mdatoms* md,
const gmx_ekindata_t* ekind,
const matrix box,
* Holian et al. Phys Rev E 52(3) : 2338, 1995
*/
- gmx::ArrayRef<const t_grp_tcstat> tcstat = ekind->tcstat;
- gmx::ArrayRef<const t_grp_acc> grpstat = ekind->grpstat;
- const unsigned short* cTC = md->cTC;
- const unsigned short* cACC = md->cACC;
+ gmx::ArrayRef<const t_grp_tcstat> tcstat = ekind->tcstat;
+ const unsigned short* cTC = md->cTC;
const rvec* gmx_restrict invMassPerDim = md->invMassPerDim;
/* Initialize group values, changed later when multiple groups are used */
- int ga = 0;
int gt = 0;
real omega_Z = 2 * static_cast<real>(M_PI) / box[ZZ][ZZ];
switch (accelerationType)
{
case AccelerationType::none: copy_rvec(v[n], vRel); break;
- case AccelerationType::group:
- if (cACC)
- {
- ga = cACC[n];
- }
- /* Avoid scaling the group velocity */
- rvec_sub(v[n], grpstat[ga].u, vRel);
- break;
case AccelerationType::cosine:
cosineZ = std::cos(x[n][ZZ] * omega_Z);
vCosine = cosineZ * ekind->cosacc.vcos;
switch (accelerationType)
{
case AccelerationType::none: break;
- case AccelerationType::group:
- /* Add back the mean velocity and apply acceleration */
- vNew += grpstat[ga].u[d] + accel[ga][d] * dt;
- break;
case AccelerationType::cosine:
if (d == XX)
{
rvec* gmx_restrict xprime,
rvec* gmx_restrict v,
const rvec* gmx_restrict f,
- const rvec* gmx_restrict accel,
const int etc,
const int epc,
const int nsttcouple,
real dtPressureCouple = (doParrinelloRahman ? nstpcouple * dt : 0);
- /* NEMD (also cosine) acceleration is applied in updateMDLeapFrogGeneral */
- bool doAcceleration = (ekind->bNEMD || ekind->cosacc.cos_accel != 0);
+ /* NEMD cosine acceleration is applied in updateMDLeapFrogGeneral */
+ bool doAcceleration = (ekind->cosacc.cos_accel != 0);
if (doNoseHoover || doPROffDiagonal || doAcceleration)
{
if (!doAcceleration)
{
updateMDLeapfrogGeneral<AccelerationType::none>(
- start, nrend, doNoseHoover, dt, dtPressureCouple, accel, md, ekind, box, x, xprime, v, f, nh_vxi, nsttcouple, stepM);
- }
- else if (ekind->bNEMD)
- {
- updateMDLeapfrogGeneral<AccelerationType::group>(
- start, nrend, doNoseHoover, dt, dtPressureCouple, accel, md, ekind, box, x, xprime, v, f, nh_vxi, nsttcouple, stepM);
+ start, nrend, doNoseHoover, dt, dtPressureCouple, md, ekind, box, x, xprime, v, f, nh_vxi, nsttcouple, stepM);
}
else
{
updateMDLeapfrogGeneral<AccelerationType::cosine>(
- start, nrend, doNoseHoover, dt, dtPressureCouple, accel, md, ekind, box, x, xprime, v, f, nh_vxi, nsttcouple, stepM);
+ start, nrend, doNoseHoover, dt, dtPressureCouple, md, ekind, box, x, xprime, v, f, nh_vxi, nsttcouple, stepM);
}
}
else
static void do_update_vv_vel(int start,
int nrend,
real dt,
- const rvec accel[],
const ivec nFreeze[],
const real invmass[],
const unsigned short ptype[],
const unsigned short cFREEZE[],
- const unsigned short cACC[],
rvec v[],
const rvec f[],
gmx_bool bExtended,
real veta,
real alpha)
{
- int gf = 0, ga = 0;
+ int gf = 0;
int n, d;
real g, mv1, mv2;
{
gf = cFREEZE[n];
}
- if (cACC)
- {
- ga = cACC[n];
- }
for (d = 0; d < DIM; d++)
{
if ((ptype[n] != eptVSite) && (ptype[n] != eptShell) && !nFreeze[gf][d])
{
- v[n][d] = mv1 * (mv1 * v[n][d] + 0.5 * (w_dt * mv2 * f[n][d])) + 0.5 * accel[ga][d] * dt;
+ v[n][d] = mv1 * (mv1 * v[n][d] + 0.5 * (w_dt * mv2 * f[n][d]));
}
else
{
int start,
int nrend,
real dt,
- const rvec accel[],
const ivec nFreeze[],
const real invmass[],
const unsigned short ptype[],
const unsigned short cFREEZE[],
- const unsigned short cACC[],
const unsigned short cTC[],
const rvec x[],
rvec xprime[],
int seed,
const int* gatindex)
{
- // cTC, cACC and cFREEZE can be nullptr any time, but various
+ // cTC and cFREEZE can be nullptr any time, but various
// instantiations do not make sense with particular pointer
// values.
if (updateType == SDUpdate::ForcesOnly)
if (updateType == SDUpdate::FrictionAndNoiseOnly)
{
GMX_ASSERT(f == nullptr, "SD update with only noise cannot handle forces");
- GMX_ASSERT(cACC == nullptr, "SD update with only noise cannot handle acceleration groups");
}
if (updateType == SDUpdate::Combined)
{
real inverseMass = invmass[n];
real invsqrtMass = std::sqrt(inverseMass);
- int freezeGroup = cFREEZE ? cFREEZE[n] : 0;
- int accelerationGroup = cACC ? cACC[n] : 0;
- int temperatureGroup = cTC ? cTC[n] : 0;
+ int freezeGroup = cFREEZE ? cFREEZE[n] : 0;
+ int temperatureGroup = cTC ? cTC[n] : 0;
for (int d = 0; d < DIM; d++)
{
{
if (updateType == SDUpdate::ForcesOnly)
{
- real vn = v[n][d] + (inverseMass * f[n][d] + accel[accelerationGroup][d]) * dt;
+ real vn = v[n][d] + inverseMass * f[n][d] * dt;
v[n][d] = vn;
// Simple position update.
xprime[n][d] = x[n][d] + v[n][d] * dt;
}
else
{
- real vn = v[n][d] + (inverseMass * f[n][d] + accel[accelerationGroup][d]) * dt;
+ real vn = v[n][d] + inverseMass * f[n][d] * dt;
v[n][d] = (vn * sd.sdc[temperatureGroup].em
+ invsqrtMass * sd.sdsig[temperatureGroup].V * dist(rng));
// Here we include half of the friction+noise
rvec* gmx_restrict xprime,
rvec* gmx_restrict v,
const rvec* gmx_restrict f,
- const rvec accel[],
const ivec nFreeze[],
const real invmass[],
const unsigned short ptype[],
const unsigned short cFREEZE[],
- const unsigned short cACC[],
const unsigned short cTC[],
int seed,
const t_commrec* cr,
{
// With constraints, the SD update is done in 2 parts
doSDUpdateGeneral<SDUpdate::ForcesOnly>(
- sd, start, nrend, dt, accel, nFreeze, invmass, ptype, cFREEZE, cACC, nullptr, x, xprime, v, f, step, seed, nullptr);
+ sd, start, nrend, dt, nFreeze, invmass, ptype, cFREEZE, nullptr, x, xprime, v, f, step, seed, nullptr);
}
else
{
start,
nrend,
dt,
- accel,
nFreeze,
invmass,
ptype,
cFREEZE,
- cACC,
cTC,
x,
xprime,
start_th,
end_th,
dt,
- inputRecord.opts.acc,
inputRecord.opts.nFreeze,
md->invmass,
md->ptype,
md->cFREEZE,
- nullptr,
md->cTC,
state->x.rvec_array(),
xp_.rvec_array(),
xp_rvec,
v_rvec,
f_rvec,
- inputRecord.opts.acc,
inputRecord.etc,
inputRecord.epc,
inputRecord.nsttcouple,
xp_rvec,
v_rvec,
f_rvec,
- inputRecord.opts.acc,
inputRecord.opts.nFreeze,
md->invmass,
md->ptype,
md->cFREEZE,
- md->cACC,
md->cTC,
inputRecord.ld_seed,
cr,
do_update_vv_vel(start_th,
end_th,
dt,
- inputRecord.opts.acc,
inputRecord.opts.nFreeze,
md->invmass,
md->ptype,
md->cFREEZE,
- md->cACC,
v_rvec,
f_rvec,
bExtended,