Use ArrayRef in special forces
authorJoe Jordan <ejjordan12@gmail.com>
Thu, 25 Mar 2021 10:48:38 +0000 (10:48 +0000)
committerPaul Bauer <paul.bauer.q@gmail.com>
Thu, 25 Mar 2021 10:48:38 +0000 (10:48 +0000)
Pass ArrayRef<RVec> instead of rvec* in edsam, pull_rotation, and
imd. Also removed unneeded header basedefinitions.

A followup change can refactor communicate_group_positions to take
ArrayRefs.

src/gromacs/essentialdynamics/edsam.cpp
src/gromacs/essentialdynamics/edsam.h
src/gromacs/imd/imd.cpp
src/gromacs/imd/imd.h
src/gromacs/mdlib/constr.cpp
src/gromacs/mdlib/sim_util.cpp
src/gromacs/mdrun/md.cpp
src/gromacs/mdrun/minimize.cpp
src/gromacs/mdrun/runner.cpp
src/gromacs/pulling/pull_rotation.cpp
src/gromacs/pulling/pull_rotation.h

index fa612d238d880bbffe15e6013c87dd036aa32475..7f33f49cbb4e5a5a99c4e39b3138826f6e77f1bc 100644 (file)
@@ -904,14 +904,14 @@ static void update_adaption(t_edpar* edi)
 }
 
 
