CarpetLib: Introduce gdata::allocated_memory_shape

Introduce a function gdata::allocated_memory_shape that determines from a given grid size how many grid points should be allocated. This allows adding padding to grid sizes, e.g. to improve cache performance.

Use this function when checking grid sizes.

18 files changed, 145 insertions, 58 deletions

diff --git a/Carpet/CarpetLib/src/accumulate_3d.cc b/Carpet/CarpetLib/src/accumulate_3d.cc
index 8fef2e616..e2b3883a7 100644
--- a/Carpet/CarpetLib/src/accumulate_3d.cc
+++ b/Carpet/CarpetLib/src/accumulate_3d.cc
@@ -7,6 +7,7 @@
 #include <cstdlib>
 #include <iostream>
 
+#include "gdata.hh"
 #include "operator_prototypes_3d.hh"
 #include "typeprops.hh"
 
@@ -74,8 +75,8 @@ namespace CarpetLib {
       }
     }
     
-    if (any (srcext != srcbbox.shape() / srcbbox.stride() or
-             dstext != dstbbox.shape() / dstbbox.stride()))
+    if (any (srcext != gdata::allocated_memory_shape(srcbbox.shape() / srcbbox.stride()) or
+             dstext != gdata::allocated_memory_shape(dstbbox.shape() / dstbbox.stride())))
     {
 #pragma omp critical
       CCTK_WARN (0, "Internal error: array sizes don't agree with bounding boxes");
diff --git a/Carpet/CarpetLib/src/copy_3d.cc b/Carpet/CarpetLib/src/copy_3d.cc
index b7fb72e40..dee9dac70 100644
--- a/Carpet/CarpetLib/src/copy_3d.cc
+++ b/Carpet/CarpetLib/src/copy_3d.cc
@@ -8,6 +8,7 @@
 #include <cstring>
 #include <iostream>
 
+#include "gdata.hh"
 #include "operator_prototypes_3d.hh"
 #include "typeprops.hh"
 
@@ -75,8 +76,8 @@ namespace CarpetLib {
       }
     }
     
