From: Jean-Paul Pelteret Date: Sat, 13 Mar 2021 18:22:31 +0000 (+0100) Subject: Add tests for MeshWorker::ScratchData::get_* X-Git-Tag: v9.3.0-rc1~326^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=bd81f3c9a824c23a80a5e2eb7514ba256d277878;p=dealii.git Add tests for MeshWorker::ScratchData::get_* --- diff --git a/tests/meshworker/scratch_data_09a.cc b/tests/meshworker/scratch_data_09a.cc new file mode 100644 index 0000000000..7a2825f5f7 --- /dev/null +++ b/tests/meshworker/scratch_data_09a.cc @@ -0,0 +1,102 @@ +/* --------------------------------------------------------------------- + * + * 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 solution values, gradients, etc. +// - Scalar valued FE + +#include +#include + +#include +#include +#include + +#include +#include + +#include + +#include + +#include "../tests.h" + + +template +void +run(const ExtractorType &extractor) +{ + LogStream::Prefix prefix("Dim " + Utilities::to_string(dim)); + std::cout << "Dim: " << dim << std::endl; + + const FE_Q fe(3); + const QGauss qf_cell(fe.degree + 1); + const QGauss qf_face(fe.degree + 1); + + Triangulation triangulation; + GridGenerator::hyper_cube(triangulation); + + DoFHandler dof_handler(triangulation); + dof_handler.distribute_dofs(fe); + + Vector solution(dof_handler.n_dofs()); + VectorTools::interpolate(dof_handler, + Functions::CosineFunction( + 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); + + const auto cell = dof_handler.begin_active(); + scratch_data.reinit(cell); + scratch_data.extract_local_dof_values("solution", solution); + + deallog << "Value: " << scratch_data.get_values("solution", extractor)[0] + << std::endl; + deallog << "Gradient: " + << scratch_data.get_gradients("solution", extractor)[0] << std::endl; + deallog << "Hessian: " << scratch_data.get_hessians("solution", extractor)[0] + << std::endl; + deallog << "Laplacian: " + << scratch_data.get_laplacians("solution", extractor)[0] << std::endl; + deallog << "Third derivative: " + << scratch_data.get_third_derivatives("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_09a.output b/tests/meshworker/scratch_data_09a.output new file mode 100644 index 0000000000..584d7eb272 --- /dev/null +++ b/tests/meshworker/scratch_data_09a.output @@ -0,0 +1,14 @@ + +DEAL:Dim 2::Value: 0.991562 +DEAL:Dim 2::Gradient: -0.158724 -0.158724 +DEAL:Dim 2::Hessian: -2.76461 0.0254076 0.0254076 -2.76461 +DEAL:Dim 2::Laplacian: -5.52923 +DEAL:Dim 2::Third derivative: 2.62952 0.442544 0.442544 0.442544 0.442544 0.442544 0.442544 2.62952 +DEAL:Dim 2::OK +DEAL:Dim 3::Value: 0.987370 +DEAL:Dim 3::Gradient: -0.158053 -0.158053 -0.158053 +DEAL:Dim 3::Hessian: -2.75292 0.0253002 0.0253002 0.0253002 -2.75292 0.0253002 0.0253002 0.0253002 -2.75292 +DEAL:Dim 3::Laplacian: -8.25877 +DEAL:Dim 3::Third derivative: 2.61841 0.440673 0.440673 0.440673 0.440673 -0.00404991 0.440673 -0.00404991 0.440673 0.440673 0.440673 -0.00404991 0.440673 2.61841 0.440673 -0.00404991 0.440673 0.440673 0.440673 -0.00404991 0.440673 -0.00404991 0.440673 0.440673 0.440673 0.440673 2.61841 +DEAL:Dim 3::OK +DEAL::OK diff --git a/tests/meshworker/scratch_data_09b.cc b/tests/meshworker/scratch_data_09b.cc new file mode 100644 index 0000000000..2315004e19 --- /dev/null +++ b/tests/meshworker/scratch_data_09b.cc @@ -0,0 +1,111 @@ +/* --------------------------------------------------------------------- + * + * 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 solution values, gradients, etc. +// - Vector valued FE + +#include +#include + +#include +#include +#include +#include + +#include +#include + +#include + +#include + +#include "../tests.h" + + +template +void +run(const ExtractorType &extractor) +{ + LogStream::Prefix prefix("Dim " + Utilities::to_string(dim)); + std::cout << "Dim: " << dim << std::endl; + + const FESystem fe(FE_Q(3), dim); + const QGauss qf_cell(fe.degree + 1); + const QGauss qf_face(fe.degree + 1); + + Triangulation triangulation; + GridGenerator::hyper_cube(triangulation); + + DoFHandler dof_handler(triangulation); + dof_handler.distribute_dofs(fe); + + Vector solution(dof_handler.n_dofs()); + VectorTools::interpolate(dof_handler, + Functions::CosineFunction( + 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); + + const auto cell = dof_handler.begin_active(); + scratch_data.reinit(cell); + scratch_data.extract_local_dof_values("solution", solution); + + deallog << "Value: " << scratch_data.get_values("solution", extractor)[0] + << std::endl; + deallog << "Gradient: " + << scratch_data.get_gradients("solution", extractor)[0] << std::endl; + deallog << "Symmetric gradient: " + << scratch_data.get_symmetric_gradients("solution", extractor)[0] + << std::endl; + deallog << "Divergence: " + << scratch_data.get_divergences("solution", extractor)[0] + << std::endl; + deallog << "Curl: " << scratch_data.get_curls("solution", extractor)[0] + << std::endl; + deallog << "Hessian: " << scratch_data.get_hessians("solution", extractor)[0] + << std::endl; + deallog << "Laplacian: " + << scratch_data.get_laplacians("solution", extractor)[0] << std::endl; + deallog << "Third derivative: " + << scratch_data.get_third_derivatives("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); + + // Just run in 3d, so we can extract the curl + run<3>(extractor); + + deallog << "OK" << std::endl; +} diff --git a/tests/meshworker/scratch_data_09b.output b/tests/meshworker/scratch_data_09b.output new file mode 100644 index 0000000000..e8d10d04de --- /dev/null +++ b/tests/meshworker/scratch_data_09b.output @@ -0,0 +1,11 @@ + +DEAL:Dim 3::Value: 0.987370 0.987370 0.987370 +DEAL:Dim 3::Gradient: -0.158053 -0.158053 -0.158053 -0.158053 -0.158053 -0.158053 -0.158053 -0.158053 -0.158053 +DEAL:Dim 3::Symmetric gradient: -0.158053 -0.158053 -0.158053 -0.158053 -0.158053 -0.158053 -0.158053 -0.158053 -0.158053 +DEAL:Dim 3::Divergence: -0.474158 +DEAL:Dim 3::Curl: 1.11890e-16 7.88432e-16 6.54858e-16 +DEAL:Dim 3::Hessian: -2.75292 0.0253002 0.0253002 0.0253002 -2.75292 0.0253002 0.0253002 0.0253002 -2.75292 -2.75292 0.0253002 0.0253002 0.0253002 -2.75292 0.0253002 0.0253002 0.0253002 -2.75292 -2.75292 0.0253002 0.0253002 0.0253002 -2.75292 0.0253002 0.0253002 0.0253002 -2.75292 +DEAL:Dim 3::Laplacian: -8.25877 -8.25877 -8.25877 +DEAL:Dim 3::Third derivative: 2.61841 0.440673 0.440673 0.440673 0.440673 -0.00404991 0.440673 -0.00404991 0.440673 0.440673 0.440673 -0.00404991 0.440673 2.61841 0.440673 -0.00404991 0.440673 0.440673 0.440673 -0.00404991 0.440673 -0.00404991 0.440673 0.440673 0.440673 0.440673 2.61841 2.61841 0.440673 0.440673 0.440673 0.440673 -0.00404991 0.440673 -0.00404991 0.440673 0.440673 0.440673 -0.00404991 0.440673 2.61841 0.440673 -0.00404991 0.440673 0.440673 0.440673 -0.00404991 0.440673 -0.00404991 0.440673 0.440673 0.440673 0.440673 2.61841 2.61841 0.440673 0.440673 0.440673 0.440673 -0.00404991 0.440673 -0.00404991 0.440673 0.440673 0.440673 -0.00404991 0.440673 2.61841 0.440673 -0.00404991 0.440673 0.440673 0.440673 -0.00404991 0.440673 -0.00404991 0.440673 0.440673 0.440673 0.440673 2.61841 +DEAL:Dim 3::OK +DEAL::OK diff --git a/tests/meshworker/scratch_data_09c.cc b/tests/meshworker/scratch_data_09c.cc new file mode 100644 index 0000000000..8f10ea7732 --- /dev/null +++ b/tests/meshworker/scratch_data_09c.cc @@ -0,0 +1,100 @@ +/* --------------------------------------------------------------------- + * + * 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 solution values, gradients, etc. +// - Tensor valued FE + +#include +#include +#include + +#include +#include +#include +#include + +#include +#include + +#include + +#include + +#include "../tests.h" + + +template +void +run(const ExtractorType &extractor) +{ + LogStream::Prefix prefix("Dim " + Utilities::to_string(dim)); + std::cout << "Dim: " << dim << std::endl; + + const FESystem fe(FE_Q(3), + Tensor<2, dim>::n_independent_components); + const QGauss qf_cell(fe.degree + 1); + const QGauss qf_face(fe.degree + 1); + + Triangulation triangulation; + GridGenerator::hyper_cube(triangulation); + + DoFHandler dof_handler(triangulation); + dof_handler.distribute_dofs(fe); + + Vector solution(dof_handler.n_dofs()); + VectorTools::interpolate(dof_handler, + Functions::CosineFunction( + fe.n_components()), + solution); + + const UpdateFlags update_flags = update_values | update_gradients; + MeshWorker::ScratchData scratch_data(fe, + qf_cell, + update_flags); + + const auto cell = dof_handler.begin_active(); + scratch_data.reinit(cell); + scratch_data.extract_local_dof_values("solution", solution); + + deallog << "Value: " << scratch_data.get_values("solution", extractor)[0] + << std::endl; + deallog << "Gradient: " + << scratch_data.get_gradients("solution", extractor)[0] << std::endl; + deallog << "Divergence: " + << scratch_data.get_divergences("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::Tensor<2> extractor(0); + + run<2>(extractor); + run<3>(extractor); + + deallog << "OK" << std::endl; +} diff --git a/tests/meshworker/scratch_data_09c.output b/tests/meshworker/scratch_data_09c.output new file mode 100644 index 0000000000..7d8e560039 --- /dev/null +++ b/tests/meshworker/scratch_data_09c.output @@ -0,0 +1,10 @@ + +DEAL:Dim 2::Value: 0.991562 0.991562 0.991562 0.991562 +DEAL:Dim 2::Gradient: -0.158724 -0.158724 -0.158724 -0.158724 -0.158724 -0.158724 -0.158724 -0.158724 +DEAL:Dim 2::Divergence: -0.317447 -0.317447 +DEAL:Dim 2::OK +DEAL:Dim 3::Value: 0.987370 0.987370 0.987370 0.987370 0.987370 0.987370 0.987370 0.987370 0.987370 +DEAL:Dim 3::Gradient: -0.158053 -0.158053 -0.158053 -0.158053 -0.158053 -0.158053 -0.158053 -0.158053 -0.158053 -0.158053 -0.158053 -0.158053 -0.158053 -0.158053 -0.158053 -0.158053 -0.158053 -0.158053 -0.158053 -0.158053 -0.158053 -0.158053 -0.158053 -0.158053 -0.158053 -0.158053 -0.158053 +DEAL:Dim 3::Divergence: -0.474158 -0.474158 -0.474158 +DEAL:Dim 3::OK +DEAL::OK diff --git a/tests/meshworker/scratch_data_09d.cc b/tests/meshworker/scratch_data_09d.cc new file mode 100644 index 0000000000..bbc6b59b4a --- /dev/null +++ b/tests/meshworker/scratch_data_09d.cc @@ -0,0 +1,98 @@ +/* --------------------------------------------------------------------- + * + * 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 solution values, gradients, etc. +// - Symmetric tensor valued FE + +#include +#include +#include + +#include +#include +#include +#include + +#include +#include + +#include + +#include + +#include "../tests.h" + + +template +void +run(const ExtractorType &extractor) +{ + LogStream::Prefix prefix("Dim " + Utilities::to_string(dim)); + std::cout << "Dim: " << dim << std::endl; + + const FESystem fe( + FE_Q(3), SymmetricTensor<2, dim>::n_independent_components); + const QGauss qf_cell(fe.degree + 1); + const QGauss qf_face(fe.degree + 1); + + Triangulation triangulation; + GridGenerator::hyper_cube(triangulation); + + DoFHandler dof_handler(triangulation); + dof_handler.distribute_dofs(fe); + + Vector solution(dof_handler.n_dofs()); + VectorTools::interpolate(dof_handler, + Functions::CosineFunction( + fe.n_components()), + solution); + + const UpdateFlags update_flags = update_values | update_gradients; + MeshWorker::ScratchData scratch_data(fe, + qf_cell, + update_flags); + + const auto cell = dof_handler.begin_active(); + scratch_data.reinit(cell); + scratch_data.extract_local_dof_values("solution", solution); + + deallog << "Value: " << scratch_data.get_values("solution", extractor)[0] + << std::endl; + deallog << "Divergence: " + << scratch_data.get_divergences("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::SymmetricTensor<2> extractor(0); + + run<2>(extractor); + run<3>(extractor); + + deallog << "OK" << std::endl; +} diff --git a/tests/meshworker/scratch_data_09d.output b/tests/meshworker/scratch_data_09d.output new file mode 100644 index 0000000000..c9b8c573ac --- /dev/null +++ b/tests/meshworker/scratch_data_09d.output @@ -0,0 +1,8 @@ + +DEAL:Dim 2::Value: 0.991562 0.991562 0.991562 0.991562 +DEAL:Dim 2::Divergence: -0.317447 -0.317447 +DEAL:Dim 2::OK +DEAL:Dim 3::Value: 0.987370 0.987370 0.987370 0.987370 0.987370 0.987370 0.987370 0.987370 0.987370 +DEAL:Dim 3::Divergence: -0.474158 -0.474158 -0.474158 +DEAL:Dim 3::OK +DEAL::OK