refactoring

src/EXTRA-FIX/fix_ttm_grid.cpp

@@ -57,6 +57,7 @@ FixTTMGrid::FixTTMGrid(LAMMPS *lmp, int narg, char **arg) :
   pergrid_flag = 1;
 
   /*
+  NOTE: uncomment this when ready to release
   if (outfile) error->all(FLERR,"Fix ttm/grid does not support outfile option - "
                           "use dump grid instead");
   */
@@ -285,7 +286,7 @@ void FixTTMGrid::read_electron_temperatures(const std::string &filename)
   int ***T_initial_set;
   memory->create3d_offset(T_initial_set, nzlo_in, nzhi_in, nylo_in, nyhi_in, nxlo_in, nxhi_in,
                           "ttm/grid:T_initial_set");
-  memset(&T_initial_set[nzlo_in][nylo_in][nxlo_in], 0, ngridmine * sizeof(int));
+  memset(&T_initial_set[nzlo_in][nylo_in][nxlo_in], 0, ngridown * sizeof(int));
 
   // proc 0 opens file
 
@@ -365,6 +366,7 @@ void FixTTMGrid::read_electron_temperatures(const std::string &filename)
 /* ----------------------------------------------------------------------
    write out current electron temperatures to user-specified file
    only written by proc 0
+   NOTE: remove this function when ready to release
 ------------------------------------------------------------------------- */
 
 void FixTTMGrid::write_electron_temperatures(const std::string &filename)
@@ -445,16 +447,10 @@ void FixTTMGrid::allocate_grid()
                     nxlo_in, nxhi_in, nylo_in, nyhi_in, nzlo_in, nzhi_in, 
                     nxlo_out, nxhi_out, nylo_out, nyhi_out, nzlo_out, nzhi_out);
 
-  // set ngridout and ngridmine and error check
+  ngridown = (nxhi_in - nxlo_in + 1) * (nyhi_in - nylo_in + 1) * 
-
+    (nzhi_in - nzlo_in + 1);
-  bigint totalmine;
+  ngridout = (nxhi_out - nxlo_out + 1) * (nyhi_out - nylo_out + 1) * 
-  totalmine =
+    (nzhi_out - nzlo_out + 1);
-      (bigint) (nxhi_out - nxlo_out + 1) * (nyhi_out - nylo_out + 1) * (nzhi_out - nzlo_out + 1);
-  if (totalmine > MAXSMALLINT) error->one(FLERR, "Too many owned+ghost grid points in fix ttm");
-  ngridout = totalmine;
-
-  totalmine = (bigint) (nxhi_in - nxlo_in + 1) * (nyhi_in - nylo_in + 1) * (nzhi_in - nzlo_in + 1);
-  ngridmine = totalmine;
 
   // setup grid communication and allocate grid data structs
 
@@ -562,6 +558,7 @@ void FixTTMGrid::restart(char *buf)
 /* ----------------------------------------------------------------------
    pack values from local grid into buf
    used by which = 0 and 1
+   NOTE: remove this function when ready to release
 ------------------------------------------------------------------------- */
 
 void FixTTMGrid::pack_gather_grid(int /*which*/, void *vbuf)
