Merge pull request #6 from pshriwise/centre_for_cylinder_spherical_meshes

Add an intermediate `PeriodicStructuredMesh` class to separate the `origin` attribute from `StructuredMesh`

include/openmc/mesh.h

@@ -75,6 +75,12 @@ public:
 
   // Methods
 
+  //! Update a position to the local coordinates of the mesh
+  virtual void local_coords(Position& r) const {};
+
+  //! Return a position in the local coordinates of the mesh
+  virtual Position local_coords(const Position& r) const { return r; };
+
   //! Determine which bins were crossed by a particle
   //
   //! \param[in] r0 Previous position of the particle
@@ -160,7 +166,7 @@ public:
     }
   };
 
-  int get_bin(Position r) const override;
+  virtual int get_bin(Position r) const;
 
   int n_bins() const override;
 
@@ -241,11 +247,27 @@ public:
   xt::xtensor<double, 1> lower_left_;  //!< Lower-left coordinates of mesh
   xt::xtensor<double, 1> upper_right_; //!< Upper-right coordinates of mesh
   std::array<int, 3> shape_; //!< Number of mesh elements in each dimension
-  Position origin_ {0.0, 0.0, 0.0};
 
 protected:
 };
 
+class PeriodicStructuredMesh : public StructuredMesh {
+
+public:
+  PeriodicStructuredMesh() = default;
+  PeriodicStructuredMesh(pugi::xml_node node) : StructuredMesh {node} {};
+
+  void local_coords(Position& r) const override { r -= origin_; };
+
+  Position local_coords(const Position& r) const override
+  {
+    return r - origin_;
+  };
+
+  // Data members
+  Position origin_ {0.0, 0.0, 0.0}; //!< Origin of the mesh
+};
+
 //==============================================================================
 //! Tessellation of n-dimensional Euclidean space by congruent squares or cubes
 //==============================================================================
@@ -336,7 +358,7 @@ public:
   int set_grid();
 };
 
-class CylindricalMesh : public StructuredMesh {
+class CylindricalMesh : public PeriodicStructuredMesh {
 public:
   // Constructors
   CylindricalMesh() = default;
@@ -390,7 +412,7 @@ private:
   }
 };
 
-class SphericalMesh : public StructuredMesh {
+class SphericalMesh : public PeriodicStructuredMesh {
 public:
   // Constructors
   SphericalMesh() = default;

src/mesh.cpp

@@ -446,8 +446,8 @@ void StructuredMesh::raytrace_mesh(
 
   // translate start and end positions,
   // this needs to come after the get_indices call because it does its own translation
-  r0 -= origin_;
-  r1 -= origin_;
+  local_coords(r0);
+  local_coords(r1);
 
   // Calculate initial distances to next surfaces in all three dimensions
   std::array<MeshDistance, 3> distances;
@@ -952,7 +952,8 @@ void RectilinearMesh::to_hdf5(hid_t group) const
 // CylindricalMesh implementation
 //==============================================================================
 
-CylindricalMesh::CylindricalMesh(pugi::xml_node node) : StructuredMesh {node}
+CylindricalMesh::CylindricalMesh(pugi::xml_node node)
+  : PeriodicStructuredMesh {node}
 {
   n_dimension_ = 3;
 
@@ -976,7 +977,7 @@ std::string CylindricalMesh::get_mesh_type() const
 StructuredMesh::MeshIndex CylindricalMesh::get_indices(
   Position r, bool& in_mesh) const
 {
-  r -= origin_;
+  local_coords(r);
 
   Position mapped_r;
   mapped_r[0] = std::hypot(r.x, r.y);
@@ -1189,7 +1190,8 @@ void CylindricalMesh::to_hdf5(hid_t group) const
 // SphericalMesh implementation
 //==============================================================================
 
-SphericalMesh::SphericalMesh(pugi::xml_node node) : StructuredMesh {node}
+SphericalMesh::SphericalMesh(pugi::xml_node node)
+  : PeriodicStructuredMesh {node}
 {
   n_dimension_ = 3;
 
@@ -1213,7 +1215,7 @@ std::string SphericalMesh::get_mesh_type() const
 StructuredMesh::MeshIndex SphericalMesh::get_indices(
   Position r, bool& in_mesh) const
 {
-  r -= origin_;
+  local_coords(r);
 
   Position mapped_r;
   mapped_r[0] = r.norm();