From: Pascal Merz Date: Thu, 5 Sep 2019 02:27:31 +0000 (-0600) Subject: Make do_constrain_first() independent of t_state X-Git-Url: http://biod.pnpi.spb.ru/gitweb/?a=commitdiff_plain;h=4b698786b3acfb0e4a4cf33d8622167f0a6eca06;p=alexxy%2Fgromacs.git Make do_constrain_first() independent of t_state This changes the do_constrain_first() function to directly take the required state data (natoms, x, v, box, lambda[efptBONDED]) as input instead of a pointer to the full t_state object. This makes subsequent changes to the t_state object easier. This also adds some trivial const qualifier for the box matrix in the constraining functions. Change-Id: I28b58df45481549cd93076334b0778e23e228154 --- diff --git a/src/gromacs/domdec/domdec.h b/src/gromacs/domdec/domdec.h index 948d6ec556..8c74c4204d 100644 --- a/src/gromacs/domdec/domdec.h +++ b/src/gromacs/domdec/domdec.h @@ -244,7 +244,7 @@ void dd_move_f_vsites(struct gmx_domdec_t *dd, rvec *f, rvec *fshift); 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, matrix box, +void dd_move_x_constraints(struct gmx_domdec_t *dd, const matrix box, rvec *x0, rvec *x1, gmx_bool bX1IsCoord); /*! \brief Communicates the coordinates involved in virtual sites */ diff --git a/src/gromacs/domdec/domdec_constraints.cpp b/src/gromacs/domdec/domdec_constraints.cpp index 5e3d382ca1..8005f7a8fd 100644 --- a/src/gromacs/domdec/domdec_constraints.cpp +++ b/src/gromacs/domdec/domdec_constraints.cpp @@ -99,7 +99,7 @@ struct gmx_domdec_constraints_t //! @endcond }; -void dd_move_x_constraints(gmx_domdec_t *dd, matrix box, +void dd_move_x_constraints(gmx_domdec_t *dd, const matrix box, rvec *x0, rvec *x1, gmx_bool bX1IsCoord) { if (dd->constraint_comm) diff --git a/src/gromacs/essentialdynamics/edsam.cpp b/src/gromacs/essentialdynamics/edsam.cpp index f4b2b2d079..2b9fd374a6 100644 --- a/src/gromacs/essentialdynamics/edsam.cpp +++ b/src/gromacs/essentialdynamics/edsam.cpp @@ -1946,7 +1946,7 @@ void dd_make_local_ed_indices(gmx_domdec_t *dd, struct gmx_edsam *ed) } -static inline void ed_unshift_single_coord(matrix box, const rvec x, const ivec is, rvec xu) +static inline void ed_unshift_single_coord(const matrix box, const rvec x, const ivec is, rvec xu) { int tx, ty, tz; @@ -3000,7 +3000,7 @@ void do_edsam(const t_inputrec *ir, const t_commrec *cr, rvec xs[], rvec v[], - matrix box, + const matrix box, gmx_edsam *ed) { int i, edinr, iupdate = 500; diff --git a/src/gromacs/essentialdynamics/edsam.h b/src/gromacs/essentialdynamics/edsam.h index 04cf6cdcd7..7fd2d28dc9 100644 --- a/src/gromacs/essentialdynamics/edsam.h +++ b/src/gromacs/essentialdynamics/edsam.h @@ -102,7 +102,7 @@ class MDLogger; * \param ed The essential dynamics data. */ void do_edsam(const t_inputrec *ir, int64_t step, - const t_commrec *cr, rvec xs[], rvec v[], matrix box, gmx_edsam *ed); + const t_commrec *cr, rvec xs[], rvec v[], const matrix box, gmx_edsam *ed); /*! \brief Initializes the essential dynamics and flooding module. diff --git a/src/gromacs/mdlib/constr.cpp b/src/gromacs/mdlib/constr.cpp index 2be5fd3c3a..6181835e4f 100644 --- a/src/gromacs/mdlib/constr.cpp +++ b/src/gromacs/mdlib/constr.cpp @@ -122,7 +122,7 @@ class Constraints::Impl rvec *x, rvec *xprime, rvec *min_proj, - matrix box, + const matrix box, real lambda, real *dvdlambda, rvec *v, @@ -238,7 +238,7 @@ void too_many_constraint_warnings(int eConstrAlg, int warncount) 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[], matrix box) + const rvec x[], const matrix box) { char fname[STRLEN]; FILE *out; @@ -295,7 +295,7 @@ static void write_constr_pdb(const char *fn, const char *title, //! 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[], matrix box) + const rvec x[], rvec xprime[], const matrix box) { char buf[STRLEN], buf2[22]; @@ -327,7 +327,7 @@ Constraints::apply(bool bLog, rvec *x, rvec *xprime, rvec *min_proj, - matrix box, + const matrix box, real lambda, real *dvdlambda, rvec *v, @@ -359,7 +359,7 @@ Constraints::Impl::apply(bool bLog, rvec *x, rvec *xprime, rvec *min_proj, - matrix box, + const matrix box, real lambda, real *dvdlambda, rvec *v, @@ -1273,7 +1273,11 @@ bool inter_charge_group_settles(const gmx_mtop_t &mtop) void do_constrain_first(FILE *fplog, gmx::Constraints *constr, const t_inputrec *ir, const t_mdatoms *md, - t_state *state) + int natoms, + ArrayRefWithPadding x, + ArrayRefWithPadding v, + const matrix box, + real lambda) { int i, m, start, end; int64_t step; @@ -1281,10 +1285,13 @@ void do_constrain_first(FILE *fplog, gmx::Constraints *constr, 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, state->natoms + 1); + snew(savex, natoms + 1); start = 0; end = md->homenr; @@ -1307,9 +1314,9 @@ void do_constrain_first(FILE *fplog, gmx::Constraints *constr, /* constrain the current position */ constr->apply(TRUE, FALSE, step, 0, 1.0, - state->x.rvec_array(), state->x.rvec_array(), nullptr, - state->box, - state->lambda[efptBONDED], &dvdl_dum, + xRvec, xRvec, nullptr, + box, + lambda, &dvdl_dum, nullptr, nullptr, gmx::ConstraintVariable::Positions); if (EI_VV(ir->eI)) { @@ -1317,24 +1324,24 @@ void do_constrain_first(FILE *fplog, gmx::Constraints *constr, /* also may be useful if we need the ekin from the halfstep for velocity verlet */ constr->apply(TRUE, FALSE, step, 0, 1.0, - state->x.rvec_array(), state->v.rvec_array(), state->v.rvec_array(), - state->box, - state->lambda[efptBONDED], &dvdl_dum, + xRvec, vRvec, vRvec, + box, + lambda, &dvdl_dum, nullptr, nullptr, gmx::ConstraintVariable::Velocities); } /* constrain the inital velocities at t-dt/2 */ if (EI_STATE_VELOCITY(ir->eI) && ir->eI != eiVV) { - auto x = makeArrayRef(state->x).subArray(start, end); - auto v = makeArrayRef(state->v).subArray(start, end); + auto subX = x.paddedArrayRef().subArray(start, end); + auto subV = v.paddedArrayRef().subArray(start, end); for (i = start; (i < end); i++) { for (m = 0; (m < DIM); m++) { /* Reverse the velocity */ - v[i][m] = -v[i][m]; + subV[i][m] = -subV[i][m]; /* Store the position at t-dt in buf */ - savex[i][m] = x[i][m] + dt*v[i][m]; + savex[i][m] = subX[i][m] + dt*subV[i][m]; } } /* Shake the positions at t=-dt with the positions at t=0 @@ -1349,17 +1356,17 @@ void do_constrain_first(FILE *fplog, gmx::Constraints *constr, dvdl_dum = 0; constr->apply(TRUE, FALSE, step, -1, 1.0, - state->x.rvec_array(), savex, nullptr, - state->box, - state->lambda[efptBONDED], &dvdl_dum, - state->v.rvec_array(), nullptr, gmx::ConstraintVariable::Positions); + xRvec, savex, nullptr, + box, + lambda, &dvdl_dum, + vRvec, nullptr, gmx::ConstraintVariable::Positions); for (i = start; i < end; i++) { for (m = 0; m < DIM; m++) { /* Re-reverse the velocities */ - v[i][m] = -v[i][m]; + subV[i][m] = -subV[i][m]; } } } diff --git a/src/gromacs/mdlib/constr.h b/src/gromacs/mdlib/constr.h index 452d2ac36d..8fd7f3003d 100644 --- a/src/gromacs/mdlib/constr.h +++ b/src/gromacs/mdlib/constr.h @@ -74,6 +74,7 @@ class t_state; namespace gmx { +template class ArrayRefWithPadding; //! Describes supported flavours of constrained updates. enum class ConstraintVariable : int @@ -172,7 +173,7 @@ class Constraints rvec *x, rvec *xprime, rvec *min_proj, - matrix box, + const matrix box, real lambda, real *dvdlambda, rvec *v, @@ -295,7 +296,11 @@ bool inter_charge_group_settles(const gmx_mtop_t &mtop); /*! \brief Constrain the initial coordinates and velocities */ void do_constrain_first(FILE *log, gmx::Constraints *constr, const t_inputrec *inputrec, const t_mdatoms *md, - t_state *state); + int natoms, + ArrayRefWithPadding x, + ArrayRefWithPadding v, + const matrix box, + real lambda); } // namespace gmx diff --git a/src/gromacs/mdlib/lincs.cpp b/src/gromacs/mdlib/lincs.cpp index 0b6243a79a..149c87590c 100644 --- a/src/gromacs/mdlib/lincs.cpp +++ b/src/gromacs/mdlib/lincs.cpp @@ -952,7 +952,7 @@ calc_dist_iter_simd(int b0, #endif // GMX_SIMD_HAVE_REAL //! Implements LINCS constraining. -static void do_lincs(const rvec *x, rvec *xp, matrix box, t_pbc *pbc, +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, @@ -2296,7 +2296,7 @@ bool constrain_lincs(bool computeRmsd, const t_commrec *cr, const gmx_multisim_t *ms, const rvec *x, rvec *xprime, rvec *min_proj, - matrix box, t_pbc *pbc, + const matrix box, t_pbc *pbc, real lambda, real *dvdlambda, real invdt, rvec *v, bool bCalcVir, tensor vir_r_m_dr, diff --git a/src/gromacs/mdlib/lincs.h b/src/gromacs/mdlib/lincs.h index c04741e45f..3405df1648 100644 --- a/src/gromacs/mdlib/lincs.h +++ b/src/gromacs/mdlib/lincs.h @@ -1,7 +1,7 @@ /* * This file is part of the GROMACS molecular simulation package. * - * Copyright (c) 2018, by the GROMACS development team, led by + * Copyright (c) 2018,2019, 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. @@ -99,7 +99,7 @@ constrain_lincs(bool computeRmsd, const t_commrec *cr, const gmx_multisim_t *ms, const rvec *x, rvec *xprime, rvec *min_proj, - matrix box, t_pbc *pbc, + const matrix box, t_pbc *pbc, real lambda, real *dvdlambda, real invdt, rvec *v, bool bCalcVir, tensor vir_r_m_dr, diff --git a/src/gromacs/mdrun/md.cpp b/src/gromacs/mdrun/md.cpp index 18dc39041d..9e33cc02dd 100644 --- a/src/gromacs/mdrun/md.cpp +++ b/src/gromacs/mdrun/md.cpp @@ -454,7 +454,11 @@ void gmx::LegacySimulator::do_md() if (constr) { /* Constrain the initial coordinates and velocities */ - do_constrain_first(fplog, constr, ir, mdatoms, state); + do_constrain_first(fplog, constr, ir, mdatoms, + state->natoms, + state->x.arrayRefWithPadding(), + state->v.arrayRefWithPadding(), + state->box, state->lambda[efptBONDED]); } if (vsite) {