Improve runtime checks on grid/domain validity

src/include/dim.hxx

@@ -3,12 +3,34 @@
 
 #include <type_traits>
 
+#include "kg/Vec3.h"
+
 template <bool IX = false, bool IY = false, bool IZ = false>
 struct Invar
 {
   using InvarX = std::integral_constant<bool, IX>;
   using InvarY = std::integral_constant<bool, IY>;
   using InvarZ = std::integral_constant<bool, IZ>;
+
+  /// @brief Gets a mask of noninvariant dimensions (i.e., 0 for invariant
+  /// dimensions, 1 for noninvariant dimensions). This is useful for multiplying
+  /// by a constant to, say, get the number of ghost cells in each dimension
+  /// (where invariant dimensions have no ghosts).
+  /// @return the mask
+  static Int3 get_noninvariant_mask()
+  {
+    return {!InvarX::value, !InvarY::value, !InvarZ::value};
+  }
+
+  static bool is_invar(int dim)
+  {
+    switch (dim) {
+      case 0: return InvarX::value;
+      case 1: return InvarY::value;
+      case 2: return InvarZ::value;
+      default: return false;
+    }
+  }
 };
 
 using dim_xyz = Invar<false, false, false>;

src/include/grid.hxx

@@ -7,7 +7,7 @@
 #include "psc_bits.h"
 #include "mrc_domain.hxx"
 #include "grid/BC.h"
-#include "grid/Domain.h"
+#include "grid/domain.hxx"
 #include <mrc_ddc.h>
 
 #include <vector>
@@ -81,10 +81,6 @@ struct Grid_
     mpi_printf(MPI_COMM_WORLD, "::: dx      = %g %g %g\n", domain.dx[0], domain.dx[1], domain.dx[2]);
 #endif
 
-    assert(domain.dx[0] > 0.);
-    assert(domain.dx[1] > 0.);
-    assert(domain.dx[2] > 0.);
-
     for (auto off : mrc_domain_.offs()) {
       patches.push_back(Patch(
         off, Vec3<double>(off) * domain.dx + domain.corner,

src/include/grid/Domain.h → src/include/grid/domain.hxx

@@ -1,9 +1,8 @@
-
-#ifndef GRID_DOMAIN_H
-#define GRID_DOMAIN_H
+#pragma once
 
 #include <mrc_common.h>
 #include <iostream>
+#include "../libpsc/vpic/psc_vpic_bits.h"
 
 namespace psc
 {
@@ -26,9 +25,23 @@ struct Domain
     : gdims(gdims), length(length), corner(corner), np(np)
   {
     for (int d = 0; d < 3; d++) {
-      assert(gdims[d] % np[d] == 0);
-      ldims[d] = gdims[d] / np[d];
+      if (gdims[d] <= 0) {
+        LOG_ERROR("dimension %d has non-positive number of cells (%d)\n", d,
+                  gdims[d]);
+      }
+
+      if (gdims[d] % np[d] != 0) {
+        LOG_ERROR("in dimension %d, number of patches (%d) doesn't divide "
+                  "number of cells (%d)\n",
+                  d, np[d], gdims[d]);
+      }
+
+      if (length[d] <= 0.0) {
+        LOG_ERROR("dimension %d has non-positive length (%f)\n", d, length[d]);
+      }
     }
+
+    ldims = gdims / np;
     dx = length / Real3(gdims);
   }
 
@@ -40,11 +53,17 @@ struct Domain
 
   bool isInvar(int d) const { return gdims[d] == 1; }
 
-  Int3 gdims;   ///< Number of grid-points in each dimension
-  Real3 length; ///< The physical size of the simulation-box
+  /// @brief Total number of grid cells in each dimension.
+  Int3 gdims;
+  /// @brief Side lengths of the domain, in physical units.
+  Real3 length;
+  /// @brief Location of lower-right corner of the domain, in physical units.
   Real3 corner;
-  Int3 np; ///< Number of patches in each dimension
+  /// @brief Number of patches in each dimension.
+  Int3 np;
+  /// @brief Number of grid cells per patch in each dimension.
   Int3 ldims;
+  /// @brief Side lengths of each grid cell, in physical units.
   Real3 dx;
 };
 
@@ -63,5 +82,3 @@ inline std::ostream& operator<<(std::ostream& os, const Domain<R>& domain)
 
 } // namespace grid
 } // namespace psc
-
-#endif

src/include/mrc_domain.hxx

@@ -2,7 +2,7 @@
 #pragma once
 
 #include "grid/BC.h"
-#include "grid/Domain.h"
+#include "grid/domain.hxx"
 
 #include <mrc_domain_private.h>
 

src/include/psc.hxx

@@ -177,9 +177,12 @@ struct Psc
   {
     time_start_ = MPI_Wtime();
 
-    assert(grid.isInvar(0) == Dim::InvarX::value);
-    assert(grid.isInvar(1) == Dim::InvarY::value);
-    assert(grid.isInvar(2) == Dim::InvarZ::value);
+    for (int d = 0; d < 3; d++) {
+      if (grid.isInvar(d) != Dim::is_invar(d)) {
+        LOG_ERROR("dimension %d is%s invariant, but gdims[%d]=%d\n", d,
+                  Dim::is_invar(d) ? "" : " not", d, grid.domain.gdims[d]);
+      }
+    }
 
     int rank;
     MPI_Comm_rank(MPI_COMM_WORLD, &rank);

src/libpsc/tests/test_boundary_injector.cxx

@@ -68,16 +68,8 @@ Grid_t* setupGrid()
   double dt = .1;
   Grid_t::Normalization norm{norm_params};
 
-  Int3 ibn = {2, 2, 2};
-  if (Dim::InvarX::value) {
-    ibn[0] = 0;
-  }
-  if (Dim::InvarY::value) {
-    ibn[1] = 0;
-  }
-  if (Dim::InvarZ::value) {
-    ibn[2] = 0;
-  }
+  int n_ghosts = 2;
+  Int3 ibn = n_ghosts * Dim::get_noninvariant_mask();
 
   return new Grid_t{domain, bc, kinds, norm, dt, -1, ibn};
 }

src/libpsc/tests/test_reflective_bcs.cxx

@@ -65,16 +65,8 @@ Grid_t* setupGrid()
   double dt = psc_params.cfl * courant_length(domain);
   Grid_t::Normalization norm{norm_params};
 
-  Int3 ibn = {2, 2, 2};
-  if (Dim::InvarX::value) {
-    ibn[0] = 0;
-  }
-  if (Dim::InvarY::value) {
-    ibn[1] = 0;
-  }
-  if (Dim::InvarZ::value) {
-    ibn[2] = 0;
-  }
+  int n_ghosts = 2;
+  Int3 ibn = n_ghosts * Dim::get_noninvariant_mask();
 
   return new Grid_t{domain, bc, kinds, norm, dt, -1, ibn};
 }

