diff --git a/include/bout/mesh.hxx b/include/bout/mesh.hxx
index 02e2a23905..4c8cb35a05 100644
--- a/include/bout/mesh.hxx
+++ b/include/bout/mesh.hxx
@@ -464,11 +464,11 @@ public:
 
   /// Is there a boundary on the lower guard cells in Y
   /// on any processor along the X direction?
-  bool hasBndryLowerY();
+  virtual bool hasBndryLowerY() const = 0;
 
   /// Is there a boundary on the upper guard cells in Y
   /// on any processor along the X direction?
-  bool hasBndryUpperY();
+  virtual bool hasBndryUpperY() const = 0;
   // Boundary regions
 
   /// Return a vector containing all the boundary regions on this processor
diff --git a/src/mesh/impls/bout/boutmesh.cxx b/src/mesh/impls/bout/boutmesh.cxx
index 9f366d9456..8917fa3753 100644
--- a/src/mesh/impls/bout/boutmesh.cxx
+++ b/src/mesh/impls/bout/boutmesh.cxx
@@ -607,6 +607,20 @@ int BoutMesh::load() {
   // Add boundary regions
   addBoundaryRegions();
 
+  // Set cached values
+  {
+    int mybndry = static_cast<int>(!(iterateBndryLowerY().isDone()));
+    int allbndry = 0;
+    mpi->MPI_Allreduce(&mybndry, &allbndry, 1, MPI_INT, MPI_BOR, getXcomm(yend));
+    has_boundary_lower_y = static_cast<bool>(allbndry);
+  }
+  {
+    int mybndry = static_cast<int>(!(iterateBndryUpperY().isDone()));
+    int allbndry = 0;
+    mpi->MPI_Allreduce(&mybndry, &allbndry, 1, MPI_INT, MPI_BOR, getXcomm(ystart));
+    has_boundary_upper_y = static_cast<bool>(allbndry);
+  }
+
   // Initialize default coordinates
   getCoordinates();
 
diff --git a/src/mesh/impls/bout/boutmesh.hxx b/src/mesh/impls/bout/boutmesh.hxx
index cc674d401a..a22e0e4473 100644
--- a/src/mesh/impls/bout/boutmesh.hxx
+++ b/src/mesh/impls/bout/boutmesh.hxx
@@ -156,6 +156,9 @@ public:
   RangeIterator iterateBndryUpperInnerY() const override;
   RangeIterator iterateBndryUpperOuterY() const override;
 
+  bool hasBndryLowerY() const override { return has_boundary_lower_y; }
+  bool hasBndryUpperY() const override { return has_boundary_upper_y; }
+
   // Boundary regions
   std::vector<BoundaryRegion*> getBoundaries() override;
   std::vector<std::shared_ptr<BoundaryRegionPar>>
@@ -400,6 +403,8 @@ private:
              static_cast<int>(BoundaryParType::SIZE)>
       par_boundary; // Vector of parallel boundary regions
 
+  bool has_boundary_lower_y{false};
+  bool has_boundary_upper_y{false};
   //////////////////////////////////////////////////
   // Communications
 
diff --git a/src/mesh/mesh.cxx b/src/mesh/mesh.cxx
index cb6cb8410c..1d3f3b0dbe 100644
--- a/src/mesh/mesh.cxx
+++ b/src/mesh/mesh.cxx
@@ -408,34 +408,6 @@ int Mesh::ySize(int jx) const {
   return all;
 }
 
-bool Mesh::hasBndryLowerY() {
-  static bool calc = false, answer;
-  if (calc) {
-    return answer; // Already calculated
-  }
-
-  int mybndry = static_cast<int>(!(iterateBndryLowerY().isDone()));
-  int allbndry;
-  mpi->MPI_Allreduce(&mybndry, &allbndry, 1, MPI_INT, MPI_BOR, getXcomm(yend));
-  answer = static_cast<bool>(allbndry);
-  calc = true;
-  return answer;
-}
-
-bool Mesh::hasBndryUpperY() {
-  static bool calc = false, answer;
-  if (calc) {
-    return answer; // Already calculated
-  }
-
-  int mybndry = static_cast<int>(!(iterateBndryUpperY().isDone()));
-  int allbndry;
-  mpi->MPI_Allreduce(&mybndry, &allbndry, 1, MPI_INT, MPI_BOR, getXcomm(ystart));
-  answer = static_cast<bool>(allbndry);
-  calc = true;
-  return answer;
-}
-
 int Mesh::localSize3D() {
   if (localNumCells3D < 0) {
     const int xs = firstX() ? xstart - 1 : xstart;
diff --git a/tests/unit/fake_mesh.hxx b/tests/unit/fake_mesh.hxx
index e6f78f8767..dc43b647ec 100644
--- a/tests/unit/fake_mesh.hxx
+++ b/tests/unit/fake_mesh.hxx
@@ -143,6 +143,8 @@ public:
   RangeIterator iterateBndryLowerInnerY() const override { return RangeIterator(); }
   RangeIterator iterateBndryUpperOuterY() const override { return RangeIterator(); }
   RangeIterator iterateBndryUpperInnerY() const override { return RangeIterator(); }
+  bool hasBndryLowerY() const override { return false; }
+  bool hasBndryUpperY() const override { return false; }
   void addBoundary(BoundaryRegion* region) override { boundaries.push_back(region); }
   std::vector<BoundaryRegion*> getBoundaries() override { return boundaries; }
   std::vector<std::shared_ptr<BoundaryRegionPar>>