From c0382ccaebfb6c333c3cca9d2caccae69267102a Mon Sep 17 00:00:00 2001 From: Jean-Paul Pelteret Date: Fri, 23 Dec 2022 17:53:01 +0100 Subject: [PATCH] Add tests --- tests/meshworker/scratch_data_10a_hp.cc | 152 +++++++++++++++++++ tests/meshworker/scratch_data_10a_hp.output | 18 +++ tests/meshworker/scratch_data_10b_hp.cc | 154 ++++++++++++++++++++ tests/meshworker/scratch_data_10b_hp.output | 18 +++ 4 files changed, 342 insertions(+) create mode 100644 tests/meshworker/scratch_data_10a_hp.cc create mode 100644 tests/meshworker/scratch_data_10a_hp.output create mode 100644 tests/meshworker/scratch_data_10b_hp.cc create mode 100644 tests/meshworker/scratch_data_10b_hp.output diff --git a/tests/meshworker/scratch_data_10a_hp.cc b/tests/meshworker/scratch_data_10a_hp.cc new file mode 100644 index 0000000000..6fc210c7f8 --- /dev/null +++ b/tests/meshworker/scratch_data_10a_hp.cc @@ -0,0 +1,152 @@ +/* --------------------------------------------------------------------- + * + * 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 +// - hp variant + +#include +#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; + + hp::FECollection fe(FE_DGQ(3), + FE_DGQ(4)); + hp::QCollection qf_cell(QGauss(fe.max_degree()), + QGauss(fe.max_degree() + 1)); + hp::QCollection qf_face(QGauss(fe.max_degree()), + QGauss(fe.max_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, + hp::QCollection(QGauss(fe.max_degree() + 2), + QGauss(fe.max_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, + cell->neighbor(face), + cell->neighbor_face_no(face)); + 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_hp.output b/tests/meshworker/scratch_data_10a_hp.output new file mode 100644 index 0000000000..090dc0c2ce --- /dev/null +++ b/tests/meshworker/scratch_data_10a_hp.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_hp.cc b/tests/meshworker/scratch_data_10b_hp.cc new file mode 100644 index 0000000000..194a8ed907 --- /dev/null +++ b/tests/meshworker/scratch_data_10b_hp.cc @@ -0,0 +1,154 @@ +/* --------------------------------------------------------------------- + * + * 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 +// - hp variant + +#include +#include + +#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; + + hp::FECollection fe( + FESystem(FE_DGQ(3), dim), + FESystem(FE_DGQ(4), dim)); + hp::QCollection qf_cell(QGauss(fe.max_degree()), + QGauss(fe.max_degree() + 1)); + hp::QCollection qf_face(QGauss(fe.max_degree()), + QGauss(fe.max_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, + hp::QCollection(QGauss(fe.max_degree() + 2), + QGauss(fe.max_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, + cell->neighbor(face), + cell->neighbor_face_no(face)); + 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_hp.output b/tests/meshworker/scratch_data_10b_hp.output new file mode 100644 index 0000000000..5773077a61 --- /dev/null +++ b/tests/meshworker/scratch_data_10b_hp.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