void dd_clear_f_vsites(struct gmx_domdec_t* dd, rvec* f);
/*! \brief Move x0 and also x1 if x1!=NULL. bX1IsCoord tells if to do PBC on x1 */
-void dd_move_x_constraints(struct gmx_domdec_t* dd, const matrix box, rvec* x0, rvec* x1, gmx_bool bX1IsCoord);
+void dd_move_x_constraints(struct gmx_domdec_t* dd,
+ const matrix box,
+ gmx::ArrayRef<gmx::RVec> x0,
+ gmx::ArrayRef<gmx::RVec> x1,
+ gmx_bool bX1IsCoord);
/*! \brief Communicates the coordinates involved in virtual sites */
void dd_move_x_vsites(struct gmx_domdec_t* dd, const matrix box, rvec* x);
//! @endcond
};
-void dd_move_x_constraints(gmx_domdec_t* dd, const matrix box, rvec* x0, rvec* x1, gmx_bool bX1IsCoord)
+void dd_move_x_constraints(gmx_domdec_t* dd,
+ const matrix box,
+ gmx::ArrayRef<gmx::RVec> x0,
+ gmx::ArrayRef<gmx::RVec> x1,
+ gmx_bool bX1IsCoord)
{
if (dd->constraint_comm)
{
- dd_move_x_specat(dd, dd->constraint_comm, box, x0, x1, bX1IsCoord);
+ dd_move_x_specat(dd, dd->constraint_comm, box, as_rvec_array(x0.data()),
+ as_rvec_array(x1.data()), bX1IsCoord);
ddReopenBalanceRegionCpu(dd);
}
#include "gromacs/utility/gmxassert.h"
#include "gromacs/utility/listoflists.h"
#include "gromacs/utility/pleasecite.h"
-#include "gromacs/utility/smalloc.h"
#include "gromacs/utility/txtdump.h"
namespace gmx
int numSettles);
~Impl();
void setConstraints(const gmx_localtop_t& top, const t_mdatoms& md);
- bool apply(bool bLog,
- bool bEner,
- int64_t step,
- int delta_step,
- real step_scaling,
- rvec* x,
- rvec* xprime,
- rvec* min_proj,
- const matrix box,
- real lambda,
- real* dvdlambda,
- rvec* v,
- tensor* vir,
- ConstraintVariable econq);
+ bool apply(bool bLog,
+ bool bEner,
+ int64_t step,
+ int delta_step,
+ real step_scaling,
+ ArrayRefWithPadding<RVec> x,
+ ArrayRefWithPadding<RVec> xprime,
+ ArrayRef<RVec> min_proj,
+ const matrix box,
+ real lambda,
+ real* dvdlambda,
+ ArrayRefWithPadding<RVec> v,
+ tensor* vir,
+ ConstraintVariable econq);
//! The total number of constraints.
int ncon_tot = 0;
//! The number of flexible constraints.
}
//! Clears constraint quantities for atoms in nonlocal region.
-static void clear_constraint_quantity_nonlocal(gmx_domdec_t* dd, rvec* q)
+static void clear_constraint_quantity_nonlocal(gmx_domdec_t* dd, ArrayRef<RVec> q)
{
int nonlocal_at_start, nonlocal_at_end, at;
}
//! Writes out coordinates.
-static void write_constr_pdb(const char* fn,
- const char* title,
- const gmx_mtop_t& mtop,
- int start,
- int homenr,
- const t_commrec* cr,
- const rvec x[],
- const matrix box)
+static void write_constr_pdb(const char* fn,
+ const char* title,
+ const gmx_mtop_t& mtop,
+ int start,
+ int homenr,
+ const t_commrec* cr,
+ ArrayRef<const RVec> x,
+ const matrix box)
{
char fname[STRLEN];
FILE* out;
}
//! Writes out domain contents to help diagnose crashes.
-static void dump_confs(FILE* log,
- int64_t step,
- const gmx_mtop_t& mtop,
- int start,
- int homenr,
- const t_commrec* cr,
- const rvec x[],
- rvec xprime[],
- const matrix box)
+static void dump_confs(FILE* log,
+ int64_t step,
+ const gmx_mtop_t& mtop,
+ int start,
+ int homenr,
+ const t_commrec* cr,
+ ArrayRef<const RVec> x,
+ ArrayRef<const RVec> xprime,
+ const matrix box)
{
char buf[STRLEN], buf2[22];
fprintf(stderr, "Wrote pdb files with previous and current coordinates\n");
}
-bool Constraints::apply(bool bLog,
- bool bEner,
- int64_t step,
- int delta_step,
- real step_scaling,
- rvec* x,
- rvec* xprime,
- rvec* min_proj,
- const matrix box,
- real lambda,
- real* dvdlambda,
- rvec* v,
- tensor* vir,
- ConstraintVariable econq)
+bool Constraints::apply(bool bLog,
+ bool bEner,
+ int64_t step,
+ int delta_step,
+ real step_scaling,
+ ArrayRefWithPadding<RVec> x,
+ ArrayRefWithPadding<RVec> xprime,
+ ArrayRef<RVec> min_proj,
+ const matrix box,
+ real lambda,
+ real* dvdlambda,
+ ArrayRefWithPadding<RVec> v,
+ tensor* vir,
+ ConstraintVariable econq)
{
- return impl_->apply(bLog, bEner, step, delta_step, step_scaling, x, xprime, min_proj, box,
- lambda, dvdlambda, v, vir, econq);
+ return impl_->apply(bLog, bEner, step, delta_step, step_scaling, std::move(x), std::move(xprime),
+ min_proj, box, lambda, dvdlambda, std::move(v), vir, econq);
}
-bool Constraints::Impl::apply(bool bLog,
- bool bEner,
- int64_t step,
- int delta_step,
- real step_scaling,
- rvec* x,
- rvec* xprime,
- rvec* min_proj,
- const matrix box,
- real lambda,
- real* dvdlambda,
- rvec* v,
- tensor* vir,
- ConstraintVariable econq)
+bool Constraints::Impl::apply(bool bLog,
+ bool bEner,
+ int64_t step,
+ int delta_step,
+ real step_scaling,
+ ArrayRefWithPadding<RVec> x,
+ ArrayRefWithPadding<RVec> xprime,
+ ArrayRef<RVec> min_proj,
+ const matrix box,
+ real lambda,
+ real* dvdlambda,
+ ArrayRefWithPadding<RVec> v,
+ tensor* vir,
+ ConstraintVariable econq)
{
bool bOK, bDump;
int start, homenr;
*/
if (cr->dd)
{
- dd_move_x_constraints(cr->dd, box, x, xprime, econq == ConstraintVariable::Positions);
+ dd_move_x_constraints(cr->dd, box, x.unpaddedArrayRef(), xprime.unpaddedArrayRef(),
+ econq == ConstraintVariable::Positions);
- if (v != nullptr)
+ if (!v.empty())
{
/* We need to initialize the non-local components of v.
* We never actually use these values, but we do increment them,
* so we should avoid uninitialized variables and overflows.
*/
- clear_constraint_quantity_nonlocal(cr->dd, v);
+ clear_constraint_quantity_nonlocal(cr->dd, v.unpaddedArrayRef());
}
}
if (lincsd != nullptr)
{
bOK = constrain_lincs(bLog || bEner, ir, step, lincsd, md, cr, ms, x, xprime, min_proj, box,
- pbc_null, lambda, dvdlambda, invdt, v, vir != nullptr, vir_r_m_dr,
- econq, nrnb, maxwarn, &warncount_lincs);
+ pbc_null, lambda, dvdlambda, invdt, v.unpaddedArrayRef(),
+ vir != nullptr, vir_r_m_dr, econq, nrnb, maxwarn, &warncount_lincs);
if (!bOK && maxwarn < INT_MAX)
{
if (log != nullptr)
if (shaked != nullptr)
{
- bOK = constrain_shake(log, shaked, md.invmass, *idef, ir, x, xprime, min_proj, nrnb, lambda,
- dvdlambda, invdt, v, vir != nullptr, vir_r_m_dr, maxwarn < INT_MAX, econq);
+ bOK = constrain_shake(
+ log, shaked, md.invmass, *idef, ir, as_rvec_array(x.unpaddedArrayRef().data()),
+ as_rvec_array(xprime.unpaddedArrayRef().data()), as_rvec_array(min_proj.data()),
+ nrnb, lambda, dvdlambda, invdt, as_rvec_array(v.unpaddedArrayRef().data()),
+ vir != nullptr, vir_r_m_dr, maxwarn < INT_MAX, econq);
if (!bOK && maxwarn < INT_MAX)
{
clear_mat(vir_r_m_dr_th[th]);
}
- csettle(settled, nth, th, pbc_null, x[0], xprime[0], invdt, v ? v[0] : nullptr,
- vir != nullptr, th == 0 ? vir_r_m_dr : vir_r_m_dr_th[th],
+ csettle(settled, nth, th, pbc_null, x, xprime, invdt, v, vir != nullptr,
+ th == 0 ? vir_r_m_dr : vir_r_m_dr_th[th],
th == 0 ? &bSettleErrorHasOccurred0 : &bSettleErrorHasOccurred[th]);
}
GMX_CATCH_ALL_AND_EXIT_WITH_FATAL_ERROR
}
inc_nrnb(nrnb, eNR_SETTLE, nsettle);
- if (v != nullptr)
+ if (!v.empty())
{
inc_nrnb(nrnb, eNR_CONSTR_V, nsettle * 3);
}
{
settle_proj(settled, econq, end_th - start_th,
settle->iatoms + start_th * (1 + NRAL(F_SETTLE)), pbc_null,
- x, xprime, min_proj, calcvir_atom_end,
- th == 0 ? vir_r_m_dr : vir_r_m_dr_th[th]);
+ x.unpaddedArrayRef(), xprime.unpaddedArrayRef(), min_proj,
+ calcvir_atom_end, th == 0 ? vir_r_m_dr : vir_r_m_dr_th[th]);
}
}
GMX_CATCH_ALL_AND_EXIT_WITH_FATAL_ERROR
if (bDump)
{
- dump_confs(log, step, mtop, start, homenr, cr, x, xprime, box);
+ dump_confs(log, step, mtop, start, homenr, cr, x.unpaddedArrayRef(),
+ xprime.unpaddedArrayRef(), box);
}
if (econq == ConstraintVariable::Positions)
t = ir.init_t;
}
set_pbc(&pbc, ir.pbcType, box);
- pull_constraint(pull_work, &md, &pbc, cr, ir.delta_t, t, x, xprime, v, *vir);
+ pull_constraint(pull_work, &md, &pbc, cr, ir.delta_t, t,
+ as_rvec_array(x.unpaddedArrayRef().data()),
+ as_rvec_array(xprime.unpaddedArrayRef().data()),
+ as_rvec_array(v.unpaddedArrayRef().data()), *vir);
}
if (ed && delta_step > 0)
{
/* apply the essential dynamics constraints here */
- do_edsam(&ir, step, cr, xprime, v, box, ed);
+ do_edsam(&ir, step, cr, as_rvec_array(xprime.unpaddedArrayRef().data()),
+ as_rvec_array(v.unpaddedArrayRef().data()), box, ed);
}
}
wallcycle_stop(wcycle, ewcCONSTR);
- if (v != nullptr && md.cFREEZE)
+ if (!v.empty() && md.cFREEZE)
{
/* Set the velocities of frozen dimensions to zero */
+ ArrayRef<RVec> vRef = v.unpaddedArrayRef();
int gmx_unused numThreads = gmx_omp_nthreads_get(emntUpdate);
{
if (ir.opts.nFreeze[freezeGroup][d])
{
- v[i][d] = 0;
+ vRef[i][d] = 0;
}
}
}
int64_t step;
real dt = ir->delta_t;
real dvdl_dum;
- rvec* savex;
- auto xRvec = as_rvec_array(x.paddedArrayRef().data());
- auto vRvec = as_rvec_array(v.paddedArrayRef().data());
-
- /* We need to allocate one element extra, since we might use
- * (unaligned) 4-wide SIMD loads to access rvec entries.
- */
- snew(savex, natoms + 1);
+ PaddedVector<RVec> savex(natoms);
start = 0;
end = md->homenr;
dvdl_dum = 0;
/* constrain the current position */
- constr->apply(TRUE, FALSE, step, 0, 1.0, xRvec, xRvec, nullptr, box, lambda, &dvdl_dum, nullptr,
- nullptr, gmx::ConstraintVariable::Positions);
+ constr->apply(TRUE, FALSE, step, 0, 1.0, x, x, ArrayRef<RVec>(), box, lambda, &dvdl_dum,
+ ArrayRefWithPadding<RVec>(), nullptr, gmx::ConstraintVariable::Positions);
if (EI_VV(ir->eI))
{
/* constrain the inital velocity, and save it */
/* also may be useful if we need the ekin from the halfstep for velocity verlet */
- constr->apply(TRUE, FALSE, step, 0, 1.0, xRvec, vRvec, vRvec, box, lambda, &dvdl_dum,
- nullptr, nullptr, gmx::ConstraintVariable::Velocities);
+ constr->apply(TRUE, FALSE, step, 0, 1.0, x, v, v.unpaddedArrayRef(), box, lambda, &dvdl_dum,
+ ArrayRefWithPadding<RVec>(), nullptr, gmx::ConstraintVariable::Velocities);
}
/* constrain the inital velocities at t-dt/2 */
if (EI_STATE_VELOCITY(ir->eI) && ir->eI != eiVV)
fprintf(fplog, "\nConstraining the coordinates at t0-dt (step %s)\n", gmx_step_str(step, buf));
}
dvdl_dum = 0;
- constr->apply(TRUE, FALSE, step, -1, 1.0, xRvec, savex, nullptr, box, lambda, &dvdl_dum,
- vRvec, nullptr, gmx::ConstraintVariable::Positions);
+ constr->apply(TRUE, FALSE, step, -1, 1.0, x, savex.arrayRefWithPadding(), ArrayRef<RVec>(),
+ box, lambda, &dvdl_dum, v, nullptr, gmx::ConstraintVariable::Positions);
for (i = start; i < end; i++)
{
}
}
}
- sfree(savex);
}
} // namespace gmx
#include <cstdio>
+#include "gromacs/math/arrayrefwithpadding.h"
#include "gromacs/math/vectypes.h"
#include "gromacs/topology/idef.h"
#include "gromacs/utility/arrayref.h"
*
* /note x is non-const, because non-local atoms need to be communicated.
*/
- bool apply(bool bLog,
- bool bEner,
- int64_t step,
- int delta_step,
- real step_scaling,
- rvec* x,
- rvec* xprime,
- rvec* min_proj,
- const matrix box,
- real lambda,
- real* dvdlambda,
- rvec* v,
- tensor* vir,
- ConstraintVariable econq);
+ bool apply(bool bLog,
+ bool bEner,
+ int64_t step,
+ int delta_step,
+ real step_scaling,
+ ArrayRefWithPadding<RVec> x,
+ ArrayRefWithPadding<RVec> xprime,
+ ArrayRef<RVec> min_proj,
+ const matrix box,
+ real lambda,
+ real* dvdlambda,
+ ArrayRefWithPadding<RVec> v,
+ tensor* vir,
+ ConstraintVariable econq);
//! Links the essentialdynamics and constraint code.
void saveEdsamPointer(gmx_edsam* ed);
//! Getter for use by domain decomposition.
#endif // GMX_SIMD_HAVE_REAL
/*! \brief LINCS projection, works on derivatives of the coordinates. */
-static void do_lincsp(const rvec* x,
- rvec* f,
- rvec* fp,
- t_pbc* pbc,
- Lincs* lincsd,
- int th,
- real* invmass,
- ConstraintVariable econq,
- bool bCalcDHDL,
- bool bCalcVir,
- tensor rmdf)
+static void do_lincsp(ArrayRefWithPadding<const RVec> xPadded,
+ ArrayRefWithPadding<RVec> fPadded,
+ ArrayRef<RVec> fp,
+ t_pbc* pbc,
+ Lincs* lincsd,
+ int th,
+ real* invmass,
+ ConstraintVariable econq,
+ bool bCalcDHDL,
+ bool bCalcVir,
+ tensor rmdf)
{
const int b0 = lincsd->task[th].b0;
const int b1 = lincsd->task[th].b1;
gmx::ArrayRef<real> rhs2 = lincsd->tmp2;
gmx::ArrayRef<real> sol = lincsd->tmp3;
+ const rvec* x = as_rvec_array(xPadded.paddedArrayRef().data());
+ rvec* f = as_rvec_array(fPadded.paddedArrayRef().data());
+
#if GMX_SIMD_HAVE_REAL
/* This SIMD code does the same as the plain-C code after the #else.
* The only difference is that we always call pbc code, as with SIMD
/* When constraining forces, we should not use mass weighting,
* so we pass invmass=NULL, which results in the use of 1 for all atoms.
*/
- lincs_update_atoms(lincsd, th, 1.0, sol, r,
- (econq != ConstraintVariable::Force) ? invmass : nullptr, fp);
+ lincs_update_atoms(lincsd, th, 1.0, sol, r, (econq != ConstraintVariable::Force) ? invmass : nullptr,
+ as_rvec_array(fp.data()));
if (bCalcDHDL)
{
#endif // GMX_SIMD_HAVE_REAL
//! Implements LINCS constraining.
-static void do_lincs(const rvec* x,
- rvec* xp,
- const matrix box,
- t_pbc* pbc,
- Lincs* lincsd,
- int th,
- const real* invmass,
- const t_commrec* cr,
- bool bCalcDHDL,
- real wangle,
- bool* bWarn,
- real invdt,
- rvec* gmx_restrict v,
- bool bCalcVir,
- tensor vir_r_m_dr)
+static void do_lincs(ArrayRefWithPadding<const RVec> xPadded,
+ ArrayRefWithPadding<RVec> xpPadded,
+ const matrix box,
+ t_pbc* pbc,
+ Lincs* lincsd,
+ int th,
+ const real* invmass,
+ const t_commrec* cr,
+ bool bCalcDHDL,
+ real wangle,
+ bool* bWarn,
+ real invdt,
+ ArrayRef<RVec> vRef,
+ bool bCalcVir,
+ tensor vir_r_m_dr)
{
+ const rvec* x = as_rvec_array(xPadded.paddedArrayRef().data());
+ rvec* xp = as_rvec_array(xpPadded.paddedArrayRef().data());
+ rvec* gmx_restrict v = as_rvec_array(vRef.data());
+
const int b0 = lincsd->task[th].b0;
const int b1 = lincsd->task[th].b1;
/* Communicate the corrected non-local coordinates */
if (DOMAINDECOMP(cr))
{
- dd_move_x_constraints(cr->dd, box, xp, nullptr, FALSE);
+ dd_move_x_constraints(cr->dd, box, xpPadded.unpaddedArrayRef(), ArrayRef<RVec>(), FALSE);
}
}
#pragma omp barrier
//! Issues a warning when LINCS constraints cannot be satisfied.
static void lincs_warning(gmx_domdec_t* dd,
- const rvec* x,
- rvec* xprime,
+ ArrayRef<const RVec> x,
+ ArrayRef<const RVec> xprime,
t_pbc* pbc,
int ncons,
gmx::ArrayRef<const AtomPair> atoms,
};
//! Determine how well the constraints have been satisfied.
-static LincsDeviations makeLincsDeviations(const Lincs& lincsd, const rvec* x, const t_pbc* pbc)
+static LincsDeviations makeLincsDeviations(const Lincs& lincsd, ArrayRef<const RVec> x, const t_pbc* pbc)
{
LincsDeviations result;
const ArrayRef<const AtomPair> atoms = lincsd.atoms;
return result;
}
-bool constrain_lincs(bool computeRmsd,
- const t_inputrec& ir,
- int64_t step,
- Lincs* lincsd,
- const t_mdatoms& md,
- const t_commrec* cr,
- const gmx_multisim_t* ms,
- const rvec* x,
- rvec* xprime,
- rvec* min_proj,
- const matrix box,
- t_pbc* pbc,
- real lambda,
- real* dvdlambda,
- real invdt,
- rvec* v,
- bool bCalcVir,
- tensor vir_r_m_dr,
- ConstraintVariable econq,
- t_nrnb* nrnb,
- int maxwarn,
- int* warncount)
+bool constrain_lincs(bool computeRmsd,
+ const t_inputrec& ir,
+ int64_t step,
+ Lincs* lincsd,
+ const t_mdatoms& md,
+ const t_commrec* cr,
+ const gmx_multisim_t* ms,
+ ArrayRefWithPadding<const RVec> xPadded,
+ ArrayRefWithPadding<RVec> xprimePadded,
+ ArrayRef<RVec> min_proj,
+ const matrix box,
+ t_pbc* pbc,
+ real lambda,
+ real* dvdlambda,
+ real invdt,
+ ArrayRef<RVec> v,
+ bool bCalcVir,
+ tensor vir_r_m_dr,
+ ConstraintVariable econq,
+ t_nrnb* nrnb,
+ int maxwarn,
+ int* warncount)
{
bool bOK = TRUE;
return bOK;
}
+ ArrayRef<const RVec> x = xPadded.unpaddedArrayRef();
+ ArrayRef<RVec> xprime = xprimePadded.unpaddedArrayRef();
+
if (econq == ConstraintVariable::Positions)
{
/* We can't use bCalcDHDL here, since NULL can be passed for dvdlambda
clear_mat(lincsd->task[th].vir_r_m_dr);
- do_lincs(x, xprime, box, pbc, lincsd, th, md.invmass, cr, bCalcDHDL, ir.LincsWarnAngle,
- &bWarn, invdt, v, bCalcVir, th == 0 ? vir_r_m_dr : lincsd->task[th].vir_r_m_dr);
+ do_lincs(xPadded, xprimePadded, box, pbc, lincsd, th, md.invmass, cr, bCalcDHDL,
+ ir.LincsWarnAngle, &bWarn, invdt, v, bCalcVir,
+ th == 0 ? vir_r_m_dr : lincsd->task[th].vir_r_m_dr);
}
GMX_CATCH_ALL_AND_EXIT_WITH_FATAL_ERROR
}
{
int th = gmx_omp_get_thread_num();
- do_lincsp(x, xprime, min_proj, pbc, lincsd, th, md.invmass, econq, bCalcDHDL,
- bCalcVir, th == 0 ? vir_r_m_dr : lincsd->task[th].vir_r_m_dr);
+ do_lincsp(xPadded, xprimePadded, min_proj, pbc, lincsd, th, md.invmass, econq,
+ bCalcDHDL, bCalcVir, th == 0 ? vir_r_m_dr : lincsd->task[th].vir_r_m_dr);
}
GMX_CATCH_ALL_AND_EXIT_WITH_FATAL_ERROR
}
{
inc_nrnb(nrnb, eNR_LINCSMAT, lincsd->nOrder * lincsd->ncc_triangle);
}
- if (v)
+ if (!v.empty())
{
inc_nrnb(nrnb, eNR_CONSTR_V, lincsd->nc_real * 2);
}
/*
* This file is part of the GROMACS molecular simulation package.
*
- * Copyright (c) 2018,2019, by the GROMACS development team, led by
+ * Copyright (c) 2018,2019,2020, 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 <cstdio>
+#include "gromacs/math/arrayrefwithpadding.h"
#include "gromacs/math/vectypes.h"
#include "gromacs/utility/arrayref.h"
#include "gromacs/utility/basedefinitions.h"
/*! \brief Applies LINCS constraints.
*
* \returns true if the constraining succeeded. */
-bool constrain_lincs(bool computeRmsd,
- const t_inputrec& ir,
- int64_t step,
- Lincs* lincsd,
- const t_mdatoms& md,
- const t_commrec* cr,
- const gmx_multisim_t* ms,
- const rvec* x,
- rvec* xprime,
- rvec* min_proj,
- const matrix box,
- t_pbc* pbc,
- real lambda,
- real* dvdlambda,
- real invdt,
- rvec* v,
- bool bCalcVir,
- tensor vir_r_m_dr,
- ConstraintVariable econq,
- t_nrnb* nrnb,
- int maxwarn,
- int* warncount);
+bool constrain_lincs(bool computeRmsd,
+ const t_inputrec& ir,
+ int64_t step,
+ Lincs* lincsd,
+ const t_mdatoms& md,
+ const t_commrec* cr,
+ const gmx_multisim_t* ms,
+ ArrayRefWithPadding<const RVec> x,
+ ArrayRefWithPadding<RVec> xprime,
+ ArrayRef<RVec> min_proj,
+ const matrix box,
+ t_pbc* pbc,
+ real lambda,
+ real* dvdlambda,
+ real invdt,
+ ArrayRef<RVec> v,
+ bool bCalcVir,
+ tensor vir_r_m_dr,
+ ConstraintVariable econq,
+ t_nrnb* nrnb,
+ int maxwarn,
+ int* warncount);
} // namespace gmx
}
}
-void settle_proj(settledata* settled,
- ConstraintVariable econq,
- int nsettle,
- const t_iatom iatoms[],
- const t_pbc* pbc,
- const rvec x[],
- rvec* der,
- rvec* derp,
- int calcvir_atom_end,
- tensor vir_r_m_dder)
+void settle_proj(settledata* settled,
+ ConstraintVariable econq,
+ int nsettle,
+ const t_iatom iatoms[],
+ const t_pbc* pbc,
+ ArrayRef<const RVec> x,
+ ArrayRef<RVec> der,
+ ArrayRef<RVec> derp,
+ int calcvir_atom_end,
+ tensor vir_r_m_dder)
{
/* Settle for projection out constraint components
* of derivatives of the coordinates.
}
}
-void csettle(settledata* settled,
- int nthread,
- int thread,
- const t_pbc* pbc,
- const real x[],
- real xprime[],
- real invdt,
- real* v,
- bool bCalcVirial,
- tensor vir_r_m_dr,
- bool* bErrorHasOccurred)
+void csettle(settledata* settled,
+ int nthread,
+ int thread,
+ const t_pbc* pbc,
+ ArrayRefWithPadding<const RVec> x,
+ ArrayRefWithPadding<RVec> xprime,
+ real invdt,
+ ArrayRefWithPadding<RVec> v,
+ bool bCalcVirial,
+ tensor vir_r_m_dr,
+ bool* bErrorHasOccurred)
{
+ const real* xPtr = as_rvec_array(x.paddedArrayRef().data())[0];
+ real* xprimePtr = as_rvec_array(xprime.paddedArrayRef().data())[0];
+ real* vPtr = as_rvec_array(v.paddedArrayRef().data())[0];
+
#if GMX_SIMD_HAVE_REAL
if (settled->bUseSimd)
{
set_pbc_simd(pbc, pbcSimd);
settleTemplateWrapper<SimdReal, SimdBool, GMX_SIMD_REAL_WIDTH, const real*>(
- settled, nthread, thread, pbcSimd, x, xprime, invdt, v, bCalcVirial, vir_r_m_dr,
- bErrorHasOccurred);
+ settled, nthread, thread, pbcSimd, xPtr, xprimePtr, invdt, vPtr, bCalcVirial,
+ vir_r_m_dr, bErrorHasOccurred);
}
else
#endif
pbcNonNull = &pbcNo;
}
- settleTemplateWrapper<real, bool, 1, const t_pbc*>(settled, nthread, thread, pbcNonNull, x,
- xprime, invdt, v, bCalcVirial,
- vir_r_m_dr, bErrorHasOccurred);
+ settleTemplateWrapper<real, bool, 1, const t_pbc*>(settled, nthread, thread, pbcNonNull,
+ &xPtr[0], &xprimePtr[0], invdt, &vPtr[0],
+ bCalcVirial, vir_r_m_dr, bErrorHasOccurred);
}
}
/*
* This file is part of the GROMACS molecular simulation package.
*
- * Copyright (c) 2018,2019, by the GROMACS development team, led by
+ * Copyright (c) 2018,2019,2020, 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.
#ifndef GMX_MDLIB_SETTLE_H
#define GMX_MDLIB_SETTLE_H
+#include "gromacs/math/arrayrefwithpadding.h"
#include "gromacs/topology/idef.h"
struct gmx_cmap_t;
/*! \brief Constrain coordinates using SETTLE.
* Can be called on any number of threads.
*/
-void csettle(settledata* settled, /* The SETTLE structure */
- int nthread, /* The number of threads used */
- int thread, /* Our thread index */
- const t_pbc* pbc, /* PBC data pointer, can be NULL */
- const real x[], /* Reference coordinates */
- real xprime[], /* New coords, to be settled */
- real invdt, /* 1/delta_t */
- real* v, /* Also constrain v if v!=NULL */
- bool bCalcVirial, /* Calculate the virial contribution */
- tensor vir_r_m_dr, /* sum r x m delta_r */
- bool* bErrorHasOccurred /* True if a settle error occurred */
+void csettle(settledata* settled, /* The SETTLE structure */
+ int nthread, /* The number of threads used */
+ int thread, /* Our thread index */
+ const t_pbc* pbc, /* PBC data pointer, can be NULL */
+ ArrayRefWithPadding<const RVec> x, /* Reference coordinates */
+ ArrayRefWithPadding<RVec> xprime, /* New coords, to be settled */
+ real invdt, /* 1/delta_t */
+ ArrayRefWithPadding<RVec> v, /* Also constrain v if v!=NULL */
+ bool bCalcVirial, /* Calculate the virial contribution */
+ tensor vir_r_m_dr, /* sum r x m delta_r */
+ bool* bErrorHasOccurred /* True if a settle error occurred */
);
/*! \brief Analytical algorithm to subtract the components of derivatives
* of coordinates working on settle type constraint.
*/
-void settle_proj(settledata* settled,
- ConstraintVariable econq,
- int nsettle,
- const t_iatom iatoms[],
- const t_pbc* pbc, /* PBC data pointer, can be NULL */
- const rvec x[],
- rvec* der,
- rvec* derp,
- int CalcVirAtomEnd,
- tensor vir_r_m_dder);
+void settle_proj(settledata* settled,
+ ConstraintVariable econq,
+ int nsettle,
+ const t_iatom iatoms[],
+ const t_pbc* pbc, /* PBC data pointer, can be NULL */
+ ArrayRef<const RVec> x,
+ ArrayRef<RVec> der,
+ ArrayRef<RVec> derp,
+ int CalcVirAtomEnd,
+ tensor vir_r_m_dder);
} // namespace gmx
// Evaluate constraints
bool success = constrain_lincs(
false, testData->ir_, 0, lincsd, testData->md_, &testData->cr_, &testData->ms_,
- as_rvec_array(testData->x_.data()), as_rvec_array(testData->xPrime_.data()),
- as_rvec_array(testData->xPrime2_.data()), pbc.box, &pbc, testData->md_.lambda,
- &testData->dHdLambda_, testData->invdt_, as_rvec_array(testData->v_.data()),
- testData->computeVirial_, testData->virialScaled_, gmx::ConstraintVariable::Positions,
- &testData->nrnb_, maxwarn, &warncount_lincs);
+ testData->x_.arrayRefWithPadding(), testData->xPrime_.arrayRefWithPadding(),
+ testData->xPrime2_.arrayRefWithPadding().unpaddedArrayRef(), pbc.box, &pbc,
+ testData->md_.lambda, &testData->dHdLambda_, testData->invdt_,
+ testData->v_.arrayRefWithPadding().unpaddedArrayRef(), testData->computeVirial_,
+ testData->virialScaled_, gmx::ConstraintVariable::Positions, &testData->nrnb_, maxwarn,
+ &warncount_lincs);
EXPECT_TRUE(success) << "Test failed with a false return value in LINCS.";
EXPECT_EQ(warncount_lincs, 0) << "There were warnings in LINCS.";
done_lincs(lincsd);
/*
* This file is part of the GROMACS molecular simulation package.
*
- * Copyright (c) 2018,2019, by the GROMACS development team, led by
+ * Copyright (c) 2018,2019,2020, 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.
bool errorOccured;
int numThreads = 1;
int threadIndex = 0;
- csettle(settled, numThreads, threadIndex, &pbc,
- reinterpret_cast<real*>(as_rvec_array(testData->x_.data())),
- reinterpret_cast<real*>(as_rvec_array(testData->xPrime_.data())), testData->reciprocalTimeStep_,
- updateVelocities ? reinterpret_cast<real*>(as_rvec_array(testData->v_.data())) : nullptr,
+ csettle(settled, numThreads, threadIndex, &pbc, testData->x_.arrayRefWithPadding(),
+ testData->xPrime_.arrayRefWithPadding(), testData->reciprocalTimeStep_,
+ updateVelocities ? testData->v_.arrayRefWithPadding() : ArrayRefWithPadding<RVec>(),
calcVirial, testData->virial_, &errorOccured);
settle_free(settled);
EXPECT_FALSE(errorOccured) << testDescription;
clear_mat(vir_part);
/* Constrain the coordinates upd->xp */
- constr->apply(do_log, do_ene, step, 1, 1.0, state->x.rvec_array(), state->v.rvec_array(),
- state->v.rvec_array(), state->box, state->lambda[efptBONDED], dvdlambda,
- nullptr, bCalcVir ? &vir_con : nullptr, ConstraintVariable::Velocities);
+ constr->apply(do_log, do_ene, step, 1, 1.0, state->x.arrayRefWithPadding(),
+ state->v.arrayRefWithPadding(), state->v.arrayRefWithPadding().unpaddedArrayRef(),
+ state->box, state->lambda[efptBONDED], dvdlambda, ArrayRefWithPadding<RVec>(),
+ bCalcVir ? &vir_con : nullptr, ConstraintVariable::Velocities);
if (bCalcVir)
{
clear_mat(vir_part);
/* Constrain the coordinates upd->xp */
- constr->apply(do_log, do_ene, step, 1, 1.0, state->x.rvec_array(), upd->xp()->rvec_array(),
- nullptr, state->box, state->lambda[efptBONDED], dvdlambda,
- as_rvec_array(state->v.data()), bCalcVir ? &vir_con : nullptr,
- ConstraintVariable::Positions);
+ constr->apply(do_log, do_ene, step, 1, 1.0, state->x.arrayRefWithPadding(),
+ upd->xp()->arrayRefWithPadding(), ArrayRef<RVec>(), state->box,
+ state->lambda[efptBONDED], dvdlambda, state->v.arrayRefWithPadding(),
+ bCalcVir ? &vir_con : nullptr, ConstraintVariable::Positions);
if (bCalcVir)
{
wallcycle_stop(wcycle, ewcUPDATE);
/* Constrain the coordinates upd->xp for half a time step */
- constr->apply(do_log, do_ene, step, 1, 0.5, state->x.rvec_array(), upd->xp()->rvec_array(),
- nullptr, state->box, state->lambda[efptBONDED], dvdlambda,
- as_rvec_array(state->v.data()), nullptr, ConstraintVariable::Positions);
+ constr->apply(do_log, do_ene, step, 1, 0.5, state->x.arrayRefWithPadding(),
+ upd->xp()->arrayRefWithPadding(), ArrayRef<RVec>(), state->box,
+ state->lambda[efptBONDED], dvdlambda, state->v.arrayRefWithPadding(), nullptr,
+ ConstraintVariable::Positions);
}
}
#include "legacysimulator.h"
#include "shellfc.h"
+using gmx::ArrayRef;
using gmx::MdrunScheduleWorkload;
+using gmx::RVec;
//! Utility structure for manipulating states during EM
typedef struct
{
/* Constrain the starting coordinates */
dvdl_constr = 0;
- constr->apply(TRUE, TRUE, -1, 0, 1.0, ems->s.x.rvec_array(), ems->s.x.rvec_array(),
- nullptr, ems->s.box, ems->s.lambda[efptFEP], &dvdl_constr, nullptr,
+ constr->apply(TRUE, TRUE, -1, 0, 1.0, ems->s.x.arrayRefWithPadding(),
+ ems->s.x.arrayRefWithPadding(), ArrayRef<RVec>(), ems->s.box,
+ ems->s.lambda[efptFEP], &dvdl_constr, gmx::ArrayRefWithPadding<RVec>(),
nullptr, gmx::ConstraintVariable::Positions);
}
}
if (constr)
{
dvdl_constr = 0;
- validStep = constr->apply(TRUE, TRUE, count, 0, 1.0, s1->x.rvec_array(), s2->x.rvec_array(),
- nullptr, s2->box, s2->lambda[efptBONDED], &dvdl_constr, nullptr,
- nullptr, gmx::ConstraintVariable::Positions);
+ validStep = constr->apply(
+ TRUE, TRUE, count, 0, 1.0, s1->x.arrayRefWithPadding(), s2->x.arrayRefWithPadding(),
+ ArrayRef<RVec>(), s2->box, s2->lambda[efptBONDED], &dvdl_constr,
+ gmx::ArrayRefWithPadding<RVec>(), nullptr, gmx::ConstraintVariable::Positions);
if (cr->nnodes > 1)
{
if (constr)
{
/* Project out the constraint components of the force */
- dvdl_constr = 0;
- rvec* f_rvec = ems->f.rvec_array();
- constr->apply(FALSE, FALSE, count, 0, 1.0, ems->s.x.rvec_array(), f_rvec, f_rvec,
- ems->s.box, ems->s.lambda[efptBONDED], &dvdl_constr, nullptr, &shake_vir,
- gmx::ConstraintVariable::ForceDispl);
+ dvdl_constr = 0;
+ auto f = ems->f.arrayRefWithPadding();
+ constr->apply(FALSE, FALSE, count, 0, 1.0, ems->s.x.arrayRefWithPadding(), f,
+ f.unpaddedArrayRef(), ems->s.box, ems->s.lambda[efptBONDED], &dvdl_constr,
+ gmx::ArrayRefWithPadding<RVec>(), &shake_vir, gmx::ConstraintVariable::ForceDispl);
enerd->term[F_DVDL_CONSTR] += dvdl_constr;
m_add(force_vir, shake_vir, vir);
}
/* Flexible constraint working data */
std::vector<RVec> acc_dir; /* Acceleration direction for flexcon */
gmx::PaddedVector<RVec> x_old; /* Old coordinates for flexcon */
- std::vector<RVec> adir_xnold; /* Work space for init_adir */
- std::vector<RVec> adir_xnew; /* Work space for init_adir */
+ gmx::PaddedVector<RVec> adir_xnold; /* Work space for init_adir */
+ gmx::PaddedVector<RVec> adir_xnew; /* Work space for init_adir */
std::int64_t numForceEvaluations; /* Total number of force evaluations */
int numConvergedIterations; /* Total number of iterations that converged */
};
{
n = end;
}
- shfc->adir_xnold.resize(n);
- shfc->adir_xnew.resize(n);
+ shfc->adir_xnold.resizeWithPadding(n);
+ shfc->adir_xnew.resizeWithPadding(n);
rvec* xnold = as_rvec_array(shfc->adir_xnold.data());
rvec* xnew = as_rvec_array(shfc->adir_xnew.data());
rvec* x_old = as_rvec_array(xOld.paddedArrayRef().data());
}
}
}
- constr->apply(FALSE, FALSE, step, 0, 1.0, x, xnold, nullptr, box, lambda[efptBONDED],
- &(dvdlambda[efptBONDED]), nullptr, nullptr, gmx::ConstraintVariable::Positions);
- constr->apply(FALSE, FALSE, step, 0, 1.0, x, xnew, nullptr, box, lambda[efptBONDED],
- &(dvdlambda[efptBONDED]), nullptr, nullptr, gmx::ConstraintVariable::Positions);
+ constr->apply(FALSE, FALSE, step, 0, 1.0, xCurrent, shfc->adir_xnold.arrayRefWithPadding(),
+ ArrayRef<RVec>(), box, lambda[efptBONDED], &(dvdlambda[efptBONDED]),
+ ArrayRefWithPadding<RVec>(), nullptr, gmx::ConstraintVariable::Positions);
+ constr->apply(FALSE, FALSE, step, 0, 1.0, xCurrent, shfc->adir_xnew.arrayRefWithPadding(),
+ ArrayRef<RVec>(), box, lambda[efptBONDED], &(dvdlambda[efptBONDED]),
+ ArrayRefWithPadding<RVec>(), nullptr, gmx::ConstraintVariable::Positions);
for (n = 0; n < end; n++)
{
}
/* Project the acceleration on the old bond directions */
- constr->apply(FALSE, FALSE, step, 0, 1.0, x_old, xnew, as_rvec_array(acc_dir.data()), box,
- lambda[efptBONDED], &(dvdlambda[efptBONDED]), nullptr, nullptr,
- gmx::ConstraintVariable::Deriv_FlexCon);
+ constr->apply(FALSE, FALSE, step, 0, 1.0, xOld, shfc->adir_xnew.arrayRefWithPadding(), acc_dir,
+ box, lambda[efptBONDED], &(dvdlambda[efptBONDED]), ArrayRefWithPadding<RVec>(),
+ nullptr, gmx::ConstraintVariable::Deriv_FlexCon);
}
void relax_shell_flexcon(FILE* fplog,
/*
* This file is part of the GROMACS molecular simulation package.
*
- * Copyright (c) 2019, by the GROMACS development team, led by
+ * Copyright (c) 2019,2020, 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.
{
tensor vir_con;
- rvec *x, *xprime, *min_proj, *v;
+ ArrayRefWithPadding<RVec> x;
+ ArrayRefWithPadding<RVec> xprime;
+ ArrayRef<RVec> min_proj;
+ ArrayRefWithPadding<RVec> v;
// disabled functionality
real lambda = 0;
switch (variable)
{
case ConstraintVariable::Positions:
- x = as_rvec_array(statePropagatorData_->previousPositionsView().paddedArrayRef().data());
- xprime = as_rvec_array(statePropagatorData_->positionsView().paddedArrayRef().data());
- min_proj = nullptr;
- v = as_rvec_array(statePropagatorData_->velocitiesView().paddedArrayRef().data());
+ x = statePropagatorData_->previousPositionsView();
+ xprime = statePropagatorData_->positionsView();
+ v = statePropagatorData_->velocitiesView();
break;
case ConstraintVariable::Velocities:
- x = as_rvec_array(statePropagatorData_->positionsView().paddedArrayRef().data());
- xprime = as_rvec_array(statePropagatorData_->velocitiesView().paddedArrayRef().data());
- min_proj = as_rvec_array(statePropagatorData_->velocitiesView().paddedArrayRef().data());
- v = nullptr;
+ x = statePropagatorData_->positionsView();
+ xprime = statePropagatorData_->velocitiesView();
+ min_proj = statePropagatorData_->velocitiesView().unpaddedArrayRef();
break;
default: gmx_fatal(FARGS, "Constraint algorithm not implemented for modular simulator.");
}