#include "gromacs/mdtypes/md_enums.h"
#include "gromacs/timing/wallcycle.h"
#include "gromacs/utility/arrayref.h"
-#include "gromacs/utility/basedefinitions.h"
#include "gromacs/utility/real.h"
class ekinstate_t;
struct t_graph;
struct t_grpopts;
struct t_inputrec;
-struct t_mdatoms;
struct t_nrnb;
class t_state;
+enum class ParticleType;
namespace gmx
{
* \returns handle to box deformation class
*/
BoxDeformation* deform() const;
- /*! \brief Resizes buffer that stores intermediate coordinates.
+ /*! \brief Sets data that changes only at domain decomposition time.
*
* \param[in] numAtoms Updated number of atoms.
+ * \param[in] cFREEZE Group index for freezing
+ * \param[in] cTC Group index for center of mass motion removal
+ * \param[in] cAcceleration Group index for constant acceleration groups
*/
- void setNumAtoms(int numAtoms);
+ void updateAfterPartition(int numAtoms,
+ gmx::ArrayRef<const unsigned short> cFREEZE,
+ gmx::ArrayRef<const unsigned short> cTC,
+ gmx::ArrayRef<const unsigned short> cAcceleration);
/*! \brief Perform numerical integration step.
*
* Selects the appropriate integrator, based on the input record and performs a numerical integration step.
*
- * \param[in] inputRecord Input record.
- * \param[in] step Current timestep.
- * \param[in] md MD atoms data.
- * \param[in] state System state object.
- * \param[in] f Buffer with atomic forces for home particles.
- * \param[in] fcdata Force calculation data to update distance and orientation restraints.
- * \param[in] ekind Kinetic energy data (for temperature coupling, energy groups, etc.).
- * \param[in] M Parrinello-Rahman velocity scaling matrix.
- * \param[in] updatePart What should be updated, coordinates or velocities. This enum only used in VV integrator.
- * \param[in] cr Comunication record (Old comment: these shouldn't be here -- need to think about it).
- * \param[in] haveConstraints If the system has constraints.
+ * \param[in] inputRecord Input record.
+ * \param[in] step Current timestep.
+ * \param[in] homenr The number of atoms on this processor.
+ * \param[in] havePartiallyFrozenAtoms Whether atoms are frozen along 1 or 2 (not 3) dimensions?
+ * \param[in] ptype The list of particle types.
+ * \param[in] invMass Inverse atomic mass per atom, 0 for vsites and shells.
+ * \param[in] invMassPerDim Inverse atomic mass per atom and dimension, 0 for vsites, shells and frozen dimensions
+ * \param[in] state System state object.
+ * \param[in] f Buffer with atomic forces for home particles.
+ * \param[in] fcdata Force calculation data to update distance and orientation restraints.
+ * \param[in] ekind Kinetic energy data (for temperature coupling, energy groups, etc.).
+ * \param[in] M Parrinello-Rahman velocity scaling matrix.
+ * \param[in] updatePart What should be updated, coordinates or velocities. This enum only used in VV integrator.
+ * \param[in] cr Comunication record (Old comment: these shouldn't be here -- need to think about it).
+ * \param[in] haveConstraints If the system has constraints.
*/
void update_coords(const t_inputrec& inputRecord,
int64_t step,
- const t_mdatoms* md,
+ int homenr,
+ bool havePartiallyFrozenAtoms,
+ gmx::ArrayRef<const ParticleType> ptype,
+ gmx::ArrayRef<const real> invMass,
+ gmx::ArrayRef<const rvec> invMassPerDim,
t_state* state,
const gmx::ArrayRefWithPadding<const gmx::RVec>& f,
- const t_fcdata& fcdata,
+ t_fcdata* fcdata,
const gmx_ekindata_t* ekind,
const matrix M,
int updatePart,
* Copy the updated coordinates to the main coordinates buffer for the atoms that are not frozen.
*
* \param[in] inputRecord Input record.
- * \param[in] md MD atoms data.
+ * \param[in] havePartiallyFrozenAtoms Whether atoms are frozen along 1 or 2 (not 3) dimensions?
+ * \param[in] homenr The number of atoms on this processor.
* \param[in] state System state object.
* \param[in] wcycle Wall-clock cycle counter.
* \param[in] haveConstraints If the system has constraints.
*/
void finish_update(const t_inputrec& inputRecord,
- const t_mdatoms* md,
+ bool havePartiallyFrozenAtoms,
+ int homenr,
t_state* state,
- gmx_wallcycle_t wcycle,
+ gmx_wallcycle* wcycle,
bool haveConstraints);
/*! \brief Secong part of the SD integrator.
* \param[in] step Current timestep.
* \param[in] dvdlambda Free energy derivative. Contribution to be added to
* the bonded interactions.
- * \param[in] md MD atoms data.
+ * \param[in] homenr The number of atoms on this processor.
+ * \param[in] ptype The list of particle types.
+ * \param[in] invMass Inverse atomic mass per atom, 0 for vsites and shells.
* \param[in] state System state object.
* \param[in] cr Comunication record.
* \param[in] nrnb Cycle counters.
* \param[in] do_log If this is logging step.
* \param[in] do_ene If this is an energy evaluation step.
*/
- void update_sd_second_half(const t_inputrec& inputRecord,
- int64_t step,
- real* dvdlambda,
- const t_mdatoms* md,
- t_state* state,
- const t_commrec* cr,
- t_nrnb* nrnb,
- gmx_wallcycle_t wcycle,
- gmx::Constraints* constr,
- bool do_log,
- bool do_ene);
+ void update_sd_second_half(const t_inputrec& inputRecord,
+ int64_t step,
+ real* dvdlambda,
+ int homenr,
+ gmx::ArrayRef<const ParticleType> ptype,
+ gmx::ArrayRef<const real> invMass,
+ t_state* state,
+ const t_commrec* cr,
+ t_nrnb* nrnb,
+ gmx_wallcycle* wcycle,
+ gmx::Constraints* constr,
+ bool do_log,
+ bool do_ene);
/*! \brief Performs a leap-frog update without updating \p state so the constrain virial
* can be computed.
*/
- void update_for_constraint_virial(const t_inputrec& inputRecord,
- const t_mdatoms& md,
- const t_state& state,
+ void update_for_constraint_virial(const t_inputrec& inputRecord,
+ int homenr,
+ bool havePartiallyFrozenAtoms,
+ gmx::ArrayRef<const real> invmass,
+ gmx::ArrayRef<const rvec> invMassPerDim,
+ const t_state& state,
const gmx::ArrayRefWithPadding<const gmx::RVec>& f,
const gmx_ekindata_t& ekind);