From: Peter Munch Date: Mon, 14 Dec 2020 09:21:40 +0000 (+0100) Subject: MatrixFree: support FE_PyramidP X-Git-Tag: v9.3.0-rc1~757^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=f4a98d4b1d82fd3ff3c06efa0d253e8b5c91517e;p=dealii.git MatrixFree: support FE_PyramidP --- diff --git a/include/deal.II/matrix_free/shape_info.templates.h b/include/deal.II/matrix_free/shape_info.templates.h index 403e25aceb..9a0bcf9b71 100644 --- a/include/deal.II/matrix_free/shape_info.templates.h +++ b/include/deal.II/matrix_free/shape_info.templates.h @@ -93,7 +93,8 @@ namespace internal if (dynamic_cast *>(&fe) != nullptr || dynamic_cast *>(&fe) != nullptr || - dynamic_cast *>(&fe) != nullptr) + dynamic_cast *>(&fe) != nullptr || + dynamic_cast *>(&fe) != nullptr) { scalar_lexicographic.resize(fe.n_dofs_per_cell()); for (unsigned int i = 0; i < scalar_lexicographic.size(); ++i) @@ -206,7 +207,8 @@ namespace internal // indicative of their support if (dynamic_cast *>(fe_poly_ptr) || dynamic_cast *>(fe_poly_ptr) || - dynamic_cast *>(fe_poly_ptr)) + dynamic_cast *>(fe_poly_ptr) || + dynamic_cast *>(fe_poly_ptr)) return true; #endif @@ -247,6 +249,8 @@ namespace internal dynamic_cast *>( &fe_in.base_element(base_element_number)) || dynamic_cast *>( + &fe_in.base_element(base_element_number)) || + dynamic_cast *>( &fe_in.base_element(base_element_number))) { // specialization for arbitrary finite elements and quadrature rules diff --git a/tests/simplex/matrix_free_01.cc b/tests/simplex/matrix_free_01.cc index 8f0df19794..daeeb53115 100644 --- a/tests/simplex/matrix_free_01.cc +++ b/tests/simplex/matrix_free_01.cc @@ -148,6 +148,14 @@ test(const unsigned int v, const unsigned int degree, const bool do_helmholtz) quad = std::make_shared>(degree + 1); fe_mapping = std::make_shared>(1); } + else if (v == 2) + { + GridGenerator::subdivided_hyper_cube_with_pyramids(tria, + dim == 2 ? 16 : 8); + fe = std::make_shared>(degree); + quad = std::make_shared>(degree + 1); + fe_mapping = std::make_shared>(1); + } else Assert(false, ExcNotImplemented()); @@ -304,16 +312,18 @@ main(int argc, char **argv) Utilities::MPI::MPI_InitFinalize mpi(argc, argv, 1); - for (unsigned int i = 0; i < 2; ++i) + for (unsigned int i = 0; i <= 2; ++i) { if (i == 0) deallog.push("SIMPLEX"); else if (i == 1) deallog.push("WEDGE "); + else if (i == 2) + deallog.push("PYRAMID"); else Assert(false, ExcNotImplemented()); - if (i == 0) + if (i == 0) // 2D makes only sense for simplex { test<2>(i, /*degree=*/1, /*do_helmholtz*/ false); test<2>(i, /*degree=*/1, /*do_helmholtz*/ true); @@ -323,8 +333,13 @@ main(int argc, char **argv) test<3>(i, /*degree=*/1, /*do_helmholtz*/ false); test<3>(i, /*degree=*/1, /*do_helmholtz*/ true); - test<3>(i, /*degree=*/2, /*do_helmholtz*/ false); - test<3>(i, /*degree=*/2, /*do_helmholtz*/ true); + + if (i != + 2) // for pyramids no quadratic elements have been implemented yet + { + test<3>(i, /*degree=*/2, /*do_helmholtz*/ false); + test<3>(i, /*degree=*/2, /*do_helmholtz*/ true); + } deallog.pop(); } diff --git a/tests/simplex/matrix_free_01.with_simplex_support=on.output b/tests/simplex/matrix_free_01.with_simplex_support=on.output index 5bcc1f1895..83ee05c055 100644 --- a/tests/simplex/matrix_free_01.with_simplex_support=on.output +++ b/tests/simplex/matrix_free_01.with_simplex_support=on.output @@ -11,3 +11,5 @@ DEAL:WEDGE ::dim=3 degree=1 Type=Possion : Convergence step 7 value 0.0002242 DEAL:WEDGE ::dim=3 degree=1 Type=Helmholtz : Convergence step 7 value 0.000211401. DEAL:WEDGE ::dim=3 degree=2 Type=Possion : Convergence step 20 value 0.000133408. DEAL:WEDGE ::dim=3 degree=2 Type=Helmholtz : Convergence step 20 value 0.000125844. +DEAL:PYRAMID::dim=3 degree=1 Type=Possion : Convergence step 10 value 0.000197898. +DEAL:PYRAMID::dim=3 degree=1 Type=Helmholtz : Convergence step 10 value 0.000186534. diff --git a/tests/simplex/simplex_grids.h b/tests/simplex/simplex_grids.h index 7c91600310..8d0e42dc1a 100644 --- a/tests/simplex/simplex_grids.h +++ b/tests/simplex/simplex_grids.h @@ -250,6 +250,31 @@ namespace dealii + template + void + subdivided_hyper_cube_with_pyramids(Triangulation &tria, + const unsigned int repetitions, + const double p1 = 0.0, + const double p2 = 1.0, + const bool colorize = false) + { + if (dim == 3) + { + subdivided_hyper_rectangle_with_pyramids( + tria, + {{repetitions, repetitions, repetitions}}, + {p1, p1, p1}, + {p2, p2, p2}, + colorize); + } + else + { + AssertThrow(false, ExcNotImplemented()) + } + } + + + template void subdivided_hyper_rectangle_with_simplices_mix(