From ad1b7e826675f8f2abdc8704bcde973eebbd1164 Mon Sep 17 00:00:00 2001 From: Martin Kronbichler Date: Fri, 14 Apr 2023 14:09:35 +0200 Subject: [PATCH] Test case --- .../compute_point_locations_try_all_02.cc | 96 +++++++++++++++++++ .../compute_point_locations_try_all_02.output | 24 +++++ 2 files changed, 120 insertions(+) create mode 100644 tests/grid/compute_point_locations_try_all_02.cc create mode 100644 tests/grid/compute_point_locations_try_all_02.output diff --git a/tests/grid/compute_point_locations_try_all_02.cc b/tests/grid/compute_point_locations_try_all_02.cc new file mode 100644 index 0000000000..3519200747 --- /dev/null +++ b/tests/grid/compute_point_locations_try_all_02.cc @@ -0,0 +1,96 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2019 - 2020 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. +// +// --------------------------------------------------------------------- + +// Test GridTools::compute_point_locations_try_all for the case where +// some points only lie outside a cell by some roundoff, and another case +// where they are farther away + +#include + +#include +#include +#include +#include + +#include "../tests.h" + +template +void +check_found(const Triangulation & tria, + const GridTools::Cache & cache, + const std::vector> &points) +{ + auto cell_qpoint_map = + GridTools::compute_point_locations_try_all(cache, points); + const auto &cells = std::get<0>(cell_qpoint_map); + const auto &qpoints = std::get<1>(cell_qpoint_map); + const auto &indices = std::get<2>(cell_qpoint_map); + const auto &other_p = std::get<3>(cell_qpoint_map); + size_t n_cells = cells.size(); + + deallog << "Points found in " << n_cells << " cells" << std::endl; + + for (unsigned int i = 0; i < n_cells; ++i) + { + deallog << "On cell " << cells[i]->id() << " found:" << std::endl; + unsigned int j = 0; + for (const Point &p : qpoints[i]) + deallog << "real " << points[indices[i][j++]] << " unit " << p + << std::endl; + } + deallog << "Points not found: "; + for (const auto i : other_p) + deallog << points[i] << " "; + deallog << std::endl << std::endl; +} + + +int +main() +{ + initlog(); + + deallog << std::setprecision(18); + + const int dim = 2; + Triangulation triangulation; + GridGenerator::hyper_cube(triangulation, -0.3, 0.7); + triangulation.refine_global(2); + + MappingQ mapping(1); + GridTools::Cache cache(triangulation, mapping); + + // Point is outside by 1e-16 + Point p0(-0.299999999999999989, -0.247168783648703205), + p1(-0.300000000000000044, -0.102831216351296786); + + check_found(triangulation, cache, {{p0, p1}}); + + // Shift point by 1e-8 + p0[0] -= 1e-8; + + check_found(triangulation, cache, {{p0, p1}}); + + // Shift point by more + p0[0] += 0.5; + p1[0] += 1e-15; + + check_found(triangulation, cache, {{p0, p1}}); + + // Shift outside + p0[0] += 1.; + + check_found(triangulation, cache, {{p0, p1}}); +} diff --git a/tests/grid/compute_point_locations_try_all_02.output b/tests/grid/compute_point_locations_try_all_02.output new file mode 100644 index 0000000000..2b11a01782 --- /dev/null +++ b/tests/grid/compute_point_locations_try_all_02.output @@ -0,0 +1,24 @@ + +DEAL::Points found in 1 cells +DEAL::On cell 0_2:00 found: +DEAL::real -0.299999999999999989 -0.247168783648703205 unit 0.00000000000000000 0.211324865405187134 +DEAL::real -0.300000000000000044 -0.102831216351296786 unit -2.22044604925031308e-16 0.788675134594812866 +DEAL::Points not found: +DEAL:: +DEAL::Points found in 1 cells +DEAL::On cell 0_2:00 found: +DEAL::real -0.300000000000000044 -0.102831216351296786 unit -2.22044604925031308e-16 0.788675134594812866 +DEAL::Points not found: -0.300000009999999984 -0.247168783648703205 +DEAL:: +DEAL::Points found in 2 cells +DEAL::On cell 0_2:01 found: +DEAL::real 0.199999990000000016 -0.247168783648703205 unit 0.999999960000000132 0.211324865405187134 +DEAL::On cell 0_2:00 found: +DEAL::real -0.299999999999999045 -0.102831216351296786 unit 3.77475828372553224e-15 0.788675134594812866 +DEAL::Points not found: +DEAL:: +DEAL::Points found in 1 cells +DEAL::On cell 0_2:00 found: +DEAL::real -0.299999999999999045 -0.102831216351296786 unit 3.77475828372553224e-15 0.788675134594812866 +DEAL::Points not found: 1.19999999000000002 -0.247168783648703205 +DEAL:: -- 2.39.5