DOLFIN-X
DOLFIN-X C++ interface
utils.h
1 // Copyright (C) 2019 Garth N. Wells
2 //
3 // This file is part of DOLFINX (https://www.fenicsproject.org)
4 //
5 // SPDX-License-Identifier: LGPL-3.0-or-later
6 
7 #pragma once
8 
9 #include "BoundingBoxTree.h"
10 #include <Eigen/Dense>
11 #include <utility>
12 #include <vector>
13 
14 namespace dolfinx
15 {
16 namespace mesh
17 {
18 class Mesh;
19 class MeshEntity;
20 } // namespace mesh
21 
22 namespace geometry
23 {
27 BoundingBoxTree create_midpoint_tree(const mesh::Mesh& mesh);
28 
29 // FIXME: document properly
31 std::pair<std::vector<int>, std::vector<int>>
32 compute_collisions(const BoundingBoxTree& tree0, const BoundingBoxTree& tree1);
33 
35 std::pair<std::vector<int>, std::vector<int>>
36 compute_entity_collisions(const BoundingBoxTree& tree0,
37  const BoundingBoxTree& tree1, const mesh::Mesh& mesh0,
38  const mesh::Mesh& mesh1);
39 
44 std::vector<int> compute_collisions(const BoundingBoxTree& tree,
45  const Eigen::Vector3d& p);
46 
52 std::vector<int> compute_entity_collisions(const BoundingBoxTree& tree,
53  const Eigen::Vector3d& p,
54  const mesh::Mesh& mesh);
55 
60 int compute_first_collision(const BoundingBoxTree& tree,
61  const Eigen::Vector3d& p);
62 
68 int compute_first_entity_collision(const BoundingBoxTree& tree,
69  const Eigen::Vector3d& p,
70  const mesh::Mesh& mesh);
71 
74 std::vector<int> compute_process_collisions(const BoundingBoxTree& tree,
75  const Eigen::Vector3d& p);
76 
78 bool bbox_in_bbox(const Eigen::Array<double, 2, 3, Eigen::RowMajor>& a,
79  const Eigen::Array<double, 2, 3, Eigen::RowMajor>& b,
80  double rtol = 1e-14);
81 
84 std::pair<int, double>
85 compute_closest_entity(const BoundingBoxTree& tree,
86  const BoundingBoxTree& tree_midpoint,
87  const Eigen::Vector3d& p, const mesh::Mesh& mesh);
88 
94 std::pair<int, double> compute_closest_point(const BoundingBoxTree& tree,
95  const Eigen::Vector3d& p);
96 
103 bool point_in_bbox(const Eigen::Array<double, 2, 3, Eigen::RowMajor>& b,
104  const Eigen::Vector3d& x, double rtol = 1e-14);
105 
109  const Eigen::Array<double, 2, 3, Eigen::RowMajor>& b,
110  const Eigen::Vector3d& x);
111 
114 double squared_distance(const mesh::MeshEntity& entity,
115  const Eigen::Vector3d& p);
116 
121 double squared_distance_triangle(const Eigen::Vector3d& point,
122  const Eigen::Vector3d& a,
123  const Eigen::Vector3d& b,
124  const Eigen::Vector3d& c);
125 
129 double squared_distance_interval(const Eigen::Vector3d& point,
130  const Eigen::Vector3d& a,
131  const Eigen::Vector3d& b);
132 
133 } // namespace geometry
134 } // namespace dolfinx
dolfinx::geometry::create_midpoint_tree
BoundingBoxTree create_midpoint_tree(const mesh::Mesh &mesh)
Create a boundary box tree for cell midpoints.
Definition: utils.cpp:342
dolfinx::geometry::squared_distance_interval
double squared_distance_interval(const Eigen::Vector3d &point, const Eigen::Vector3d &a, const Eigen::Vector3d &b)
Compute squared distance to given point. This version takes the two vertex coordinates as 3D points....
Definition: utils.cpp:749
dolfinx::geometry::squared_distance
double squared_distance(const mesh::MeshEntity &entity, const Eigen::Vector3d &p)
Compute squared distance from a given point to the nearest point on a cell (only simplex cells are su...
Definition: utils.cpp:545
dolfinx::geometry::compute_closest_point
std::pair< int, double > compute_closest_point(const BoundingBoxTree &tree, const Eigen::Vector3d &p)
Compute closest point and distance to a given point.
Definition: utils.cpp:510
dolfinx::geometry::compute_entity_collisions
std::pair< std::vector< int >, std::vector< int > > compute_entity_collisions(const BoundingBoxTree &tree0, const BoundingBoxTree &tree1, const mesh::Mesh &mesh0, const mesh::Mesh &mesh1)
Compute all collisions between entities and BoundingBoxTree.
Definition: utils.cpp:381
dolfinx::geometry::compute_closest_entity
std::pair< int, double > compute_closest_entity(const BoundingBoxTree &tree, const BoundingBoxTree &tree_midpoint, const Eigen::Vector3d &p, const mesh::Mesh &mesh)
Compute closest mesh entity and distance to the point. The tree must have been initialised with topol...
Definition: utils.cpp:479
dolfinx::geometry::compute_first_collision
int compute_first_collision(const BoundingBoxTree &tree, const Eigen::Vector3d &p)
Compute first collision between bounding boxes and point.
Definition: utils.cpp:423
dolfinx::geometry::point_in_bbox
bool point_in_bbox(const Eigen::Array< double, 2, 3, Eigen::RowMajor > &b, const Eigen::Vector3d &x, double rtol=1e-14)
Check whether point (x) is in bounding box.
Definition: utils.cpp:536
dolfinx::geometry::compute_squared_distance_bbox
double compute_squared_distance_bbox(const Eigen::Array< double, 2, 3, Eigen::RowMajor > &b, const Eigen::Vector3d &x)
Compute squared distance between point and bounding box wih index "node". Returns zero if point is in...
Definition: utils.cpp:469
dolfinx::geometry::bbox_in_bbox
bool bbox_in_bbox(const Eigen::Array< double, 2, 3, Eigen::RowMajor > &a, const Eigen::Array< double, 2, 3, Eigen::RowMajor > &b, double rtol=1e-14)
Check whether bounding box a collides with bounding box (b)
Definition: utils.cpp:460
dolfinx::geometry::compute_process_collisions
std::vector< int > compute_process_collisions(const BoundingBoxTree &tree, const Eigen::Vector3d &p)
Compute all collisions between processes and Point returning a list of process ranks.
Definition: utils.cpp:446
dolfinx::geometry::compute_collisions
std::pair< std::vector< int >, std::vector< int > > compute_collisions(const BoundingBoxTree &tree0, const BoundingBoxTree &tree1)
Compute all collisions between bounding boxes and BoundingBoxTree.
Definition: utils.cpp:365
dolfinx::geometry::compute_first_entity_collision
int compute_first_entity_collision(const BoundingBoxTree &tree, const Eigen::Vector3d &p, const mesh::Mesh &mesh)
Compute first collision between entities and point.
Definition: utils.cpp:429
dolfinx::geometry::squared_distance_triangle
double squared_distance_triangle(const Eigen::Vector3d &point, const Eigen::Vector3d &a, const Eigen::Vector3d &b, const Eigen::Vector3d &c)
Compute squared distance to given point. This version takes the three vertex coordinates as 3D points...
Definition: utils.cpp:677