From 8d9ad233cd78f151a01e20242d8e56fc1751c113 Mon Sep 17 00:00:00 2001 From: David Wells Date: Fri, 19 May 2023 08:59:33 -0400 Subject: [PATCH] Relax the check for points inside a pyramid. --- include/deal.II/grid/reference_cell.h | 2 +- tests/grid/reference_cell_type_04.cc | 6 ++++++ 2 files changed, 7 insertions(+), 1 deletion(-) diff --git a/include/deal.II/grid/reference_cell.h b/include/deal.II/grid/reference_cell.h index b769c476d6..db61250178 100644 --- a/include/deal.II/grid/reference_cell.h +++ b/include/deal.II/grid/reference_cell.h @@ -2622,7 +2622,7 @@ ReferenceCell::contains_point(const Point &p, const double tolerance) const // We are inside the pyramid if the distance from the axis is less // than (1-z) - return (distance_from_axis < 1 + tolerance - p[z_coordinate]); + return (distance_from_axis <= 1 + tolerance - p[z_coordinate]); } case ReferenceCells::Wedge: { diff --git a/tests/grid/reference_cell_type_04.cc b/tests/grid/reference_cell_type_04.cc index 07c69f4b44..cbb8a9a2d6 100644 --- a/tests/grid/reference_cell_type_04.cc +++ b/tests/grid/reference_cell_type_04.cc @@ -35,6 +35,12 @@ test(const ReferenceCell &reference_cell) unsigned int n_samples_inside = 0; const unsigned int n_samples = 200000; + // sanity check: does the reference cell contain its own nodes? + for (unsigned int vertex_no : reference_cell.vertex_indices()) + AssertThrow(reference_cell.contains_point( + reference_cell.template vertex(vertex_no)), + ExcInternalError()); + for (unsigned int n = 0; n < n_samples; ++n) { // Choose a random point in the box [-1,1]^d that contains all -- 2.39.5