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 <Eigen/Dense>
10 #include <utility>
11 #include <vector>
12 
13 namespace dolfinx
14 {
15 namespace mesh
16 {
17 class Mesh;
18 } // namespace mesh
19 
20 namespace geometry
21 {
22 class BoundingBoxTree;
23 
27 BoundingBoxTree create_midpoint_tree(const mesh::Mesh& mesh);
28 
33 std::vector<std::array<int, 2>>
34 compute_collisions(const BoundingBoxTree& tree0, const BoundingBoxTree& tree1);
35 
40 std::vector<int> compute_collisions(const BoundingBoxTree& tree,
41  const Eigen::Vector3d& p);
42 
45 std::vector<int> compute_process_collisions(const BoundingBoxTree& tree,
46  const Eigen::Vector3d& p);
47 
50 std::pair<int, double>
51 compute_closest_entity(const BoundingBoxTree& tree,
52  const BoundingBoxTree& tree_midpoint,
53  const Eigen::Vector3d& p, const mesh::Mesh& mesh);
54 
60 std::pair<int, double> compute_closest_point(const BoundingBoxTree& tree,
61  const Eigen::Vector3d& p);
62 
66  const Eigen::Array<double, 2, 3, Eigen::RowMajor>& b,
67  const Eigen::Vector3d& x);
68 
82 double squared_distance(const mesh::Mesh& mesh, int dim, std::int32_t index,
83  const Eigen::Vector3d& p);
84 
93 std::vector<int> select_colliding_cells(const dolfinx::mesh::Mesh& mesh,
94  const std::vector<int>& candidate_cells,
95  const Eigen::Vector3d& point, int n);
96 } // namespace geometry
97 } // namespace dolfinx
A Mesh consists of a set of connected and numbered mesh topological entities, and geometry data.
Definition: Mesh.h:47
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:316
std::vector< std::array< int, 2 > > compute_collisions(const BoundingBoxTree &tree0, const BoundingBoxTree &tree1)
Compute all collisions between two BoundingBoxTrees.
Definition: utils.cpp:246
std::vector< int > select_colliding_cells(const dolfinx::mesh::Mesh &mesh, const std::vector< int > &candidate_cells, const Eigen::Vector3d &point, int n)
From the given Mesh, select up to n cells from the list which actually collide with point p....
Definition: utils.cpp:394
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:265
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:289
double squared_distance(const mesh::Mesh &mesh, int dim, std::int32_t index, const Eigen::Vector3d &p)
Compute squared distance from a given point to the nearest point on a cell (only first order convex c...
Definition: utils.cpp:342
BoundingBoxTree create_midpoint_tree(const mesh::Mesh &mesh)
Create a boundary box tree for cell midpoints.
Definition: utils.cpp:223
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:279