}
real Awh::applyBiasForcesAndUpdateBias(PbcType pbcType,
- const real* masses,
+ ArrayRef<const real> masses,
ArrayRef<const double> neighborLambdaEnergies,
ArrayRef<const double> neighborLambdaDhdl,
const matrix box,
* \returns the potential energy for the bias.
*/
real applyBiasForcesAndUpdateBias(PbcType pbcType,
- const real* masses,
+ ArrayRef<const real> masses,
ArrayRef<const double> neighborLambdaEnergies,
ArrayRef<const double> neighborLambdaDhdl,
const matrix box,
#include "gromacs/pbcutil/pbc.h"
#include "gromacs/topology/mtop_util.h"
#include "gromacs/topology/topology.h"
+#include "gromacs/utility/arrayref.h"
#include "gromacs/utility/gmxassert.h"
#include "gromacs/utility/smalloc.h"
numHomeAtoms,
mdAtoms);
- auto mdatoms = mdAtoms->mdatoms();
+ t_mdatoms* mdatoms = mdAtoms->mdatoms();
if (usingDomDec)
{
dd_sort_local_top(*cr->dd, mdatoms, top);
constr->setConstraints(top,
mdatoms->nr,
mdatoms->homenr,
- mdatoms->massT,
- mdatoms->invmass,
+ gmx::arrayRefFromArray(mdatoms->massT, mdatoms->nr),
+ gmx::arrayRefFromArray(mdatoms->invmass, mdatoms->nr),
mdatoms->nMassPerturbed != 0,
mdatoms->lambda,
- mdatoms->cFREEZE);
+ mdatoms->cFREEZE ? gmx::arrayRefFromArray(mdatoms->cFREEZE, mdatoms->nr)
+ : gmx::ArrayRef<const unsigned short>());
}
}
if (ir->bPull)
{
- pull = set_pull_init(ir,
- sys,
- state.x.rvec_array(),
- state.box,
- state.lambda[FreeEnergyPerturbationCouplingType::Mass],
- wi);
+ pull = set_pull_init(
+ ir, sys, state.x, state.box, state.lambda[FreeEnergyPerturbationCouplingType::Mass], wi);
}
/* Modules that supply external potential for pull coordinates
gmx::ArrayRef<const t_pull_coord> pullCoords);
/* Process the pull coordinates after reading the pull groups */
-pull_t* set_pull_init(t_inputrec* ir, const gmx_mtop_t& mtop, rvec* x, matrix box, real lambda, warninp_t wi);
+pull_t* set_pull_init(t_inputrec* ir,
+ const gmx_mtop_t& mtop,
+ gmx::ArrayRef<const gmx::RVec> x,
+ matrix box,
+ real lambda,
+ warninp_t wi);
/* Prints the initial pull group distances in x.
* If requested, adds the current distance to the initial reference location.
* Returns the pull_t pull work struct. This should be passed to finish_pull()
}
}
-pull_t* set_pull_init(t_inputrec* ir, const gmx_mtop_t& mtop, rvec* x, matrix box, real lambda, warninp_t wi)
+pull_t* set_pull_init(t_inputrec* ir,
+ const gmx_mtop_t& mtop,
+ gmx::ArrayRef<const gmx::RVec> x,
+ matrix box,
+ real lambda,
+ warninp_t wi)
{
pull_t* pull_work;
t_pbc pbc;
if (pull->bSetPbcRefToPrevStepCOM)
{
- initPullComFromPrevStep(nullptr, pull_work, md->massT, &pbc, x);
+ initPullComFromPrevStep(nullptr, pull_work, gmx::arrayRefFromArray(md->massT, md->nr), &pbc, x);
}
- pull_calc_coms(nullptr, pull_work, md->massT, &pbc, t_start, x, nullptr);
+ pull_calc_coms(nullptr, pull_work, gmx::arrayRefFromArray(md->massT, md->nr), &pbc, t_start, x, {});
for (int g = 0; g < pull->ngroup; g++)
{
- bool groupObeysPbc = pullCheckPbcWithinGroup(
- *pull_work,
- gmx::arrayRefFromArray(reinterpret_cast<gmx::RVec*>(x), mtop.natoms),
- pbc,
- g,
- c_pullGroupSmallGroupThreshold);
+ bool groupObeysPbc =
+ pullCheckPbcWithinGroup(*pull_work, x, pbc, g, c_pullGroupSmallGroupThreshold);
if (!groupObeysPbc)
{
char buf[STRLEN];
}
if (groupObeysPbc)
{
- groupObeysPbc = pullCheckPbcWithinGroup(
- *pull_work,
- gmx::arrayRefFromArray(reinterpret_cast<gmx::RVec*>(x), mtop.natoms),
- pbc,
- g,
- c_pullGroupPbcMargin);
+ groupObeysPbc = pullCheckPbcWithinGroup(*pull_work, x, pbc, g, c_pullGroupPbcMargin);
if (!groupObeysPbc)
{
char buf[STRLEN];
#include "gromacs/fileio/gmxfio.h"
#include "gromacs/fileio/pdbio.h"
#include "gromacs/gmxlib/nrnb.h"
+#include "gromacs/math/arrayrefwithpadding.h"
#include "gromacs/math/utilities.h"
#include "gromacs/math/vec.h"
#include "gromacs/mdlib/gmx_omp_nthreads.h"
#include "gromacs/topology/ifunc.h"
#include "gromacs/topology/mtop_lookup.h"
#include "gromacs/topology/mtop_util.h"
+#include "gromacs/utility/arrayref.h"
#include "gromacs/utility/exceptions.h"
#include "gromacs/utility/fatalerror.h"
#include "gromacs/utility/gmxassert.h"
int numConstraints,
int numSettles);
~Impl();
- void setConstraints(gmx_localtop_t* top,
- int numAtoms,
- int numHomeAtoms,
- real* masses,
- real* inverseMasses,
- bool hasMassPerturbedAtoms,
- real lambda,
- unsigned short* cFREEZE);
+ void setConstraints(gmx_localtop_t* top,
+ int numAtoms,
+ int numHomeAtoms,
+ gmx::ArrayRef<const real> masses,
+ gmx::ArrayRef<const real> inverseMasses,
+ bool hasMassPerturbedAtoms,
+ real lambda,
+ gmx::ArrayRef<const unsigned short> cFREEZE);
bool apply(bool bLog,
bool bEner,
int64_t step,
//! Number of local atoms.
int numHomeAtoms_ = 0;
//! Atoms masses.
- real* masses_;
+ gmx::ArrayRef<const real> masses_;
//! Inverse masses.
- real* inverseMasses_;
+ gmx::ArrayRef<const real> inverseMasses_;
//! If there are atoms with perturbed mass (for FEP).
bool hasMassPerturbedAtoms_;
//! FEP lambda value.
real lambda_;
//! Freeze groups data
- unsigned short* cFREEZE_;
+ gmx::ArrayRef<const unsigned short> cFREEZE_;
//! Whether we need to do pbc for handling bonds.
bool pbcHandlingRequired_ = false;
cr,
ir.delta_t,
t,
- as_rvec_array(x.unpaddedArrayRef().data()),
- as_rvec_array(xprime.unpaddedArrayRef().data()),
- as_rvec_array(v.unpaddedArrayRef().data()),
+ x.unpaddedArrayRef(),
+ xprime.unpaddedArrayRef(),
+ v.unpaddedArrayRef(),
constraintsVirial);
}
if (ed && delta_step > 0)
wallcycle_stop(wcycle, ewcCONSTR);
const bool haveVelocities = (!v.empty() || econq == ConstraintVariable::Velocities);
- if (haveVelocities && cFREEZE_)
+ if (haveVelocities && !cFREEZE_.empty())
{
/* Set the velocities of frozen dimensions to zero */
ArrayRef<RVec> vRef;
return at2s;
}
-void Constraints::Impl::setConstraints(gmx_localtop_t* top,
- int numAtoms,
- int numHomeAtoms,
- real* masses,
- real* inverseMasses,
- const bool hasMassPerturbedAtoms,
- const real lambda,
- unsigned short* cFREEZE)
+void Constraints::Impl::setConstraints(gmx_localtop_t* top,
+ int numAtoms,
+ int numHomeAtoms,
+ gmx::ArrayRef<const real> masses,
+ gmx::ArrayRef<const real> inverseMasses,
+ const bool hasMassPerturbedAtoms,
+ const real lambda,
+ gmx::ArrayRef<const unsigned short> cFREEZE)
{
numAtoms_ = numAtoms;
numHomeAtoms_ = numHomeAtoms;
}
}
-void Constraints::setConstraints(gmx_localtop_t* top,
- const int numAtoms,
- const int numHomeAtoms,
- real* masses,
- real* inverseMasses,
- const bool hasMassPerturbedAtoms,
- const real lambda,
- unsigned short* cFREEZE)
+void Constraints::setConstraints(gmx_localtop_t* top,
+ const int numAtoms,
+ const int numHomeAtoms,
+ gmx::ArrayRef<const real> masses,
+ gmx::ArrayRef<const real> inverseMasses,
+ const bool hasMassPerturbedAtoms,
+ const real lambda,
+ gmx::ArrayRef<const unsigned short> cFREEZE)
{
impl_->setConstraints(
top, numAtoms, numHomeAtoms, masses, inverseMasses, hasMassPerturbedAtoms, lambda, cFREEZE);
*
* \todo Make this a callback that is called automatically
* once a new domain has been made. */
- void setConstraints(gmx_localtop_t* top,
- int numAtoms,
- int numHomeAtoms,
- real* masses,
- real* inverseMasses,
- bool hasMassPerturbedAtoms,
- real lambda,
- unsigned short* cFREEZE);
+ void setConstraints(gmx_localtop_t* top,
+ int numAtoms,
+ int numHomeAtoms,
+ gmx::ArrayRef<const real> masses,
+ gmx::ArrayRef<const real> inverseMasses,
+ bool hasMassPerturbedAtoms,
+ real lambda,
+ gmx::ArrayRef<const unsigned short> cFREEZE);
/*! \brief Applies constraints to coordinates.
*
real preFactor,
gmx::ArrayRef<const real> fac,
gmx::ArrayRef<const gmx::RVec> r,
- const real* invmass,
+ gmx::ArrayRef<const real> invmass,
rvec* x)
{
- if (invmass != nullptr)
+ if (!invmass.empty())
{
for (int b = 0; b < ncons; b++)
{
real preFactor,
gmx::ArrayRef<const real> fac,
gmx::ArrayRef<const gmx::RVec> r,
- const real* invmass,
+ gmx::ArrayRef<const real> invmass,
rvec* x)
{
- if (invmass != nullptr)
+ if (!invmass.empty())
{
for (int b : ind)
{
real preFactor,
gmx::ArrayRef<const real> fac,
gmx::ArrayRef<const gmx::RVec> r,
- const real* invmass,
+ gmx::ArrayRef<const real> invmass,
rvec* x)
{
if (li->ntask == 1)
t_pbc* pbc,
Lincs* lincsd,
int th,
- const real* invmass,
+ ArrayRef<const real> invmass,
ConstraintVariable econq,
bool bCalcDHDL,
bool bCalcVir,
1.0,
sol,
r,
- (econq != ConstraintVariable::Force) ? invmass : nullptr,
+ (econq != ConstraintVariable::Force) ? invmass : gmx::ArrayRef<real>(),
as_rvec_array(fp.data()));
if (bCalcDHDL)
t_pbc* pbc,
Lincs* lincsd,
int th,
- const real* invmass,
+ ArrayRef<const real> invmass,
const t_commrec* cr,
bool bCalcDHDL,
real wangle,
}
/*! \brief Sets the elements in the LINCS matrix for task task. */
-static void set_lincs_matrix_task(Lincs* li, Task* li_task, const real* invmass, int* ncc_triangle, int* nCrossTaskTriangles)
+static void set_lincs_matrix_task(Lincs* li,
+ Task* li_task,
+ ArrayRef<const real> invmass,
+ int* ncc_triangle,
+ int* nCrossTaskTriangles)
{
/* Construct the coupling coefficient matrix blmf */
li_task->ntriangle = 0;
}
/*! \brief Sets the elements in the LINCS matrix. */
-static void set_lincs_matrix(Lincs* li, const real* invmass, real lambda)
+static void set_lincs_matrix(Lincs* li, ArrayRef<const real> invmass, real lambda)
{
const real invsqrt2 = 0.7071067811865475244;
void set_lincs(const InteractionDefinitions& idef,
const int numAtoms,
- const real* invmass,
+ ArrayRef<const real> invmass,
const real lambda,
bool bDynamics,
const t_commrec* cr,
const t_inputrec& ir,
int64_t step,
Lincs* lincsd,
- const real* invmass,
+ ArrayRef<const real> invmass,
const t_commrec* cr,
const gmx_multisim_t* ms,
ArrayRefWithPadding<const RVec> xPadded,
/*
* This file is part of the GROMACS molecular simulation package.
*
- * 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.
namespace gmx
{
-
template<typename>
class ArrayRefWithPadding;
+template<typename>
+class ArrayRef;
enum class ConstraintVariable : int;
class Lincs;
template<typename>
/*! \brief Initialize lincs stuff */
void set_lincs(const InteractionDefinitions& idef,
int numAtoms,
- const real* invmass,
+ ArrayRef<const real> invmass,
real lambda,
bool bDynamics,
const t_commrec* cr,
const t_inputrec& ir,
int64_t step,
Lincs* lincsd,
- const real* invmass,
+ ArrayRef<const real> invmass,
const t_commrec* cr,
const gmx_multisim_t* ms,
ArrayRefWithPadding<const RVec> x,
* 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.
#include "gromacs/math/invertmatrix.h"
#include "gromacs/math/vec.h"
#include "gromacs/mdlib/constr.h"
-#include "gromacs/mdtypes/mdatom.h"
#include "gromacs/pbcutil/ishift.h"
#include "gromacs/pbcutil/pbc.h"
#include "gromacs/pbcutil/pbc_simd.h"
#include "gromacs/topology/idef.h"
#include "gromacs/topology/ifunc.h"
#include "gromacs/topology/mtop_util.h"
+#include "gromacs/utility/arrayref.h"
#include "gromacs/utility/fatalerror.h"
#include "gromacs/utility/gmxassert.h"
parametersAllMasses1_ = settleParameters(1.0, 1.0, 1.0, 1.0, dOH, dHH);
}
-void SettleData::setConstraints(const InteractionList& il_settle,
- const int numHomeAtoms,
- const real* masses,
- const real* inverseMasses)
+void SettleData::setConstraints(const InteractionList& il_settle,
+ const int numHomeAtoms,
+ gmx::ArrayRef<const real> masses,
+ gmx::ArrayRef<const real> inverseMasses)
{
#if GMX_SIMD_HAVE_REAL
const int pack_size = GMX_SIMD_REAL_WIDTH;
/*
* This file is part of the GROMACS molecular simulation package.
*
- * 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.
SettleData(const gmx_mtop_t& mtop);
//! Sets the constraints from the interaction list and the masses
- void setConstraints(const InteractionList& il_settle,
- int numHomeAtoms,
- const real* masses,
- const real* inverseMasses);
+ void setConstraints(const InteractionList& il_settle,
+ int numHomeAtoms,
+ gmx::ArrayRef<const real> masses,
+ gmx::ArrayRef<const real> inverseMasses);
//! Returns settle parameters for constraining coordinates and forces
const SettleParameters& parametersMassWeighted() const { return parametersMassWeighted_; }
#include "gromacs/mdtypes/md_enums.h"
#include "gromacs/pbcutil/pbc.h"
#include "gromacs/topology/invblock.h"
+#include "gromacs/utility/arrayref.h"
#include "gromacs/utility/fatalerror.h"
#include "gromacs/utility/smalloc.h"
ArrayRef<const RVec> initial_displacements,
ArrayRef<const real> half_of_reduced_mass,
real omega,
- const real invmass[],
+ ArrayRef<const real> invmass,
ArrayRef<const real> distance_squared_tolerance,
ArrayRef<real> scaled_lagrange_multiplier,
int* nerror)
ArrayRef<const RVec> rij,
ArrayRef<const real> m2,
real omega,
- const real invmass[],
+ ArrayRef<const real> invmass,
ArrayRef<const real> distance_squared_tolerance,
ArrayRef<real> scaled_lagrange_multiplier,
int* nerror,
//! Applies SHAKE
static int vec_shakef(FILE* fplog,
shakedata* shaked,
- const real invmass[],
+ ArrayRef<const real> invmass,
int ncon,
ArrayRef<const t_iparams> ip,
const int* iatom,
const t_pbc* pbc,
ArrayRef<const t_iparams> ip,
const int* iatom,
- const real invmass[],
+ ArrayRef<const real> invmass,
ConstraintVariable econq)
{
int ai, aj;
//! Applies SHAKE.
static bool bshakef(FILE* log,
shakedata* shaked,
- const real invmass[],
+ ArrayRef<const real> invmass,
const InteractionDefinitions& idef,
const t_inputrec& ir,
ArrayRef<const RVec> x_s,
bool constrain_shake(FILE* log,
shakedata* shaked,
- const real invmass[],
+ ArrayRef<const real> invmass,
const InteractionDefinitions& idef,
const t_inputrec& ir,
ArrayRef<const RVec> x_s,
/*
* This file is part of the GROMACS molecular simulation package.
*
- * 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.
* sblock[n] to sblock[n+1]. Array sblock should be large enough.
* Return TRUE when OK, FALSE when shake-error
*/
-bool constrain_shake(FILE* log, /* Log file */
- shakedata* shaked, /* Total number of atoms */
- const real invmass[], /* Atomic masses */
- const InteractionDefinitions& idef, /* The interaction def */
- const t_inputrec& ir, /* Input record */
- ArrayRef<const RVec> x_s, /* Coords before update */
- ArrayRef<RVec> xprime, /* Output coords when constraining x */
- ArrayRef<RVec> vprime, /* Output coords when constraining v */
- const t_pbc* pbc, /* PBC information */
- t_nrnb* nrnb, /* Performance measure */
- real lambda, /* FEP lambda */
+bool constrain_shake(FILE* log, /* Log file */
+ shakedata* shaked, /* Total number of atoms */
+ gmx::ArrayRef<const real> invmass, /* Atomic masses */
+ const InteractionDefinitions& idef, /* The interaction def */
+ const t_inputrec& ir, /* Input record */
+ ArrayRef<const RVec> x_s, /* Coords before update */
+ ArrayRef<RVec> xprime, /* Output coords when constraining x */
+ ArrayRef<RVec> vprime, /* Output coords when constraining v */
+ const t_pbc* pbc, /* PBC information */
+ t_nrnb* nrnb, /* Performance measure */
+ real lambda, /* FEP lambda */
real* dvdlambda, /* FEP force */
real invdt, /* 1/delta_t */
ArrayRef<RVec> v, /* Also constrain v if not empty */
ArrayRef<const RVec> rij,
ArrayRef<const real> half_of_reduced_mass,
real omega,
- const real invmass[],
+ ArrayRef<const real> invmass,
ArrayRef<const real> distance_squared_tolerance,
ArrayRef<real> scaled_lagrange_multiplier,
int* nerror);
dvdl = 0;
enerd->term[F_COM_PULL] +=
pull_potential(pull_work,
- mdatoms->massT,
+ gmx::arrayRefFromArray(mdatoms->massT, mdatoms->nr),
&pbc,
cr,
t,
lambda[static_cast<int>(FreeEnergyPerturbationCouplingType::Restraint)],
- as_rvec_array(x.data()),
+ x,
force,
&dvdl);
enerd->dvdl_lin[FreeEnergyPerturbationCouplingType::Restraint] += dvdl;
}
auto& forceWithVirial = (mtsLevel == 0) ? forceWithVirialMtsLevel0 : forceWithVirialMtsLevel1;
- enerd->term[F_COM_PULL] += awh->applyBiasForcesAndUpdateBias(inputrec.pbcType,
- mdatoms->massT,
- foreignLambdaDeltaH,
- foreignLambdaDhDl,
- box,
- forceWithVirial,
- t,
- step,
- wcycle,
- fplog);
+ enerd->term[F_COM_PULL] += awh->applyBiasForcesAndUpdateBias(
+ inputrec.pbcType,
+ gmx::arrayRefFromArray(mdatoms->massT, mdatoms->nr),
+ foreignLambdaDeltaH,
+ foreignLambdaDhDl,
+ box,
+ forceWithVirial,
+ t,
+ step,
+ wcycle,
+ fplog);
}
}
/*
* This file is part of the GROMACS molecular simulation package.
*
- * 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.
make_shake_sblock_serial(&shaked, testData->idef_.get(), testData->numAtoms_);
bool success = constrain_shake(nullptr,
&shaked,
- testData->invmass_.data(),
+ testData->invmass_,
*testData->idef_,
testData->ir_,
testData->x_,
testData->ir_.nProjOrder);
set_lincs(*testData->idef_,
testData->numAtoms_,
- testData->invmass_.data(),
+ testData->invmass_,
testData->lambda_,
EI_DYNAMICS(testData->ir_.eI),
&cr,
testData->ir_,
0,
lincsd,
- testData->invmass_.data(),
+ testData->invmass_,
&cr,
&ms,
testData->x_.arrayRefWithPadding(),
/*
* This file is part of the GROMACS molecular simulation package.
*
- * 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.
{
SettleData settled(testData->mtop_);
- settled.setConstraints(testData->idef_->il[F_SETTLE],
- testData->numAtoms_,
- testData->masses_.data(),
- testData->inverseMasses_.data());
+ settled.setConstraints(
+ testData->idef_->il[F_SETTLE], testData->numAtoms_, testData->masses_, testData->inverseMasses_);
bool errorOccured;
int numThreads = 1;
/*
* This file is part of the GROMACS molecular simulation package.
*
- * Copyright (c) 2014,2015,2017,2018,2019,2020, by the GROMACS development team, led by
+ * Copyright (c) 2014,2015,2017,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.
initialDisplacements,
halfOfReducedMasses,
omega_,
- inverseMasses.data(),
+ inverseMasses,
distanceSquaredTolerances,
lagrangianValues,
&numErrors);
EnergyData::initializeEnergyHistory(startingBehavior, observablesHistory, &energyOutput);
}
- preparePrevStepPullCom(
- ir, pull_work, mdatoms->massT, state, state_global, cr, startingBehavior != StartingBehavior::NewSimulation);
+ preparePrevStepPullCom(ir,
+ pull_work,
+ gmx::arrayRefFromArray(mdatoms->massT, mdatoms->nr),
+ state,
+ state_global,
+ cr,
+ startingBehavior != StartingBehavior::NewSimulation);
// TODO: Remove this by converting AWH into a ForceProvider
auto awh = prepareAwhModule(fplog,
*/
#include "gmxpre.h"
-#include "gromacs/utility/stringutil.h"
#include "pull.h"
#include "config.h"
#include "gromacs/pbcutil/pbc.h"
#include "gromacs/topology/mtop_lookup.h"
#include "gromacs/topology/topology.h"
+#include "gromacs/utility/arrayref.h"
#include "gromacs/utility/cstringutil.h"
#include "gromacs/utility/exceptions.h"
#include "gromacs/utility/fatalerror.h"
#include "gromacs/utility/real.h"
#include "gromacs/utility/smalloc.h"
#include "gromacs/utility/strconvert.h"
+#include "gromacs/utility/stringutil.h"
#include "pull_internal.h"
}
/* Apply forces in a mass weighted fashion for part of the pull group */
-static void apply_forces_grp_part(const pull_group_work_t* pgrp,
- int ind_start,
- int ind_end,
- const real* masses,
- const dvec f_pull,
- int sign,
- rvec* f)
+static void apply_forces_grp_part(const pull_group_work_t* pgrp,
+ int ind_start,
+ int ind_end,
+ gmx::ArrayRef<const real> masses,
+ const dvec f_pull,
+ int sign,
+ rvec* f)
{
double inv_wm = pgrp->mwscale;
}
/* Apply forces in a mass weighted fashion */
-static void apply_forces_grp(const pull_group_work_t* pgrp,
- const real* masses,
- const dvec f_pull,
- int sign,
- rvec* f,
- int nthreads)
+static void apply_forces_grp(const pull_group_work_t* pgrp,
+ gmx::ArrayRef<const real> masses,
+ const dvec f_pull,
+ int sign,
+ rvec* f,
+ int nthreads)
{
auto localAtomIndices = pgrp->atomSet.localIndex();
}
/* Apply forces in a mass weighted fashion to a cylinder group */
-static void apply_forces_cyl_grp(const pull_group_work_t* pgrp,
- const double dv_corr,
- const real* masses,
- const dvec f_pull,
- double f_scal,
- int sign,
- rvec* f,
+static void apply_forces_cyl_grp(const pull_group_work_t* pgrp,
+ const double dv_corr,
+ gmx::ArrayRef<const real> masses,
+ const dvec f_pull,
+ double f_scal,
+ int sign,
+ rvec* f,
int gmx_unused nthreads)
{
double inv_wm = pgrp->mwscale;
/* Apply torque forces in a mass weighted fashion to the groups that define
* the pull vector direction for pull coordinate pcrd.
*/
-static void apply_forces_vec_torque(const struct pull_t* pull,
- const pull_coord_work_t* pcrd,
- const real* masses,
- rvec* f)
+static void apply_forces_vec_torque(const struct pull_t* pull,
+ const pull_coord_work_t* pcrd,
+ gmx::ArrayRef<const real> masses,
+ rvec* f)
{
const PullCoordSpatialData& spatialData = pcrd->spatialData;
static void apply_forces_coord(struct pull_t* pull,
int coord,
const PullCoordVectorForces& forces,
- const real* masses,
+ gmx::ArrayRef<const real> masses,
rvec* f)
{
/* Here it would be more efficient to use one large thread-parallel
}
/* Apply constraint using SHAKE */
-static void
-do_constraint(struct pull_t* pull, t_pbc* pbc, rvec* x, rvec* v, gmx_bool bMaster, tensor vir, double dt, double t)
+static void do_constraint(struct pull_t* pull,
+ t_pbc* pbc,
+ gmx::ArrayRef<gmx::RVec> x,
+ gmx::ArrayRef<gmx::RVec> v,
+ gmx_bool bMaster,
+ tensor vir,
+ double dt,
+ double t)
{
dvec* r_ij; /* x[i] com of i in prev. step. Obeys constr. -> r_ij[i] */
/* DONE ITERATING, NOW UPDATE COORDINATES AND CALC. CONSTRAINT FORCES */
- if (v)
+ if (!v.empty())
{
invdt = 1 / dt;
}
{
x[ii][m] += tmp[m];
}
- if (v)
+ if (!v.empty())
{
for (m = 0; m < DIM; m++)
{
* potential energy is added either to the pull term or to a term
* specific to the external module.
*/
-void apply_external_pull_coord_force(struct pull_t* pull,
- int coord_index,
- double coord_force,
- const real* masses,
- gmx::ForceWithVirial* forceWithVirial)
+void apply_external_pull_coord_force(struct pull_t* pull,
+ int coord_index,
+ double coord_force,
+ gmx::ArrayRef<const real> masses,
+ gmx::ForceWithVirial* forceWithVirial)
{
pull_coord_work_t* pcrd;
return pullCoordForces;
}
-real pull_potential(struct pull_t* pull,
- const real* masses,
- t_pbc* pbc,
- const t_commrec* cr,
- double t,
- real lambda,
- const rvec* x,
- gmx::ForceWithVirial* force,
- real* dvdlambda)
+real pull_potential(struct pull_t* pull,
+ gmx::ArrayRef<const real> masses,
+ t_pbc* pbc,
+ const t_commrec* cr,
+ double t,
+ real lambda,
+ gmx::ArrayRef<const gmx::RVec> x,
+ gmx::ForceWithVirial* force,
+ real* dvdlambda)
{
real V = 0;
{
real dVdl = 0;
- pull_calc_coms(cr, pull, masses, pbc, t, x, nullptr);
+ pull_calc_coms(cr, pull, masses, pbc, t, x, {});
rvec* f = as_rvec_array(force->force_.data());
matrix virial = { { 0 } };
return (MASTER(cr) ? V : 0.0);
}
-void pull_constraint(struct pull_t* pull,
- const real* masses,
- t_pbc* pbc,
- const t_commrec* cr,
- double dt,
- double t,
- rvec* x,
- rvec* xp,
- rvec* v,
- tensor vir)
+void pull_constraint(struct pull_t* pull,
+ gmx::ArrayRef<const real> masses,
+ t_pbc* pbc,
+ const t_commrec* cr,
+ double dt,
+ double t,
+ gmx::ArrayRef<gmx::RVec> x,
+ gmx::ArrayRef<gmx::RVec> xp,
+ gmx::ArrayRef<gmx::RVec> v,
+ tensor vir)
{
assert(pull != nullptr);
delete pull;
}
-void preparePrevStepPullCom(const t_inputrec* ir,
- pull_t* pull_work,
- const real* masses,
- t_state* state,
- const t_state* state_global,
- const t_commrec* cr,
- bool startingFromCheckpoint)
+void preparePrevStepPullCom(const t_inputrec* ir,
+ pull_t* pull_work,
+ gmx::ArrayRef<const real> masses,
+ t_state* state,
+ const t_state* state_global,
+ const t_commrec* cr,
+ bool startingFromCheckpoint)
{
if (!ir->pull || !ir->pull->bSetPbcRefToPrevStepCOM)
{
{
t_pbc pbc;
set_pbc(&pbc, ir->pbcType, state->box);
- initPullComFromPrevStep(cr, pull_work, masses, &pbc, state->x.rvec_array());
+ initPullComFromPrevStep(
+ cr, pull_work, masses, &pbc, state->x.arrayRefWithPadding().unpaddedArrayRef());
updatePrevStepPullCom(pull_work, state);
}
}
namespace gmx
{
+template<typename>
+class ArrayRef;
class ForceWithVirial;
class LocalAtomSetManager;
} // namespace gmx
* \param[in] masses Atoms masses.
* \param[in,out] forceWithVirial Force and virial buffers.
*/
-void apply_external_pull_coord_force(struct pull_t* pull,
- int coord_index,
- double coord_force,
- const real* masses,
- gmx::ForceWithVirial* forceWithVirial);
-
+void apply_external_pull_coord_force(struct pull_t* pull,
+ int coord_index,
+ double coord_force,
+ gmx::ArrayRef<const real> masses,
+ gmx::ForceWithVirial* forceWithVirial);
/*! \brief Set the all the pull forces to zero.
*
*
* \returns The pull potential energy.
*/
-real pull_potential(struct pull_t* pull,
- const real* masses,
- struct t_pbc* pbc,
- const t_commrec* cr,
- double t,
- real lambda,
- const rvec* x,
- gmx::ForceWithVirial* force,
- real* dvdlambda);
+real pull_potential(struct pull_t* pull,
+ gmx::ArrayRef<const real> masses,
+ struct t_pbc* pbc,
+ const t_commrec* cr,
+ double t,
+ real lambda,
+ gmx::ArrayRef<const gmx::RVec> x,
+ gmx::ForceWithVirial* force,
+ real* dvdlambda);
/*! \brief Constrain the coordinates xp in the directions in x
* \param[in,out] v Velocities, which may get a pull correction.
* \param[in,out] vir The virial, which, if != NULL, gets a pull correction.
*/
-void pull_constraint(struct pull_t* pull,
- const real* masses,
- struct t_pbc* pbc,
- const t_commrec* cr,
- double dt,
- double t,
- rvec* x,
- rvec* xp,
- rvec* v,
- tensor vir);
+void pull_constraint(struct pull_t* pull,
+ gmx::ArrayRef<const real> masses,
+ struct t_pbc* pbc,
+ const t_commrec* cr,
+ double dt,
+ double t,
+ gmx::ArrayRef<gmx::RVec> x,
+ gmx::ArrayRef<gmx::RVec> xp,
+ gmx::ArrayRef<gmx::RVec> v,
+ tensor vir);
/*! \brief Make a selection of the home atoms for all pull groups.
* \param[in,out] xp Updated x, can be NULL.
*
*/
-void pull_calc_coms(const t_commrec* cr,
- pull_t* pull,
- const real* masses,
- t_pbc* pbc,
- double t,
- const rvec x[],
- rvec* xp);
+void pull_calc_coms(const t_commrec* cr,
+ pull_t* pull,
+ gmx::ArrayRef<const real> masses,
+ t_pbc* pbc,
+ double t,
+ gmx::ArrayRef<const gmx::RVec> x,
+ gmx::ArrayRef<gmx::RVec> xp);
/*! \brief Margin for checking pull group PBC distances compared to half the box size */
static constexpr real c_pullGroupPbcMargin = 0.9;
* \param[in] pbcMargin The minimum margin (as a fraction) to half the box size
* \returns -1 when all groups obey PBC or the first group index that fails PBC
*/
-int pullCheckPbcWithinGroups(const pull_t& pull, const rvec* x, const t_pbc& pbc, real pbcMargin);
+int pullCheckPbcWithinGroups(const pull_t& pull, gmx::ArrayRef<const gmx::RVec> x, const t_pbc& pbc, real pbcMargin);
/*! \brief Checks whether a specific group that uses a reference atom is within PBC restrictions
*
* \param[in] cr Struct for communication info.
* \param[in] startingFromCheckpoint Is the simulation starting from a checkpoint?
*/
-void preparePrevStepPullCom(const t_inputrec* ir,
- pull_t* pull_work,
- const real* masses,
- t_state* state,
- const t_state* state_global,
- const t_commrec* cr,
- bool startingFromCheckpoint);
+void preparePrevStepPullCom(const t_inputrec* ir,
+ pull_t* pull_work,
+ gmx::ArrayRef<const real> masses,
+ t_state* state,
+ const t_state* state_global,
+ const t_commrec* cr,
+ bool startingFromCheckpoint);
/*! \brief Initializes the COM of the previous step (set to initial COM)
*
* \param[in] pbc Information struct about periodicity.
* \param[in] x The local positions.
*/
-void initPullComFromPrevStep(const t_commrec* cr, pull_t* pull, const real* masses, t_pbc* pbc, const rvec x[]);
+void initPullComFromPrevStep(const t_commrec* cr,
+ pull_t* pull,
+ gmx::ArrayRef<const real> masses,
+ t_pbc* pbc,
+ gmx::ArrayRef<const gmx::RVec> x);
#endif
#include "gromacs/math/units.h"
#include "gromacs/math/utilities.h"
#include "gromacs/math/vec.h"
+#include "gromacs/math/vectypes.h"
#include "gromacs/mdtypes/commrec.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/pbc.h"
#include "gromacs/pulling/pull.h"
+#include "gromacs/utility/arrayref.h"
#include "gromacs/utility/fatalerror.h"
#include "gromacs/utility/futil.h"
#include "gromacs/utility/gmxassert.h"
/* Copies the coordinates of the PBC atom of pgrp to x_pbc.
* When those coordinates are not available on this rank, clears x_pbc.
*/
-static void setPbcAtomCoords(const pull_group_work_t& pgrp, const rvec* x, rvec x_pbc)
+static void setPbcAtomCoords(const pull_group_work_t& pgrp, gmx::ArrayRef<const gmx::RVec> x, rvec x_pbc)
{
if (pgrp.pbcAtomSet != nullptr)
{
}
}
-static void pull_set_pbcatoms(const t_commrec* cr, struct pull_t* pull, const rvec* x, gmx::ArrayRef<gmx::RVec> x_pbc)
+static void pull_set_pbcatoms(const t_commrec* cr,
+ struct pull_t* pull,
+ gmx::ArrayRef<const gmx::RVec> x,
+ gmx::ArrayRef<gmx::RVec> x_pbc)
{
int numPbcAtoms = 0;
for (size_t g = 0; g < pull->group.size(); g++)
}
}
-static void
-make_cyl_refgrps(const t_commrec* cr, pull_t* pull, const real* masses, t_pbc* pbc, double t, const rvec* x)
+static void make_cyl_refgrps(const t_commrec* cr,
+ pull_t* pull,
+ gmx::ArrayRef<const real> masses,
+ t_pbc* pbc,
+ double t,
+ gmx::ArrayRef<const gmx::RVec> x)
{
pull_comm_t* comm = &pull->comm;
return a;
}
-static void sum_com_part(const pull_group_work_t* pgrp,
- int ind_start,
- int ind_end,
- const rvec* x,
- const rvec* xp,
- const real* mass,
- const t_pbc* pbc,
- const rvec x_pbc,
- ComSums* sum_com)
+static void sum_com_part(const pull_group_work_t* pgrp,
+ int ind_start,
+ int ind_end,
+ gmx::ArrayRef<const gmx::RVec> x,
+ gmx::ArrayRef<const gmx::RVec> xp,
+ gmx::ArrayRef<const real> mass,
+ const t_pbc* pbc,
+ const rvec x_pbc,
+ ComSums* sum_com)
{
double sum_wm = 0;
double sum_wwm = 0;
{
sum_wmx[d] += wm * x[ii][d];
}
- if (xp)
+ if (!xp.empty())
{
for (int d = 0; d < DIM; d++)
{
{
sum_wmx[d] += wm * dx[d];
}
- if (xp)
+ if (!xp.empty())
{
/* For xp add the difference between xp and x to dx,
* such that we use the same periodic image,
sum_com->sum_wm = sum_wm;
sum_com->sum_wwm = sum_wwm;
copy_dvec(sum_wmx, sum_com->sum_wmx);
- if (xp)
+ if (!xp.empty())
{
copy_dvec(sum_wmxp, sum_com->sum_wmxp);
}
}
-static void sum_com_part_cosweight(const pull_group_work_t* pgrp,
- int ind_start,
- int ind_end,
- int cosdim,
- real twopi_box,
- const rvec* x,
- const rvec* xp,
- const real* mass,
- ComSums* sum_com)
+static void sum_com_part_cosweight(const pull_group_work_t* pgrp,
+ int ind_start,
+ int ind_end,
+ int cosdim,
+ real twopi_box,
+ gmx::ArrayRef<const gmx::RVec> x,
+ gmx::ArrayRef<const gmx::RVec> xp,
+ gmx::ArrayRef<const real> mass,
+ ComSums* sum_com)
{
/* Cosine weighting geometry */
double sum_cm = 0;
sum_csm += static_cast<double>(cw * sw * m);
sum_ssm += static_cast<double>(sw * sw * m);
- if (xp != nullptr)
+ if (!xp.empty())
{
real cw = std::cos(xp[ii][cosdim] * twopi_box);
real sw = std::sin(xp[ii][cosdim] * twopi_box);
}
/* calculates center of mass of selection index from all coordinates x */
-void pull_calc_coms(const t_commrec* cr, pull_t* pull, const real* masses, t_pbc* pbc, double t, const rvec x[], rvec* xp)
+void pull_calc_coms(const t_commrec* cr,
+ pull_t* pull,
+ gmx::ArrayRef<const real> masses,
+ t_pbc* pbc,
+ double t,
+ gmx::ArrayRef<const gmx::RVec> x,
+ gmx::ArrayRef<gmx::RVec> xp)
{
real twopi_box = 0;
pull_comm_t* comm;
if (pgrp->params.ind.size() == 1 && pgrp->atomSet.numAtomsLocal() == 1
&& masses[pgrp->atomSet.localIndex()[0]] == 0)
{
- GMX_ASSERT(xp == nullptr,
+ GMX_ASSERT(xp.empty(),
"We should not have groups with zero mass with constraints, i.e. "
- "xp!=NULL");
+ "xp not empty");
/* Copy the single atom coordinate */
for (int d = 0; d < DIM; d++)
for (m = 0; m < DIM; m++)
{
pgrp->x[m] = comBuffer[0][m] * pgrp->mwscale;
- if (xp)
+ if (!xp.empty())
{
pgrp->xp[m] = comBuffer[1][m] * pgrp->mwscale;
}
if (pgrp->epgrppbc == epgrppbcREFAT || pgrp->epgrppbc == epgrppbcPREVSTEPCOM)
{
pgrp->x[m] += comm->pbcAtomBuffer[g][m];
- if (xp)
+ if (!xp.empty())
{
pgrp->xp[m] += comm->pbcAtomBuffer[g][m];
}
pgrp->localWeights[i] = csw * std::cos(twopi_box * x[ii][pull->cosdim])
+ snw * std::sin(twopi_box * x[ii][pull->cosdim]);
}
- if (xp)
+ if (!xp.empty())
{
csw = comBuffer[2][0];
snw = comBuffer[2][1];
using BoolVec = gmx::BasicVector<bool>;
/* Returns whether the pull group obeys the PBC restrictions */
-static bool pullGroupObeysPbcRestrictions(const pull_group_work_t& group,
- const BoolVec& dimUsed,
- const rvec* x,
- const t_pbc& pbc,
- const gmx::RVec& x_pbc,
- const real pbcMargin)
+static bool pullGroupObeysPbcRestrictions(const pull_group_work_t& group,
+ const BoolVec& dimUsed,
+ gmx::ArrayRef<const gmx::RVec> x,
+ const t_pbc& pbc,
+ const gmx::RVec& x_pbc,
+ const real pbcMargin)
{
/* Determine which dimensions are relevant for PBC */
BoolVec dimUsesPbc = { false, false, false };
return true;
}
-int pullCheckPbcWithinGroups(const pull_t& pull, const rvec* x, const t_pbc& pbc, real pbcMargin)
+int pullCheckPbcWithinGroups(const pull_t& pull, gmx::ArrayRef<const gmx::RVec> x, const t_pbc& pbc, real pbcMargin)
{
if (pbc.pbcType == PbcType::No)
{
}
return (pullGroupObeysPbcRestrictions(
- group, dimUsed, as_rvec_array(x.data()), pbc, pull.comm.pbcAtomBuffer[groupNr], pbcMargin));
+ group, dimUsed, x, pbc, pull.comm.pbcAtomBuffer[groupNr], pbcMargin));
}
void setPrevStepPullComFromState(struct pull_t* pull, const t_state* state)
}
}
-void initPullComFromPrevStep(const t_commrec* cr, pull_t* pull, const real* masses, t_pbc* pbc, const rvec x[])
+void initPullComFromPrevStep(const t_commrec* cr,
+ pull_t* pull,
+ gmx::ArrayRef<const real> masses,
+ t_pbc* pbc,
+ gmx::ArrayRef<const gmx::RVec> x)
{
pull_comm_t* comm = &pull->comm;
size_t ngroup = pull->group.size();
"Groups with no atoms, or only one atom, should not "
"use the COM from the previous step as reference.");
- rvec x_pbc = { 0, 0, 0 };
+ gmx::RVec x_pbc = { 0, 0, 0 };
copy_rvec(comm->pbcAtomBuffer[g], x_pbc);
if (debug)
if (pgrp->atomSet.numAtomsLocal() <= c_pullMaxNumLocalAtomsSingleThreaded)
{
- sum_com_part(pgrp, 0, pgrp->atomSet.numAtomsLocal(), x, nullptr, masses, pbc, x_pbc, &comSumsTotal);
+ sum_com_part(pgrp, 0, pgrp->atomSet.numAtomsLocal(), x, {}, masses, pbc, x_pbc, &comSumsTotal);
}
else
{
{
int ind_start = (pgrp->atomSet.numAtomsLocal() * (t + 0)) / pull->nthreads;
int ind_end = (pgrp->atomSet.numAtomsLocal() * (t + 1)) / pull->nthreads;
- sum_com_part(
- pgrp, ind_start, ind_end, x, nullptr, masses, pbc, x_pbc, &pull->comSums[t]);
+ sum_com_part(pgrp, ind_start, ind_end, x, {}, masses, pbc, x_pbc, &pull->comSums[t]);
}
/* Reduce the thread contributions to sum_com[0] */