21#ifdef DEAL_II_WITH_CGAL
31# include <CGAL/Boolean_set_operations_2.h>
34# include <CGAL/Cartesian.h>
35# include <CGAL/Circular_kernel_intersections.h>
36# include <CGAL/Constrained_Delaunay_triangulation_2.h>
37# include <CGAL/Delaunay_mesh_face_base_2.h>
38# include <CGAL/Delaunay_mesh_size_criteria_2.h>
39# include <CGAL/Delaunay_mesher_2.h>
40# include <CGAL/Delaunay_triangulation_2.h>
41# include <CGAL/Exact_predicates_exact_constructions_kernel_with_sqrt.h>
42# include <CGAL/Kernel_traits.h>
43# include <CGAL/Polygon_2.h>
44# include <CGAL/Polygon_with_holes_2.h>
45# include <CGAL/Projection_traits_xy_3.h>
46# include <CGAL/Segment_3.h>
47# include <CGAL/Simple_cartesian.h>
48# include <CGAL/Surface_mesh/Surface_mesh.h>
49# include <CGAL/Tetrahedron_3.h>
50# include <CGAL/Triangle_2.h>
51# include <CGAL/Triangle_3.h>
52# include <CGAL/Triangulation_2.h>
53# include <CGAL/Triangulation_3.h>
54# include <CGAL/Triangulation_face_base_with_id_2.h>
55# include <CGAL/Triangulation_face_base_with_info_2.h>
107 template <
class T,
class... Types>
111 return std::get_if<T>(v);
114 template <
class T,
class... Types>
118 return boost::get<T>(v);
129 template <
typename TargetVariant>
130 struct Repackage : boost::static_visitor<TargetVariant>
132 template <
typename T>
134 operator()(
const T &t)
const
146 template <
typename... Types>
147 std::optional<std::variant<Types...>>
168 template <
typename... Types>
169 const std::optional<std::variant<Types...>> &
182 std::list<CDT::Edge> &
border)
184 if (start->info().nesting_level != -1)
188 std::list<Face_handle> queue;
189 queue.push_back(start);
190 while (!queue.empty())
194 if (
fh->info().nesting_level == -1)
196 fh->info().nesting_level = index;
197 for (
int i = 0; i < 3; i++)
201 if (n->info().nesting_level == -1)
203 if (
ct.is_constrained(e))
218 for (CDT::Face_handle f :
cdt.all_face_handles())
220 f->info().nesting_level = -1;
222 std::list<CDT::Edge>
border;
226 CDT::Edge e =
border.front();
229 if (n->info().nesting_level == -1)
248 std::vector<CGALPoint2>>>
256 std::array<CGALPoint2, 3>
pts0,
pts1;
261 &CGALWrappers::dealii_point_to_cgal_point<CGALPoint2, 2>);
266 &CGALWrappers::dealii_point_to_cgal_point<CGALPoint2, 2>);
275 std::optional<std::variant<CGALPoint2, CGALSegment2>>
283 std::array<CGALPoint2, 3>
pts0;
284 std::array<CGALPoint2, 2>
pts1;
289 &CGALWrappers::dealii_point_to_cgal_point<CGALPoint2, 2>);
291 std::transform(
segment.begin(),
294 &CGALWrappers::dealii_point_to_cgal_point<CGALPoint2, 2>);
305 std::vector<Polygon_with_holes_2>
312 std::array<CGALPoint2, 4>
pts0,
pts1;
317 &CGALWrappers::dealii_point_to_cgal_point<CGALPoint2, 2>);
322 &CGALWrappers::dealii_point_to_cgal_point<CGALPoint2, 2>);
327 std::vector<Polygon_with_holes_2>
poly_list;
336 std::optional<std::variant<CGALPoint3, CGALSegment3>>
341# if DEAL_II_CGAL_VERSION_GTE(5, 5, 0)
346 std::array<CGALPoint3, 4>
pts0;
347 std::array<CGALPoint3, 2>
pts1;
352 &CGALWrappers::dealii_point_to_cgal_point<CGALPoint3, 3>);
354 std::transform(
segment.begin(),
357 &CGALWrappers::dealii_point_to_cgal_point<CGALPoint3, 3>);
367 "This function requires a version of CGAL greater or equal than 5.5."));
379 std::vector<CGALPoint3>>>
384# if DEAL_II_CGAL_VERSION_GTE(5, 5, 0)
389 std::array<CGALPoint3, 4>
pts0;
390 std::array<CGALPoint3, 3>
pts1;
395 &CGALWrappers::dealii_point_to_cgal_point<CGALPoint3, 3>);
400 &CGALWrappers::dealii_point_to_cgal_point<CGALPoint3, 3>);
411 "This function requires a version of CGAL greater or equal than 5.5."));
419 std::vector<std::array<Point<2>, 3>>
433 const unsigned int size_poly = poly.size();
439 {{CGALWrappers::cgal_point_to_dealii_point<2>(poly.vertex(0)),
440 CGALWrappers::cgal_point_to_dealii_point<2>(poly.vertex(1)),
441 CGALWrappers::cgal_point_to_dealii_point<2>(
447 std::vector<std::array<Point<2>, 3>>
collection;
450 cdt.insert_constraint(poly.vertices_begin(),
458 if (f->info().in_domain() &&
462 {{CGALWrappers::cgal_point_to_dealii_point<2>(
463 cdt.triangle(f).vertex(0)),
464 CGALWrappers::cgal_point_to_dealii_point<2>(
465 cdt.triangle(f).vertex(1)),
466 CGALWrappers::cgal_point_to_dealii_point<2>(
467 cdt.triangle(f).vertex(2))}});
485 std::vector<std::array<Point<2>, 2>>
493 std::array<CGALPoint2, 4>
pts;
495 std::transform(quad.begin(),
498 &CGALWrappers::dealii_point_to_cgal_point<CGALPoint2, 2>);
503 CGALWrappers::dealii_point_to_cgal_point<CGALPoint2>(line[0]),
504 CGALWrappers::dealii_point_to_cgal_point<CGALPoint2>(line[1]));
506 cdt.insert_constraint(poly.vertices_begin(), poly.vertices_end(),
true);
507 std::vector<std::array<Point<2>, 2>> vertices;
511 if (f->info().in_domain() &&
520 {{CGALWrappers::cgal_point_to_dealii_point<2>((*s)[0]),
521 CGALWrappers::cgal_point_to_dealii_point<2>((*s)[1])}});
530 std::vector<std::array<Point<3>, 2>>
535# if DEAL_II_CGAL_VERSION_GTE(5, 5, 0)
540 std::array<CGALPoint3_exact, 8>
pts;
546 &CGALWrappers::dealii_point_to_cgal_point<CGALPoint3_exact, 3>);
549 CGALWrappers::dealii_point_to_cgal_point<CGALPoint3_exact>(line[0]),
550 CGALWrappers::dealii_point_to_cgal_point<CGALPoint3_exact>(line[1]));
554 std::vector<std::array<Point<3>, 2>> vertices;
567 if (s->squared_length() > tol * tol)
570 {{CGALWrappers::cgal_point_to_dealii_point<3>(
572 CGALWrappers::cgal_point_to_dealii_point<3>(
583 "This function requires a version of CGAL greater or equal than 5.5."));
591 std::vector<std::array<Point<3>, 3>>
596# if DEAL_II_CGAL_VERSION_GTE(5, 5, 0)
601 std::array<CGALPoint3_exact, 8>
pts_hex;
602 std::array<CGALPoint3_exact, 4>
pts_quad;
608 &CGALWrappers::dealii_point_to_cgal_point<CGALPoint3_exact, 3>);
614 &CGALWrappers::dealii_point_to_cgal_point<CGALPoint3_exact, 3>);
617 std::vector<std::array<Point<3>, 3>> vertices;
648 if (
const std::vector<CGALPoint3_exact> *
vps =
662 std::array<Point<3>, 3>
verts = {
663 {CGALWrappers::cgal_point_to_dealii_point<3>(
665 CGALWrappers::cgal_point_to_dealii_point<3>(
667 CGALWrappers::cgal_point_to_dealii_point<3>(
670 vertices.push_back(
verts);
683 "This function requires a version of CGAL greater or equal than 5.5."));
691 std::vector<std::array<Point<3>, 4>>
699 std::array<CGALPoint3_exact, 8>
pts_hex0;
700 std::array<CGALPoint3_exact, 8>
pts_hex1;
706 &CGALWrappers::dealii_point_to_cgal_point<CGALPoint3_exact, 3>);
712 &CGALWrappers::dealii_point_to_cgal_point<CGALPoint3_exact, 3>);
716 std::vector<std::array<Point<3>, 4>> vertices;
722 for (
const auto &
c0 :
tria0.finite_cell_handles())
731 for (
const auto &
c1 :
tria1.finite_cell_handles())
753 {{CGALWrappers::cgal_point_to_dealii_point<3>(
755 CGALWrappers::cgal_point_to_dealii_point<3>(
757 CGALWrappers::cgal_point_to_dealii_point<3>(
759 CGALWrappers::cgal_point_to_dealii_point<3>(
774 template <
int structdim0,
int structdim1,
int spacedim>
775 std::vector<std::array<Point<spacedim>,
structdim1 + 1>>
787 "The intersection cannot be computed as at least one of the two cells has no vertices."));
844 template <
int structdim0,
int structdim1,
int spacedim>
845 std::vector<std::array<Point<spacedim>,
structdim1 + 1>>
854 ReferenceCells::get_hypercube<structdim0>().n_vertices(),
857 ReferenceCells::get_hypercube<structdim1>().n_vertices(),
876# include "cgal/intersections.inst"
892std::vector<std::array<Point<spacedim>,
structdim1 + 1>>
904template <
int structdim0,
int structdim1,
int spacedim>
905std::vector<std::array<Point<spacedim>,
structdim1 + 1>>
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_DISABLE_EXTRA_DIAGNOSTICS
#define DEAL_II_NAMESPACE_CLOSE
#define DEAL_II_ENABLE_EXTRA_DIAGNOSTICS
static ::ExceptionBase & ExcNotImplemented()
#define Assert(cond, exc)
#define AssertDimension(dim1, dim2)
static ::ExceptionBase & ExcNeedsCGAL()
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
#define DEAL_II_NOT_IMPLEMENTED()
std::vector< std::array< Point< 2 >, 2 > > compute_intersection_quad_line(const ArrayView< const Point< 2 > > &quad, const ArrayView< const Point< 2 > > &line, const double tol)
std::vector< std::array< Point< 3 >, 3 > > compute_intersection_hexa_quad(const ArrayView< const Point< 3 > > &hexa, const ArrayView< const Point< 3 > > &quad, const double tol)
void mark_domains(CDT &ct, Face_handle start, int index, std::list< CDT::Edge > &border)
std::optional< std::variant< CGALPoint3, CGALSegment3 > > compute_intersection_tetra_segment(const ArrayView< const Point< 3 > > &tetrahedron, const ArrayView< const Point< 3 > > &segment)
std::vector< Polygon_with_holes_2 > compute_intersection_rect_rect(const ArrayView< const Point< 2 > > &rectangle0, const ArrayView< const Point< 2 > > &rectangle1)
std::vector< std::array< Point< 3 >, 4 > > compute_intersection_hexa_hexa(const ArrayView< const Point< 3 > > &hexa0, const ArrayView< const Point< 3 > > &hexa1, const double tol)
std::vector< std::array< Point< 3 >, 2 > > compute_intersection_hexa_line(const ArrayView< const Point< 3 > > &hexa, const ArrayView< const Point< 3 > > &line, const double tol)
std::optional< std::variant< CGALPoint3, CGALSegment3, CGALTriangle3, std::vector< CGALPoint3 > > > compute_intersection_tetra_triangle(const ArrayView< const Point< 3 > > &tetrahedron, const ArrayView< const Point< 3 > > &triangle)
std::optional< std::variant< CGALPoint2, CGALSegment2 > > compute_intersection_triangle_segment(const ArrayView< const Point< 2 > > &triangle, const ArrayView< const Point< 2 > > &segment)
std::optional< std::variant< CGALPoint2, CGALSegment2, CGALTriangle2, std::vector< CGALPoint2 > > > compute_intersection_triangle_triangle(const ArrayView< const Point< 2 > > &triangle0, const ArrayView< const Point< 2 > > &triangle1)
std::vector< std::array< Point< 2 >, 3 > > compute_intersection_quad_quad(const ArrayView< const Point< 2 > > &quad0, const ArrayView< const Point< 2 > > &quad1, const double tol)
K::Tetrahedron_3 CGALTetra
CGAL::Surface_mesh< K_exact::Point_3 > Surface_mesh
K::Segment_2 CGALSegment2
const T * get_if_(const std::variant< Types... > *v)
CGAL::Exact_predicates_tag Itag
K_exact::Tetrahedron_3 CGALTetra_exact
K_exact::Triangle_3 CGALTriangle3_exact
CDT::Vertex_handle Vertex_handle
K::Triangle_3 CGALTriangle3
K::Segment_3 CGALSegment3
CDT::Face_handle Face_handle
CGAL::Triangulation_3< K_exact > Triangulation3_exact
CGAL::Constrained_Delaunay_triangulation_2< K, Tds, Itag > CDT
std::vector< std::array< Point< spacedim >, structdim1+1 > > compute_intersection_of_cells(const typename Triangulation< structdim0, spacedim >::cell_iterator &cell0, const typename Triangulation< structdim1, spacedim >::cell_iterator &cell1, const Mapping< structdim0, spacedim > &mapping0, const Mapping< structdim1, spacedim > &mapping1, const double tol=1e-9)
CGAL::Delaunay_mesh_size_criteria_2< CDT > Criteria
CGAL::Triangulation_vertex_base_2< K > Vb
CGAL::Polygon_with_holes_2< K > Polygon_with_holes_2
CGAL::Triangulation_data_structure_2< Vb, Fb > Tds
CGAL::Polygon_2< K > CGALPolygon
CGAL::Triangulation_2< K > Triangulation2
CGAL::Exact_predicates_exact_constructions_kernel K_exact
CGAL::Triangulation_3< K > Triangulation3
CGAL::Constrained_triangulation_face_base_2< K, Fbb > CFb
CGAL::Exact_predicates_exact_constructions_kernel_with_sqrt K
K::Triangle_2 CGALTriangle2
CGAL::Delaunay_mesh_face_base_2< K, CFb > Fb
K_exact::Segment_3 CGALSegment3_exact
K_exact::Point_3 CGALPoint3_exact
CGAL::Triangulation_face_base_with_info_2< FaceInfo2, K > Fbb