-    if (any (srcext != srcbbox.shape() / srcbbox.stride() or
-             dstext != dstbbox.shape() / dstbbox.stride()))
+    if (any (srcext != gdata::allocated_memory_shape(srcbbox.shape() / srcbbox.stride()) or
+             dstext != gdata::allocated_memory_shape(dstbbox.shape() / dstbbox.stride())))
     {
 #pragma omp critical
       CCTK_WARN (0, "Internal error: array sizes don't agree with bounding boxes");
diff --git a/Carpet/CarpetLib/src/copy_4d.cc b/Carpet/CarpetLib/src/copy_4d.cc
index dfebb5acc..376b81428 100644
--- a/Carpet/CarpetLib/src/copy_4d.cc
+++ b/Carpet/CarpetLib/src/copy_4d.cc
@@ -7,6 +7,7 @@
 #include <cstdlib>
 #include <iostream>
 
+#include "gdata.hh"
 #include "operator_prototypes_4d.hh"
 #include "typeprops.hh"
 
@@ -63,8 +64,8 @@ namespace CarpetLib {
       CCTK_WARN (0, "Internal error: region extent is not contained in array extent");
     }
     
-    if (any (srcext != srcbbox.shape() / srcbbox.stride() or
-             dstext != dstbbox.shape() / dstbbox.stride()))
+    if (any (srcext != gdata::allocated_memory_shape(srcbbox.shape() / srcbbox.stride()) or
+             dstext != gdata::allocated_memory_shape(dstbbox.shape() / dstbbox.stride())))
     {
       CCTK_WARN (0, "Internal error: array sizes don't agree with bounding boxes");
     }
diff --git a/Carpet/CarpetLib/src/gdata.cc b/Carpet/CarpetLib/src/gdata.cc
index e6c76050c..958153784 100644
--- a/Carpet/CarpetLib/src/gdata.cc
+++ b/Carpet/CarpetLib/src/gdata.cc
@@ -3,6 +3,8 @@
 #include <util_ErrorCodes.h>
 #include <util_Table.h>
 
+#include <vectors.h>
+
 #include <cassert>
 #include <cstdlib>
 #include <iomanip>
@@ -69,22 +71,67 @@ gdata::~gdata ()
 
 // Storage management
 
-ivect
+template<int D>
+vect<int,D>
 gdata::
-allocated_memory_shape (ibbox const& extent)
+allocated_memory_shape (bbox<int,D> const& extent)
+{
+  vect<int,D> const shape =
+    max (vect<int,D>(0), extent.shape() / extent.stride());
+  return allocated_memory_shape (shape);
+}
+
+template<int D>
+vect<int,D>
+gdata::
+allocated_memory_shape (vect<int,D> shape)
 {
   DECLARE_CCTK_PARAMETERS;
-  ivect shape = max (ivect(0), extent.shape() / extent.stride());
-  // Enlarge shape to avoid multiples of cache line colours
-  if (avoid_arraysize_bytes > 0) {
-    for (int d=0; d<dim; ++d) {
-      if (shape[d] > 0 and
-          shape[d] * sizeof(CCTK_REAL) % avoid_arraysize_bytes == 0)
-      {
-        ++shape[d];
+  bool did_change;
+  for (int count=0; count<10; ++count) {
+    did_change = false;
+    int magic_size = -1;
+    // Enlarge shape to avoid multiples of cache line colours
+    if (avoid_arraysize_bytes > 0) {
+      magic_size = avoid_arraysize_bytes;
+      if (avoid_arraysize_bytes == sizeof(CCTK_REAL)) {
+        CCTK_VWarn (CCTK_WARN_ABORT, __LINE__, __FILE__, CCTK_THORNSTRING,
+                    "The run-time parameter avoid_arraysize_bytes=%d requires that the size of a grid variable in the x direction is not a multiple of %d bytes, but the size of CCTK_REAL is %d. This is inconsistent -- aborting",
+                    int(avoid_arraysize_bytes),
+                    int(avoid_arraysize_bytes),
+                    int(sizeof(CCTK_REAL)));
+      }
+      for (int d=0; d<D; ++d) {
+        if (shape[d] > 0 and
+            shape[d] * sizeof(CCTK_REAL) % avoid_arraysize_bytes == 0)
+        {
+          ++shape[d];
+          did_change = true;
+        }
       }
     }
+#if defined(VECTORISE_ALIGNED_ARRAYS)
+    // Enlarge shape in the x direction to ensure it is a multiple of
+    // the vector size
+#warning "TODO: Support other datatypes as well, don't target only CCTK_REAL"
+    if (sizeof(CCTK_REAL) * CCTK_REAL_VEC_SIZE == magic_size) {
+      CCTK_VWarn (CCTK_WARN_ABORT, __LINE__, __FILE__, CCTK_THORNSTRING,
+                  "The build-time option KRANC_ALIGNED_ARRAYS requires that the size of a grid variable in the x direction is a multiple of %d bytes, while the run-time parameter avoid_arraysize_bytes=%d requests the opposite. This is inconsistent -- aborting",
+                  int(avoid_arraysize_bytes),
+                  int(avoid_arraysize_bytes));
+    }
+    if (shape[0] % CCTK_REAL_VEC_SIZE != 0) {
+      shape[0] =
+        (shape[0] + CCTK_REAL_VEC_SIZE - 1) /
+        CCTK_REAL_VEC_SIZE *
+        CCTK_REAL_VEC_SIZE;
+      did_change = true;
+    }
+#endif
+    if (not did_change) goto done;
   }
+  CCTK_WARN (CCTK_WARN_ABORT, "endless loop");
+ done:;
   return shape;
 }
 
@@ -168,8 +215,8 @@ transfer_from (comm_state & state,
   find_source_timelevel
     (times, time, order_time, my_transport_operator, timelevel0, ntimelevels);
   if (is_src) assert (int (srcs.size()) >= ntimelevels);
-  int const dstpoints = dstbox.size();
-  int const srcpoints = srcbox.size() * ntimelevels;
+  int const dstpoints = prod(allocated_memory_shape(dstbox));
+  int const srcpoints = prod(allocated_memory_shape(srcbox)) * ntimelevels;
   bool const interp_on_src = dstpoints <= srcpoints;
   int const npoints = interp_on_src ? dstpoints : srcpoints;
   
@@ -202,29 +249,36 @@ transfer_from (comm_state & state,
             ioffset = 1 - is_centered;
           }
           ibbox const bufbox = dstbox.shift(ioffset, 2);
-          assert (bufbox.size() == dstbox.size());
-          size_t const sendbufsize = src->c_datatype_size() * dstbox.size();
+          assert (prod(allocated_memory_shape(bufbox)) ==
+                  prod(allocated_memory_shape(dstbox)));
+          size_t const sendbufsize =
+            src->c_datatype_size() * prod(allocated_memory_shape(dstbox));
           void * const sendbuf =
-            state.send_buffer (src->c_datatype(), dstproc, dstbox.size());
+            state.send_buffer (src->c_datatype(), dstproc,
+                               prod(allocated_memory_shape(dstbox)));
           gdata * const buf =
             src->make_typed (src->varindex, src->cent, src->transport_operator);
           buf->allocate (bufbox, srcproc, sendbuf, sendbufsize);
           buf->transfer_from_innerloop
             (srcs, times, dstbox, time, order_space, order_time);
           delete buf;
-          state.commit_send_space (src->c_datatype(), dstproc, dstbox.size());
+          state.commit_send_space
+            (src->c_datatype(), dstproc, prod(allocated_memory_shape(dstbox)));
         } else {
           for (int tl = timelevel0; tl < timelevel0 + ntimelevels; ++ tl) {
-            size_t const sendbufsize = src->c_datatype_size() * srcbox.size();
+            size_t const sendbufsize =
+              src->c_datatype_size() * prod(allocated_memory_shape(srcbox));
             void * const sendbuf =
-              state.send_buffer (src->c_datatype(), dstproc, srcbox.size());
+              state.send_buffer (src->c_datatype(), dstproc,
+                                 prod(allocated_memory_shape(srcbox)));
             gdata * const buf =
               src->make_typed (src->varindex, src->cent,
                                src->transport_operator);
             buf->allocate (srcbox, srcproc, sendbuf, sendbufsize);
             buf->copy_from_innerloop (srcs.AT(tl), srcbox);
             delete buf;
-            state.commit_send_space (src->c_datatype(), dstproc, srcbox.size());
+            state.commit_send_space (src->c_datatype(), dstproc,
+                                     prod(allocated_memory_shape(srcbox)));
           }
         }
       }
@@ -244,12 +298,15 @@ transfer_from (comm_state & state,
       if (is_dst) {
         // copy from the recv buffer
         if (interp_on_src) {
-          size_t const recvbufsize = c_datatype_size() * dstbox.size();
+          size_t const recvbufsize =
+            c_datatype_size() * prod(allocated_memory_shape(dstbox));
           void * const recvbuf =
-            state.recv_buffer (c_datatype(), srcproc, dstbox.size());
+            state.recv_buffer (c_datatype(), srcproc,
+                               prod(allocated_memory_shape(dstbox)));
           gdata * const buf = make_typed (varindex, cent, transport_operator);
           buf->allocate (dstbox, dstproc, recvbuf, recvbufsize);
-          state.commit_recv_space (c_datatype(), srcproc, dstbox.size());
+          state.commit_recv_space (c_datatype(), srcproc,
+                                   prod(allocated_memory_shape(dstbox)));
           copy_from_innerloop (buf, dstbox);
           delete buf;
         } else {
@@ -264,12 +321,15 @@ transfer_from (comm_state & state,
           vector <gdata const *> bufs (ntimelevels, null);
           vector <CCTK_REAL> timebuf (ntimelevels);
           for (int tl = 0; tl < ntimelevels; ++ tl) {
-            size_t const recvbufsize = c_datatype_size() * srcbox.size();
+            size_t const recvbufsize =
+              c_datatype_size() * prod(allocated_memory_shape(srcbox));
             void * const recvbuf =
-              state.recv_buffer (c_datatype(), srcproc, srcbox.size());
+              state.recv_buffer (c_datatype(), srcproc,
+                                 prod(allocated_memory_shape(srcbox)));
             gdata * const buf = make_typed (varindex, cent, transport_operator);
             buf->allocate (srcbox, dstproc, recvbuf, recvbufsize);
-            state.commit_recv_space (c_datatype(), srcproc, srcbox.size());
+            state.commit_recv_space
+              (c_datatype(), srcproc, prod(allocated_memory_shape(srcbox)));
             bufs.AT(tl) = buf;
             timebuf.AT(tl) = times.AT(timelevel0 + tl);
           }
@@ -377,3 +437,11 @@ allmemory ()
   }
   return mem;
 }
+
+
+
+template vect<int,3> gdata::allocated_memory_shape (bbox<int,3> const& extent);
+template vect<int,3> gdata::allocated_memory_shape (vect<int,3> shape);
+
+template vect<int,4> gdata::allocated_memory_shape (bbox<int,4> const& extent);
+template vect<int,4> gdata::allocated_memory_shape (vect<int,4> shape);
diff --git a/Carpet/CarpetLib/src/gdata.hh b/Carpet/CarpetLib/src/gdata.hh
index 634f5a917..a06509bd7 100644
--- a/Carpet/CarpetLib/src/gdata.hh
+++ b/Carpet/CarpetLib/src/gdata.hh
@@ -75,10 +75,14 @@ public:
   
   // Storage management
   virtual void allocate (const ibbox& extent, const int proc,
-			 void* const memptr = NULL, size_t const memsize = 0) = 0;
+			 void* const memptr = NULL, size_t const memsize = 0) =
+    0;
   virtual void free () = 0;
   virtual size_t allocsize (const ibbox& extent, const int proc) const = 0;
-  static ivect allocated_memory_shape (ibbox const& extent);
+  template<int D>
+  static vect<int,D> allocated_memory_shape (bbox<int,D> const& extent);
+  template<int D>
+  static vect<int,D> allocated_memory_shape (vect<int,D> shape);
   
   // Accessors
   bool has_storage () const {
diff --git a/Carpet/CarpetLib/src/interpolate_3d_2tl.cc b/Carpet/CarpetLib/src/interpolate_3d_2tl.cc
index 3627dee7d..754b3a48c 100644
--- a/Carpet/CarpetLib/src/interpolate_3d_2tl.cc
+++ b/Carpet/CarpetLib/src/interpolate_3d_2tl.cc
@@ -8,6 +8,7 @@
 
 #include <loopcontrol.h>
 
+#include "gdata.hh"
 #include "operator_prototypes_3d.hh"
 #include "typeprops.hh"
 
@@ -68,8 +69,8 @@ namespace CarpetLib {
       CCTK_WARN (0, "Internal error: region extent is not contained in array extent");
     }
     
-    if (any (srcext != srcbbox.shape() / srcbbox.stride() or
-             dstext != dstbbox.shape() / dstbbox.stride()))
+    if (any (srcext != gdata::allocated_memory_shape(srcbbox.shape() / srcbbox.stride()) or
+             dstext != gdata::allocated_memory_shape(dstbbox.shape() / dstbbox.stride())))
     {
       CCTK_WARN (0, "Internal error: array sizes don't agree with bounding boxes");
     }
diff --git a/Carpet/CarpetLib/src/interpolate_3d_3tl.cc b/Carpet/CarpetLib/src/interpolate_3d_3tl.cc
index 28527ce3e..f0b2a78a2 100644
--- a/Carpet/CarpetLib/src/interpolate_3d_3tl.cc
+++ b/Carpet/CarpetLib/src/interpolate_3d_3tl.cc
@@ -8,6 +8,7 @@
 
 #include <loopcontrol.h>
 
+#include "gdata.hh"
 #include "operator_prototypes_3d.hh"
 #include "typeprops.hh"
 
@@ -70,8 +71,8 @@ namespace CarpetLib {
       CCTK_WARN (0, "Internal error: region extent is not contained in array extent");
     }
     
-    if (any (srcext != srcbbox.shape() / srcbbox.stride() or
-             dstext != dstbbox.shape() / dstbbox.stride()))
+    if (any (srcext != gdata::allocated_memory_shape(srcbbox.shape() / srcbbox.stride()) or
+             dstext != gdata::allocated_memory_shape(dstbbox.shape() / dstbbox.stride())))
     {
       CCTK_WARN (0, "Internal error: array sizes don't agree with bounding boxes");
     }
diff --git a/Carpet/CarpetLib/src/interpolate_3d_4tl.cc b/Carpet/CarpetLib/src/interpolate_3d_4tl.cc
index c5d44ca82..0e5439893 100644
--- a/Carpet/CarpetLib/src/interpolate_3d_4tl.cc
+++ b/Carpet/CarpetLib/src/interpolate_3d_4tl.cc
@@ -8,6 +8,7 @@
 
 #include <loopcontrol.h>
 
+#include "gdata.hh"
 #include "operator_prototypes_3d.hh"
 #include "typeprops.hh"
 
@@ -72,8 +73,8 @@ namespace CarpetLib {
       CCTK_WARN (0, "Internal error: region extent is not contained in array extent");
     }
     
-    if (any (srcext != srcbbox.shape() / srcbbox.stride() or
-             dstext != dstbbox.shape() / dstbbox.stride()))
+    if (any (srcext != gdata::allocated_memory_shape(srcbbox.shape() / srcbbox.stride()) or
+             dstext != gdata::allocated_memory_shape(dstbbox.shape() / dstbbox.stride())))
     {
       CCTK_WARN (0, "Internal error: array sizes don't agree with bounding boxes");
     }
diff --git a/Carpet/CarpetLib/src/interpolate_3d_5tl.cc b/Carpet/CarpetLib/src/interpolate_3d_5tl.cc
index 431351b37..2cb033f37 100644
--- a/Carpet/CarpetLib/src/interpolate_3d_5tl.cc
+++ b/Carpet/CarpetLib/src/interpolate_3d_5tl.cc
@@ -8,6 +8,7 @@
 
 #include <loopcontrol.h>
 
+#include "gdata.hh"
 #include "operator_prototypes_3d.hh"
 #include "typeprops.hh"
 
@@ -74,8 +75,8 @@ namespace CarpetLib {
       CCTK_WARN (0, "Internal error: region extent is not contained in array extent");
     }
     
-    if (any (srcext != srcbbox.shape() / srcbbox.stride() or
-             dstext != dstbbox.shape() / dstbbox.stride()))
+    if (any (srcext != gdata::allocated_memory_shape(srcbbox.shape() / srcbbox.stride()) or
+             dstext != gdata::allocated_memory_shape(dstbbox.shape() / dstbbox.stride())))
     {
       CCTK_WARN (0, "Internal error: array sizes don't agree with bounding boxes");
     }
diff --git a/Carpet/CarpetLib/src/interpolate_eno_3d_3tl.cc b/Carpet/CarpetLib/src/interpolate_eno_3d_3tl.cc
index bd5f53129..8682eb1a0 100644
--- a/Carpet/CarpetLib/src/interpolate_eno_3d_3tl.cc
+++ b/Carpet/CarpetLib/src/interpolate_eno_3d_3tl.cc
@@ -8,6 +8,7 @@
 
 #include <loopcontrol.h>
 
+#include "gdata.hh"
 #include "operator_prototypes_3d.hh"
 #include "typeprops.hh"
 
@@ -88,8 +89,8 @@ namespace CarpetLib {
       CCTK_WARN (0, "Internal error: region extent is not contained in array extent");
     }
     
-    if (any (srcext != srcbbox.shape() / srcbbox.stride() or
-             dstext != dstbbox.shape() / dstbbox.stride()))
+    if (any (srcext != gdata::allocated_memory_shape(srcbbox.shape() / srcbbox.stride()) or
+             dstext != gdata::allocated_memory_shape(dstbbox.shape() / dstbbox.stride())))
     {
       CCTK_WARN (0, "Internal error: array sizes don't agree with bounding boxes");
     }
diff --git a/Carpet/CarpetLib/src/mem.hh b/Carpet/CarpetLib/src/mem.hh
index 5fdc70018..04748f436 100644
--- a/Carpet/CarpetLib/src/mem.hh
+++ b/Carpet/CarpetLib/src/mem.hh
@@ -33,7 +33,7 @@ public:
 
 
 template<typename T>
-class mem: gmem
+class mem: public gmem
 {
   T * storage_;
   size_t nelems_;
diff --git a/Carpet/CarpetLib/src/prolongate_3d_o5_monotone_rf2.cc b/Carpet/CarpetLib/src/prolongate_3d_o5_monotone_rf2.cc
index dac0870d2..b6c9b8d90 100644
--- a/Carpet/CarpetLib/src/prolongate_3d_o5_monotone_rf2.cc
+++ b/Carpet/CarpetLib/src/prolongate_3d_o5_monotone_rf2.cc
@@ -23,6 +23,7 @@
 #include <cmath>
 #include <cstdlib>
 
+#include "gdata.hh"
 #include "operator_prototypes_3d.hh"
 #include "typeprops.hh"
 
@@ -132,8 +133,8 @@ namespace CarpetLib {
       CCTK_WARN (0, "Internal error: region extent is not contained in array extent");
     }
     
-    if (any (srcext != srcbbox.shape() / srcbbox.stride() or
-             dstext != dstbbox.shape() / dstbbox.stride()))
+    if (any (srcext != gdata::allocated_memory_shape(srcbbox.shape() / srcbbox.stride()) or
+             dstext != gdata::allocated_memory_shape(dstbbox.shape() / dstbbox.stride())))
     {
       CCTK_WARN (0, "Internal error: array sizes don't agree with bounding boxes");
     }
diff --git a/Carpet/CarpetLib/src/prolongate_3d_rf2.cc b/Carpet/CarpetLib/src/prolongate_3d_rf2.cc
index 3b632eb6c..4f669d6ea 100644
--- a/Carpet/CarpetLib/src/prolongate_3d_rf2.cc
+++ b/Carpet/CarpetLib/src/prolongate_3d_rf2.cc
@@ -8,6 +8,7 @@
 #include <cstdlib>
 #include <iostream>
 
+#include "gdata.hh"
 #include "vectors.h"
 #include "operator_prototypes_3d.hh"
 #include "typeprops.hh"
@@ -419,8 +420,8 @@ namespace CarpetLib {
       CCTK_WARN (0, "Internal error: region extent is not contained in array extent");
     }
     
-    if (any (srcext != srcbbox.shape() / srcbbox.stride() or
-             dstext != dstbbox.shape() / dstbbox.stride()))
+    if (any (srcext != gdata::allocated_memory_shape(srcbbox.shape() / srcbbox.stride()) or
+             dstext != gdata::allocated_memory_shape(dstbbox.shape() / dstbbox.stride())))
     {
       CCTK_WARN (0, "Internal error: array sizes don't agree with bounding boxes");
     }
diff --git a/Carpet/CarpetLib/src/prolongate_4d_o1_rf2.cc b/Carpet/CarpetLib/src/prolongate_4d_o1_rf2.cc
index 366bb88ef..af4b1696b 100644
--- a/Carpet/CarpetLib/src/prolongate_4d_o1_rf2.cc
+++ b/Carpet/CarpetLib/src/prolongate_4d_o1_rf2.cc
@@ -6,6 +6,7 @@
 #include <cmath>
 #include <cstdlib>
 
+#include "gdata.hh"
 #include "operator_prototypes_4d.hh"
 #include "typeprops.hh"
 
@@ -79,8 +80,8 @@ namespace CarpetLib {
       CCTK_WARN (0, "Internal error: region extent is not contained in array extent");
     }
     
-    if (any (srcext != srcbbox.shape() / srcbbox.stride() or
-             dstext != dstbbox.shape() / dstbbox.stride()))
+    if (any (srcext != gdata::allocated_memory_shape(srcbbox.shape() / srcbbox.stride()) or
+             dstext != gdata::allocated_memory_shape(dstbbox.shape() / dstbbox.stride())))
     {
       CCTK_WARN (0, "Internal error: array sizes don't agree with bounding boxes");
     }
diff --git a/Carpet/CarpetLib/src/restrict_3d_cc_rf2.cc b/Carpet/CarpetLib/src/restrict_3d_cc_rf2.cc
index 5e9ff5e84..29b574354 100644
--- a/Carpet/CarpetLib/src/restrict_3d_cc_rf2.cc
+++ b/Carpet/CarpetLib/src/restrict_3d_cc_rf2.cc
@@ -5,6 +5,7 @@
 #include <cassert>
 #include <cmath>
 
+#include "gdata.hh"
 #include "operator_prototypes_3d.hh"
 #include "typeprops.hh"
 
@@ -63,8 +64,8 @@ namespace CarpetLib {
       CCTK_WARN (0, "Internal error: region extent is not contained in array extent");
     }
     
-    if (any (srcext != srcbbox.shape() / srcbbox.stride() or
-             dstext != dstbbox.shape() / dstbbox.stride()))
+    if (any (srcext != gdata::allocated_memory_shape(srcbbox.shape() / srcbbox.stride()) or
+             dstext != gdata::allocated_memory_shape(dstbbox.shape() / dstbbox.stride())))
     {
       CCTK_WARN (0, "Internal error: array sizes don't agree with bounding boxes");
     }
diff --git a/Carpet/CarpetLib/src/restrict_3d_rf2.cc b/Carpet/CarpetLib/src/restrict_3d_rf2.cc
index a7df73e62..4286cd2ef 100644
--- a/Carpet/CarpetLib/src/restrict_3d_rf2.cc
+++ b/Carpet/CarpetLib/src/restrict_3d_rf2.cc
@@ -7,6 +7,7 @@
 #include <cstdlib>
 #include <iostream>
 
+#include "gdata.hh"
 #include "operator_prototypes_3d.hh"
 #include "typeprops.hh"
 
@@ -62,8 +63,8 @@ namespace CarpetLib {
       CCTK_WARN (0, "Internal error: region extent is not contained in array extent");
     }
     
-    if (any (srcext != srcbbox.shape() / srcbbox.stride() or
-             dstext != dstbbox.shape() / dstbbox.stride()))
+    if (any (srcext != gdata::allocated_memory_shape(srcbbox.shape() / srcbbox.stride()) or
+             dstext != gdata::allocated_memory_shape(dstbbox.shape() / dstbbox.stride())))
     {
       CCTK_WARN (0, "Internal error: array sizes don't agree with bounding boxes");
     }
diff --git a/Carpet/CarpetLib/src/restrict_3d_vc_rf2.cc b/Carpet/CarpetLib/src/restrict_3d_vc_rf2.cc
index a17462de0..8d80b40e0 100644
--- a/Carpet/CarpetLib/src/restrict_3d_vc_rf2.cc
+++ b/Carpet/CarpetLib/src/restrict_3d_vc_rf2.cc
@@ -5,6 +5,7 @@
 #include <cassert>
 #include <cmath>
 
+#include "gdata.hh"
 #include "operator_prototypes_3d.hh"
 #include "typeprops.hh"
 
@@ -169,8 +170,8 @@ namespace CarpetLib {
       CCTK_WARN (0, "Internal error: region extent is not contained in array extent");
     }
     
-    if (any (srcext != srcbbox.shape() / srcbbox.stride() or
-             dstext != dstbbox.shape() / dstbbox.stride()))
+    if (any (srcext != gdata::allocated_memory_shape(srcbbox.shape() / srcbbox.stride()) or
+             dstext != gdata::allocated_memory_shape(dstbbox.shape() / dstbbox.stride())))
     {
       CCTK_WARN (0, "Internal error: array sizes don't agree with bounding boxes");
     }
diff --git a/Carpet/CarpetLib/src/restrict_4d_rf2.cc b/Carpet/CarpetLib/src/restrict_4d_rf2.cc
index ba1370b93..8d4574278 100644
--- a/Carpet/CarpetLib/src/restrict_4d_rf2.cc
+++ b/Carpet/CarpetLib/src/restrict_4d_rf2.cc
@@ -7,6 +7,7 @@
 #include <cstdlib>
 #include <iostream>
 
+#include "gdata.hh"
 #include "operator_prototypes_4d.hh"
 #include "typeprops.hh"
 
@@ -62,8 +63,8 @@ namespace CarpetLib {
       CCTK_WARN (0, "Internal error: region extent is not contained in array extent");
     }
     
-    if (any (srcext != srcbbox.shape() / srcbbox.stride() or
-             dstext != dstbbox.shape() / dstbbox.stride()))
+    if (any (srcext != gdata::allocated_memory_shape(srcbbox.shape() / srcbbox.stride()) or
+             dstext != gdata::allocated_memory_shape(dstbbox.shape() / dstbbox.stride())))
     {
       CCTK_WARN (0, "Internal error: array sizes don't agree with bounding boxes");
     }