src/libpsc/tests/test_reflective_bcs_integration.cxx

@@ -65,16 +65,8 @@ Grid_t* setupGrid()
   double dt = psc_params.cfl * courant_length(domain);
   Grid_t::Normalization norm{norm_params};
 
-  Int3 ibn = {2, 2, 2};
-  if (Dim::InvarX::value) {
-    ibn[0] = 0;
-  }
-  if (Dim::InvarY::value) {
-    ibn[1] = 0;
-  }
-  if (Dim::InvarZ::value) {
-    ibn[2] = 0;
-  }
+  int n_ghosts = 2;
+  Int3 ibn = n_ghosts * Dim::get_noninvariant_mask();
 
   return new Grid_t{domain, bc, kinds, norm, dt, -1, ibn};
 }

src/psc_2d_shock.cxx

@@ -320,16 +320,8 @@ Grid_t* setupGrid()
   double dt = psc_params.cfl * courant_length(domain);
   Grid_t::Normalization norm{norm_params};
 
-  Int3 ibn = {2, 2, 2};
-  if (Dim::InvarX::value) {
-    ibn[0] = 0;
-  }
-  if (Dim::InvarY::value) {
-    ibn[1] = 0;
-  }
-  if (Dim::InvarZ::value) {
-    ibn[2] = 0;
-  }
+  int n_ghosts = 2;
+  Int3 ibn = n_ghosts * Dim::get_noninvariant_mask();
 
   return new Grid_t{domain, bc, kinds, norm, dt, -1, ibn};
 }

src/psc_bgk.cxx

@@ -129,16 +129,8 @@ Grid_t* setupGrid()
   double dt = psc_params.cfl * courant_length(domain);
   Grid_t::Normalization norm{norm_params};
 
-  Int3 ibn = {2, 2, 2};
-  if (Dim::InvarX::value) {
-    ibn[0] = 0;
-  }
-  if (Dim::InvarY::value) {
-    ibn[1] = 0;
-  }
-  if (Dim::InvarZ::value) {
-    ibn[2] = 0;
-  }
+  int n_ghosts = 2;
+  Int3 ibn = n_ghosts * Dim::get_noninvariant_mask();
 
   return new Grid_t{domain, bc, kinds, norm, dt, -1, ibn};
 }

src/psc_bubble_yz.cxx

@@ -141,16 +141,8 @@ Grid_t* setupGrid()
   double dt = psc_params.cfl * courant_length(domain);
   Grid_t::Normalization norm{norm_params};
 
