Use RVec and std::array in gmx_domdec_comm_t
authorJoe Jordan <ejjordan12@gmail.com>
Tue, 16 Mar 2021 18:03:19 +0000 (18:03 +0000)
committerArtem Zhmurov <zhmurov@gmail.com>
Tue, 16 Mar 2021 18:03:19 +0000 (18:03 +0000)
Where possible, rvecs have been changed to RVecs and pointers
changed to std::arrays.

src/gromacs/domdec/domdec_internal.h
src/gromacs/domdec/partition.cpp
src/gromacs/mdlib/nsgrid.cpp
src/gromacs/mdlib/nsgrid.h

index 680515174d538bdde84ec7a920fb282be541a8e5..41730c5c1b346842e10585ea02717727cac5d758 100644 (file)
@@ -583,35 +583,35 @@ struct gmx_domdec_comm_t // NOLINT (clang-analyzer-optin.performance.Padding)
     /**< Cut-off for multi-body interactions, also 2-body bonded when \p cutoff_mody > \p cutoff */
     real cutoff_mbody = 0;
     /**< The minimum guaranteed cell-size, Cartesian indexing */
-    rvec cellsize_min = {};
+    gmx::RVec cellsize_min = { 0, 0, 0 };
     /**< The minimum guaranteed cell-size with dlb=auto */
-    rvec cellsize_min_dlb = {};
+    gmx::RVec cellsize_min_dlb = { 0, 0, 0 };
     /**< The lower limit for the DD cell size with DLB */
     real cellsize_limit = 0;
     /**< Effectively no NB cut-off limit with DLB for systems without PBC? */
-    gmx_bool bVacDLBNoLimit = false;
+    bool bVacDLBNoLimit = false;
 
     /** With PME load balancing we set limits on DLB */
-    gmx_bool bPMELoadBalDLBLimits = false;
+    bool bPMELoadBalDLBLimits = false;
     /** DLB needs to take into account that we want to allow this maximum
      *  cut-off (for PME load balancing), this could limit cell boundaries.
      */
     real PMELoadBal_max_cutoff = 0;
 
     /**< box lower corner, required with dim's without pbc and -gcom */
-    rvec box0 = {};
+    gmx::RVec box0 = { 0, 0, 0 };
     /**< box size, required with dim's without pbc and -gcom */
-    rvec box_size = {};
+    gmx::RVec box_size = { 0, 0, 0 };
 
     /**< The DD cell lower corner, in triclinic space */
-    rvec cell_x0 = {};
+    gmx::RVec cell_x0 = { 0, 0, 0 };
     /**< The DD cell upper corner, in triclinic space */
-    rvec cell_x1 = {};
+    gmx::RVec cell_x1 = { 0, 0, 0 };
 
     /**< The old \p cell_x0, to check cg displacements */
-    rvec old_cell_x0 = {};
+    gmx::RVec old_cell_x0 = { 0, 0, 0 };
     /**< The old \p cell_x1, to check cg displacements */
-    rvec old_cell_x1 = {};
+    gmx::RVec old_cell_x1 = { 0, 0, 0 };
 
     /** The communication setup and charge group boundaries for the zones */
     gmx_domdec_zones_t zones;
@@ -621,12 +621,12 @@ struct gmx_domdec_comm_t // NOLINT (clang-analyzer-optin.performance.Padding)
      * dynamic load balancing.
      */
     /**< Zone limits for dim 1 with staggered grids */
-    gmx_ddzone_t zone_d1[2];
+    std::array<gmx_ddzone_t, 2> zone_d1;
     /**< Zone limits for dim 2 with staggered grids */
     gmx_ddzone_t zone_d2[2][2];
 
     /** The coordinate/force communication setup and indices */
-    gmx_domdec_comm_dim_t cd[DIM];
+    std::array<gmx_domdec_comm_dim_t, DIM> cd;
     /** Restricts the maximum number of cells to communicate with in one dimension
      *
      * Dynamic load balancing is not permitted to change sizes if it
@@ -639,7 +639,7 @@ struct gmx_domdec_comm_t // NOLINT (clang-analyzer-optin.performance.Padding)
     int64_t master_cg_ddp_count = 0;
 
     /** The number of cg's received from the direct neighbors */
