diff --git a/include/openmc/mesh.h b/include/openmc/mesh.h index 9554a51436d..f19a22c09f2 100644 --- a/include/openmc/mesh.h +++ b/include/openmc/mesh.h @@ -75,6 +75,12 @@ class Mesh { // 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 @@ class StructuredMesh : public Mesh { } }; - int get_bin(Position r) const override; + virtual int get_bin(Position r) const; int n_bins() const override; @@ -241,11 +247,27 @@ class StructuredMesh : public Mesh { xt::xtensor lower_left_; //!< Lower-left coordinates of mesh xt::xtensor upper_right_; //!< Upper-right coordinates of mesh std::array 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 @@ class RectilinearMesh : public StructuredMesh { int set_grid(); }; -class CylindricalMesh : public StructuredMesh { +class CylindricalMesh : public PeriodicStructuredMesh { public: // Constructors CylindricalMesh() = default; @@ -390,7 +412,7 @@ class CylindricalMesh : public StructuredMesh { } }; -class SphericalMesh : public StructuredMesh { +class SphericalMesh : public PeriodicStructuredMesh { public: // Constructors SphericalMesh() = default; diff --git a/src/mesh.cpp b/src/mesh.cpp index e0fa305edd4..9020a95aafc 100644 --- a/src/mesh.cpp +++ b/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 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();