From: Daniel Garcia-Sanchez Date: Fri, 21 May 2021 19:30:24 +0000 (+0200) Subject: Add point_value_03 test X-Git-Tag: v9.4.0-rc1~72^2~2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=076095d618a59f7a0fb2a6d314745c9951182ce3;p=dealii.git Add point_value_03 test Parallel test for find_active_cell_around_point() --- diff --git a/tests/bits/point_value_03.cc b/tests/bits/point_value_03.cc new file mode 100644 index 0000000000..f02afcad25 --- /dev/null +++ b/tests/bits/point_value_03.cc @@ -0,0 +1,213 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2004 - 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. +// +// --------------------------------------------------------------------- + + +// check VectorTools::point_value +// This test is the parallel version of point_value_01 + + + +#include +#include +#include + +#include + +#include + +#include + +#include +#include +#include +#include +#include + +#include +#include + +#include + +#include "../tests.h" + + +using VectorType = LinearAlgebraPETSc::MPI::Vector; + +template +class MySquareFunction : public Function +{ +public: + MySquareFunction() + : Function() + {} + + virtual VectorType::value_type + value(const Point &p, const unsigned int component = 0) const + { + return (component + 1) * p.square() + 1; + } + + virtual void + vector_value(const Point & p, + Vector &values) const + { + values(0) = value(p, 0); + } +}; + + +template +class MyExpFunction : public Function +{ +public: + MyExpFunction() + : Function() + {} + + virtual double + value(const Point &p, const unsigned int component) const + { + return std::exp(p(0)); + } + + virtual void + vector_value(const Point &p, Vector &values) const + { + values(0) = value(p, 0); + } +}; + + + +template +void +make_mesh(Triangulation &tria) +{ + GridGenerator::hyper_cube(tria, -1, 1); + + // refine the mesh in a random way so as to + // generate as many cells with + // hanging nodes as possible + tria.refine_global(4 - dim); + const double steps[4] = {/*d=0*/ 0, 7, 3, 3}; + for (unsigned int i = 0; i < steps[dim]; ++i) + { + typename Triangulation::active_cell_iterator cell = + tria.begin_active(); + for (unsigned int index = 0; cell != tria.end(); ++cell, ++index) + if (index % (3 * dim) == 0) + cell->set_refine_flag(); + tria.execute_coarsening_and_refinement(); + } +} + + + +template +void +check() +{ + parallel::distributed::Triangulation tria(MPI_COMM_WORLD); + make_mesh(tria); + + FE_Q element(QIterated<1>(QTrapez<1>(), 3)); + DoFHandler dof(tria); + dof.distribute_dofs(element); + + static const MySquareFunction function; + + + IndexSet locally_owned_dofs; + + VectorType v; + + locally_owned_dofs = dof.locally_owned_dofs(); + + v.reinit(locally_owned_dofs, MPI_COMM_WORLD); + VectorTools::interpolate(dof, function, v); + + // v = 1.; + + // for the following points, check the + // function value, output it, and + // verify that the value retrieved from + // the interpolated function is close + // enough to that of the real function + // + // also verify that the actual value is + // roughly correct + Point p[3]; + for (unsigned int d = 0; d < dim; ++d) + { + p[0][d] = 0; + p[1][d] = 0.5; + p[2][d] = 1. / 3.; + } + Vector value(1); + const auto & mapping = get_default_linear_mapping(tria); + GridTools::Cache cache(tria, mapping); + typename Triangulation::active_cell_iterator cell_hint{}; + std::vector marked_vertices = {}; + const double tolerance = 1.e-10; + for (unsigned int i = 0; i < 3; ++i) + { + { + bool point_in_locally_owned_cell = false; + value(0) = 0; + { + auto cell_and_ref_point = GridTools::find_active_cell_around_point( + cache, p[i], cell_hint, marked_vertices, tolerance); + if (cell_and_ref_point.first.state() == IteratorState::valid) + { + cell_hint = cell_and_ref_point.first; + point_in_locally_owned_cell = + cell_and_ref_point.first->is_locally_owned(); + } + } + if (point_in_locally_owned_cell) + { + VectorTools::point_value(dof, v, p[i], value); + } + } + value(0) = Utilities::MPI::sum(value(0), MPI_COMM_WORLD); + deallog << -value(0) << std::endl; + + Assert(std::abs(value(0) - function.value(p[i])) < 2e-4, + ExcInternalError()); + } + + deallog << "OK" << std::endl; +} + + +int +main(int argc, char **argv) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); + + MPI_Comm mpi_communicator = MPI_COMM_WORLD; + + + MPILogInitAll log; + + deallog << std::setprecision(4); + + deallog.push("2d"); + check<2>(); + deallog.pop(); + deallog.push("3d"); + check<3>(); + deallog.pop(); +} diff --git a/tests/bits/point_value_03.with_petsc=true.with_petsc_with_complex=false.mpirun=1.output b/tests/bits/point_value_03.with_petsc=true.with_petsc_with_complex=false.mpirun=1.output new file mode 100644 index 0000000000..2a7da639d8 --- /dev/null +++ b/tests/bits/point_value_03.with_petsc=true.with_petsc_with_complex=false.mpirun=1.output @@ -0,0 +1,9 @@ + +DEAL:0:2d::-1.000 +DEAL:0:2d::-1.500 +DEAL:0:2d::-1.222 +DEAL:0:2d::OK +DEAL:0:3d::-1.000 +DEAL:0:3d::-1.750 +DEAL:0:3d::-1.333 +DEAL:0:3d::OK diff --git a/tests/bits/point_value_03.with_petsc=true.with_petsc_with_complex=false.mpirun=7.output b/tests/bits/point_value_03.with_petsc=true.with_petsc_with_complex=false.mpirun=7.output new file mode 100644 index 0000000000..2a7da639d8 --- /dev/null +++ b/tests/bits/point_value_03.with_petsc=true.with_petsc_with_complex=false.mpirun=7.output @@ -0,0 +1,9 @@ + +DEAL:0:2d::-1.000 +DEAL:0:2d::-1.500 +DEAL:0:2d::-1.222 +DEAL:0:2d::OK +DEAL:0:3d::-1.000 +DEAL:0:3d::-1.750 +DEAL:0:3d::-1.333 +DEAL:0:3d::OK diff --git a/tests/bits/point_value_03.with_petsc=true.with_petsc_with_complex=true.mpirun=1.output b/tests/bits/point_value_03.with_petsc=true.with_petsc_with_complex=true.mpirun=1.output new file mode 100644 index 0000000000..435c47dfb7 --- /dev/null +++ b/tests/bits/point_value_03.with_petsc=true.with_petsc_with_complex=true.mpirun=1.output @@ -0,0 +1,9 @@ + +DEAL:0:2d::(-1.000,0.000) +DEAL:0:2d::(-1.500,0.000) +DEAL:0:2d::(-1.222,0.000) +DEAL:0:2d::OK +DEAL:0:3d::(-1.000,0.000) +DEAL:0:3d::(-1.750,0.000) +DEAL:0:3d::(-1.333,0.000) +DEAL:0:3d::OK diff --git a/tests/bits/point_value_03.with_petsc=true.with_petsc_with_complex=true.mpirun=7.output b/tests/bits/point_value_03.with_petsc=true.with_petsc_with_complex=true.mpirun=7.output new file mode 100644 index 0000000000..435c47dfb7 --- /dev/null +++ b/tests/bits/point_value_03.with_petsc=true.with_petsc_with_complex=true.mpirun=7.output @@ -0,0 +1,9 @@ + +DEAL:0:2d::(-1.000,0.000) +DEAL:0:2d::(-1.500,0.000) +DEAL:0:2d::(-1.222,0.000) +DEAL:0:2d::OK +DEAL:0:3d::(-1.000,0.000) +DEAL:0:3d::(-1.750,0.000) +DEAL:0:3d::(-1.333,0.000) +DEAL:0:3d::OK