-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 */
@@ -932,7 +932,7 @@ static void do_single_flood(FILE*            edo,
                                 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,
@@ -948,7 +948,7 @@ static void do_single_flood(FILE*            edo,
                                     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,
@@ -1028,14 +1028,14 @@ static void do_single_flood(FILE*            edo,
 
 
 /* Main flooding routine, called from do_force */
-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,
-              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. */
@@ -1054,7 +1054,7 @@ void do_flood(const t_commrec*  cr,
         /* 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);
         }
     }
 }
@@ -3142,7 +3142,13 @@ std::unique_ptr<gmx::EssentialDynamics> init_edsam(const gmx::MDLogger&        m
 }
 
 
-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 */
@@ -3191,7 +3197,7 @@ void do_edsam(const t_inputrec* ir, int64_t step, const t_commrec* cr, rvec xs[]
                                         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,
@@ -3207,7 +3213,7 @@ void do_edsam(const t_inputrec* ir, int64_t step, const t_commrec* cr, rvec xs[]
                                             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,
@@ -3313,17 +3319,17 @@ void do_edsam(const t_inputrec* ir, int64_t step, const t_commrec* cr, rvec xs[]
                             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) ) */
index ab26dfaf668c8860ffb8707c76935195d3a88c55..957f8c8f0297d73e20db96c39b6f4863bfe65b8b 100644 (file)
@@ -52,7 +52,6 @@
 #include <memory>
 
 #include "gromacs/math/vectypes.h"
-#include "gromacs/utility/basedefinitions.h"
 
 /*! \brief Abstract type for essential dynamics
  *
@@ -72,6 +71,9 @@ namespace gmx
 {
 enum class StartingBehavior;
 class Constraints;
+template<typename>
+class ArrayRef;
+
 class EssentialDynamics
 {
 public:
@@ -97,18 +99,18 @@ class MDLogger;
  * \param ir                MD input parameter record.
  * \param step              Number of the time step.
  * \param cr                Data needed for MPI communication.
- * \param xs                The local positions on this processor.
- * \param v                 The local velocities.
+ * \param coords            The local positions on this processor.
+ * \param velocities        The local velocities.
  * \param box               The simulation box.
  * \param ed                The essential dynamics data.
  */
-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);
 
 
 /*! \brief Initializes the essential dynamics and flooding module.
@@ -153,20 +155,20 @@ void dd_make_local_ed_indices(gmx_domdec_t* dd, gmx_edsam* ed);
  *
  * \param cr                Data needed for MPI communication.
  * \param ir                MD input parameter record.
- * \param x                 Positions on the local processor.
+ * \param coords            Positions on the local processor.
  * \param force             Forcefield forces to which the flooding forces are added.
  * \param ed                The essential dynamics data.
  * \param box               The simulation box.
  * \param step              Number of the time step.
  * \param bNS               Are we in a neighbor searching step?
  */
-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,
-              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);
 
 #endif
index edc6a2fed9196013da8c9034c92fa01442916e5a..c75f8d7ce0503b43978f3c789e4970af9084b02f 100644 (file)
@@ -203,9 +203,9 @@ public:
     /*! \brief Removes shifts of molecules diffused outside of the box. */
     void removeMolecularShifts(const matrix box);
     /*! \brief Initialize arrays used to assemble the positions from the other nodes. */
-    void prepareForPositionAssembly(const t_commrec* cr, const rvec x[]);
+    void prepareForPositionAssembly(const t_commrec* cr, gmx::ArrayRef<const gmx::RVec> coords);
     /*! \brief Interact with any connected VMD session */
-    bool run(int64_t step, bool bNS, const matrix box, const rvec x[], double t);
+    bool run(int64_t step, bool bNS, const matrix box, gmx::ArrayRef<const gmx::RVec> coords, double t);
 
     // TODO rename all the data members to have underscore suffixes
 
@@ -1242,7 +1242,7 @@ void ImdSession::Impl::removeMolecularShifts(const matrix box)
 }
 
 
-void ImdSession::Impl::prepareForPositionAssembly(const t_commrec* cr, const rvec x[])
+void ImdSession::Impl::prepareForPositionAssembly(const t_commrec* cr, gmx::ArrayRef<const gmx::RVec> coords)
 {
     snew(xa, nat);
     snew(xa_ind, nat);
@@ -1257,7 +1257,7 @@ void ImdSession::Impl::prepareForPositionAssembly(const t_commrec* cr, const rve
         for (int i = 0; i < nat; i++)
         {
             int ii = ind[i];
-            copy_rvec(x[ii], xa_old[i]);
+            copy_rvec(coords[ii], xa_old[i]);
         }
     }
 
@@ -1297,19 +1297,19 @@ static void imd_check_integrator_parallel(const t_inputrec* ir, const t_commrec*
     }
 }
 
-std::unique_ptr<ImdSession> makeImdSession(const t_inputrec*           ir,
-                                           const t_commrec*            cr,
-                                           gmx_wallcycle*              wcycle,
-                                           gmx_enerdata_t*             enerd,
-                                           const gmx_multisim_t*       ms,
-                                           const gmx_mtop_t&           top_global,
-                                           const MDLogger&             mdlog,
-                                           const rvec                  x[],
-                                           int                         nfile,
-                                           const t_filenm              fnm[],
-                                           const gmx_output_env_t*     oenv,
-                                           const ImdOptions&           options,
-                                           const gmx::StartingBehavior startingBehavior)
+std::unique_ptr<ImdSession> makeImdSession(const t_inputrec*              ir,
+                                           const t_commrec*               cr,
+                                           gmx_wallcycle*                 wcycle,
+                                           gmx_enerdata_t*                enerd,
+                                           const gmx_multisim_t*          ms,
+                                           const gmx_mtop_t&              top_global,
+                                           const MDLogger&                mdlog,
+                                           gmx::ArrayRef<const gmx::RVec> coords,
+                                           int                            nfile,
+                                           const t_filenm                 fnm[],
+                                           const gmx_output_env_t*        oenv,
+                                           const ImdOptions&              options,
+                                           const gmx::StartingBehavior    startingBehavior)
 {
     std::unique_ptr<ImdSession> session(new ImdSession(mdlog));
     auto                        impl = session->impl_.get();
@@ -1505,7 +1505,7 @@ std::unique_ptr<ImdSession> makeImdSession(const t_inputrec*           ir,
     impl->syncNodes(cr, 0);
 
     /* Initialize arrays used to assemble the positions from the other nodes */
-    impl->prepareForPositionAssembly(cr, x);
+    impl->prepareForPositionAssembly(cr, coords);
 
     /* Initialize molecule blocks to make them whole later...*/
     if (MASTER(cr))
@@ -1517,7 +1517,7 @@ std::unique_ptr<ImdSession> makeImdSession(const t_inputrec*           ir,
 }
 
 
-bool ImdSession::Impl::run(int64_t step, bool bNS, const matrix box, const rvec x[], double t)
+bool ImdSession::Impl::run(int64_t step, bool bNS, const matrix box, gmx::ArrayRef<const gmx::RVec> coords, double t)
 {
     /* IMD at all? */
     if (!sessionPossible)
@@ -1567,7 +1567,7 @@ bool ImdSession::Impl::run(int64_t step, bool bNS, const matrix box, const rvec
         /* Transfer the IMD positions to the master node. Every node contributes
          * its local positions x and stores them in the assembled xa array. */
         communicate_group_positions(
-                cr, xa, xa_shifts, xa_eshifts, true, x, nat, nat_loc, ind_loc, xa_ind, xa_old, box);
+                cr, xa, xa_shifts, xa_eshifts, true, as_rvec_array(coords.data()), nat, nat_loc, ind_loc, xa_ind, xa_old, box);
 
         /* If connected and master -> remove shifts */
         if ((imdstep && bConnected) && MASTER(cr))
@@ -1581,9 +1581,9 @@ bool ImdSession::Impl::run(int64_t step, bool bNS, const matrix box, const rvec
     return imdstep;
 }
 
-bool ImdSession::run(int64_t step, bool bNS, const matrix box, const rvec x[], double t)
+bool ImdSession::run(int64_t step, bool bNS, const matrix box, gmx::ArrayRef<const gmx::RVec> coords, double t)
 {
-    return impl_->run(step, bNS, box, x, t);
+    return impl_->run(step, bNS, box, coords, t);
 }
 
 void ImdSession::fillEnergyRecord(int64_t step, bool bHaveNewEnergies)
@@ -1655,7 +1655,7 @@ void ImdSession::updateEnergyRecordAndSendPositionsAndEnergies(bool bIMDstep, in
     wallcycle_stop(impl_->wcycle, ewcIMD);
 }
 
-void ImdSession::applyForces(rvec* f)
+void ImdSession::applyForces(gmx::ArrayRef<gmx::RVec> force)
 {
     if (!impl_->sessionPossible || !impl_->bForceActivated)
     {
@@ -1677,7 +1677,7 @@ void ImdSession::applyForces(rvec* f)
             j = *locndx;
         }
 
-        rvec_inc(f[j], impl_->f[i]);
+        rvec_inc(force[j], impl_->f[i]);
     }
 
     wallcycle_stop(impl_->wcycle, ewcIMD);
index 6c6517717c0c2f85a10c3818a7c49813669c9974..1d5c4fa56d60bb5e62d3a1574bec68ecf6c113f3 100644 (file)
@@ -67,7 +67,6 @@
 #include <memory>
 
 #include "gromacs/math/vectypes.h"
-#include "gromacs/utility/basedefinitions.h"
 
 struct gmx_domdec_t;
 struct gmx_enerdata_t;
@@ -90,6 +89,8 @@ class InteractiveMolecularDynamics;
 class MDLogger;
 struct ImdOptions;
 struct MdrunOptions;
+template<typename>
+class ArrayRef;
 
 /*! \brief
  * Creates a module for interactive molecular dynamics.
@@ -130,26 +131,26 @@ void write_IMDgroup_to_file(bool              bIMD,
  * \param ms           Handler for multi-simulations.
  * \param top_global   The topology of the whole system.
  * \param mdlog        Logger
- * \param x            The starting positions of the atoms.
+ * \param coords       The starting positions of the atoms.
  * \param nfile        Number of files.
  * \param fnm          Struct containing file names etc.
  * \param oenv         Output options.
  * \param options      Options for interactive MD.
  * \param startingBehavior  Describes whether this is a restart appending to output files
  */
-std::unique_ptr<ImdSession> makeImdSession(const t_inputrec*       ir,
-                                           const t_commrec*        cr,
-                                           gmx_wallcycle*          wcycle,
-                                           gmx_enerdata_t*         enerd,
-                                           const gmx_multisim_t*   ms,
-                                           const gmx_mtop_t&       top_global,
-                                           const MDLogger&         mdlog,
-                                           const rvec              x[],
-                                           int                     nfile,
-                                           const t_filenm          fnm[],
-                                           const gmx_output_env_t* oenv,
-                                           const ImdOptions&       options,
-                                           StartingBehavior        startingBehavior);
+std::unique_ptr<ImdSession> makeImdSession(const t_inputrec*              ir,
+                                           const t_commrec*               cr,
+                                           gmx_wallcycle*                 wcycle,
+                                           gmx_enerdata_t*                enerd,
+                                           const gmx_multisim_t*          ms,
+                                           const gmx_mtop_t&              top_global,
+                                           const MDLogger&                mdlog,
+                                           gmx::ArrayRef<const gmx::RVec> coords,
+                                           int                            nfile,
+                                           const t_filenm                 fnm[],
+                                           const gmx_output_env_t*        oenv,
+                                           const ImdOptions&              options,
+                                           StartingBehavior               startingBehavior);
 
 class ImdSession
 {
@@ -179,18 +180,18 @@ public:
      * \param step         The time step.
      * \param bNS          Is this a neighbor searching step?
      * \param box          The simulation box.
-     * \param x            The local atomic positions on this node.
+     * \param coords       The local atomic positions on this node.
      * \param t            The time.
      *
      * \returns            Whether or not we have to do IMD communication at this step.
      */
-    bool run(int64_t step, bool bNS, const matrix box, const rvec x[], double t);
+    bool run(int64_t step, bool bNS, const matrix box, gmx::ArrayRef<const gmx::RVec> coords, double t);
 
     /*! \brief Add external forces from a running interactive molecular dynamics session.
      *
-     * \param f            The forces.
+     * \param force The forces.
      */
-    void applyForces(rvec* f);
+    void applyForces(gmx::ArrayRef<gmx::RVec> force);
 
     /*! \brief Copy energies and convert to float from enerdata to the IMD energy record.
      *
@@ -220,19 +221,19 @@ private:
 
 public:
     // Befriend the factory function.
-    friend std::unique_ptr<ImdSession> makeImdSession(const t_inputrec*       ir,
-                                                      const t_commrec*        cr,
-                                                      gmx_wallcycle*          wcycle,
-                                                      gmx_enerdata_t*         enerd,
-                                                      const gmx_multisim_t*   ms,
-                                                      const gmx_mtop_t&       top_global,
-                                                      const MDLogger&         mdlog,
-                                                      const rvec              x[],
-                                                      int                     nfile,
-                                                      const t_filenm          fnm[],
-                                                      const gmx_output_env_t* oenv,
-                                                      const ImdOptions&       options,
-                                                      StartingBehavior        startingBehavior);
+    friend std::unique_ptr<ImdSession> makeImdSession(const t_inputrec*              ir,
+                                                      const t_commrec*               cr,
+                                                      gmx_wallcycle*                 wcycle,
+                                                      gmx_enerdata_t*                enerd,
+                                                      const gmx_multisim_t*          ms,
+                                                      const gmx_mtop_t&              top_global,
+                                                      const MDLogger&                mdlog,
+                                                      gmx::ArrayRef<const gmx::RVec> coords,
+                                                      int                            nfile,
+                                                      const t_filenm                 fnm[],
+                                                      const gmx_output_env_t*        oenv,
+                                                      const ImdOptions&              options,
+                                                      StartingBehavior startingBehavior);
 };
 
 } // namespace gmx
index d8cddba26ed7dede05b2f0fe9f7a0136b3e364e7..398841179cf4d26e5744ada8e2fc9aeaec2bd921 100644 (file)
@@ -778,13 +778,7 @@ bool Constraints::Impl::apply(bool                      bLog,
         if (ed && delta_step > 0)
         {
             /* apply the essential dynamics constraints here */
-            do_edsam(&ir,
-                     step,
-                     cr,
-                     as_rvec_array(xprime.unpaddedArrayRef().data()),
-                     as_rvec_array(v.unpaddedArrayRef().data()),
-                     box,
-                     ed);
+            do_edsam(&ir, step, cr, xprime.unpaddedArrayRef(), v.unpaddedArrayRef(), box, ed);
         }
     }
     wallcycle_stop(wcycle, ewcCONSTR);
index 424e9844fc4e0e8a138b64e07ac1d978aa7b6699..362cbaccf4b90f517577d22ea9506cbba3c0442d 100644 (file)
@@ -693,14 +693,12 @@ static void computeSpecialForces(FILE*                          fplog,
                     fplog);
         }
     }
-
-    rvec* f = as_rvec_array(forceWithVirialMtsLevel0->force_.data());
-
     /* Add the forces from enforced rotation potentials (if any) */
     if (inputrec.bRot)
     {
         wallcycle_start(wcycle, ewcROTadd);
-        enerd->term[F_COM_PULL] += add_rot_forces(enforcedRotation, f, cr, step, t);
+        enerd->term[F_COM_PULL] +=
+                add_rot_forces(enforcedRotation, forceWithVirialMtsLevel0->force_, cr, step, t);
         wallcycle_stop(wcycle, ewcROTadd);
     }
 
@@ -711,13 +709,13 @@ static void computeSpecialForces(FILE*                          fplog,
          * Thus if no other algorithm (e.g. PME) requires it, the forces
          * here will contribute to the virial.
          */
-        do_flood(cr, inputrec, as_rvec_array(x.data()), f, ed, box, step, didNeighborSearch);
+        do_flood(cr, inputrec, x, forceWithVirialMtsLevel0->force_, ed, box, step, didNeighborSearch);
     }
 
     /* Add forces from interactive molecular dynamics (IMD), if any */
     if (inputrec.bIMD && stepWork.computeForces)
     {
-        imdSession->applyForces(f);
+        imdSession->applyForces(forceWithVirialMtsLevel0->force_);
     }
 }
 
@@ -1697,7 +1695,7 @@ void do_force(FILE*                               fplog,
     if (inputrec.bRot)
     {
         wallcycle_start(wcycle, ewcROT);
-        do_rotation(cr, enforcedRotation, box, as_rvec_array(x.unpaddedArrayRef().data()), t, step, stepWork.doNeighborSearch);
+        do_rotation(cr, enforcedRotation, box, x.unpaddedConstArrayRef(), t, step, stepWork.doNeighborSearch);
         wallcycle_stop(wcycle, ewcROT);
     }
 
index 5584b60232d7e89090ea58cacc193c2144edd0e5..8c58a5bbc034cfbf39fa632f69c49c720972ba5a 100644 (file)
@@ -1359,7 +1359,7 @@ void gmx::LegacySimulator::do_md()
                                  mdrunOptions.writeConfout,
                                  bSumEkinhOld);
         /* Check if IMD step and do IMD communication, if bIMD is TRUE. */
-        bInteractiveMDstep = imdSession->run(step, bNS, state->box, state->x.rvec_array(), t);
+        bInteractiveMDstep = imdSession->run(step, bNS, state->box, state->x, t);
 
         /* kludge -- virial is lost with restart for MTTK NPT control. Must reload (saved earlier). */
         if (startingBehavior != StartingBehavior::NewSimulation && bFirstStep
index 7c59a560058da0cd06e7bfce55eaadbdb59eae46..f1bf1197940868dc03e1c35e4301a8695e849fbd 100644 (file)
@@ -1777,7 +1777,7 @@ void LegacySimulator::do_cg()
         }
 
         /* Send energies and positions to the IMD client if bIMD is TRUE. */
-        if (MASTER(cr) && imdSession->run(step, TRUE, state_global->box, state_global->x.rvec_array(), 0))
+        if (MASTER(cr) && imdSession->run(step, TRUE, state_global->box, state_global->x, 0))
         {
             imdSession->sendPositionsAndEnergies();
         }
@@ -2591,7 +2591,7 @@ void LegacySimulator::do_lbfgs()
         }
 
         /* Send x and E to IMD client, if bIMD is TRUE. */
