* 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/gmxlib/nrnb.h"
#include "gromacs/linearalgebra/nrjac.h"
#include "gromacs/math/functions.h"
+#include "gromacs/math/units.h"
#include "gromacs/math/utilities.h"
#include "gromacs/math/vec.h"
#include "gromacs/math/vectypes.h"
rad += gmx::square((vec.refproj[i] - vec.xproj[i]));
}
- return rad = sqrt(rad);
+ return sqrt(rad);
}
} // namespace
}
-static void do_single_flood(FILE* edo,
- const rvec x[],
- rvec force[],
- t_edpar* edi,
- int64_t step,
- const matrix box,
- const t_commrec* cr,
- gmx_bool bNS) /* Are we in a neighbor searching step? */
+static void do_single_flood(FILE* edo,
+ gmx::ArrayRef<const gmx::RVec> coords,
+ gmx::ArrayRef<gmx::RVec> force,
+ t_edpar* edi,
+ int64_t step,
+ const matrix box,
+ const t_commrec* cr,
+ gmx_bool bNS) /* Are we in a neighbor searching step? */
{
int i;
matrix rotmat; /* rotation matrix */
buf->shifts_xcoll,
buf->extra_shifts_xcoll,
bNS,
- x,
+ as_rvec_array(coords.data()),
edi->sav.nr,
edi->sav.nr_loc,
edi->sav.anrs_loc,
buf->shifts_xc_ref,
buf->extra_shifts_xc_ref,
bNS,
- x,
+ as_rvec_array(coords.data()),
edi->sref.nr,
edi->sref.nr_loc,
edi->sref.anrs_loc,
/* Main flooding routine, called from do_force */
-extern void do_flood(const t_commrec* cr,
- const t_inputrec* ir,
- const rvec x[],
- rvec force[],
- gmx_edsam* ed,
- const matrix box,
- int64_t step,
- gmx_bool bNS)
+void do_flood(const t_commrec* cr,
+ const t_inputrec& ir,
+ gmx::ArrayRef<const gmx::RVec> coords,
+ gmx::ArrayRef<gmx::RVec> force,
+ gmx_edsam* ed,
+ const matrix box,
+ int64_t step,
+ bool bNS)
{
/* Write time to edo, when required. Output the time anyhow since we need
* it in the output file for ED constraints. */
if (MASTER(cr) && do_per_step(step, ed->edpar.begin()->outfrq))
{
- fprintf(ed->edo, "\n%12f", ir->init_t + step * ir->delta_t);
+ fprintf(ed->edo, "\n%12f", ir.init_t + step * ir.delta_t);
}
if (ed->eEDtype != EssentialDynamicsType::Flooding)
/* Call flooding for one matrix */
if (edi.flood.vecs.neig)
{
- do_single_flood(ed->edo, x, force, &edi, step, box, cr, bNS);
+ do_single_flood(ed->edo, coords, force, &edi, step, box, cr, bNS);
}
}
}
const gmx_output_env_t* oenv,
const t_commrec* cr)
{
- auto edHandle = std::make_unique<gmx::EssentialDynamics>();
- auto ed = edHandle->getLegacyED();
+ auto edHandle = std::make_unique<gmx::EssentialDynamics>();
+ auto* ed = edHandle->getLegacyED();
/* We want to perform ED (this switch might later be upgraded to EssentialDynamicsType::Flooding) */
ed->eEDtype = EssentialDynamicsType::EDSampling;
/* init-routine called for every *.edi-cycle, initialises t_edpar structure */
-static void init_edi(const gmx_mtop_t* mtop, t_edpar* edi)
+static void init_edi(const gmx_mtop_t& mtop, t_edpar* edi)
{
int i;
real totalmass = 0.0;
/* Gets the rms deviation of the positions to the structure s */
/* fit_to_structure has to be called before calling this routine! */
-static real rmsd_from_structure(rvec* x, /* The positions under consideration */
+static real rmsd_from_structure(rvec* x, /* The positions under consideration */
struct gmx_edx* s) /* The structure from which the rmsd shall be computed */
{
real rmsd = 0.0;
std::unique_ptr<gmx::EssentialDynamics> init_edsam(const gmx::MDLogger& mdlog,
const char* ediFileName,
const char* edoFileName,
- const gmx_mtop_t* mtop,
- const t_inputrec* ir,
+ const gmx_mtop_t& mtop,
+ const t_inputrec& ir,
const t_commrec* cr,
gmx::Constraints* constr,
const t_state* globalState,
"gmx grompp and the related .mdp options may change also.");
/* Open input and output files, allocate space for ED data structure */
- auto edHandle = ed_open(mtop->natoms, oh, ediFileName, edoFileName, startingBehavior, oenv, cr);
- auto ed = edHandle->getLegacyED();
+ auto edHandle = ed_open(mtop.natoms, oh, ediFileName, edoFileName, startingBehavior, oenv, cr);
+ auto* ed = edHandle->getLegacyED();
GMX_RELEASE_ASSERT(constr != nullptr, "Must have valid constraints object");
constr->saveEdsamPointer(ed);
for (auto& edi : ed->edpar)
{
init_edi(mtop, &edi);
- init_flood(&edi, ed, ir->delta_t);
+ init_flood(&edi, ed, ir.delta_t);
}
}
if (!EDstate->bFromCpt)
{
/* Remove PBC, make molecule(s) subject to ED whole. */
- snew(x_pbc, mtop->natoms);
- copy_rvecn(globalState->x.rvec_array(), x_pbc, 0, mtop->natoms);
- do_pbc_first_mtop(nullptr, ir->pbcType, globalState->box, mtop, x_pbc);
+ snew(x_pbc, mtop.natoms);
+ copy_rvecn(globalState->x.rvec_array(), x_pbc, 0, mtop.natoms);
+ do_pbc_first_mtop(nullptr, ir.pbcType, globalState->box, &mtop, x_pbc);
}
/* Reset pointer to first ED data set which contains the actual ED data */
auto edi = ed->edpar.begin();
} /* end of MASTER only section */
- if (PAR(cr))
+ if (haveDDAtomOrdering(*cr))
{
- /* Broadcast the essential dynamics / flooding data to all nodes */
+ /* Broadcast the essential dynamics / flooding data to all nodes.
+ * In a single-rank case, only the necessary memory allocation is done. */
broadcast_ed_data(cr, ed);
}
else
{
- /* In the single-CPU case, point the local atom numbers pointers to the global
+ /* In the non-DD case, point the local atom numbers pointers to the global
* one, so that we can use the same notation in serial and parallel case: */
/* Loop over all ED data sets (usually only one, though) */
for (auto edi = ed->edpar.begin(); edi != ed->edpar.end(); ++edi)
}
-void do_edsam(const t_inputrec* ir, int64_t step, const t_commrec* cr, rvec xs[], rvec v[], const matrix box, gmx_edsam* ed)
+void do_edsam(const t_inputrec* ir,
+ int64_t step,
+ const t_commrec* cr,
+ gmx::ArrayRef<gmx::RVec> coords,
+ gmx::ArrayRef<gmx::RVec> velocities,
+ const matrix box,
+ gmx_edsam* ed)
{
int i, edinr, iupdate = 500;
matrix rotmat; /* rotation matrix */
buf->shifts_xcoll,
buf->extra_shifts_xcoll,
PAR(cr) ? buf->bUpdateShifts : TRUE,
- xs,
+ as_rvec_array(coords.data()),
edi.sav.nr,
edi.sav.nr_loc,
edi.sav.anrs_loc,
buf->shifts_xc_ref,
buf->extra_shifts_xc_ref,
PAR(cr) ? buf->bUpdateShifts : TRUE,
- xs,
+ as_rvec_array(coords.data()),
edi.sref.nr,
edi.sref.nr_loc,
edi.sref.anrs_loc,
box, buf->xcoll[edi.sav.c_ind[i]], buf->shifts_xcoll[edi.sav.c_ind[i]], x_unsh);
/* dx is the ED correction to the positions: */
- rvec_sub(x_unsh, xs[edi.sav.anrs_loc[i]], dx);
+ rvec_sub(x_unsh, coords[edi.sav.anrs_loc[i]], dx);
- if (v != nullptr)
+ if (!velocities.empty())
{
/* dv is the ED correction to the velocity: */
svmul(dt_1, dx, dv);
/* apply the velocity correction: */
- rvec_inc(v[edi.sav.anrs_loc[i]], dv);
+ rvec_inc(velocities[edi.sav.anrs_loc[i]], dv);
}
/* Finally apply the position correction due to ED: */
- copy_rvec(x_unsh, xs[edi.sav.anrs_loc[i]]);
+ copy_rvec(x_unsh, coords[edi.sav.anrs_loc[i]]);
}
}
} /* END of if ( bNeedDoEdsam(edi) ) */
biod.pnpi.spb.ru / alexxy / gromacs.git / blobdiff