-  Int3 ibn = {2, 2, 2};
-  if (Dim::InvarX::value) {
-    ibn[0] = 0;
-  }
-  if (Dim::InvarY::value) {
-    ibn[1] = 0;
-  }
-  if (Dim::InvarZ::value) {
-    ibn[2] = 0;
-  }
+  int n_ghosts = 2;
+  Int3 ibn = n_ghosts * Dim::get_noninvariant_mask();
 
   return new Grid_t{domain, bc, kinds, norm, dt, -1, ibn};
 }

src/psc_flatfoil_yz.cxx

@@ -339,16 +339,8 @@ Grid_t* setupGrid()
 
   mpi_printf(MPI_COMM_WORLD, "dt = %g\n", dt);
 
-  Int3 ibn = {2, 2, 2};
-  if (Dim::InvarX::value) {
-    ibn[0] = 0;
-  }
-  if (Dim::InvarY::value) {
-    ibn[1] = 0;
-  }
-  if (Dim::InvarZ::value) {
-    ibn[2] = 0;
-  }
+  int n_ghosts = 2;
+  Int3 ibn = n_ghosts * Dim::get_noninvariant_mask();
 
   return new Grid_t{domain, bc, kinds, norm, dt, -1, ibn};
 }

src/psc_harris_xz.cxx

@@ -368,16 +368,8 @@ Grid_t* setupGrid()
 #ifdef VPIC
   Int3 ibn = {1, 1, 1};
 #else
-  Int3 ibn = {2, 2, 2};
-  if (Dim::InvarX::value) {
-    ibn[0] = 0;
-  }
-  if (Dim::InvarY::value) {
-    ibn[1] = 0;
-  }
-  if (Dim::InvarZ::value) {
-    ibn[2] = 0;
-  }
+  int n_ghosts = 2;
+  Int3 ibn = n_ghosts * Dim::get_noninvariant_mask();
 #endif
 
   auto grid_ptr = new Grid_t{domain, bc, kinds, norm, dt, -1, ibn};

src/psc_harris_yz.cxx

@@ -247,16 +247,8 @@ Grid_t* setupGrid()
 
   mpi_printf(MPI_COMM_WORLD, "dt = %g\n", dt);
 
-  Int3 ibn = {2, 2, 2};
-  if (Dim::InvarX::value) {
-    ibn[0] = 0;
-  }
-  if (Dim::InvarY::value) {
-    ibn[1] = 0;
-  }
-  if (Dim::InvarZ::value) {
-    ibn[2] = 0;
-  }
+  int n_ghosts = 2;
+  Int3 ibn = n_ghosts * Dim::get_noninvariant_mask();
 
   return new Grid_t{domain, bc, kinds, norm, dt, -1, ibn};
 }

src/psc_radiation.cxx

@@ -150,16 +150,8 @@ Grid_t* setupGrid()
 
   mpi_printf(MPI_COMM_WORLD, "dt = %g\n", dt);
 
-  Int3 ibn = {2, 2, 2};
-  if (Dim::InvarX::value) {
-    ibn[0] = 0;
-  }
-  if (Dim::InvarY::value) {
-    ibn[1] = 0;
-  }
-  if (Dim::InvarZ::value) {
-    ibn[2] = 0;
-  }
+  int n_ghosts = 2;
+  Int3 ibn = n_ghosts * Dim::get_noninvariant_mask();
 
   return new Grid_t{domain, bc, kinds, norm, dt, -1, ibn};
 }

src/psc_shock.cxx

@@ -175,16 +175,8 @@ Grid_t* setupGrid()
   double dt = psc_params.cfl * courant_length(domain);
   Grid_t::Normalization norm{norm_params};
 
-  Int3 ibn = {2, 2, 2};
-  if (Dim::InvarX::value) {
-    ibn[0] = 0;
-  }
-  if (Dim::InvarY::value) {
-    ibn[1] = 0;
-  }
-  if (Dim::InvarZ::value) {
-    ibn[2] = 0;
-  }
+  int n_ghosts = 2;
+  Int3 ibn = n_ghosts * Dim::get_noninvariant_mask();
 
   return new Grid_t{domain, bc, kinds, norm, dt, -1, ibn};
 }

src/psc_whistler.cxx

@@ -163,16 +163,8 @@ Grid_t* setupGrid()
   double dt = psc_params.cfl * courant_length(domain);
   Grid_t::Normalization norm{norm_params};
 
-  Int3 ibn = {2, 2, 2};
-  if (Dim::InvarX::value) {
-    ibn[0] = 0;
-  }
-  if (Dim::InvarY::value) {
-    ibn[1] = 0;
-  }
-  if (Dim::InvarZ::value) {
-    ibn[2] = 0;
-  }
+  int n_ghosts = 2;
+  Int3 ibn = n_ghosts * Dim::get_noninvariant_mask();
 
   return new Grid_t{domain, bc, kinds, norm, dt, -1, ibn};
 }