Line data    Source code

       1              : // Copyright Amazon.com, Inc. or its affiliates. All Rights Reserved.
       2              : // SPDX-License-Identifier: Apache-2.0
       3              : 
       4              : #ifndef PALACE_UTILS_GEODATA_IMPL_HPP
       5              : #define PALACE_UTILS_GEODATA_IMPL_HPP
       6              : 
       7              : #include <array>
       8              : #include <memory>
       9              : #include <vector>
      10              : #include <Eigen/Dense>
      11              : #include <mfem.hpp>
      12              : #include "utils/iodata.hpp"
      13              : 
      14              : namespace palace
      15              : {
      16              : 
      17              : struct BoundingBox;
      18              : 
      19              : namespace mesh
      20              : {
      21              : 
      22              : //
      23              : // Implementations Functions for mesh related functionality. Isolated to avoid Eigen
      24              : // propagation.
      25              : //
      26              : 
      27              : // For the public interface, we can use a BoundingBox as a generalization of a BoundingBall.
      28              : // Internally, however, it's nice to work with a specific ball data type.
      29            0 : struct BoundingBall
      30              : {
      31              :   Eigen::Vector3d origin;
      32              :   double radius;
      33              :   bool planar;
      34              : };
      35              : 
      36              : // Helper for collecting a point cloud from a mesh, used in calculating bounding boxes and
      37              : // bounding balls. Returns the dominant rank, for which the vertices argument will be
      38              : // filled, while all other ranks will have an empty vector. Vertices are de-duplicated to a
      39              : // certain floating point precision.
      40              : int CollectPointCloudOnRoot(const mfem::ParMesh &mesh, const mfem::Array<int> &marker,
      41              :                             bool bdr, std::vector<Eigen::Vector3d> &vertices);
      42              : 
      43              : // Calculates a bounding box from a point cloud, result is broadcast across all processes.
      44              : BoundingBox BoundingBoxFromPointCloud(MPI_Comm comm,
      45              :                                       const std::vector<Eigen::Vector3d> &vertices,
      46              :                                       int dominant_rank);
      47              : 
      48              : // Compute the distance from a point orthogonal to the list of normal axes, relative to
      49              : // the given origin.
      50              : auto PerpendicularDistance(const std::initializer_list<Eigen::Vector3d> &normals,
      51              :                            const Eigen::Vector3d &origin, const Eigen::Vector3d &v);
      52              : 
      53              : // Use 4 points to define a sphere in 3D. If the points are coplanar, 3 of them are used to
      54              : // define a circle which is interpreted as the equator of the sphere. We assume the points
      55              : // are unique and not collinear.
      56              : BoundingBall SphereFromPoints(const std::vector<std::size_t> &indices,
      57              :                               const std::vector<Eigen::Vector3d> &vertices);
      58              : 
      59              : // Implementation of the recursive Welzl algorithm kernel.
      60              : BoundingBall Welzl(std::vector<std::size_t> P, std::vector<std::size_t> R,
      61              :                    const std::vector<Eigen::Vector3d> &vertices);
      62              : 
      63              : // Calculates a bounding ball from a point cloud using Welzl's algorithm, result is
      64              : // broadcast across all processes. We don't operate on the convex hull, since the number of
      65              : // points should be small enough that operating on the full set should be OK. If only three
      66              : // points are provided, the bounding circle is computed (likewise for if the points are
      67              : // coplanar).
      68              : BoundingBox BoundingBallFromPointCloud(MPI_Comm comm,
      69              :                                        const std::vector<Eigen::Vector3d> &vertices,
      70              :                                        int dominant_rank);
      71              : 
      72              : // Compute a normal vector from an element transformation, optionally ensure aligned
      73              : // (| normal ⋅ align | > 0)
      74              : void Normal(mfem::ElementTransformation &T, mfem::Vector &normal,
      75              :             const mfem::Vector *const align);
      76              : 
      77              : // Determine the vertex mapping between donor and receiver boundary attributes.
      78              : // Uses the translation vector or affine transformation matrix specified in the
      79              : // configuration file. If not provided, attempts to automatically detect the
      80              : // affine transformation between donor and receiver boundary vertices.
      81              : std::vector<int>
      82              : DeterminePeriodicVertexMapping(std::unique_ptr<mfem::Mesh> &mesh,
      83              :                                const struct palace::config::PeriodicData &data,
      84              :                                const double tol = 1e-8);
      85              : 
      86              : }  // namespace mesh
      87              : 
      88              : }  // namespace palace
      89              : 
      90              : #endif  // PALACE_UTILS_GEODATA_IMPL_HPP