From 6eb34acc99175437f4c335edf0b2dbf163ca6ba1 Mon Sep 17 00:00:00 2001 From: Johannes Heinz Date: Wed, 24 May 2023 14:45:43 +0200 Subject: [PATCH] drop inexact CGAL kernels Co-authored-by: Marco Feder --- source/cgal/intersections.cc | 54 +++++------ ...c => cgal_intersect_simplices_3d_3d_01.cc} | 0 ... cgal_intersect_simplices_3d_3d_01.output} | 0 .../cgal/cgal_intersect_simplices_3d_3d_02.cc | 94 +++++++++++++++++++ .../cgal_intersect_simplices_3d_3d_02.output | 2 + 5 files changed, 121 insertions(+), 29 deletions(-) rename tests/cgal/{cgal_intersect_simplices_3d_3d.cc => cgal_intersect_simplices_3d_3d_01.cc} (100%) rename tests/cgal/{cgal_intersect_simplices_3d_3d.output => cgal_intersect_simplices_3d_3d_01.output} (100%) create mode 100644 tests/cgal/cgal_intersect_simplices_3d_3d_02.cc create mode 100644 tests/cgal/cgal_intersect_simplices_3d_3d_02.output diff --git a/source/cgal/intersections.cc b/source/cgal/intersections.cc index 245103a748..f2cab839da 100644 --- a/source/cgal/intersections.cc +++ b/source/cgal/intersections.cc @@ -39,7 +39,6 @@ # include # include # include -# include # include # include # include @@ -62,28 +61,25 @@ DEAL_II_NAMESPACE_OPEN namespace CGALWrappers { - using K = CGAL::Exact_predicates_exact_constructions_kernel_with_sqrt; - using K_exact = CGAL::Exact_predicates_exact_constructions_kernel; - using K_inexact = CGAL::Exact_predicates_inexact_constructions_kernel; - using CGALPolygon = CGAL::Polygon_2; - using Polygon_with_holes_2 = CGAL::Polygon_with_holes_2; - using CGALTriangle2 = K::Triangle_2; - using CGALTriangle3 = K::Triangle_3; - using CGALTriangle3_exact = K_exact::Triangle_3; - using CGALPoint2 = K::Point_2; - using CGALPoint3 = K::Point_3; - using CGALPoint3_exact = K_exact::Point_3; - using CGALPoint3_inexact = K_inexact::Point_3; - using CGALSegment2 = K::Segment_2; - using Surface_mesh = CGAL::Surface_mesh; - using CGALSegment3 = K::Segment_3; - using CGALSegment3_exact = K_exact::Segment_3; - using CGALTetra = K::Tetrahedron_3; - using CGALTetra_exact = K_exact::Tetrahedron_3; - using Triangulation2 = CGAL::Triangulation_2; - using Triangulation3 = CGAL::Triangulation_3; - using Triangulation3_exact = CGAL::Triangulation_3; - using Triangulation3_inexact = CGAL::Triangulation_3; + using K = CGAL::Exact_predicates_exact_constructions_kernel_with_sqrt; + using K_exact = CGAL::Exact_predicates_exact_constructions_kernel; + using CGALPolygon = CGAL::Polygon_2; + using Polygon_with_holes_2 = CGAL::Polygon_with_holes_2; + using CGALTriangle2 = K::Triangle_2; + using CGALTriangle3 = K::Triangle_3; + using CGALTriangle3_exact = K_exact::Triangle_3; + using CGALPoint2 = K::Point_2; + using CGALPoint3 = K::Point_3; + using CGALPoint3_exact = K_exact::Point_3; + using CGALSegment2 = K::Segment_2; + using Surface_mesh = CGAL::Surface_mesh; + using CGALSegment3 = K::Segment_3; + using CGALSegment3_exact = K_exact::Segment_3; + using CGALTetra = K::Tetrahedron_3; + using CGALTetra_exact = K_exact::Tetrahedron_3; + using Triangulation2 = CGAL::Triangulation_2; + using Triangulation3 = CGAL::Triangulation_3; + using Triangulation3_exact = CGAL::Triangulation_3; struct FaceInfo2 { @@ -682,25 +678,25 @@ namespace CGALWrappers AssertDimension(hexa0.size(), 8); AssertDimension(hexa0.size(), hexa1.size()); - std::array pts_hex0; - std::array pts_hex1; + std::array pts_hex0; + std::array pts_hex1; std::transform( hexa0.begin(), hexa0.end(), pts_hex0.begin(), - &CGALWrappers::dealii_point_to_cgal_point); + &CGALWrappers::dealii_point_to_cgal_point); std::transform( hexa1.begin(), hexa1.end(), pts_hex1.begin(), - &CGALWrappers::dealii_point_to_cgal_point); + &CGALWrappers::dealii_point_to_cgal_point); Surface_mesh surf0, surf1, sm; // Subdivide hex into tetrahedrons std::vector, 4>> vertices; - Triangulation3_inexact tria0, tria1; + Triangulation3_exact tria0, tria1; tria0.insert(pts_hex0.begin(), pts_hex0.end()); tria1.insert(pts_hex1.begin(), pts_hex1.end()); @@ -729,7 +725,7 @@ namespace CGALWrappers if (PMP::volume(sm) > tol && test_intersection) { // Collect tetrahedrons - Triangulation3_inexact triangulation_hexa; + Triangulation3_exact triangulation_hexa; triangulation_hexa.insert(sm.points().begin(), sm.points().end()); for (const auto &c : triangulation_hexa.finite_cell_handles()) diff --git a/tests/cgal/cgal_intersect_simplices_3d_3d.cc b/tests/cgal/cgal_intersect_simplices_3d_3d_01.cc similarity index 100% rename from tests/cgal/cgal_intersect_simplices_3d_3d.cc rename to tests/cgal/cgal_intersect_simplices_3d_3d_01.cc diff --git a/tests/cgal/cgal_intersect_simplices_3d_3d.output b/tests/cgal/cgal_intersect_simplices_3d_3d_01.output similarity index 100% rename from tests/cgal/cgal_intersect_simplices_3d_3d.output rename to tests/cgal/cgal_intersect_simplices_3d_3d_01.output diff --git a/tests/cgal/cgal_intersect_simplices_3d_3d_02.cc b/tests/cgal/cgal_intersect_simplices_3d_3d_02.cc new file mode 100644 index 0000000000..d594b6add1 --- /dev/null +++ b/tests/cgal/cgal_intersect_simplices_3d_3d_02.cc @@ -0,0 +1,94 @@ +// --------------------------------------------------------------------- + +// Copyright (C) 2022 by the deal.II authors + +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE.md at +// the top level directory of deal.II. +// +// --------------------------------------------------------------------- + +// Compute intersection of two 3D cells. This additional test is added becase +// intersections are not found with inexact kernels. + +#include + +#include +#include +#include + +#include +#include + +#include "../tests.h" + +void +test_intersection() +{ + MappingQ1<3> mapping; + + Triangulation<3> tria0; + GridGenerator::hyper_cube(tria0); + const auto cell0 = tria0.begin_active(); + { + cell0->vertex(0) = + Point<3>(2.600000000000e-01, 0.000000000000e+00, 0.000000000000e+00); + cell0->vertex(1) = + Point<3>(2.600000000000e-01, 0.000000000000e+00, 6.955875000000e-03); + cell0->vertex(2) = + Point<3>(2.600000000000e-01, 6.955875000000e-03, 0.000000000000e+00); + cell0->vertex(3) = + Point<3>(2.600000000000e-01, 6.158125000000e-03, 6.158125000000e-03); + cell0->vertex(4) = + Point<3>(2.550000000000e-01, 0.000000000000e+00, 0.000000000000e+00); + cell0->vertex(5) = + Point<3>(2.550000000000e-01, 0.000000000000e+00, 6.955875000000e-03); + cell0->vertex(6) = + Point<3>(2.550000000000e-01, 6.955875000000e-03, 0.000000000000e+00); + cell0->vertex(7) = + Point<3>(2.550000000000e-01, 6.158125000000e-03, 6.158125000000e-03); + } + + Triangulation<3> tria1; + GridGenerator::hyper_cube(tria1); + const auto cell1 = tria1.begin_active(); + { + cell1->vertex(0) = + Point<3>(2.600000000000e-01, 0.000000000000e+00, 0.000000000000e+00); + cell1->vertex(1) = + Point<3>(2.600000000000e-01, 0.000000000000e+00, 1.391175000000e-02); + cell1->vertex(2) = + Point<3>(2.600000000000e-01, 1.391175000000e-02, 0.000000000000e+00); + cell1->vertex(3) = + Point<3>(2.600000000000e-01, 1.072075000000e-02, 1.072075000000e-02); + cell1->vertex(4) = + Point<3>(2.500000000000e-01, 0.000000000000e+00, 0.000000000000e+00); + cell1->vertex(5) = + Point<3>(2.500000000000e-01, 0.000000000000e+00, 1.391175000000e-02); + cell1->vertex(6) = + Point<3>(2.500000000000e-01, 1.391175000000e-02, 0.000000000000e+00); + cell1->vertex(7) = + Point<3>(2.500000000000e-01, 1.072075000000e-02, 1.072075000000e-02); + } + + const auto &vec_of_arrays = CGALWrappers::compute_intersection_of_cells( + cell1, cell0, mapping, mapping, 1e-9); + + assert(vec_of_arrays.size() == 18); + + deallog << "OK" << std::endl; +} + +int +main() +{ + initlog(); + test_intersection(); + + return 0; +} diff --git a/tests/cgal/cgal_intersect_simplices_3d_3d_02.output b/tests/cgal/cgal_intersect_simplices_3d_3d_02.output new file mode 100644 index 0000000000..0fd8fc12f0 --- /dev/null +++ b/tests/cgal/cgal_intersect_simplices_3d_3d_02.output @@ -0,0 +1,2 @@ + +DEAL::OK -- 2.39.5