DOLFINx
DOLFINx C++ interface
utils.h
1// Copyright (C) 2019-2021 Garth N. Wells and Jørgen S. Dokken
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 <array>
10#include <dolfinx/graph/AdjacencyList.h>
11#include <vector>
12#include <xtensor/xfixed.hpp>
13#include <xtensor/xshape.hpp>
14#include <xtensor/xtensor.hpp>
15#include <xtl/xspan.hpp>
16
17namespace dolfinx::mesh
18{
19class Mesh;
20}
21
22namespace dolfinx::geometry
23{
24class BoundingBoxTree;
25
32BoundingBoxTree
33create_midpoint_tree(const mesh::Mesh& mesh, int tdim,
34 const xtl::span<const std::int32_t>& entity_indices);
35
41std::vector<std::array<int, 2>>
42compute_collisions(const BoundingBoxTree& tree0, const BoundingBoxTree& tree1);
43
50graph::AdjacencyList<std::int32_t>
51compute_collisions(const BoundingBoxTree& tree,
52 const xt::xtensor<double, 2>& points);
53
63std::vector<std::int32_t> compute_closest_entity(
64 const BoundingBoxTree& tree, const BoundingBoxTree& midpoint_tree,
65 const mesh::Mesh& mesh, const xt::xtensor<double, 2>& points);
66
73 const xt::xtensor_fixed<double, xt::xshape<2, 3>>& b,
74 const xt::xtensor_fixed<double, xt::xshape<3>>& x);
75
83xt::xtensor<double, 2>
84shortest_vector(const mesh::Mesh& mesh, int dim,
85 const xtl::span<const std::int32_t>& entities,
86 const xt::xtensor<double, 2>& points);
87
100xt::xtensor<double, 1>
101squared_distance(const mesh::Mesh& mesh, int dim,
102 const xtl::span<const std::int32_t>& entities,
103 const xt::xtensor<double, 2>& points);
104
116graph::AdjacencyList<int> compute_colliding_cells(
117 const mesh::Mesh& mesh,
118 const graph::AdjacencyList<std::int32_t>& candidate_cells,
119 const xt::xtensor<double, 2>& points);
120} // namespace dolfinx::geometry
Geometry data structures and algorithms.
Definition: BoundingBoxTree.h:21
std::vector< std::int32_t > compute_closest_entity(const BoundingBoxTree &tree, const BoundingBoxTree &midpoint_tree, const mesh::Mesh &mesh, const xt::xtensor< double, 2 > &points)
Compute closest mesh entity to a point.
Definition: utils.cpp:307
std::vector< std::array< int, 2 > > compute_collisions(const BoundingBoxTree &tree0, const BoundingBoxTree &tree1)
Compute all collisions between two BoundingBoxTrees (local to process)
Definition: utils.cpp:268
BoundingBoxTree create_midpoint_tree(const mesh::Mesh &mesh, int tdim, const xtl::span< const std::int32_t > &entity_indices)
Create a bounding box tree for a subset of entities (local to process) based on the entity midpoints.
Definition: utils.cpp:247
xt::xtensor< double, 2 > shortest_vector(const mesh::Mesh &mesh, int dim, const xtl::span< const std::int32_t > &entities, const xt::xtensor< double, 2 > &points)
Compute the shortest vector from a mesh entity to a point.
Definition: utils.cpp:371
double compute_squared_distance_bbox(const xt::xtensor_fixed< double, xt::xshape< 2, 3 > > &b, const xt::xtensor_fixed< double, xt::xshape< 3 > > &x)
Compute squared distance between point and bounding box.
Definition: utils.cpp:359
xt::xtensor< double, 1 > squared_distance(const mesh::Mesh &mesh, int dim, const xtl::span< const std::int32_t > &entities, const xt::xtensor< double, 2 > &points)
Compute the squared distance between a point and a mesh entity. The distance is computed between the ...
Definition: utils.cpp:442
graph::AdjacencyList< int > compute_colliding_cells(const mesh::Mesh &mesh, const graph::AdjacencyList< std::int32_t > &candidate_cells, const xt::xtensor< double, 2 > &points)
From a Mesh, find which cells collide with a set of points.
Definition: utils.cpp:449
Mesh data structures and algorithms on meshes.
Definition: DofMap.h:30