@@ -579,6 +576,7 @@ void FixTTMGrid::pack_gather_grid(int /*which*/, void *vbuf)
 /* ----------------------------------------------------------------------
    which = 0: unpack values from buf into global gbuf based on their indices
    which = 1: print values from buf to FPout file
+   NOTE: remove this function when ready to release
 ------------------------------------------------------------------------- */
 
 void FixTTMGrid::unpack_gather_grid(int which, void *vbuf, void *vgbuf, int xlo, int xhi, int ylo,
@@ -617,6 +615,7 @@ void FixTTMGrid::unpack_gather_grid(int which, void *vbuf, void *vgbuf, int xlo,
    return index of grid associated with name
    this class can store M named grids, indexed 0 to M-1
    also set dim for 2d vs 3d grid
+   return -1 if grid name not found
 ------------------------------------------------------------------------- */
 
 int FixTTMGrid::get_grid_by_name(char *name, int &dim)
@@ -632,6 +631,7 @@ int FixTTMGrid::get_grid_by_name(char *name, int &dim)
 /* ----------------------------------------------------------------------
    return ptr to Grid data struct for grid with index
    this class can store M named grids, indexed 0 to M-1
+   return nullptr if index is invalid
 ------------------------------------------------------------------------- */
 
 void *FixTTMGrid::get_grid_by_index(int index)
@@ -645,8 +645,9 @@ void *FixTTMGrid::get_grid_by_index(int index)
    this class can store M named grids, indexed 0 to M-1
    each grid can store G named data sets, indexed 0 to G-1
      a data set name can be associated with multiple grids
-   also set ncol for data set, 0 = vector, 1-N for array with N columns
+   set ncol for data set, 0 = vector, 1-N for array with N columns
      vector = single value per grid pt, array = N values per grid pt
+   return -1 if data name not found
 ------------------------------------------------------------------------- */
 
 int FixTTMGrid::get_griddata_by_name(int igrid, char *name, int &ncol)
@@ -662,6 +663,7 @@ int FixTTMGrid::get_griddata_by_name(int igrid, char *name, int &ncol)
 /* ----------------------------------------------------------------------
    return ptr to multidim data array associated with index
    this class can store G named data sets, indexed 0 to M-1
+   return nullptr if index is invalid
 ------------------------------------------------------------------------- */
 
 void *FixTTMGrid::get_griddata_by_index(int index)

src/EXTRA-FIX/fix_ttm_grid.h

@@ -54,7 +54,7 @@ class FixTTMGrid : public FixTTM {
   double memory_usage() override;
 
  private:
-  int ngridmine, ngridout;
+  int ngridown, ngridout;
   int nxlo_in, nxhi_in, nylo_in, nyhi_in, nzlo_in, nzhi_in;
   int nxlo_out, nxhi_out, nylo_out, nyhi_out, nzlo_out, nzhi_out;
   double delxinv, delyinv, delzinv;

src/compute_property_grid.cpp

@@ -109,9 +109,18 @@ ComputePropertyGrid::ComputePropertyGrid(LAMMPS *lmp, int narg, char **arg) :
   }
 
   // instantiate the Grid class and allocate per-grid memory
-  // NOTE: need new memory create methods for 2d
 
   if (dimension == 2) {
+    grid2d = new Grid2d(lmp, world, nx, ny, 0, 0.0, 0.0,
+                        nxlo_in, nxhi_in, nylo_in, nyhi_in,
+                        nxlo_out, nxhi_out, nylo_out, nyhi_out);
+    if (nvalues == 1)
+      memory->create2d_offset(vec2d, nylo_out, nyhi_out, nxlo_out, nxhi_out, 
+                              "property/grid:vec2d");
+    else
+      memory->create3d_offset_last(array2d, nylo_out, nyhi_out, nxlo_out,
+                                   nxhi_out, nvalues, "property/grid:array2d");
+    ngridout = (nxhi_out - nxlo_out + 1) * (nyhi_out - nylo_out + 1);
     
   } else {
     grid3d = new Grid3d(lmp, world, nx, ny, nz, 0, 0.0, 0.0,
@@ -126,6 +135,8 @@ ComputePropertyGrid::ComputePropertyGrid(LAMMPS *lmp, int narg, char **arg) :
       memory->create4d_offset_last(array3d, nzlo_out, nzhi_out, nylo_out, 
                                    nyhi_out, nxlo_out,
                                    nxhi_out, nvalues, "property/grid:array3d");
+    ngridout = (nxhi_out - nxlo_out + 1) * (nyhi_out - nylo_out + 1) * 
+      (nzhi_out - nzlo_out + 1);
   }
 }
 
@@ -137,8 +148,8 @@ ComputePropertyGrid::~ComputePropertyGrid()
 
   delete grid2d;
   delete grid3d;
-  //memory->destroy2d_offset(vec2d);
+  memory->destroy2d_offset(vec2d,nylo_out,nxlo_out);
-  //memory->destroy2d_offset(array2d);
+  memory->destroy2d_offset(array2d,nylo_out,nxlo_out);
   memory->destroy3d_offset(vec3d,nzlo_out,nylo_out,nxlo_out);
   memory->destroy4d_offset_last(array3d,nzlo_out,nylo_out,nxlo_out);
 }
@@ -153,12 +164,12 @@ void ComputePropertyGrid::compute_pergrid()
   // may change between compute invocations due to load balancing
   
   if (dimension == 2)
-    grid2d->query_bounds(nxlo_in,nxhi_in,nylo_in,nyhi_in);
+    grid2d->get_bounds(nxlo_in,nxhi_in,nylo_in,nyhi_in);
   else
-    grid3d->query_bounds(nxlo_in,nxhi_in,nylo_in,nyhi_in,nzlo_in,nzhi_in);
+    grid3d->get_bounds(nxlo_in,nxhi_in,nylo_in,nyhi_in,nzlo_in,nzhi_in);
 
   // reallocate data vector or array if changed
-
+  // NOTE: still need to implement
 
 
 
@@ -175,6 +186,7 @@ void ComputePropertyGrid::compute_pergrid()
    return index of grid associated with name
    this class can store M named grids, indexed 0 to M-1
    also set dim for 2d vs 3d grid
+   return -1 if grid name not found
 ------------------------------------------------------------------------- */
 
 int ComputePropertyGrid::get_grid_by_name(char *name, int &dim)
@@ -190,6 +202,7 @@ int ComputePropertyGrid::get_grid_by_name(char *name, int &dim)
 /* ----------------------------------------------------------------------
    return ptr to Grid data struct for grid with index
    this class can store M named grids, indexed 0 to M-1
+   return nullptr if index is invalid
 ------------------------------------------------------------------------- */
 
 void *ComputePropertyGrid::get_grid_by_index(int index)
@@ -206,8 +219,9 @@ void *ComputePropertyGrid::get_grid_by_index(int index)
    this class can store M named grids, indexed 0 to M-1
    each grid can store G named data sets, indexed 0 to G-1
      a data set name can be associated with multiple grids
-   also set ncol for data set, 0 = vector, 1-N for array with N columns
+   set ncol for data set, 0 = vector, 1-N for array with N columns
      vector = single value per grid pt, array = N values per grid pt
+   return -1 if data name not found
 ------------------------------------------------------------------------- */
 
 int ComputePropertyGrid::get_griddata_by_name(int igrid, char *name, int &ncol)
@@ -224,6 +238,7 @@ int ComputePropertyGrid::get_griddata_by_name(int igrid, char *name, int &ncol)
 /* ----------------------------------------------------------------------
    return ptr to multidim data array associated with index
    this class can store G named data sets, indexed 0 to M-1
+   return nullptr if index is invalid
 ------------------------------------------------------------------------- */
 
 void *ComputePropertyGrid::get_griddata_by_index(int index)
@@ -246,9 +261,7 @@ void *ComputePropertyGrid::get_griddata_by_index(int index)
 
 double ComputePropertyGrid::memory_usage()
 {
-  double bytes = 0.0;
+  double bytes = (double) ngridout * nvalues * sizeof(double);
-  //double bytes = (double) nmax * nvalues * sizeof(double);
-  //bytes += (double) nmax * 2 * sizeof(int);
   return bytes;
 }
 
@@ -365,7 +378,7 @@ void ComputePropertyGrid::pack_x(int n)
   double boxlo,dx;
 
   if (dimension == 2) {
-    grid2d->query_box(0,boxlo,dx);
+    grid2d->get_box(0,boxlo,dx);
     if (nvalues == 0) {
       for (int iy = nylo_in; iy <= nyhi_in; iy++)
         for (int ix = nxlo_in; ix <= nxhi_in; ix++)
@@ -376,7 +389,7 @@ void ComputePropertyGrid::pack_x(int n)
           array2d[iy][ix][n] = boxlo + ix*dx;
     }
   } else if (dimension == 3) {
-    grid3d->query_box(0,boxlo,dx);
+    grid3d->get_box(0,boxlo,dx);
     if (nvalues == 0) {
       for (int iz = nzlo_in; iz <= nzhi_in; iz++)
         for (int iy = nylo_in; iy <= nyhi_in; iy++)
@@ -398,7 +411,7 @@ void ComputePropertyGrid::pack_y(int n)
   double boxlo,dy;
 
   if (dimension == 2) {
-    grid2d->query_box(1,boxlo,dy);
+    grid2d->get_box(1,boxlo,dy);
     if (nvalues == 0) {
       for (int iy = nylo_in; iy <= nyhi_in; iy++)
         for (int ix = nxlo_in; ix <= nxhi_in; ix++)
@@ -409,7 +422,7 @@ void ComputePropertyGrid::pack_y(int n)
           array2d[iy][ix][n] = boxlo + iy*dy;
     }
   } else if (dimension == 3) {
-    grid3d->query_box(1,boxlo,dy);
+    grid3d->get_box(1,boxlo,dy);
     if (nvalues == 0) {
       for (int iz = nzlo_in; iz <= nzhi_in; iz++)
         for (int iy = nylo_in; iy <= nyhi_in; iy++)
@@ -429,7 +442,7 @@ void ComputePropertyGrid::pack_y(int n)
 void ComputePropertyGrid::pack_z(int n)
 {
   double boxlo,dz;
-  grid3d->query_box(2,boxlo,dz);
+  grid3d->get_box(2,boxlo,dz);
 
   if (nvalues == 0) {
     for (int iz = nzlo_in; iz <= nzhi_in; iz++)
@@ -532,7 +545,7 @@ void ComputePropertyGrid::pack_xc(int n)
   double boxlo,dx;
 
   if (dimension == 2) {
-    grid2d->query_box(0,boxlo,dx);
+    grid2d->get_box(0,boxlo,dx);
     if (nvalues == 0) {
       for (int iy = nylo_in; iy <= nyhi_in; iy++)
         for (int ix = nxlo_in; ix <= nxhi_in; ix++)
@@ -543,7 +556,7 @@ void ComputePropertyGrid::pack_xc(int n)
           array2d[iy][ix][n] = boxlo + (ix+0.5)*dx;
     }
   } else if (dimension == 3) {
-    grid3d->query_box(0,boxlo,dx);
+    grid3d->get_box(0,boxlo,dx);
     if (nvalues == 0) {
       for (int iz = nzlo_in; iz <= nzhi_in; iz++)
         for (int iy = nylo_in; iy <= nyhi_in; iy++)
@@ -565,7 +578,7 @@ void ComputePropertyGrid::pack_yc(int n)
   double boxlo,dy;
 
   if (dimension == 2) {
-    grid2d->query_box(1,boxlo,dy);
+    grid2d->get_box(1,boxlo,dy);
     if (nvalues == 0) {
       for (int iy = nylo_in; iy <= nyhi_in; iy++)
         for (int ix = nxlo_in; ix <= nxhi_in; ix++)
@@ -576,7 +589,7 @@ void ComputePropertyGrid::pack_yc(int n)
           array2d[iy][ix][n] = boxlo + (iy+0.5)*dy;
     }
   } else if (dimension == 3) {
-    grid3d->query_box(1,boxlo,dy);
+    grid3d->get_box(1,boxlo,dy);
     if (nvalues == 0) {
       for (int iz = nzlo_in; iz <= nzhi_in; iz++)
         for (int iy = nylo_in; iy <= nyhi_in; iy++)
@@ -596,7 +609,7 @@ void ComputePropertyGrid::pack_yc(int n)
 void ComputePropertyGrid::pack_zc(int n)
 {
   double boxlo,dz;
-  grid3d->query_box(2,boxlo,dz);
+  grid3d->get_box(2,boxlo,dz);
 
   if (nvalues == 0) {
     for (int iz = nzlo_in; iz <= nzhi_in; iz++)

src/compute_property_grid.h

@@ -48,6 +48,7 @@ class ComputePropertyGrid : public Compute {
 
   int nxlo_in,nxhi_in,nylo_in,nyhi_in,nzlo_in,nzhi_in;
   int nxlo_out,nxhi_out,nylo_out,nyhi_out,nzlo_out,nzhi_out;
+  int ngridout;
 
   double **vec2d,***vec3d;
   double ***array2d,****array3d;

src/dump_grid.cpp

@@ -283,9 +283,9 @@ void DumpGrid::init_style()
         ifix = fix[field2index[i]];
         grid2d = (Grid2d *) ifix->get_grid_by_index(field2grid[i]);
       }
-      if (i == 0) grid2d->query_size(nx,ny);
+      if (i == 0) grid2d->get_size(nx,ny);
       else {
-        grid2d->query_size(nxtmp,nytmp);
+        grid2d->get_size(nxtmp,nytmp);
         if (nxtmp != nx || nytmp != ny) 
           error->all(FLERR,"Dump grid field grid sizes do not match");
       }
@@ -298,9 +298,9 @@ void DumpGrid::init_style()
         ifix = fix[field2index[i]];
         grid3d = (Grid3d *) ifix->get_grid_by_index(field2grid[i]);
       }
-      if (i == 0) grid3d->query_size(nx,ny,nz);
+      if (i == 0) grid3d->get_size(nx,ny,nz);
       else {
-        grid3d->query_size(nxtmp,nytmp,nztmp);
+        grid3d->get_size(nxtmp,nytmp,nztmp);
         if (nxtmp != nx || nytmp != ny || nztmp != nz)
           error->all(FLERR,"Dump grid field grid sizes do not match");
       }
@@ -535,7 +535,7 @@ int DumpGrid::count()
     else if (field2source[0] == FIX)
       grid2d = (Grid2d *) 
         fix[field2index[0]]->get_grid_by_index(field2grid[0]);
-    grid2d->query_bounds(nxlo_in,nxhi_in,nylo_in,nyhi_in);
+    grid2d->get_bounds(nxlo_in,nxhi_in,nylo_in,nyhi_in);
   } else {
     if (field2source[0] == COMPUTE)
       grid3d = (Grid3d *) 
@@ -543,7 +543,7 @@ int DumpGrid::count()
     else if (field2source[0] == FIX)
       grid3d = (Grid3d *) 
         fix[field2index[0]]->get_grid_by_index(field2grid[0]);
-    grid3d->query_bounds(nxlo_in,nxhi_in,nylo_in,nyhi_in,nzlo_in,nzhi_in);
+    grid3d->get_bounds(nxlo_in,nxhi_in,nylo_in,nyhi_in,nzlo_in,nzhi_in);
   }
 
   // invoke Computes for per-grid quantities
@@ -614,6 +614,7 @@ int DumpGrid::count()
 void DumpGrid::pack(tagint *ids)
 {
   for (int n = 0; n < size_one; n++) (this->*pack_choice[n])(n);
+
   // NOTE: this needs to be grid IDs ?
   /*
   if (ids) {
@@ -731,29 +732,10 @@ int DumpGrid::parse_fields(int narg, char **arg)
         vtype[iarg] = Dump::DOUBLE;
 	field2source[iarg] = COMPUTE;
         
-        // name = idcompute:gname:fname, split into 3 strings
+        // split name = idcompute:gname:dname into 3 strings
 
-        char *ptr1 = strchr((char *) name,':');
+        char *idcompute,*gname,*dname;
-        if (!ptr1) 
+        utils::grid_parse(FLERR,name,idcompute,gname,dname,error);
-          error->all(FLERR,"Dump grid fix {} does not contain 2 ':' chars");
-        *ptr1 = '\0';
-        char *ptr2 = strchr(ptr1+1,':');
-        if (!ptr2) 
-          error->all(FLERR,"Dump grid fix {} does not contain 2 ':' chars");
-        *ptr2 = '\0';
-
-        int n = strlen(name) + 1;
-        char *idcompute = new char[n];
-        strcpy(idcompute,name);
-        n = strlen(ptr1+1) + 1;
-        char *gname = new char[n];
-        strcpy(gname,ptr1+1);
-        n = strlen(ptr2+1) + 1;
-        char *dname = new char[n];
-        strcpy(dname,ptr2+1);
-
-        *ptr1 = ':';
-        *ptr2 = ':';
 
         icompute = modify->get_compute_by_id(idcompute);
         if (!icompute) 
@@ -806,29 +788,10 @@ int DumpGrid::parse_fields(int narg, char **arg)
         vtype[iarg] = Dump::DOUBLE;
 	field2source[iarg] = FIX;
 
-        // name = idfix:gname:fname, split into 3 strings
+        // split name = idfix:gname:dname into 3 strings
 
-        char *ptr1 = strchr((char *) name,':');
+        char *idfix,*gname,*dname;
-        if (!ptr1) 
+        utils::grid_parse(FLERR,name,idfix,gname,dname,error);
-          error->all(FLERR,"Dump grid fix {} does not contain 2 ':' chars");
-        *ptr1 = '\0';
-        char *ptr2 = strchr(ptr1+1,':');
-        if (!ptr2) 
-          error->all(FLERR,"Dump grid fix {} does not contain 2 ':' chars");
-        *ptr2 = '\0';
-
-        int n = strlen(name) + 1;
-        char *idfix = new char[n];
-        strcpy(idfix,name);
-        n = strlen(ptr1+1) + 1;
-        char *gname = new char[n];
-        strcpy(gname,ptr1+1);
-        n = strlen(ptr2+1) + 1;
-        char *dname = new char[n];
-        strcpy(dname,ptr2+1);
-
-        *ptr1 = ':';
-        *ptr2 = ':';
 
         ifix = modify->get_fix_by_id(idfix);
         if (!ifix) error->all(FLERR,"Could not find dump grid fix ID: {}",idfix);

src/grid2d.cpp

@@ -122,6 +122,12 @@ Grid2d::Grid2d(LAMMPS *lmp, MPI_Comm gcomm,
     oyhi = MIN(gnx-1,oyhi);
   }
 
+  // error check on size of grid stored by this proc
+
+  bigint total = (bigint) (oxhi - oxlo + 1) * (oyhi - oylo + 1);
+  if (total > MAXSMALLINT) 
+    error->one(FLERR, "Too many owned+ghost grid2d points");
+
   // store grid bounds and proc neighs
 
   if (layout == REGULAR) {
@@ -336,7 +342,7 @@ void Grid2d::store(int ixlo, int ixhi, int iylo, int iyhi,
 
 /* ---------------------------------------------------------------------- */
 
-void Grid2d::query_size(int &nxgrid, int &nygrid)
+void Grid2d::get_size(int &nxgrid, int &nygrid)
 {
   nxgrid = nx;
   nygrid = ny;
@@ -344,7 +350,7 @@ void Grid2d::query_size(int &nxgrid, int &nygrid)
 
 /* ---------------------------------------------------------------------- */
 
-void Grid2d::query_bounds(int &xlo, int &xhi, int &ylo, int &yhi)
+void Grid2d::get_bounds(int &xlo, int &xhi, int &ylo, int &yhi)
 {
   xlo = inxlo;
   xhi = inxhi;
@@ -354,7 +360,7 @@ void Grid2d::query_bounds(int &xlo, int &xhi, int &ylo, int &yhi)
 
 /* ---------------------------------------------------------------------- */
 
-void Grid2d::query_box(int dim, double &lo, double &delta)
+void Grid2d::get_box(int dim, double &lo, double &delta)
 {
   lo = boxlo[dim];
   delta = prd[dim] / ngrid[dim];

src/grid2d.h

@@ -29,9 +29,9 @@ class Grid2d : protected Pointers {
   Grid2d(class LAMMPS *, MPI_Comm, int, int, int, int, int, int, int, int, int, int,
 	 int, int, int, int, int);
   ~Grid2d() override;
-  void query_size(int &, int &);
+  void get_size(int &, int &);
-  void query_bounds(int &, int &, int &, int &);
+  void get_bounds(int &, int &, int &, int &);
-  void query_box(int, double &, double &);
+  void get_box(int, double &, double &);
   void setup(int &, int &);
   int ghost_adjacent();
   void forward_comm(int, void *, int, int, int, void *, void *, MPI_Datatype);

src/grid3d.cpp

@@ -133,6 +133,13 @@ Grid3d::Grid3d(LAMMPS *lmp, MPI_Comm gcomm,
     ozhi = MIN(gnx-1,ozhi);
   }
   
+  // error check on size of grid stored by this proc
+
+  bigint total = (bigint) 
+    (oxhi - oxlo + 1) * (oyhi - oylo + 1) * (ozhi - ozlo + 1);
+  if (total > MAXSMALLINT) 
+    error->one(FLERR, "Too many owned+ghost grid3d points");
+
   // store grid bounds and proc neighs
 
   if (layout == REGULAR) {
@@ -364,7 +371,7 @@ void Grid3d::store(int ixlo, int ixhi, int iylo, int iyhi,
 
 /* ---------------------------------------------------------------------- */
 
-void Grid3d::query_size(int &nxgrid, int &nygrid, int &nzgrid)
+void Grid3d::get_size(int &nxgrid, int &nygrid, int &nzgrid)
 {
   nxgrid = nx;
   nygrid = ny;
@@ -373,7 +380,7 @@ void Grid3d::query_size(int &nxgrid, int &nygrid, int &nzgrid)
 
 /* ---------------------------------------------------------------------- */
 
-void Grid3d::query_bounds(int &xlo, int &xhi, int &ylo, int &yhi, 
+void Grid3d::get_bounds(int &xlo, int &xhi, int &ylo, int &yhi, 
                         int &zlo, int &zhi)
 {
   xlo = inxlo;
@@ -386,7 +393,7 @@ void Grid3d::query_bounds(int &xlo, int &xhi, int &ylo, int &yhi,
 
 /* ---------------------------------------------------------------------- */
 
-void Grid3d::query_box(int dim, double &lo, double &delta)
+void Grid3d::get_box(int dim, double &lo, double &delta)
 {
   lo = boxlo[dim];
   delta = prd[dim] / ngrid[dim];

src/grid3d.h

@@ -31,9 +31,9 @@ class Grid3d : protected Pointers {
          int, int, int, int, int, int, int, int, int, int, int, int, 
          int, int, int, int, int, int);
   ~Grid3d() override;
-  void query_size(int &, int &, int &);
+  void get_size(int &, int &, int &);
-  void query_bounds(int &, int &, int &, int &, int &, int &);
+  void get_bounds(int &, int &, int &, int &, int &, int &);
-  void query_box(int, double &, double &);
+  void get_box(int, double &, double &);
   void setup(int &, int &);
   int ghost_adjacent();
   void forward_comm(int, void *, int, int, int, void *, void *, MPI_Datatype);

src/memory.h

@@ -36,6 +36,10 @@ class Memory : protected Pointers {
    for these other cases, use smalloc/srealloc/sfree directly
 ------------------------------------------------------------------------- */
 
+// -------------------------------------------------------------------------
+// -------------------------------------------------------------------------
+// -------------------------------------------------------------------------
+
 /* ----------------------------------------------------------------------
    create a 1d array
 ------------------------------------------------------------------------- */
@@ -91,6 +95,10 @@ class Memory : protected Pointers {
     array = nullptr;
   }
 
+// -------------------------------------------------------------------------
+// -------------------------------------------------------------------------
+// -------------------------------------------------------------------------
+
 /* ----------------------------------------------------------------------
    create a 1d array with index from nlo to nhi inclusive
    cannot grow it
@@ -124,6 +132,10 @@ class Memory : protected Pointers {
     array = nullptr;
   }
 
+// -------------------------------------------------------------------------
+// -------------------------------------------------------------------------
+// -------------------------------------------------------------------------
+
 /* ----------------------------------------------------------------------
    create a 2d array
 ------------------------------------------------------------------------- */
@@ -200,6 +212,10 @@ class Memory : protected Pointers {
     array = nullptr;
   }
 
+// -------------------------------------------------------------------------
+// -------------------------------------------------------------------------
+// -------------------------------------------------------------------------
+
 /* ----------------------------------------------------------------------
    create a 2d array with a ragged 2nd dimension
 ------------------------------------------------------------------------- */
@@ -259,6 +275,10 @@ class Memory : protected Pointers {
     return array;
   }
 
+// -------------------------------------------------------------------------
+// -------------------------------------------------------------------------
+// -------------------------------------------------------------------------
+
   /* ----------------------------------------------------------------------
    create a 2d array with 2nd index from n2lo to n2hi inclusive
    cannot grow it
@@ -296,6 +316,57 @@ class Memory : protected Pointers {
     array = nullptr;
   }
 
+// -------------------------------------------------------------------------
+// -------------------------------------------------------------------------
+// -------------------------------------------------------------------------
+
+/* ----------------------------------------------------------------------
+   create a 2d array with both indices offset
+   1st index from n1lo to n1hi inclusive
+   2nd index from n2lo to n2hi inclusive
+   cannot grow it
+------------------------------------------------------------------------- */
+
+  template <typename TYPE>
+  TYPE **create2d_offset(TYPE **&array, int n1lo, int n1hi, int n2lo, int n2hi,
+                         const char *name)
+  {
+    if (n1lo > n1hi || n2lo > n2hi) return nullptr;
+
+    int n1 = n1hi - n1lo + 1;
+    int n2 = n2hi - n2lo + 1;
+    create(array, n1, n2, name);
+
+    for (int i = 0; i < n1; i++) array[i] -= n2lo;
+    array -= n1lo;
+    return array;
+  }
+
+  template <typename TYPE>
+  TYPE ***create2d_offset(TYPE ***& /*array*/, int /*n1lo*/, int /*n1hi*/, int /*n2lo*/,
+                           int /*n2hi*/, const char *name)
+  {
+    fail(name);
+    return nullptr;
+  }
+
+  /* ----------------------------------------------------------------------
+   free a 2d array with both indices offset
+------------------------------------------------------------------------- */
+
+  template <typename TYPE>
+  void destroy2d_offset(TYPE **&array, int n1_offset, int n2_offset)
+  {
+    if (array == nullptr) return;
+    sfree(&array[n1_offset][n2_offset]);
+    sfree(&array[n1_offset]);
+    array = nullptr;
+  }
+
+// -------------------------------------------------------------------------
+// -------------------------------------------------------------------------
+// -------------------------------------------------------------------------
+
   /* ----------------------------------------------------------------------
    create a 3d array
 ------------------------------------------------------------------------- */
@@ -389,6 +460,10 @@ class Memory : protected Pointers {
     array = nullptr;
   }
 
+// -------------------------------------------------------------------------
+// -------------------------------------------------------------------------
+// -------------------------------------------------------------------------
+
   /* ----------------------------------------------------------------------
    create a 3d array with 1st index from n1lo to n1hi inclusive
    cannot grow it
@@ -426,6 +501,10 @@ class Memory : protected Pointers {
     array = nullptr;
   }
 
+// -------------------------------------------------------------------------
+// -------------------------------------------------------------------------
+// -------------------------------------------------------------------------
+
 /* ----------------------------------------------------------------------
    create a 3d array with all 3 indices offset
    1st index from n1lo to n1hi inclusive
@@ -474,6 +553,57 @@ class Memory : protected Pointers {
     array = nullptr;
   }
 
+// -------------------------------------------------------------------------
+// -------------------------------------------------------------------------
+// -------------------------------------------------------------------------
+
+/* ----------------------------------------------------------------------
+   create a 3d array with indices 1,2 offset, but not last
+   1st index from n1lo to n1hi inclusive
+   2nd index from n2lo to n2hi inclusive
+   cannot grow it
+------------------------------------------------------------------------- */
+
+  template <typename TYPE>
+  TYPE ***create3d_offset_last(TYPE ***&array, int n1lo, int n1hi,
+                               int n2lo, int n2hi, int n3,
+                               const char *name)
+  {
+    if (n1lo > n1hi || n2lo > n2hi || n3 < 0) return nullptr;
+
+    int n1 = n1hi - n1lo + 1;
+    int n2 = n2hi - n2lo + 1;
+    create(array,n1,n2,n3,name);
+
+    for (int i = 0; i < n1; i++) array[i] -= n2lo;
+    array -= n1lo;
+    return array;
+  }
+
+  template <typename TYPE>
+  TYPE ***create3d_offset_last(TYPE ****& /*array*/, int /*n1lo*/, int /*n1hi*/,
+                               int /*n2lo*/, int /*n2hi*/, int /*n3*/,
+                               const char *name)
+  {fail(name); return nullptr;}
+
+/* ----------------------------------------------------------------------
+   free a 3d array with indices 1,2 offset, but not last
+------------------------------------------------------------------------- */
+
+  template <typename TYPE>
+  void destroy3d_offset_last(TYPE ***&array, int n1_offset, int n2_offset)
+  {
+    if (array == nullptr) return;
+    sfree(&array[n1_offset][n2_offset][0]);
+    sfree(&array[n1_offset][n2_offset]);
+    sfree(&array[n1_offset]);
+    array = nullptr;
+  }
+
+// -------------------------------------------------------------------------
+// -------------------------------------------------------------------------
+// -------------------------------------------------------------------------
+
   /* ----------------------------------------------------------------------
    create a 4d array
 ------------------------------------------------------------------------- */
@@ -588,8 +718,12 @@ class Memory : protected Pointers {
     array = nullptr;
   }
 
+// -------------------------------------------------------------------------
+// -------------------------------------------------------------------------
+// -------------------------------------------------------------------------
+
 /* ----------------------------------------------------------------------
-   create a 4d array with indices 2,3,4 offset
+   create a 4d array with indices 2,3,4 offset, but not first
    2nd index from n2lo to n2hi inclusive
    3rd index from n3lo to n3hi inclusive
    4th index from n4lo to n4hi inclusive
@@ -624,7 +758,7 @@ class Memory : protected Pointers {
   }
 
   /* ----------------------------------------------------------------------
-   free a 4d array with indices 2,3,4 offset
+   free a 4d array with indices 2,3,4 offset, but not first
 ------------------------------------------------------------------------- */
 
   template <typename TYPE>
@@ -638,11 +772,15 @@ class Memory : protected Pointers {
     array = nullptr;
   }
 
+// -------------------------------------------------------------------------
+// -------------------------------------------------------------------------
+// -------------------------------------------------------------------------
+
 /* ----------------------------------------------------------------------
-   create a 4d array with indices 1,2,3 offset
+   create a 4d array with indices 1,2,3 offset, but not last
-   1st index from n1lo to n2hi inclusive
+   1st index from n1lo to n1hi inclusive
-   2nd index from n3lo to n3hi inclusive
+   2nd index from n2lo to n2hi inclusive
-   3rd index from n4lo to n4hi inclusive
+   3rd index from n3lo to n3hi inclusive
    cannot grow it
 ------------------------------------------------------------------------- */
 
@@ -673,7 +811,7 @@ class Memory : protected Pointers {
   {fail(name); return nullptr;}
 
 /* ----------------------------------------------------------------------
-   free a 4d array with indices 1,2,3 offset
+   free a 4d array with indices 1,2,3 offset, but not last
 ------------------------------------------------------------------------- */
 
   template <typename TYPE>
@@ -688,6 +826,10 @@ class Memory : protected Pointers {
     array = nullptr;
   }
 
+// -------------------------------------------------------------------------
+// -------------------------------------------------------------------------
+// -------------------------------------------------------------------------
+
 /* ----------------------------------------------------------------------
    create a 5d array
 ------------------------------------------------------------------------- */
@@ -757,6 +899,10 @@ class Memory : protected Pointers {
     array = nullptr;
   }
 
+// -------------------------------------------------------------------------
+// -------------------------------------------------------------------------
+// -------------------------------------------------------------------------
+
 /* ----------------------------------------------------------------------
    memory usage of arrays, including pointers
 ------------------------------------------------------------------------- */

src/utils.cpp

@@ -766,6 +766,30 @@ int utils::expand_args(const char *file, int line, int narg, char **arg, int mod
   return newarg;
 }
 
+/* ----------------------------------------------------------------------
+   Parse grid reference into id:gridname:dataname
+   return ptrs to 3 substrings
+------------------------------------------------------------------------- */
+
+void utils::grid_parse(const char *file, int line, const std::string &name, 
+                       char *id, char *gridname, char *dataname, Error *error)
+{
+  char *copy = strdup(name);
+
+  char *ptr1 = strchr(copy,':');
+  if (!ptr1) 
+    error->all(FLERR,"Grid reference {} does not contain 2 ':' chars",name);
+  *ptr1 = '\0';
+  char *ptr2 = strchr(ptr1+1,':');
+  if (!ptr2) 
+    error->all(FLERR,"Grid reference {} does not contain 2 ':' chars",name);
+  *ptr2 = '\0';
+
+  id = strdup(copy);
+  gridname = strdup(ptr1+1);
+  dataname = strdup(ptr2+1);
+}
+
 /* ----------------------------------------------------------------------
    Make copy of string in new storage. Works like the (non-portable)
    C-style strdup() but also accepts a C++ string as argument.
@@ -1087,6 +1111,7 @@ std::vector<std::string> utils::split_words(const std::string &text)
 /* ----------------------------------------------------------------------
    Convert multi-line string into lines
 ------------------------------------------------------------------------- */
+
 std::vector<std::string> utils::split_lines(const std::string &text)
 {
   return Tokenizer(text, "\r\n").as_vector();
@@ -1204,6 +1229,7 @@ std::string utils::get_potential_units(const std::string &path, const std::strin
 /* ----------------------------------------------------------------------
    return bitmask of supported conversions for a given property
 ------------------------------------------------------------------------- */
+
 int utils::get_supported_conversions(const int property)
 {
   if (property == ENERGY)

src/utils.h

@@ -368,6 +368,19 @@ namespace utils {
    * \param text  string that should be copied
    * \return new buffer with copy of string */
 
+  void grid_parse(const char *file, int line, const std::string &name,
+                  char *id, char *gridname, char *dataname, Error *error);
+
+  /*! Parse grid reference name into 3 sub-strings
+   *
+   * Format of grid reference name = id:gridname:dataname
+   * Return ptrs to the 3 sub-strings
+   *
+   * \param name = grid reference
+   * \param id = ptr to 1st substring
+   * \param gridname = ptr to 2nd substring
+   * \param dataname = ptr to 3rd substring */
+
   char *strdup(const std::string &text);
 
   /*! Convert string to lowercase