From eba520c87ae9aefa2f37b4cbf499e238e18789b7 Mon Sep 17 00:00:00 2001 From: Jean-Paul Pelteret Date: Sat, 14 Aug 2021 21:28:32 +0200 Subject: [PATCH] Add some tests --- tests/meshworker/scratch_data_10a.cc | 145 ++++++++++++++++++++++ tests/meshworker/scratch_data_10a.output | 18 +++ tests/meshworker/scratch_data_10b.cc | 146 +++++++++++++++++++++++ tests/meshworker/scratch_data_10b.output | 18 +++ 4 files changed, 327 insertions(+) create mode 100644 tests/meshworker/scratch_data_10a.cc create mode 100644 tests/meshworker/scratch_data_10a.output create mode 100644 tests/meshworker/scratch_data_10b.cc create mode 100644 tests/meshworker/scratch_data_10b.output diff --git a/tests/meshworker/scratch_data_10a.cc b/tests/meshworker/scratch_data_10a.cc new file mode 100644 index 0000000000..360c5e6c7e --- /dev/null +++ b/tests/meshworker/scratch_data_10a.cc @@ -0,0 +1,145 @@ +/* --------------------------------------------------------------------- + * + * Copyright (C) 2021 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 that ScratchData returns the correct jumps in, and averages of, +// solution values, gradients, etc. at interfaces +// - Scalar valued FE + +#include +#include + +#include +#include +#include + +#include +#include + +#include + +#include + +#include "../tests.h" + +template +class DiscontinuousFunction : public Function +{ +public: + DiscontinuousFunction(const unsigned int n_components = 1) + : Function(n_components) + {} + + virtual double + value(const Point &p, const unsigned int component = 0) const override + { + const double value = (p[0] * p[1] > 0 ? +3 : -1); // Non-trivial average + return value * (1.0 + std::sin(p[0]) * std::cos(p[1])); + } +}; + + +template +void +run(const ExtractorType &extractor) +{ + LogStream::Prefix prefix("Dim " + Utilities::to_string(dim)); + std::cout << "Dim: " << dim << std::endl; + + const FE_DGQ fe(3); + const QGauss qf_cell(fe.degree + 1); + const QGauss qf_face(fe.degree + 1); + + Triangulation triangulation; + GridGenerator::hyper_cube(triangulation, -1, 1); + triangulation.refine_global(1); + + DoFHandler dof_handler(triangulation); + dof_handler.distribute_dofs(fe); + + Vector solution(dof_handler.n_dofs()); + // Don't interpolate, as points precisely on the interfaces get evaluated and + // we don't end up with a jump in the values across it. + AffineConstraints constraints; + constraints.close(); + VectorTools::project(dof_handler, + constraints, + QGauss(fe.degree + 3), + DiscontinuousFunction(fe.n_components()), + solution); + + const UpdateFlags update_flags = + update_values | update_gradients | update_hessians | update_3rd_derivatives; + MeshWorker::ScratchData scratch_data( + fe, qf_cell, update_flags, qf_face, update_flags); + + const auto cell = dof_handler.begin_active(); + unsigned int face = 0; + while (cell->face(face)->at_boundary()) + ++face; + + scratch_data.reinit(cell, + face, + numbers::invalid_unsigned_int, + cell->neighbor(face), + cell->neighbor_face_no(face), + numbers::invalid_unsigned_int); + scratch_data.extract_local_dof_values("solution", solution); + + deallog << "Jumps in values: " + << scratch_data.get_jumps_in_values("solution", extractor)[0] + << std::endl; + deallog << "Jumps in gradients: " + << scratch_data.get_jumps_in_gradients("solution", extractor)[0] + << std::endl; + deallog << "Jumps in Hessians: " + << scratch_data.get_jumps_in_hessians("solution", extractor)[0] + << std::endl; + deallog << "Jumps in third derivatives: " + << scratch_data.get_jumps_in_third_derivatives("solution", + extractor)[0] + << std::endl; + + deallog << "Averages of values: " + << scratch_data.get_averages_of_values("solution", extractor)[0] + << std::endl; + deallog << "Averages of gradients: " + << scratch_data.get_averages_of_gradients("solution", extractor)[0] + << std::endl; + deallog << "Averages of Hessians: " + << scratch_data.get_averages_of_hessians("solution", extractor)[0] + << std::endl; + + deallog << "OK" << std::endl; +} + + +int +main(int argc, char *argv[]) +{ + initlog(); + Utilities::MPI::MPI_InitFinalize mpi_initialization( + argc, argv, testing_max_num_threads()); + + FEValuesExtractors::Scalar extractor(0); + + run<2>(extractor); + run<3>(extractor); + + deallog << "OK" << std::endl; +} diff --git a/tests/meshworker/scratch_data_10a.output b/tests/meshworker/scratch_data_10a.output new file mode 100644 index 0000000000..090dc0c2ce --- /dev/null +++ b/tests/meshworker/scratch_data_10a.output @@ -0,0 +1,18 @@ + +DEAL:Dim 2::Jumps in values: 4.00030 +DEAL:Dim 2::Jumps in gradients: 2.40085 0.000407768 +DEAL:Dim 2::Jumps in Hessians: 0.0456300 3.24213 3.24213 -0.000332829 +DEAL:Dim 2::Jumps in third derivatives: -2.06802 0.0616191 0.0616191 -2.64630 0.0616191 -2.64630 -2.64630 -0.000238992 +DEAL:Dim 2::Averages of values: 1.00030 +DEAL:Dim 2::Averages of gradients: 0.600214 0.000407768 +DEAL:Dim 2::Averages of Hessians: 0.0456300 0.810533 0.810533 -0.000332829 +DEAL:Dim 2::OK +DEAL:Dim 3::Jumps in values: 4.00030 +DEAL:Dim 3::Jumps in gradients: 2.40085 0.000407768 -1.31631e-14 +DEAL:Dim 3::Jumps in Hessians: 0.0456300 3.24213 2.93071e-13 3.24213 -0.000332829 -6.39488e-13 2.93071e-13 -6.39488e-13 -2.08722e-13 +DEAL:Dim 3::Jumps in third derivatives: -2.06802 0.0616191 -4.40314e-13 0.0616191 -2.64630 -7.67264e-12 -4.40314e-13 -7.67264e-12 -3.20455e-12 0.0616191 -2.64630 -7.67264e-12 -2.64630 -0.000238992 4.35918e-12 -7.67264e-12 4.35918e-12 4.07763e-12 -4.40314e-13 -7.67264e-12 -3.20455e-12 -7.67264e-12 4.35918e-12 4.07763e-12 -3.20455e-12 4.07763e-12 7.54952e-13 +DEAL:Dim 3::Averages of values: 1.00030 +DEAL:Dim 3::Averages of gradients: 0.600214 0.000407768 5.69336e-15 +DEAL:Dim 3::Averages of Hessians: 0.0456300 0.810533 3.11001e-14 0.810533 -0.000332829 -4.40536e-13 3.11001e-14 -4.40536e-13 -1.72751e-13 +DEAL:Dim 3::OK +DEAL::OK diff --git a/tests/meshworker/scratch_data_10b.cc b/tests/meshworker/scratch_data_10b.cc new file mode 100644 index 0000000000..be8664645a --- /dev/null +++ b/tests/meshworker/scratch_data_10b.cc @@ -0,0 +1,146 @@ +/* --------------------------------------------------------------------- + * + * Copyright (C) 2021 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 that ScratchData returns the correct jumps in, and averages of, +// solution values, gradients, etc. at interfaces +// - Vector valued FE + +#include +#include + +#include +#include +#include +#include + +#include +#include + +#include + +#include + +#include "../tests.h" + +template +class DiscontinuousFunction : public Function +{ +public: + DiscontinuousFunction(const unsigned int n_components = 1) + : Function(n_components) + {} + + virtual double + value(const Point &p, const unsigned int component = 0) const override + { + const double value = (p[0] * p[1] > 0 ? +3 : -1); // Non-trivial average + return value * (1.0 + std::sin(p[0]) * std::cos(p[1])); + } +}; + + +template +void +run(const ExtractorType &extractor) +{ + LogStream::Prefix prefix("Dim " + Utilities::to_string(dim)); + std::cout << "Dim: " << dim << std::endl; + + const FESystem fe(FE_DGQ(3), dim); + const QGauss qf_cell(fe.degree + 1); + const QGauss qf_face(fe.degree + 1); + + Triangulation triangulation; + GridGenerator::hyper_cube(triangulation, -1, 1); + triangulation.refine_global(1); + + DoFHandler dof_handler(triangulation); + dof_handler.distribute_dofs(fe); + + Vector solution(dof_handler.n_dofs()); + // Don't interpolate, as points precisely on the interfaces get evaluated and + // we don't end up with a jump in the values across it. + AffineConstraints constraints; + constraints.close(); + VectorTools::project(dof_handler, + constraints, + QGauss(fe.degree + 3), + DiscontinuousFunction(fe.n_components()), + solution); + + const UpdateFlags update_flags = + update_values | update_gradients | update_hessians | update_3rd_derivatives; + MeshWorker::ScratchData scratch_data( + fe, qf_cell, update_flags, qf_face, update_flags); + + const auto cell = dof_handler.begin_active(); + unsigned int face = 0; + while (cell->face(face)->at_boundary()) + ++face; + + scratch_data.reinit(cell, + face, + numbers::invalid_unsigned_int, + cell->neighbor(face), + cell->neighbor_face_no(face), + numbers::invalid_unsigned_int); + scratch_data.extract_local_dof_values("solution", solution); + + deallog << "Jumps in values: " + << scratch_data.get_jumps_in_values("solution", extractor)[0] + << std::endl; + deallog << "Jumps in gradients: " + << scratch_data.get_jumps_in_gradients("solution", extractor)[0] + << std::endl; + deallog << "Jumps in Hessians: " + << scratch_data.get_jumps_in_hessians("solution", extractor)[0] + << std::endl; + deallog << "Jumps in third derivatives: " + << scratch_data.get_jumps_in_third_derivatives("solution", + extractor)[0] + << std::endl; + + deallog << "Averages of values: " + << scratch_data.get_averages_of_values("solution", extractor)[0] + << std::endl; + deallog << "Averages of gradients: " + << scratch_data.get_averages_of_gradients("solution", extractor)[0] + << std::endl; + deallog << "Averages of Hessians: " + << scratch_data.get_averages_of_hessians("solution", extractor)[0] + << std::endl; + + deallog << "OK" << std::endl; +} + + +int +main(int argc, char *argv[]) +{ + initlog(); + Utilities::MPI::MPI_InitFinalize mpi_initialization( + argc, argv, testing_max_num_threads()); + + FEValuesExtractors::Vector extractor(0); + + run<2>(extractor); + run<3>(extractor); + + deallog << "OK" << std::endl; +} diff --git a/tests/meshworker/scratch_data_10b.output b/tests/meshworker/scratch_data_10b.output new file mode 100644 index 0000000000..5773077a61 --- /dev/null +++ b/tests/meshworker/scratch_data_10b.output @@ -0,0 +1,18 @@ + +DEAL:Dim 2::Jumps in values: 4.00030 4.00030 +DEAL:Dim 2::Jumps in gradients: 2.40085 0.000407768 2.40085 0.000407768 +DEAL:Dim 2::Jumps in Hessians: 0.0456300 3.24213 3.24213 -0.000332829 0.0456300 3.24213 3.24213 -0.000332829 +DEAL:Dim 2::Jumps in third derivatives: -2.06802 0.0616191 0.0616191 -2.64630 0.0616191 -2.64630 -2.64630 -0.000238992 -2.06802 0.0616191 0.0616191 -2.64630 0.0616191 -2.64630 -2.64630 -0.000238992 +DEAL:Dim 2::Averages of values: 1.00030 1.00030 +DEAL:Dim 2::Averages of gradients: 0.600214 0.000407768 0.600214 0.000407768 +DEAL:Dim 2::Averages of Hessians: 0.0456300 0.810533 0.810533 -0.000332829 0.0456300 0.810533 0.810533 -0.000332829 +DEAL:Dim 2::OK +DEAL:Dim 3::Jumps in values: 4.00030 4.00030 4.00030 +DEAL:Dim 3::Jumps in gradients: 2.40085 0.000407768 -9.48894e-15 2.40085 0.000407768 -9.48894e-15 2.40085 0.000407768 -9.48894e-15 +DEAL:Dim 3::Jumps in Hessians: 0.0456300 3.24213 2.80921e-13 3.24213 -0.000332829 -6.32244e-13 2.80921e-13 -6.32244e-13 -2.22489e-13 0.0456300 3.24213 2.80921e-13 3.24213 -0.000332829 -6.32244e-13 2.80921e-13 -6.32244e-13 -2.22489e-13 0.0456300 3.24213 2.80921e-13 3.24213 -0.000332829 -6.32244e-13 2.80921e-13 -6.32244e-13 -2.22489e-13 +DEAL:Dim 3::Jumps in third derivatives: -2.06802 0.0616191 -3.77309e-13 0.0616191 -2.64630 -7.53642e-12 -3.77309e-13 -7.53642e-12 -3.14782e-12 0.0616191 -2.64630 -7.53642e-12 -2.64630 -0.000238992 4.33609e-12 -7.53642e-12 4.33609e-12 4.07230e-12 -3.77309e-13 -7.53642e-12 -3.14782e-12 -7.53642e-12 4.33609e-12 4.07230e-12 -3.14782e-12 4.07230e-12 7.48956e-13 -2.06802 0.0616191 -3.77309e-13 0.0616191 -2.64630 -7.53642e-12 -3.77309e-13 -7.53642e-12 -3.14782e-12 0.0616191 -2.64630 -7.53642e-12 -2.64630 -0.000238992 4.33609e-12 -7.53642e-12 4.33609e-12 4.07230e-12 -3.77309e-13 -7.53642e-12 -3.14782e-12 -7.53642e-12 4.33609e-12 4.07230e-12 -3.14782e-12 4.07230e-12 7.48956e-13 -2.06802 0.0616191 -3.77309e-13 0.0616191 -2.64630 -7.53642e-12 -3.77309e-13 -7.53642e-12 -3.14782e-12 0.0616191 -2.64630 -7.53642e-12 -2.64630 -0.000238992 4.33609e-12 -7.53642e-12 4.33609e-12 4.07230e-12 -3.77309e-13 -7.53642e-12 -3.14782e-12 -7.53642e-12 4.33609e-12 4.07230e-12 -3.14782e-12 4.07230e-12 7.48956e-13 +DEAL:Dim 3::Averages of values: 1.00030 1.00030 1.00030 +DEAL:Dim 3::Averages of gradients: 0.600214 0.000407768 4.74447e-15 0.600214 0.000407768 4.74447e-15 0.600214 0.000407768 4.74447e-15 +DEAL:Dim 3::Averages of Hessians: 0.0456300 0.810533 3.71751e-14 0.810533 -0.000332829 -4.44159e-13 3.71751e-14 -4.44159e-13 -1.74749e-13 0.0456300 0.810533 3.71751e-14 0.810533 -0.000332829 -4.44159e-13 3.71751e-14 -4.44159e-13 -1.74749e-13 0.0456300 0.810533 3.71751e-14 0.810533 -0.000332829 -4.44159e-13 3.71751e-14 -4.44159e-13 -1.74749e-13 +DEAL:Dim 3::OK +DEAL::OK -- 2.39.5