-    int zone_ncg1[DD_MAXZONE] = { 0 };
+    std::array<int, DD_MAXZONE> zone_ncg1 = { 0 };
 
     /** The atom ranges in the local state */
     DDAtomRanges atomRanges;
@@ -687,11 +687,11 @@ struct gmx_domdec_comm_t // NOLINT (clang-analyzer-optin.performance.Padding)
 
     /* Cycle counters over nstlist steps */
     /**< Total cycles counted */
-    float cycl[ddCyclNr] = {};
+    std::array<float, ddCyclNr> cycl = { 0 };
     /**< The number of cycle recordings */
-    int cycl_n[ddCyclNr] = {};
+    std::array<int, ddCyclNr> cycl_n = { 0 };
     /**< The maximum cycle count */
-    float cycl_max[ddCyclNr] = {};
+    std::array<float, ddCyclNr> cycl_max = { 0 };
     /**< Total flops counted */
     double flop = 0.0;
     /**< The number of flop recordings */
@@ -727,7 +727,7 @@ struct gmx_domdec_comm_t // NOLINT (clang-analyzer-optin.performance.Padding)
     /**< Max \p load_sum over the ranks */
     double load_max = 0.0;
     /**< Was load balancing limited, per DD dim */
-    ivec load_lim = {};
+    gmx::IVec load_lim = { 0, 0, 0 };
     /**< Total time on PP done during PME overlap time */
     double load_mdf = 0.0;
     /**< Total time on our PME-only rank */
index cbf2ff5c303b0a3607dab886d5bf07e768d21470..a4e994e490c723a0e8c1705ea31c4a2972f464e5 100644 (file)
@@ -491,8 +491,8 @@ static void dd_set_cginfo(gmx::ArrayRef<const int> index_gl, int cg0, int cg1, t
 static void make_dd_indices(gmx_domdec_t* dd, const int atomStart)
 {
     const int                numZones               = dd->comm->zones.n;
-    const int*               zone2cg                = dd->comm->zones.cg_range.data();
-    const int*               zone_ncg1              = dd->comm->zone_ncg1;
+    gmx::ArrayRef<const int> zone2cg                = dd->comm->zones.cg_range;
+    gmx::ArrayRef<const int> zone_ncg1              = dd->comm->zone_ncg1;
     gmx::ArrayRef<const int> globalAtomGroupIndices = dd->globalAtomGroupIndices;
 
     std::vector<int>& globalAtomIndices = dd->globalAtomIndices;
index 1bfe85acf07edf6d063a62dc9747cdb9d85384fc..4b821c6484e66d4dec2bfb5376047da4a1847d91 100644 (file)
@@ -112,8 +112,8 @@ void get_nsgrid_boundaries(int           nboundeddim,
                            matrix        box,
                            gmx_domdec_t* dd,
                            gmx_ddbox_t*  ddbox,
-                           rvec*         gr0,
-                           rvec*         gr1,
+                           gmx::RVec*    gr0,
+                           gmx::RVec*    gr1,
                            int           ncg,
                            rvec*         cgcm,
                            rvec          grid_x0,
index 7de1184bb2ad9a8c113f01fa1a9f9aaa39f31ab5..421b917b309f68e9a4bf73446a55bdfee31670d0 100644 (file)
@@ -73,8 +73,8 @@ void get_nsgrid_boundaries(int                  nboundeddim,
                            matrix               box,
                            struct gmx_domdec_t* dd,
                            gmx_ddbox_t*         ddbox,
-                           rvec*                gr0,
-                           rvec*                gr1,
+                           gmx::RVec*           gr0,
+                           gmx::RVec*           gr1,
                            int                  ncg,
                            rvec*                cgcm,
                            rvec                 grid_x0,