#include "gromacs/mdtypes/group.h"
#include "gromacs/mdtypes/inputrec.h"
#include "gromacs/mdtypes/md_enums.h"
-#include "gromacs/mdtypes/mdatom.h"
#include "gromacs/mdtypes/state.h"
#include "gromacs/pbcutil/boxutilities.h"
#include "gromacs/pbcutil/pbc.h"
sy_const_3, nullptr, sy_const_5 };
-void update_tcouple(int64_t step,
- const t_inputrec* inputrec,
- t_state* state,
- gmx_ekindata_t* ekind,
- const t_extmass* MassQ,
- const t_mdatoms* md)
+void update_tcouple(int64_t step,
+ const t_inputrec* inputrec,
+ t_state* state,
+ gmx_ekindata_t* ekind,
+ const t_extmass* MassQ,
+ int homenr,
+ gmx::ArrayRef<const unsigned short> cTC)
{
// This condition was explicitly checked in previous version, but should have never been satisfied
berendsen_tcoupl(inputrec, ekind, dttc, state->therm_integral);
break;
case TemperatureCoupling::NoseHoover:
- nosehoover_tcoupl(&(inputrec->opts),
- ekind,
- dttc,
- state->nosehoover_xi.data(),
- state->nosehoover_vxi.data(),
- MassQ);
+ nosehoover_tcoupl(
+ &(inputrec->opts), ekind, dttc, state->nosehoover_xi, state->nosehoover_vxi, MassQ);
break;
case TemperatureCoupling::VRescale:
- vrescale_tcoupl(inputrec, step, ekind, dttc, state->therm_integral.data());
+ vrescale_tcoupl(inputrec, step, ekind, dttc, state->therm_integral);
break;
default: gmx_fatal(FARGS, "Unknown temperature coupling algorithm");
}
/* rescale in place here */
if (EI_VV(inputrec->eI))
{
- rescale_velocities(ekind, md, 0, md->homenr, state->v.rvec_array());
+ rescale_velocities(ekind, cTC, 0, homenr, state->v);
}
}
else
t_state* state,
matrix parrinellorahmanMu,
matrix M,
- gmx_bool bInitStep)
+ bool bInitStep)
{
/* Berendsen P-coupling is completely handled after the coordinate update.
* Trotter P-coupling is handled by separate calls to trotter_update().
}
}
-void update_pcouple_after_coordinates(FILE* fplog,
- int64_t step,
- const t_inputrec* inputrec,
- const t_mdatoms* md,
- const matrix pressure,
- const matrix forceVirial,
- const matrix constraintVirial,
- matrix pressureCouplingMu,
- t_state* state,
- t_nrnb* nrnb,
- gmx::BoxDeformation* boxDeformation,
- const bool scaleCoordinates)
+void update_pcouple_after_coordinates(FILE* fplog,
+ int64_t step,
+ const t_inputrec* inputrec,
+ const int homenr,
+ gmx::ArrayRef<const unsigned short> cFREEZE,
+ const matrix pressure,
+ const matrix forceVirial,
+ const matrix constraintVirial,
+ matrix pressureCouplingMu,
+ t_state* state,
+ t_nrnb* nrnb,
+ gmx::BoxDeformation* boxDeformation,
+ const bool scaleCoordinates)
{
- int start = 0;
- int homenr = md->homenr;
+ int start = 0;
/* Cast to real for faster code, no loss in precision (see comment above) */
real dt = inputrec->delta_t;
constraintVirial,
pressureCouplingMu,
&state->baros_integral);
- pressureCouplingScaleBoxAndCoordinates<PressureCoupling::Berendsen>(inputrec,
- pressureCouplingMu,
- state->box,
- state->box_rel,
- start,
- homenr,
- state->x.rvec_array(),
- nullptr,
- md->cFREEZE,
- nrnb,
- scaleCoordinates);
+ pressureCouplingScaleBoxAndCoordinates<PressureCoupling::Berendsen>(
+ inputrec,
+ pressureCouplingMu,
+ state->box,
+ state->box_rel,
+ start,
+ homenr,
+ state->x,
+ gmx::ArrayRef<gmx::RVec>(),
+ cFREEZE,
+ nrnb,
+ scaleCoordinates);
}
break;
case (PressureCoupling::CRescale):
constraintVirial,
pressureCouplingMu,
&state->baros_integral);
- pressureCouplingScaleBoxAndCoordinates<PressureCoupling::CRescale>(
- inputrec,
- pressureCouplingMu,
- state->box,
- state->box_rel,
- start,
- homenr,
- state->x.rvec_array(),
- state->v.rvec_array(),
- md->cFREEZE,
- nrnb,
- scaleCoordinates);
+ pressureCouplingScaleBoxAndCoordinates<PressureCoupling::CRescale>(inputrec,
+ pressureCouplingMu,
+ state->box,
+ state->box_rel,
+ start,
+ homenr,
+ state->x,
+ state->v,
+ cFREEZE,
+ nrnb,
+ scaleCoordinates);
}
break;
case (PressureCoupling::ParrinelloRahman):
}
}
-extern gmx_bool update_randomize_velocities(const t_inputrec* ir,
- int64_t step,
- const t_commrec* cr,
- const t_mdatoms* md,
- gmx::ArrayRef<gmx::RVec> v,
- const gmx::Update* upd,
- const gmx::Constraints* constr)
+extern bool update_randomize_velocities(const t_inputrec* ir,
+ int64_t step,
+ const t_commrec* cr,
+ int homenr,
+ gmx::ArrayRef<const unsigned short> cTC,
+ gmx::ArrayRef<const real> invMass,
+ gmx::ArrayRef<gmx::RVec> v,
+ const gmx::Update* upd,
+ const gmx::Constraints* constr)
{
real rate = (ir->delta_t) / ir->opts.tau_t[0];
if ((ir->etc == TemperatureCoupling::Andersen) || do_per_step(step, gmx::roundToInt(1.0 / rate)))
{
andersen_tcoupl(
- ir, step, cr, md, v, rate, upd->getAndersenRandomizeGroup(), upd->getBoltzmanFactor());
+ ir, step, cr, homenr, cTC, invMass, v, rate, upd->getAndersenRandomizeGroup(), upd->getBoltzmanFactor());
return TRUE;
}
return FALSE;
}
-/*
- static const double sy_const[MAX_SUZUKI_YOSHIDA_NUM+1][MAX_SUZUKI_YOSHIDA_NUM+1] = {
- {},
- {1},
- {},
- {0.828981543588751,-0.657963087177502,0.828981543588751},
- {},
- {0.2967324292201065,0.2967324292201065,-0.186929716880426,0.2967324292201065,0.2967324292201065}
- };*/
-
/* these integration routines are only referenced inside this file */
static void NHC_trotter(const t_grpopts* opts,
int nvar,
double scalefac[],
real* veta,
const t_extmass* MassQ,
- gmx_bool bEkinAveVel)
+ bool bEkinAveVel)
{
/* general routine for both barostat and thermostat nose hoover chains */
- int i, j, mi, mj;
- double Ekin, Efac, reft, kT, nd;
- double dt;
- double * ivxi, *ixi;
- double* GQ;
- gmx_bool bBarostat;
- int mstepsi, mstepsj;
- int ns = SUZUKI_YOSHIDA_NUM; /* set the degree of integration in the types/state.h file */
- int nh = opts->nhchainlength;
+ int i, j, mi, mj;
+ double Ekin, Efac, reft, kT, nd;
+ double dt;
+ double *ivxi, *ixi;
+ double* GQ;
+ bool bBarostat;
+ int mstepsi, mstepsj;
+ int ns = SUZUKI_YOSHIDA_NUM; /* set the degree of integration in the types/state.h file */
+ int nh = opts->nhchainlength;
snew(GQ, nh);
mstepsi = mstepsj = ns;
tensor boxv,
tensor M,
matrix mu,
- gmx_bool bFirstStep)
+ bool bFirstStep)
{
/* This doesn't do any coordinate updating. It just
* integrates the box vector equations from the calculated
}
template<PressureCoupling pressureCouplingType>
-void pressureCouplingScaleBoxAndCoordinates(const t_inputrec* ir,
- const matrix mu,
- matrix box,
- matrix box_rel,
- int start,
- int nr_atoms,
- rvec x[],
- rvec v[],
- const unsigned short cFREEZE[],
- t_nrnb* nrnb,
- const bool scaleCoordinates)
+void pressureCouplingScaleBoxAndCoordinates(const t_inputrec* ir,
+ const matrix mu,
+ matrix box,
+ matrix box_rel,
+ int start,
+ int nr_atoms,
+ gmx::ArrayRef<gmx::RVec> x,
+ gmx::ArrayRef<gmx::RVec> v,
+ gmx::ArrayRef<const unsigned short> cFREEZE,
+ t_nrnb* nrnb,
+ const bool scaleCoordinates)
{
static_assert(pressureCouplingType == PressureCoupling::Berendsen
|| pressureCouplingType == PressureCoupling::CRescale,
{
// Trivial OpenMP region that does not throw
int g = 0;
- if (cFREEZE != nullptr)
+ if (!cFREEZE.empty())
{
g = cFREEZE[n];
}
}
}
-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 std::vector<bool>& randomize,
- gmx::ArrayRef<const real> boltzfac)
+void andersen_tcoupl(const t_inputrec* ir,
+ int64_t step,
+ const t_commrec* cr,
+ const int homenr,
+ gmx::ArrayRef<const unsigned short> cTC,
+ gmx::ArrayRef<const real> invMass,
+ gmx::ArrayRef<gmx::RVec> v,
+ real rate,
+ const std::vector<bool>& randomize,
+ gmx::ArrayRef<const real> boltzfac)
{
const int* gatindex = (DOMAINDECOMP(cr) ? cr->dd->globalAtomIndices.data() : nullptr);
int i;
/* randomize the velocities of the selected particles */
- for (i = 0; i < md->homenr; i++) /* now loop over the list of atoms */
+ for (i = 0; i < homenr; i++) /* now loop over the list of atoms */
{
- int ng = gatindex ? gatindex[i] : i;
- gmx_bool bRandomize;
+ int ng = gatindex ? gatindex[i] : i;
+ bool bRandomize;
rng.restart(step, ng);
- if (md->cTC)
+ if (!cTC.empty())
{
- gc = md->cTC[i]; /* assign the atom to a temperature group if there are more than one */
+ gc = cTC[i]; /* assign the atom to a temperature group if there are more than one */
}
if (randomize[gc])
{
real scal;
int d;
- scal = std::sqrt(boltzfac[gc] * md->invmass[i]);
+ scal = std::sqrt(boltzfac[gc] * invMass[i]);
normalDist.reset();
void nosehoover_tcoupl(const t_grpopts* opts,
const gmx_ekindata_t* ekind,
real dt,
- double xi[],
- double vxi[],
+ gmx::ArrayRef<double> xi,
+ gmx::ArrayRef<double> vxi,
const t_extmass* MassQ)
{
int i;
}
}
-void trotter_update(const t_inputrec* ir,
- int64_t step,
- gmx_ekindata_t* ekind,
- const gmx_enerdata_t* enerd,
- t_state* state,
- const tensor vir,
- const t_mdatoms* md,
- const t_extmass* MassQ,
- gmx::ArrayRef<std::vector<int>> trotter_seqlist,
- int trotter_seqno)
+void trotter_update(const t_inputrec* ir,
+ int64_t step,
+ gmx_ekindata_t* ekind,
+ const gmx_enerdata_t* enerd,
+ t_state* state,
+ const tensor vir,
+ int homenr,
+ gmx::ArrayRef<const unsigned short> cTC,
+ gmx::ArrayRef<const real> invMass,
+ const t_extmass* MassQ,
+ gmx::ArrayRef<std::vector<int>> trotter_seqlist,
+ int trotter_seqno)
{
int n, i, d, ngtc, gc = 0, t;
real dt;
double * scalefac, dtc;
rvec sumv = { 0, 0, 0 };
- gmx_bool bCouple;
+ bool bCouple;
if (trotter_seqno <= ettTSEQ2)
{
/* but do we actually need the total? */
/* modify the velocities as well */
- for (n = 0; n < md->homenr; n++)
+ for (n = 0; n < homenr; n++)
{
- if (md->cTC) /* does this conditional need to be here? is this always true?*/
+ if (!cTC.empty()) /* does this conditional need to be here? is this always true?*/
{
- gc = md->cTC[n];
+ gc = cTC[n];
}
for (d = 0; d < DIM; d++)
{
{
for (d = 0; d < DIM; d++)
{
- sumv[d] += (v[n][d]) / md->invmass[n];
+ sumv[d] += (v[n][d]) / invMass[n];
}
}
}
}
-extern void init_npt_masses(const t_inputrec* ir, t_state* state, t_extmass* MassQ, gmx_bool bInit)
+extern void init_npt_masses(const t_inputrec* ir, t_state* state, t_extmass* MassQ, bool bInit)
{
int n, i, j, d, ngtc, nh;
const t_grpopts* opts;
}
std::array<std::vector<int>, ettTSEQMAX>
-init_npt_vars(const t_inputrec* ir, t_state* state, t_extmass* MassQ, gmx_bool bTrotter)
+init_npt_vars(const t_inputrec* ir, t_state* state, t_extmass* MassQ, bool bTrotter)
{
int i, j, nnhpres, nh;
const t_grpopts* opts;
return ekin_new;
}
-void vrescale_tcoupl(const t_inputrec* ir, int64_t step, gmx_ekindata_t* ekind, real dt, double therm_integral[])
+void vrescale_tcoupl(const t_inputrec* ir, int64_t step, gmx_ekindata_t* ekind, real dt, gmx::ArrayRef<double> therm_integral)
{
const t_grpopts* opts;
int i;
}
}
-void rescale_velocities(const gmx_ekindata_t* ekind, const t_mdatoms* mdatoms, int start, int end, rvec v[])
+void rescale_velocities(const gmx_ekindata_t* ekind,
+ gmx::ArrayRef<const unsigned short> cTC,
+ int start,
+ int end,
+ gmx::ArrayRef<gmx::RVec> v)
{
- const unsigned short* cTC = mdatoms->cTC;
gmx::ArrayRef<const t_grp_tcstat> tcstat = ekind->tcstat;
for (int n = start; n < end; n++)
{
int gt = 0;
- if (cTC)
+ if (!cTC.empty())
{
gt = cTC[n];
}
#ifndef GMX_MDLIB_COUPLING_H
#define GMX_MDLIB_COUPLING_H
-#include "gromacs/math/paddedvector.h"
+#include <cstdio>
+
+#include <array>
+#include <vector>
+
#include "gromacs/math/vectypes.h"
#include "gromacs/mdtypes/md_enums.h"
-#include "gromacs/utility/arrayref.h"
-#include "gromacs/utility/basedefinitions.h"
#include "gromacs/utility/real.h"
class gmx_ekindata_t;
struct t_extmass;
struct t_grpopts;
struct t_inputrec;
-struct t_mdatoms;
struct t_nrnb;
class t_state;
class BoxDeformation;
class Constraints;
class Update;
+template<typename>
+class ArrayRef;
}; // namespace gmx
/* Update the size of per-atom arrays (e.g. after DD re-partitioning,
which might increase the number of home atoms). */
-void update_tcouple(int64_t step,
- const t_inputrec* inputrec,
- t_state* state,
- gmx_ekindata_t* ekind,
- const t_extmass* MassQ,
- const t_mdatoms* md);
+void update_tcouple(int64_t step,
+ const t_inputrec* inputrec,
+ t_state* state,
+ gmx_ekindata_t* ekind,
+ const t_extmass* MassQ,
+ int homenr,
+ gmx::ArrayRef<const unsigned short> cTC);
/* Update Parrinello-Rahman, to be called before the coordinate update */
void update_pcouple_before_coordinates(FILE* fplog,
t_state* state,
matrix parrinellorahmanMu,
matrix M,
- gmx_bool bInitStep);
+ bool bInitStep);
/* Update the box, to be called after the coordinate update.
* For Berendsen P-coupling, also calculates the scaling factor
* and scales the coordinates.
* When the deform option is used, scales coordinates and box here.
*/
-void update_pcouple_after_coordinates(FILE* fplog,
- int64_t step,
- const t_inputrec* inputrec,
- const t_mdatoms* md,
- const matrix pressure,
- const matrix forceVirial,
- const matrix constraintVirial,
- matrix pressureCouplingMu,
- t_state* state,
- t_nrnb* nrnb,
- gmx::BoxDeformation* boxDeformation,
- bool scaleCoordinates);
+void update_pcouple_after_coordinates(FILE* fplog,
+ int64_t step,
+ const t_inputrec* inputrec,
+ int homenr,
+ gmx::ArrayRef<const unsigned short> cFREEZE,
+ const matrix pressure,
+ const matrix forceVirial,
+ const matrix constraintVirial,
+ matrix pressureCouplingMu,
+ t_state* state,
+ t_nrnb* nrnb,
+ gmx::BoxDeformation* boxDeformation,
+ bool scaleCoordinates);
/* Return TRUE if OK, FALSE in case of Shake Error */
-extern gmx_bool update_randomize_velocities(const t_inputrec* ir,
- int64_t step,
- const t_commrec* cr,
- const t_mdatoms* md,
- gmx::ArrayRef<gmx::RVec> v,
- const gmx::Update* upd,
- const gmx::Constraints* constr);
+extern bool update_randomize_velocities(const t_inputrec* ir,
+ int64_t step,
+ const t_commrec* cr,
+ int homenr,
+ gmx::ArrayRef<const unsigned short> cTC,
+ gmx::ArrayRef<const real> invMass,
+ gmx::ArrayRef<gmx::RVec> v,
+ const gmx::Update* upd,
+ const gmx::Constraints* constr);
void berendsen_tcoupl(const t_inputrec* ir,
gmx_ekindata_t* ekind,
real dt,
std::vector<double>& therm_integral); //NOLINT(google-runtime-references)
-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 std::vector<bool>& randomize,
- gmx::ArrayRef<const real> boltzfac);
+void andersen_tcoupl(const t_inputrec* ir,
+ int64_t step,
+ const t_commrec* cr,
+ int homenr,
+ gmx::ArrayRef<const unsigned short> cTC,
+ gmx::ArrayRef<const real> invMass,
+ gmx::ArrayRef<gmx::RVec> v,
+ 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[],
+ gmx::ArrayRef<double> xi,
+ gmx::ArrayRef<double> vxi,
const t_extmass* MassQ);
-void trotter_update(const t_inputrec* ir,
- int64_t step,
- gmx_ekindata_t* ekind,
- const gmx_enerdata_t* enerd,
- t_state* state,
- const tensor vir,
- const t_mdatoms* md,
- const t_extmass* MassQ,
- gmx::ArrayRef<std::vector<int>> trotter_seqlist,
- int trotter_seqno);
+void trotter_update(const t_inputrec* ir,
+ int64_t step,
+ gmx_ekindata_t* ekind,
+ const gmx_enerdata_t* enerd,
+ t_state* state,
+ const tensor vir,
+ int homenr,
+ gmx::ArrayRef<const unsigned short> cTC,
+ gmx::ArrayRef<const real> invMass,
+ const t_extmass* MassQ,
+ gmx::ArrayRef<std::vector<int>> trotter_seqlist,
+ int trotter_seqno);
std::array<std::vector<int>, ettTSEQMAX>
-init_npt_vars(const t_inputrec* ir, t_state* state, t_extmass* Mass, gmx_bool bTrotter);
+init_npt_vars(const t_inputrec* ir, t_state* state, t_extmass* Mass, bool bTrotter);
real NPT_energy(const t_inputrec* ir, const t_state* state, const t_extmass* MassQ);
/* computes all the pressure/tempertature control energy terms to get a conserved energy */
-void vrescale_tcoupl(const t_inputrec* ir, int64_t step, gmx_ekindata_t* ekind, real dt, double therm_integral[]);
+void vrescale_tcoupl(const t_inputrec* ir,
+ int64_t step,
+ gmx_ekindata_t* ekind,
+ real dt,
+ gmx::ArrayRef<double> therm_integral);
/* Compute temperature scaling. For V-rescale it is done in update. */
-void rescale_velocities(const gmx_ekindata_t* ekind, const t_mdatoms* mdatoms, int start, int end, rvec v[]);
+void rescale_velocities(const gmx_ekindata_t* ekind,
+ gmx::ArrayRef<const unsigned short> cTC,
+ int start,
+ int end,
+ gmx::ArrayRef<gmx::RVec> v);
/* Rescale the velocities with the scaling factor in ekind */
//! Check whether we do simulated annealing.
tensor boxv,
tensor M,
matrix mu,
- gmx_bool bFirstStep);
+ bool bFirstStep);
/*! \brief Calculate the pressure coupling scaling matrix
*
* coupling algorithm.
*/
template<PressureCoupling pressureCouplingType>
-void pressureCouplingScaleBoxAndCoordinates(const t_inputrec* ir,
- const matrix mu,
- matrix box,
- matrix box_rel,
- int start,
- int nr_atoms,
- rvec x[],
- rvec v[],
- const unsigned short cFREEZE[],
- t_nrnb* nrnb,
- bool scaleCoordinates);
+void pressureCouplingScaleBoxAndCoordinates(const t_inputrec* ir,
+ const matrix mu,
+ matrix box,
+ matrix box_rel,
+ int start,
+ int nr_atoms,
+ gmx::ArrayRef<gmx::RVec> x,
+ gmx::ArrayRef<gmx::RVec> v,
+ gmx::ArrayRef<const unsigned short> cFREEZE,
+ t_nrnb* nrnb,
+ bool scaleCoordinates);
void pleaseCiteCouplingAlgorithms(FILE* fplog, const t_inputrec& ir);
#include "gromacs/mdtypes/inputrec.h"
-real sum_ekin(const t_grpopts* opts, gmx_ekindata_t* ekind, real* dekindlambda, gmx_bool bEkinAveVel, gmx_bool bScaleEkin)
+real sum_ekin(const t_grpopts* opts, gmx_ekindata_t* ekind, real* dekindlambda, bool bEkinAveVel, bool bScaleEkin)
{
int i, j, m, ngtc;
real T;
else
{
/* this is for NHC in the Ekin(t+dt/2) version of vv */
- trotter_update(ir, step, ekind, enerd, state, total_vir, mdatoms, MassQ, trotter_seq, ettTSEQ1);
+ trotter_update(ir,
+ step,
+ ekind,
+ enerd,
+ state,
+ total_vir,
+ mdatoms->homenr,
+ mdatoms->cTC ? gmx::arrayRefFromArray(mdatoms->cTC, mdatoms->nr)
+ : gmx::ArrayRef<const unsigned short>(),
+ gmx::arrayRefFromArray(mdatoms->invmass, mdatoms->nr),
+ MassQ,
+ trotter_seq,
+ ettTSEQ1);
}
upd->update_coords(*ir,
if (bTrotter)
{
m_add(force_vir, shake_vir, total_vir); /* we need the un-dispersion corrected total vir here */
- trotter_update(ir, step, ekind, enerd, state, total_vir, mdatoms, MassQ, trotter_seq, ettTSEQ2);
+ trotter_update(ir,
+ step,
+ ekind,
+ enerd,
+ state,
+ total_vir,
+ mdatoms->homenr,
+ mdatoms->cTC ? gmx::arrayRefFromArray(mdatoms->cTC, mdatoms->nr)
+ : gmx::ArrayRef<const unsigned short>(),
+ gmx::arrayRefFromArray(mdatoms->invmass, mdatoms->nr),
+ MassQ,
+ trotter_seq,
+ ettTSEQ2);
/* TODO This is only needed when we're about to write
* a checkpoint, because we use it after the restart
bSumEkinhOld,
(bGStat ? CGLO_GSTAT : 0) | CGLO_TEMPERATURE);
wallcycle_start(wcycle, WallCycleCounter::Update);
- trotter_update(ir, step, ekind, enerd, state, total_vir, mdatoms, MassQ, trotter_seq, ettTSEQ4);
+ trotter_update(ir,
+ step,
+ ekind,
+ enerd,
+ state,
+ total_vir,
+ mdatoms->homenr,
+ mdatoms->cTC ? gmx::arrayRefFromArray(mdatoms->cTC, mdatoms->nr)
+ : gmx::ArrayRef<const unsigned short>(),
+ gmx::arrayRefFromArray(mdatoms->invmass, mdatoms->nr),
+ MassQ,
+ trotter_seq,
+ ettTSEQ4);
/* now we know the scaling, we can compute the positions again */
std::copy(cbuf->begin(), cbuf->end(), state->x.begin());
if (ETC_ANDERSEN(ir->etc)) /* keep this outside of update_tcouple because of the extra info required to pass */
{
gmx_bool bIfRandomize;
- bIfRandomize = update_randomize_velocities(ir, step, cr, md, state->v, &upd, constr);
+ bIfRandomize = update_randomize_velocities(ir,
+ step,
+ cr,
+ md->homenr,
+ md->cTC ? gmx::arrayRefFromArray(md->cTC, md->nr)
+ : gmx::ArrayRef<const unsigned short>(),
+ gmx::arrayRefFromArray(md->invmass, md->nr),
+ state->v,
+ &upd,
+ constr);
/* if we have constraints, we have to remove the kinetic energy parallel to the bonds */
if (constr && bIfRandomize)
{
/* UPDATE PRESSURE VARIABLES IN TROTTER FORMULATION WITH CONSTRAINTS */
if (bTrotter)
{
- trotter_update(ir, step, ekind, enerd, state, total_vir, md, &MassQ, trotter_seq, ettTSEQ3);
+ trotter_update(ir,
+ step,
+ ekind,
+ enerd,
+ state,
+ total_vir,
+ md->homenr,
+ md->cTC ? gmx::arrayRefFromArray(md->cTC, md->nr)
+ : gmx::ArrayRef<const unsigned short>(),
+ gmx::arrayRefFromArray(md->invmass, md->nr),
+ &MassQ,
+ trotter_seq,
+ ettTSEQ3);
/* We can only do Berendsen coupling after we have summed
* the kinetic energy or virial. Since the happens
* in global_state after update, we should only do it at
}
else
{
- update_tcouple(step, ir, state, ekind, &MassQ, md);
+ update_tcouple(step,
+ ir,
+ state,
+ ekind,
+ &MassQ,
+ md->homenr,
+ md->cTC ? gmx::arrayRefFromArray(md->cTC, md->nr)
+ : gmx::ArrayRef<const unsigned short>());
update_pcouple_before_coordinates(fplog, step, ir, state, pressureCouplingMu, M, bInitStep);
}
accumulateKineticLambdaComponents(enerd, state->lambda, *ir->fepvals);
}
- update_pcouple_after_coordinates(
- fplog, step, ir, md, pres, force_vir, shake_vir, pressureCouplingMu, state, nrnb, upd.deform(), !useGpuForUpdate);
+ update_pcouple_after_coordinates(fplog,
+ step,
+ ir,
+ md->homenr,
+ md->cFREEZE ? gmx::arrayRefFromArray(md->cFREEZE, md->nr)
+ : gmx::ArrayRef<const unsigned short>(),
+ pres,
+ force_vir,
+ shake_vir,
+ pressureCouplingMu,
+ state,
+ nrnb,
+ upd.deform(),
+ !useGpuForUpdate);
const bool doBerendsenPressureCoupling = (inputrec->epc == PressureCoupling::Berendsen
&& do_per_step(step, inputrec->nstpcouple));
biod.pnpi.spb.ru / alexxy / gromacs.git / commitdiff