Line data    Source code

       1             : //* This file is part of the MOOSE framework
       2             : //* https://mooseframework.inl.gov
       3             : //*
       4             : //* All rights reserved, see COPYRIGHT for full restrictions
       5             : //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
       6             : //*
       7             : //* Licensed under LGPL 2.1, please see LICENSE for details
       8             : //* https://www.gnu.org/licenses/lgpl-2.1.html
       9             : 
      10             : #pragma once
      11             : 
      12             : #include "MooseTypes.h"
      13             : 
      14             : #include "libmesh/bounding_box.h"
      15             : #include "libmesh/face_tri3.h"
      16             : #include "libmesh/mesh_tet_interface.h"
      17             : 
      18             : #include <iosfwd>
      19             : #include <memory>
      20             : #include <cstdint>
      21             : #include <string>
      22             : #include <unordered_map>
      23             : #include <vector>
      24             : 
      25             : /**
      26             :  * Utility for querying point containment against a closed triangulated surface mesh.
      27             :  *
      28             :  * The class takes ownership of a reference to a prepared 2D surface mesh, validates that it forms
      29             :  * a closed 2-manifold (every triangle has exactly three neighbors, all elements are Tri3, and the
      30             :  * mesh is consistently oriented), builds a lightweight yz-plane acceleration grid, and exposes
      31             :  * `contains()` for point-in-solid queries.
      32             :  *
      33             :  * Containment is resolved in three stages:
      34             :  *
      35             :  * 1. Cheap global bounding-box rejection.
      36             :  * 2. Near-surface detection via a pre-built point locator (points within `surface_tolerance` of
      37             :  *    the mesh surface are treated as inside).
      38             :  * 3. Odd/even parity counting on a fixed +x ray, with automatic fallback to a solid-angle
      39             :  *    accumulation test when the ray grazes a triangle edge or vertex.
      40             :  *
      41             :  * The referenced mesh must outlive this object. Any geometric transforms (scale, rotation,
      42             :  * translation) should be applied to the mesh before constructing a `TriangleManifold`.
      43             :  */
      44             : class TriangleManifold
      45             : {
      46             : public:
      47             :   /**
      48             :    * Build a manifold classifier from a prepared surface mesh.
      49             :    *
      50             :    * @param mesh Serialized 2D surface mesh that defines the closed manifold. Must outlive this
      51             :    * object. Any desired transforms should already be applied to the mesh.
      52             :    * @param surface_tolerance Absolute tolerance used for manifold validation and near-surface
      53             :    * classification. Choose this relative to the mesh length scale and expected coordinate noise
      54             :    * from the export pipeline.
      55             :    */
      56             :   TriangleManifold(MeshBase & mesh, const Real surface_tolerance);
      57             : 
      58             :   /**
      59             :    * @return True if the point is inside the manifold or lies on the surface.
      60             :    */
      61             :   bool contains(const Point & point) const;
      62             : 
      63             :   /**
      64             :    * @return The manifold bounding box.
      65             :    */
      66      192438 :   const libMesh::BoundingBox & boundingBox() const { return _bounding_box; };
      67             : 
      68             :   /**
      69             :    * @return The number of triangles in the loaded manifold.
      70             :    */
      71          83 :   std::size_t numTriangles() const { return _mesh.n_active_elem(); }
      72             : 
      73             : private:
      74             :   /**
      75             :    * Result of intersecting the positive x-direction ray with a triangle.
      76             :    *
      77             :    * `Ambiguous` is returned when the hit is too close to an edge, vertex, or the ray origin. In
      78             :    * that case we abandon parity counting and fall back to the more expensive but robust solid-angle
      79             :    * test.
      80             :    */
      81             :   enum class RayIntersection
      82             :   {
      83             :     Miss,
      84             :     Hit,
      85             :     Ambiguous
      86             :   };
      87             : 
      88             :   /// Complete post-parse validation and acceleration-structure setup.
      89             :   void finalize();
      90             : 
      91             :   /// Build the yz-plane lookup grid used to accelerate +x ray queries.
      92             :   void buildCandidateGrid();
      93             : 
      94             :   /// Cheap global bounding-box rejection for containment queries.
      95             :   bool pointInsideBoundingBox(const Point & point) const;
      96             : 
      97             :   /// Detect whether a query point lies on or extremely near the manifold surface.
      98             :   bool pointOnSurface(const Point & point) const;
      99             : 
     100             :   /// Intersect a positive x-direction ray with a single triangle.
     101             :   RayIntersection rayIntersectsTriangle(const Point & point, const libMesh::Elem & tri) const;
     102             : 
     103             :   /// Robust fallback containment query based on accumulated solid angle.
     104             :   bool containsBySolidAngle(const Point & point) const;
     105             : 
     106             :   /// Get the subset of triangles whose yz extents may intersect the query ray.
     107             :   std::vector<dof_id_type> rayCandidates(const Point & point) const;
     108             : 
     109             :   /// Number of yz-grid bins in the y direction.
     110             :   std::size_t _num_y_cells = 1;
     111             : 
     112             :   /// Number of yz-grid bins in the z direction.
     113             :   std::size_t _num_z_cells = 1;
     114             : 
     115             :   /// Minimum global y coordinate used to map query points into yz-grid bins.
     116             :   Real _y_min = 0.0;
     117             : 
     118             :   /// Minimum global z coordinate used to map query points into yz-grid bins.
     119             :   Real _z_min = 0.0;
     120             : 
     121             :   /// Width of one yz-grid cell in the y direction.
     122             :   Real _y_cell_size = 1.0;
     123             : 
     124             :   /// Width of one yz-grid cell in the z direction.
     125             :   Real _z_cell_size = 1.0;
     126             : 
     127             :   /// Lookup from packed yz-grid cell index to triangles that could intersect the +x query ray.
     128             :   std::unordered_map<dof_id_type, std::vector<dof_id_type>> _ray_grid;
     129             : 
     130             :   MeshBase & _mesh;
     131             : 
     132             :   /// Absolute tolerance used throughout validation and geometric classification.
     133             :   const Real _surface_tolerance;
     134             : 
     135             :   /// Global bounding box of the transformed manifold.
     136             :   const libMesh::BoundingBox _bounding_box;
     137             : 
     138             :   /// Pre-built point locator for fast proximity-to-surface detection.
     139             :   const std::unique_ptr<libMesh::PointLocatorBase> _point_locator;
     140             : };