DOLFIN-X
DOLFIN-X C++ interface
CollisionPredicates.h
1 // Copyright (C) 2014-2016 Anders Logg, August Johansson and Benjamin Kehlet
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 <cstddef>
11 
12 namespace dolfinx
13 {
14 namespace mesh
15 {
16 class MeshEntity;
17 }
18 
19 namespace geometry
20 {
21 
24 
26 {
27 public:
28  //--- High-level collision detection predicates ---
29 
35  static bool collides(const mesh::MeshEntity& entity,
36  const Eigen::Vector3d& point);
37 
43  static bool collides(const mesh::MeshEntity& entity_0,
44  const mesh::MeshEntity& entity_1);
45 
46  //--- Low-level collision detection predicates ---
47 
49  static bool collides_segment_point(const Eigen::Vector3d& p0,
50  const Eigen::Vector3d& p1,
51  const Eigen::Vector3d& point,
52  std::size_t gdim);
53 
55  static bool collides_segment_point_1d(double p0, double p1, double point);
56 
58  static bool collides_segment_point_2d(const Eigen::Vector3d& p0,
59  const Eigen::Vector3d& p1,
60  const Eigen::Vector3d& point);
61 
63  static bool collides_segment_point_3d(const Eigen::Vector3d& p0,
64  const Eigen::Vector3d& p1,
65  const Eigen::Vector3d& point);
66 
68  static bool collides_segment_segment(const Eigen::Vector3d& p0,
69  const Eigen::Vector3d& p1,
70  const Eigen::Vector3d& q0,
71  const Eigen::Vector3d& q1,
72  std::size_t gdim);
73 
76  static bool collides_segment_segment_1d(double p0, double p1, double q0,
77  double q1);
78 
81  static bool collides_segment_segment_2d(const Eigen::Vector3d& p0,
82  const Eigen::Vector3d& p1,
83  const Eigen::Vector3d& q0,
84  const Eigen::Vector3d& q1);
85 
88  static bool collides_segment_segment_3d(const Eigen::Vector3d& p0,
89  const Eigen::Vector3d& p1,
90  const Eigen::Vector3d& q0,
91  const Eigen::Vector3d& q1);
92 
94  static bool collides_triangle_point(const Eigen::Vector3d& p0,
95  const Eigen::Vector3d& p1,
96  const Eigen::Vector3d& p2,
97  const Eigen::Vector3d& point,
98  std::size_t gdim);
99 
101  static bool collides_triangle_point_2d(const Eigen::Vector3d& p0,
102  const Eigen::Vector3d& p1,
103  const Eigen::Vector3d& p2,
104  const Eigen::Vector3d& point);
105 
107  static bool collides_quad_point_2d(const Eigen::Vector3d& p0,
108  const Eigen::Vector3d& p1,
109  const Eigen::Vector3d& p2,
110  const Eigen::Vector3d& p3,
111  const Eigen::Vector3d& point);
112 
114  static bool collides_triangle_point_3d(const Eigen::Vector3d& p0,
115  const Eigen::Vector3d& p1,
116  const Eigen::Vector3d& p2,
117  const Eigen::Vector3d& point);
118 
120  static bool collides_triangle_segment(const Eigen::Vector3d& p0,
121  const Eigen::Vector3d& p1,
122  const Eigen::Vector3d& p2,
123  const Eigen::Vector3d& q0,
124  const Eigen::Vector3d& q1,
125  std::size_t gdim);
126 
129  static bool collides_triangle_segment_2d(const Eigen::Vector3d& p0,
130  const Eigen::Vector3d& p1,
131  const Eigen::Vector3d& p2,
132  const Eigen::Vector3d& q0,
133  const Eigen::Vector3d& q1);
134 
137  static bool collides_triangle_segment_3d(const Eigen::Vector3d& p0,
138  const Eigen::Vector3d& p1,
139  const Eigen::Vector3d& p2,
140  const Eigen::Vector3d& q0,
141  const Eigen::Vector3d& q1);
142 
144  static bool collides_triangle_triangle(
145  const Eigen::Vector3d& p0, const Eigen::Vector3d& p1,
146  const Eigen::Vector3d& p2, const Eigen::Vector3d& q0,
147  const Eigen::Vector3d& q1, const Eigen::Vector3d& q2, std::size_t gdim);
148 
151  static bool collides_triangle_triangle_2d(const Eigen::Vector3d& p0,
152  const Eigen::Vector3d& p1,
153  const Eigen::Vector3d& p2,
154  const Eigen::Vector3d& q0,
155  const Eigen::Vector3d& q1,
156  const Eigen::Vector3d& q2);
157 
160  static bool collides_triangle_triangle_3d(const Eigen::Vector3d& p0,
161  const Eigen::Vector3d& p1,
162  const Eigen::Vector3d& p2,
163  const Eigen::Vector3d& q0,
164  const Eigen::Vector3d& q1,
165  const Eigen::Vector3d& q2);
166 
168  static bool collides_tetrahedron_point_3d(const Eigen::Vector3d& p0,
169  const Eigen::Vector3d& p1,
170  const Eigen::Vector3d& p2,
171  const Eigen::Vector3d& p3,
172  const Eigen::Vector3d& point);
173 
175  static bool collides_tetrahedron_segment_3d(const Eigen::Vector3d& p0,
176  const Eigen::Vector3d& p1,
177  const Eigen::Vector3d& p2,
178  const Eigen::Vector3d& p3,
179  const Eigen::Vector3d& q0,
180  const Eigen::Vector3d& q1);
181 
184  static bool collides_tetrahedron_triangle_3d(const Eigen::Vector3d& p0,
185  const Eigen::Vector3d& p1,
186  const Eigen::Vector3d& p2,
187  const Eigen::Vector3d& p3,
188  const Eigen::Vector3d& q0,
189  const Eigen::Vector3d& q1,
190  const Eigen::Vector3d& q2);
191 
195  const Eigen::Vector3d& p0, const Eigen::Vector3d& p1,
196  const Eigen::Vector3d& p2, const Eigen::Vector3d& p3,
197  const Eigen::Vector3d& q0, const Eigen::Vector3d& q1,
198  const Eigen::Vector3d& q2, const Eigen::Vector3d& q3);
199 };
200 } // namespace geometry
201 } // namespace dolfinx
dolfinx::geometry::CollisionPredicates::collides_triangle_point_2d
static bool collides_triangle_point_2d(const Eigen::Vector3d &p0, const Eigen::Vector3d &p1, const Eigen::Vector3d &p2, const Eigen::Vector3d &point)
Check whether triangle p0-p1-p2 collides with point (2D version)
Definition: CollisionPredicates.cpp:541
dolfinx::geometry::CollisionPredicates::collides_triangle_point
static bool collides_triangle_point(const Eigen::Vector3d &p0, const Eigen::Vector3d &p1, const Eigen::Vector3d &p2, const Eigen::Vector3d &point, std::size_t gdim)
Check whether triangle p0-p1-p2 collides with point.
Definition: CollisionPredicates.cpp:305
dolfinx::geometry::CollisionPredicates::collides_triangle_point_3d
static bool collides_triangle_point_3d(const Eigen::Vector3d &p0, const Eigen::Vector3d &p1, const Eigen::Vector3d &p2, const Eigen::Vector3d &point)
Check whether triangle p0-p1-p2 collides with point (3D version)
Definition: CollisionPredicates.cpp:597
dolfinx::geometry::CollisionPredicates::collides_tetrahedron_point_3d
static bool collides_tetrahedron_point_3d(const Eigen::Vector3d &p0, const Eigen::Vector3d &p1, const Eigen::Vector3d &p2, const Eigen::Vector3d &p3, const Eigen::Vector3d &point)
Check whether tetrahedron p0-p1-p2-p3 collides with point.
Definition: CollisionPredicates.cpp:817
dolfinx::geometry::CollisionPredicates::collides_tetrahedron_segment_3d
static bool collides_tetrahedron_segment_3d(const Eigen::Vector3d &p0, const Eigen::Vector3d &p1, const Eigen::Vector3d &p2, const Eigen::Vector3d &p3, const Eigen::Vector3d &q0, const Eigen::Vector3d &q1)
Check whether tetrahedron p0-p1-p2-p3 collides with segment q0-q1.
Definition: CollisionPredicates.cpp:847
dolfinx::geometry::CollisionPredicates::collides_triangle_segment_2d
static bool collides_triangle_segment_2d(const Eigen::Vector3d &p0, const Eigen::Vector3d &p1, const Eigen::Vector3d &p2, const Eigen::Vector3d &q0, const Eigen::Vector3d &q1)
Check whether triangle p0-p1-p2 collides with segment q0-q1 (2D version)
Definition: CollisionPredicates.cpp:615
dolfinx::geometry::CollisionPredicates::collides_triangle_triangle_3d
static bool collides_triangle_triangle_3d(const Eigen::Vector3d &p0, const Eigen::Vector3d &p1, const Eigen::Vector3d &p2, const Eigen::Vector3d &q0, const Eigen::Vector3d &q1, const Eigen::Vector3d &q2)
Check whether triangle p0-p1-p2 collides with triangle q0-q1-q2 (3D version)
Definition: CollisionPredicates.cpp:764
dolfinx::geometry::CollisionPredicates::collides_segment_point_2d
static bool collides_segment_point_2d(const Eigen::Vector3d &p0, const Eigen::Vector3d &p1, const Eigen::Vector3d &point)
Check whether segment p0-p1 collides with point (2D version)
Definition: CollisionPredicates.cpp:373
dolfinx::geometry::CollisionPredicates::collides_segment_point
static bool collides_segment_point(const Eigen::Vector3d &p0, const Eigen::Vector3d &p1, const Eigen::Vector3d &point, std::size_t gdim)
Check whether segment p0-p1 collides with point.
Definition: CollisionPredicates.cpp:262
dolfinx::geometry::CollisionPredicates::collides_segment_point_3d
static bool collides_segment_point_3d(const Eigen::Vector3d &p0, const Eigen::Vector3d &p1, const Eigen::Vector3d &point)
Check whether segment p0-p1 collides with point (3D version)
Definition: CollisionPredicates.cpp:385
dolfinx::geometry::CollisionPredicates::collides
static bool collides(const mesh::MeshEntity &entity, const Eigen::Vector3d &point)
Check whether entity collides with point.
Definition: CollisionPredicates.cpp:43
dolfinx::mesh::MeshEntity
A MeshEntity represents a mesh entity associated with a specific topological dimension of some Mesh....
Definition: MeshEntity.h:21
dolfinx::geometry::CollisionPredicates::collides_segment_segment_2d
static bool collides_segment_segment_2d(const Eigen::Vector3d &p0, const Eigen::Vector3d &p1, const Eigen::Vector3d &q0, const Eigen::Vector3d &q1)
Check whether segment p0-p1 collides with segment q0-q1 (2D version)
Definition: CollisionPredicates.cpp:433
dolfinx::geometry::CollisionPredicates::collides_triangle_triangle
static bool collides_triangle_triangle(const Eigen::Vector3d &p0, const Eigen::Vector3d &p1, const Eigen::Vector3d &p2, const Eigen::Vector3d &q0, const Eigen::Vector3d &q1, const Eigen::Vector3d &q2, std::size_t gdim)
Check whether triangle p0-p1-p2 collides with triangle q0-q1-q2.
Definition: CollisionPredicates.cpp:346
dolfinx::geometry::CollisionPredicates::collides_triangle_segment
static bool collides_triangle_segment(const Eigen::Vector3d &p0, const Eigen::Vector3d &p1, const Eigen::Vector3d &p2, const Eigen::Vector3d &q0, const Eigen::Vector3d &q1, std::size_t gdim)
Check whether triangle p0-p1-p2 collides with segment q0-q1.
Definition: CollisionPredicates.cpp:325
dolfinx::geometry::CollisionPredicates::collides_segment_point_1d
static bool collides_segment_point_1d(double p0, double p1, double point)
Check whether segment p0-p1 collides with point (1D version)
Definition: CollisionPredicates.cpp:365
dolfinx::geometry::CollisionPredicates
This class implements algorithms for detecting pairwise collisions between mesh entities of varying d...
Definition: CollisionPredicates.h:25
dolfinx::geometry::CollisionPredicates::collides_tetrahedron_triangle_3d
static bool collides_tetrahedron_triangle_3d(const Eigen::Vector3d &p0, const Eigen::Vector3d &p1, const Eigen::Vector3d &p2, const Eigen::Vector3d &p3, const Eigen::Vector3d &q0, const Eigen::Vector3d &q1, const Eigen::Vector3d &q2)
Check whether tetrahedron p0-p1-p2-p3 collides with triangle q0-q1-q2.
Definition: CollisionPredicates.cpp:873
dolfinx::geometry::CollisionPredicates::collides_quad_point_2d
static bool collides_quad_point_2d(const Eigen::Vector3d &p0, const Eigen::Vector3d &p1, const Eigen::Vector3d &p2, const Eigen::Vector3d &p3, const Eigen::Vector3d &point)
Check whether quad p0-p1-p2-p3 collides with point (2D version)
Definition: CollisionPredicates.cpp:571
dolfinx::geometry::CollisionPredicates::collides_segment_segment_3d
static bool collides_segment_segment_3d(const Eigen::Vector3d &p0, const Eigen::Vector3d &p1, const Eigen::Vector3d &q0, const Eigen::Vector3d &q1)
Check whether segment p0-p1 collides with segment q0-q1 (3D version)
Definition: CollisionPredicates.cpp:461
dolfinx::geometry::CollisionPredicates::collides_triangle_triangle_2d
static bool collides_triangle_triangle_2d(const Eigen::Vector3d &p0, const Eigen::Vector3d &p1, const Eigen::Vector3d &p2, const Eigen::Vector3d &q0, const Eigen::Vector3d &q1, const Eigen::Vector3d &q2)
Check whether triangle p0-p1-p2 collides with triangle q0-q1-q2 (2D version)
Definition: CollisionPredicates.cpp:708
dolfinx::geometry::CollisionPredicates::collides_tetrahedron_tetrahedron_3d
static bool collides_tetrahedron_tetrahedron_3d(const Eigen::Vector3d &p0, const Eigen::Vector3d &p1, const Eigen::Vector3d &p2, const Eigen::Vector3d &p3, const Eigen::Vector3d &q0, const Eigen::Vector3d &q1, const Eigen::Vector3d &q2, const Eigen::Vector3d &q3)
Check whether tetrahedron p0-p1-p2-p3 collides with tetrahedron q0-q1-q2.
Definition: CollisionPredicates.cpp:902
dolfinx::geometry::CollisionPredicates::collides_segment_segment
static bool collides_segment_segment(const Eigen::Vector3d &p0, const Eigen::Vector3d &p1, const Eigen::Vector3d &q0, const Eigen::Vector3d &q1, std::size_t gdim)
Check whether segment p0-p1 collides with segment q0-q1.
Definition: CollisionPredicates.cpp:283
dolfinx::geometry::CollisionPredicates::collides_triangle_segment_3d
static bool collides_triangle_segment_3d(const Eigen::Vector3d &p0, const Eigen::Vector3d &p1, const Eigen::Vector3d &p2, const Eigen::Vector3d &q0, const Eigen::Vector3d &q1)
Check whether triangle p0-p1-p2 collides with segment q0-q1 (3D version)
Definition: CollisionPredicates.cpp:639
dolfinx::geometry::CollisionPredicates::collides_segment_segment_1d
static bool collides_segment_segment_1d(double p0, double p1, double q0, double q1)
Check whether segment p0-p1 collides with segment q0-q1 (1D version)
Definition: CollisionPredicates.cpp:419