-        if (imdSession->run(step, TRUE, state_global->box, state_global->x.rvec_array(), 0) && MASTER(cr))
+        if (imdSession->run(step, TRUE, state_global->box, state_global->x, 0) && MASTER(cr))
         {
             imdSession->sendPositionsAndEnergies();
         }
@@ -2940,7 +2940,7 @@ void LegacySimulator::do_steep()
         if (imdSession->run(count,
                             TRUE,
                             MASTER(cr) ? state_global->box : nullptr,
-                            MASTER(cr) ? state_global->x.rvec_array() : nullptr,
+                            MASTER(cr) ? state_global->x : gmx::ArrayRef<gmx::RVec>(),
                             0)
             && MASTER(cr))
         {
index 488951c568a3aecc545cf149caee8e6831ff7af0..f732c60ede37c85914da1450daac5590bc5d7db7 100644 (file)
@@ -1879,7 +1879,7 @@ int Mdrunner::mdrunner()
                                          ms,
                                          mtop,
                                          mdlog,
-                                         MASTER(cr) ? globalState->x.rvec_array() : nullptr,
+                                         MASTER(cr) ? globalState->x : gmx::ArrayRef<gmx::RVec>(),
                                          filenames.size(),
                                          filenames.data(),
                                          oenv,
index ab10a7a8a8740def7f291adfd2688b6b4afd2916..27adb6c158a5a59d2cd2c525a98ccd36f5ec7cba 100644 (file)
@@ -609,7 +609,7 @@ static void reduce_output(const t_commrec* cr, gmx_enfrot* er, real t, int64_t s
 
 /* Add the forces from enforced rotation potential to the local forces.
  * Should be called after the SR forces have been evaluated */
-real add_rot_forces(gmx_enfrot* er, rvec f[], const t_commrec* cr, int64_t step, real t)
+real add_rot_forces(gmx_enfrot* er, gmx::ArrayRef<gmx::RVec> force, const t_commrec* cr, int64_t step, real t)
 {
     real Vrot = 0.0; /* If more than one rotation group is present, Vrot
                         assembles the local parts from all groups         */
@@ -626,7 +626,7 @@ real add_rot_forces(gmx_enfrot* er, rvec f[], const t_commrec* cr, int64_t step,
             /* Get the right index of the local force */
             int ii = localRotationGroupIndex[l];
             /* Add */
-            rvec_inc(f[ii], erg->f_rot_loc[l]);
+            rvec_inc(force[ii], erg->f_rot_loc[l]);
         }
     }
 
@@ -2010,12 +2010,12 @@ static void flex_precalc_inner_sum(const gmx_enfrotgrp* erg)
 }
 
 
-static real do_flex2_lowlevel(gmx_enfrotgrp* erg,
-                              real           sigma, /* The Gaussian width sigma */
-                              rvec           x[],
-                              gmx_bool       bOutstepRot,
-                              gmx_bool       bOutstepSlab,
-                              const matrix   box)
+static real do_flex2_lowlevel(gmx_enfrotgrp*                 erg,
+                              real                           sigma, /* The Gaussian width sigma */
+                              gmx::ArrayRef<const gmx::RVec> coords,
+                              gmx_bool                       bOutstepRot,
+                              gmx_bool                       bOutstepSlab,
+                              const matrix                   box)
 {
     int  count, ii, iigrp;
     rvec xj;          /* position in the i-sum                         */
@@ -2076,7 +2076,7 @@ static real do_flex2_lowlevel(gmx_enfrotgrp* erg,
          * Note that erg->xc_center contains the center of mass in case the flex2-t
          * potential was chosen. For the flex2 potential erg->xc_center must be
          * zero. */
-        rvec_sub(x[ii], erg->xc_center, xj);
+        rvec_sub(coords[ii], erg->xc_center, xj);
 
         /* Shift this atom such that it is near its reference */
         shift_single_coord(box, xj, erg->xc_shifts[iigrp]);
@@ -2278,11 +2278,11 @@ static real do_flex2_lowlevel(gmx_enfrotgrp* erg,
 
 
 static real do_flex_lowlevel(gmx_enfrotgrp* erg,
-                             real         sigma, /* The Gaussian width sigma                      */
-                             rvec         x[],
-                             gmx_bool     bOutstepRot,
-                             gmx_bool     bOutstepSlab,
-                             const matrix box)
+                             real sigma, /* The Gaussian width sigma                      */
+                             gmx::ArrayRef<const gmx::RVec> coords,
+                             gmx_bool                       bOutstepRot,
+                             gmx_bool                       bOutstepSlab,
+                             const matrix                   box)
 {
     int      count, iigrp;
     rvec     xj, yj0;        /* current and reference position                */
@@ -2335,7 +2335,7 @@ static real do_flex_lowlevel(gmx_enfrotgrp* erg,
          * Note that erg->xc_center contains the center of mass in case the flex-t
          * potential was chosen. For the flex potential erg->xc_center must be
          * zero. */
-        rvec_sub(x[ii], erg->xc_center, xj);
+        rvec_sub(coords[ii], erg->xc_center, xj);
 
         /* Shift this atom such that it is near its reference */
         shift_single_coord(box, xj, erg->xc_shifts[iigrp]);
@@ -2608,11 +2608,11 @@ static void get_firstlast_slab_check(gmx_enfrotgrp* erg, /* The rotation group (
 static void do_flexible(gmx_bool       bMaster,
                         gmx_enfrot*    enfrot, /* Other rotation data                        */
                         gmx_enfrotgrp* erg,
-                        rvec           x[], /* The local positions                        */
-                        const matrix   box,
-                        double         t,           /* Time in picoseconds                        */
-                        gmx_bool       bOutstepRot, /* Output to main rotation output file        */
-                        gmx_bool bOutstepSlab)      /* Output per-slab data                       */
+                        gmx::ArrayRef<const gmx::RVec> coords, /* The local positions */
+                        const matrix                   box,
+                        double   t,            /* Time in picoseconds                        */
+                        gmx_bool bOutstepRot,  /* Output to main rotation output file        */
+                        gmx_bool bOutstepSlab) /* Output per-slab data                       */
 {
     int  nslabs;
     real sigma; /* The Gaussian width sigma */
@@ -2646,12 +2646,12 @@ static void do_flexible(gmx_bool       bMaster,
     if (erg->rotg->eType == EnforcedRotationGroupType::Flex
         || erg->rotg->eType == EnforcedRotationGroupType::Flext)
     {
-        erg->V = do_flex_lowlevel(erg, sigma, x, bOutstepRot, bOutstepSlab, box);
+        erg->V = do_flex_lowlevel(erg, sigma, coords, bOutstepRot, bOutstepSlab, box);
     }
     else if (erg->rotg->eType == EnforcedRotationGroupType::Flex2
              || erg->rotg->eType == EnforcedRotationGroupType::Flex2t)
     {
-        erg->V = do_flex2_lowlevel(erg, sigma, x, bOutstepRot, bOutstepSlab, box);
+        erg->V = do_flex2_lowlevel(erg, sigma, coords, bOutstepRot, bOutstepSlab, box);
     }
     else
     {
@@ -2936,11 +2936,11 @@ static void do_radial_motion(gmx_enfrotgrp* erg,
 
 
 /* Calculate the radial motion pivot-free potential and forces */
-static void do_radial_motion_pf(gmx_enfrotgrp* erg,
-                                rvec           x[], /* The positions                              */
-                                const matrix   box, /* The simulation box                         */
-                                gmx_bool bOutstepRot,  /* Output to main rotation output file  */
-                                gmx_bool bOutstepSlab) /* Output per-slab data */
+static void do_radial_motion_pf(gmx_enfrotgrp*                 erg,
+                                gmx::ArrayRef<const gmx::RVec> coords, /* The positions */
+                                const matrix box, /* The simulation box                         */
+                                gmx_bool     bOutstepRot, /* Output to main rotation output file  */
+                                gmx_bool     bOutstepSlab)    /* Output per-slab data */
 {
     rvec     xj;      /* Current position */
     rvec     xj_xc;   /* xj  - xc  */
@@ -3005,7 +3005,7 @@ static void do_radial_motion_pf(gmx_enfrotgrp* erg,
         wj = N_M * mj;
 
         /* Current position of this atom: x[ii][XX/YY/ZZ] */
-        copy_rvec(x[ii], xj);
+        copy_rvec(coords[ii], xj);
 
         /* Shift this atom such that it is near its reference */
         shift_single_coord(box, xj, erg->xc_shifts[iigrp]);
@@ -3138,9 +3138,9 @@ static void radial_motion2_precalc_inner_sum(const gmx_enfrotgrp* erg, rvec inne
 
 
 /* Calculate the radial motion 2 potential and forces */
-static void do_radial_motion2(gmx_enfrotgrp* erg,
-                              rvec           x[],   /* The positions                              */
-                              const matrix   box,   /* The simulation box                         */
+static void do_radial_motion2(gmx_enfrotgrp*                 erg,
+                              gmx::ArrayRef<const gmx::RVec> coords, /* The positions */
+                              const matrix box,     /* The simulation box                         */
                               gmx_bool bOutstepRot, /* Output to main rotation output file        */
                               gmx_bool bOutstepSlab) /* Output per-slab data */
 {
@@ -3195,7 +3195,7 @@ static void do_radial_motion2(gmx_enfrotgrp* erg,
             mj = erg->mc[iigrp];
 
             /* Current position of this atom: x[ii] */
-            copy_rvec(x[ii], xj);
+            copy_rvec(coords[ii], xj);
 
             /* Shift this atom such that it is near its reference */
             shift_single_coord(box, xj, erg->xc_shifts[iigrp]);
@@ -3856,7 +3856,7 @@ static void rotate_local_reference(gmx_enfrotgrp* erg)
 /* Select the PBC representation for each local x position and store that
  * for later usage. We assume the right PBC image of an x is the one nearest to
  * its rotated reference */
-static void choose_pbc_image(rvec x[], gmx_enfrotgrp* erg, const matrix box, int npbcdim)
+static void choose_pbc_image(gmx::ArrayRef<const gmx::RVec> coords, gmx_enfrotgrp* erg, const matrix box, int npbcdim)
 {
     const auto& localRotationGroupIndex = erg->atomSet->localIndex();
     for (gmx::index i = 0; i < localRotationGroupIndex.ssize(); i++)
@@ -3873,12 +3873,18 @@ static void choose_pbc_image(rvec x[], gmx_enfrotgrp* erg, const matrix box, int
         copy_rvec(erg->xr_loc[i], xref);
         rvec_inc(xref, erg->xc_ref_center);
 
-        copy_correct_pbc_image(x[ii], erg->x_loc_pbc[i], xref, box, npbcdim);
+        copy_correct_pbc_image(coords[ii], erg->x_loc_pbc[i], xref, box, npbcdim);
     }
 }
 
 
-void do_rotation(const t_commrec* cr, gmx_enfrot* er, const matrix box, rvec x[], real t, int64_t step, gmx_bool bNS)
+void do_rotation(const t_commrec*               cr,
+                 gmx_enfrot*                    er,
+                 const matrix                   box,
+                 gmx::ArrayRef<const gmx::RVec> coords,
+                 real                           t,
+                 int64_t                        step,
+                 gmx_bool                       bNS)
 {
     gmx_bool    outstep_slab, outstep_rot;
     gmx_bool    bColl;
@@ -3925,7 +3931,7 @@ void do_rotation(const t_commrec* cr, gmx_enfrot* er, const matrix box, rvec x[]
                                         erg->xc_shifts,
                                         erg->xc_eshifts,
                                         bNS,
-                                        x,
+                                        as_rvec_array(coords.data()),
                                         rotg->nat,
                                         erg->atomSet->numAtomsLocal(),
                                         erg->atomSet->localIndex().data(),
@@ -3953,7 +3959,7 @@ void do_rotation(const t_commrec* cr, gmx_enfrot* er, const matrix box, rvec x[]
 
             /* Choose the nearest PBC images of the group atoms with respect
              * to the rotated reference positions */
-            choose_pbc_image(x, erg, box, 3);
+            choose_pbc_image(coords, erg, box, 3);
 
             /* Get the center of the rotation group */
             if ((rotg->eType == EnforcedRotationGroupType::Isopf)
@@ -4016,11 +4022,11 @@ void do_rotation(const t_commrec* cr, gmx_enfrot* er, const matrix box, rvec x[]
                 do_radial_motion(erg, outstep_rot, outstep_slab);
                 break;
             case EnforcedRotationGroupType::Rmpf:
-                do_radial_motion_pf(erg, x, box, outstep_rot, outstep_slab);
+                do_radial_motion_pf(erg, coords, box, outstep_rot, outstep_slab);
                 break;
             case EnforcedRotationGroupType::Rm2:
             case EnforcedRotationGroupType::Rm2pf:
-                do_radial_motion2(erg, x, box, outstep_rot, outstep_slab);
+                do_radial_motion2(erg, coords, box, outstep_rot, outstep_slab);
                 break;
             case EnforcedRotationGroupType::Flext:
             case EnforcedRotationGroupType::Flex2t:
@@ -4030,13 +4036,13 @@ void do_rotation(const t_commrec* cr, gmx_enfrot* er, const matrix box, rvec x[]
                 get_center(erg->xc, erg->mc, rotg->nat, erg->xc_center);
                 svmul(-1.0, erg->xc_center, transvec);
                 translate_x(erg->xc, rotg->nat, transvec);
-                do_flexible(MASTER(cr), er, erg, x, box, t, outstep_rot, outstep_slab);
+                do_flexible(MASTER(cr), er, erg, coords, box, t, outstep_rot, outstep_slab);
                 break;
             case EnforcedRotationGroupType::Flex:
             case EnforcedRotationGroupType::Flex2:
                 /* Do NOT subtract the center of mass in the low level routines! */
                 clear_rvec(erg->xc_center);
-                do_flexible(MASTER(cr), er, erg, x, box, t, outstep_rot, outstep_slab);
+                do_flexible(MASTER(cr), er, erg, coords, box, t, outstep_rot, outstep_slab);
                 break;
             default: gmx_fatal(FARGS, "No such rotation potential.");
         }
index de5f1388a06590ab7fcfa205c01c57cc72320350..d367317dc755e5330c155c1b5ac0f24046fe726c 100644 (file)
@@ -54,7 +54,6 @@
 #include <memory>
 
 #include "gromacs/math/vectypes.h"
-#include "gromacs/utility/basedefinitions.h"
 
 struct gmx_domdec_t;
 struct gmx_enfrot;
@@ -71,6 +70,8 @@ namespace gmx
 enum class StartingBehavior;
 class LocalAtomSetManager;
 struct MdrunOptions;
+template<typename>
+class ArrayRef;
 
 class EnforcedRotation
 {
@@ -133,13 +134,19 @@ std::unique_ptr<gmx::EnforcedRotation> init_rot(FILE*                     fplog,
  * \param cr      Pointer to MPI communication data.
  * \param er      Pointer to the enforced rotation working data.
  * \param box     Simulation box, needed to make group whole.
- * \param x       The positions of all the local particles.
+ * \param coords  The positions of all the local particles.
  * \param t       Time.
  * \param step    The time step.
  * \param bNS     After domain decomposition / neighbor searching several
  *                local arrays have to be updated (masses, shifts)
  */
-void do_rotation(const t_commrec* cr, gmx_enfrot* er, const matrix box, rvec x[], real t, int64_t step, gmx_bool bNS);
+void do_rotation(const t_commrec*               cr,
+                 gmx_enfrot*                    er,
+                 const matrix                   box,
+                 gmx::ArrayRef<const gmx::RVec> coords,
+                 real                           t,
+                 int64_t                        step,
+                 bool                           bNS);
 
 
 /*! \brief Add the enforced rotation forces to the official force array.
@@ -152,14 +159,14 @@ void do_rotation(const t_commrec* cr, gmx_enfrot* er, const matrix box, rvec x[]
  * the potential, the angle of the group(s), and torques).
  *
  * \param er      Pointer to the enforced rotation working data.
- * \param f       The local forces to which the rotational forces have
+ * \param force   The local forces to which the rotational forces have
  *                to be added.
  * \param cr      Pointer to MPI communication data.
  * \param step    The time step, used for output.
  * \param t       Time, used for output.
  * \returns       The potential energy of the rotation potentials.
  */
-real add_rot_forces(gmx_enfrot* er, rvec f[], const t_commrec* cr, int64_t step, real t);
+real add_rot_forces(gmx_enfrot* er, gmx::ArrayRef<gmx::RVec> force, const t_commrec* cr, int64_t step, real t);
 
 
 #endif