From e73354e26b14dd8f222abc0666327a8e140a51f3 Mon Sep 17 00:00:00 2001 From: Peter Munch Date: Sun, 25 Oct 2020 19:05:00 +0100 Subject: [PATCH] Pass n_q_points_1D to Simplex::QGauss --- include/deal.II/simplex/quadrature_lib.h | 14 +++++++++----- source/simplex/quadrature_lib.cc | 16 ++++++++-------- tests/simplex/fe_lib_02.cc | 16 ++++++++-------- tests/simplex/poisson_01.cc | 6 ++---- tests/simplex/poisson_02.cc | 6 ++---- tests/simplex/q_projection_01.cc | 12 ++++++------ tests/simplex/quadrature_lib_01.cc | 10 +++++----- tests/simplex/quadrature_lib_01.output | 5 +++-- tests/simplex/step-12.cc | 6 ++---- tests/simplex/step-18.cc | 2 +- 10 files changed, 46 insertions(+), 47 deletions(-) diff --git a/include/deal.II/simplex/quadrature_lib.h b/include/deal.II/simplex/quadrature_lib.h index 41c02dc02f..09e2e97abd 100644 --- a/include/deal.II/simplex/quadrature_lib.h +++ b/include/deal.II/simplex/quadrature_lib.h @@ -28,12 +28,15 @@ namespace Simplex /** * Integration rule for simplex entities. * - * Following number of quadrature points are currently supported for 2D and - * 3D: + * Users specify a number `n_points_1D` as an indication of what polynomial + * degree to be integrated exactly, similarly to the number of points in a + * QGauss quadrature object, even though the present quadrature formula is not + * a tensor product. The given value is translated for n_points_1D=1,2,3 to + * following number of quadrature points for 2D and 3D: * - 2D: 1, 3, 7 * - 3D: 1, 4, 10 * - * For 1D, the quadrature rule degenerates to a `QGauss<1>(n_points)`. + * For 1D, the quadrature rule degenerates to a `QGauss<1>(n_points_1D)`. * * @ingroup simplex */ @@ -42,9 +45,10 @@ namespace Simplex { public: /** - * Constructor taking the number of quadrature points @p n_points. + * Constructor taking the number of quadrature points in 1D direction + * @p n_points_1D. */ - explicit QGauss(const unsigned int n_points); + explicit QGauss(const unsigned int n_points_1D); }; } // namespace Simplex diff --git a/source/simplex/quadrature_lib.cc b/source/simplex/quadrature_lib.cc index 683c9ab7cd..bb60716e6c 100644 --- a/source/simplex/quadrature_lib.cc +++ b/source/simplex/quadrature_lib.cc @@ -29,26 +29,26 @@ DEAL_II_NAMESPACE_OPEN namespace Simplex { template - QGauss::QGauss(const unsigned int n_points) + QGauss::QGauss(const unsigned int n_points_1D) : QSimplex(Quadrature()) { // fill quadrature points and quadrature weights if (dim == 1) { - const dealii::QGauss quad(n_points); + const dealii::QGauss quad(n_points_1D); this->quadrature_points = quad.get_points(); this->weights = quad.get_weights(); } else if (dim == 2) { - if (n_points == 1) + if (n_points_1D == 1) { const double p = 1.0 / 3.0; this->quadrature_points.emplace_back(p, p); this->weights.emplace_back(1.0); } - else if (n_points == 3) + else if (n_points_1D == 2) { const double Q23 = 2.0 / 3.0; const double Q16 = 1.0 / 6.0; @@ -60,7 +60,7 @@ namespace Simplex this->weights.emplace_back(Q16); this->weights.emplace_back(Q16); } - else if (n_points == 7) + else if (n_points_1D == 3) { const double q12 = 0.5; @@ -85,7 +85,7 @@ namespace Simplex } else if (dim == 3) { - if (n_points == 1) + if (n_points_1D == 1) { const double Q14 = 1.0 / 4.0; const double Q16 = 1.0 / 6.0; @@ -93,7 +93,7 @@ namespace Simplex this->quadrature_points.emplace_back(Q14, Q14, Q14); this->weights.emplace_back(Q16); } - else if (n_points == 4) + else if (n_points_1D == 2) { const double Q124 = 1.0 / 6.0 / 4.0; @@ -108,7 +108,7 @@ namespace Simplex this->weights.emplace_back(Q124); this->weights.emplace_back(Q124); } - else if (n_points == 10) + else if (n_points_1D == 3) { const double Q16 = 1.0 / 6.0; diff --git a/tests/simplex/fe_lib_02.cc b/tests/simplex/fe_lib_02.cc index 8e4c141a94..feee438b32 100644 --- a/tests/simplex/fe_lib_02.cc +++ b/tests/simplex/fe_lib_02.cc @@ -49,12 +49,12 @@ main() { initlog(); - test(Simplex::FE_P<2>(1), Simplex::QGauss<2>(3)); - test(Simplex::FE_P<2>(2), Simplex::QGauss<2>(7)); - test(Simplex::FE_P<3>(1), Simplex::QGauss<3>(4)); - test(Simplex::FE_P<3>(2), Simplex::QGauss<3>(10)); - test(Simplex::FE_DGP<2>(1), Simplex::QGauss<2>(3)); - test(Simplex::FE_DGP<2>(2), Simplex::QGauss<2>(7)); - test(Simplex::FE_DGP<3>(1), Simplex::QGauss<3>(4)); - test(Simplex::FE_DGP<3>(2), Simplex::QGauss<3>(10)); + test(Simplex::FE_P<2>(1), Simplex::QGauss<2>(2)); + test(Simplex::FE_P<2>(2), Simplex::QGauss<2>(3)); + test(Simplex::FE_P<3>(1), Simplex::QGauss<3>(2)); + test(Simplex::FE_P<3>(2), Simplex::QGauss<3>(3)); + test(Simplex::FE_DGP<2>(1), Simplex::QGauss<2>(2)); + test(Simplex::FE_DGP<2>(2), Simplex::QGauss<2>(3)); + test(Simplex::FE_DGP<3>(1), Simplex::QGauss<3>(2)); + test(Simplex::FE_DGP<3>(2), Simplex::QGauss<3>(3)); } diff --git a/tests/simplex/poisson_01.cc b/tests/simplex/poisson_01.cc index e920c6e539..e020592b7d 100644 --- a/tests/simplex/poisson_01.cc +++ b/tests/simplex/poisson_01.cc @@ -336,11 +336,9 @@ test_tet(const MPI_Comm &comm, const Parameters ¶ms) // 3) Select components Simplex::FE_P fe(params.degree); - Simplex::QGauss quad(dim == 2 ? (params.degree == 1 ? 3 : 7) : - (params.degree == 1 ? 4 : 10)); + Simplex::QGauss quad(params.degree + 1); - Simplex::QGauss face_quad(dim == 2 ? (params.degree == 1 ? 2 : 3) : - (params.degree == 1 ? 3 : 7)); + Simplex::QGauss face_quad(params.degree + 1); Simplex::FE_P fe_mapping(1); MappingFE mapping(fe_mapping); diff --git a/tests/simplex/poisson_02.cc b/tests/simplex/poisson_02.cc index 494dee8715..c117671b3d 100644 --- a/tests/simplex/poisson_02.cc +++ b/tests/simplex/poisson_02.cc @@ -189,10 +189,8 @@ public: false, new Simplex::FE_DGP(degree), new MappingFE(Simplex::FE_P(1)), - new Simplex::QGauss(dim == 2 ? (degree == 1 ? 3 : 7) : - (degree == 1 ? 4 : 10)), - new Simplex::QGauss(dim == 2 ? (degree == 1 ? 2 : 3) : - (degree == 1 ? 3 : 7)), + new Simplex::QGauss(degree + 1), + new Simplex::QGauss(degree + 1), initial_refinement, number_refinement); } diff --git a/tests/simplex/q_projection_01.cc b/tests/simplex/q_projection_01.cc index 324cdc6080..7639b3a0c7 100644 --- a/tests/simplex/q_projection_01.cc +++ b/tests/simplex/q_projection_01.cc @@ -56,8 +56,8 @@ test<2>(const unsigned int n_points) face_orientation, false, false, - n_points); - q < n_points; + quad_ref.size()); + q < quad_ref.size(); ++q, ++i) { deallog << quad.point(i) << " "; @@ -102,8 +102,8 @@ test<3>(const unsigned int n_points) face_orientation, face_flip, face_rotation, - n_points); - q < n_points; + quad_ref.size()); + q < quad_ref.size(); ++q, ++i) { deallog << quad.point(i) << " "; @@ -152,12 +152,12 @@ main() } { deallog.push("3d-4"); - test<3>(3); + test<3>(2); deallog.pop(); } { deallog.push("3d-10"); - test<3>(7); + test<3>(3); deallog.pop(); } } diff --git a/tests/simplex/quadrature_lib_01.cc b/tests/simplex/quadrature_lib_01.cc index 85d6b4a0aa..fa05a42fd0 100644 --- a/tests/simplex/quadrature_lib_01.cc +++ b/tests/simplex/quadrature_lib_01.cc @@ -54,7 +54,7 @@ main() } { deallog.push("1d-3"); - test<1>(2); + test<1>(3); deallog.pop(); } { @@ -64,12 +64,12 @@ main() } { deallog.push("2d-3"); - test<2>(3); + test<2>(2); deallog.pop(); } { deallog.push("2d-7"); - test<2>(7); + test<2>(3); deallog.pop(); } @@ -80,12 +80,12 @@ main() } { deallog.push("3d-4"); - test<3>(4); + test<3>(2); deallog.pop(); } { deallog.push("3d-10"); - test<3>(10); + test<3>(3); deallog.pop(); } } diff --git a/tests/simplex/quadrature_lib_01.output b/tests/simplex/quadrature_lib_01.output index 24a8420bbd..94888c477f 100644 --- a/tests/simplex/quadrature_lib_01.output +++ b/tests/simplex/quadrature_lib_01.output @@ -2,8 +2,9 @@ DEAL:1d-1::0.500000 1.00000 DEAL:1d-2::0.211325 0.500000 DEAL:1d-2::0.788675 0.500000 -DEAL:1d-3::0.211325 0.500000 -DEAL:1d-3::0.788675 0.500000 +DEAL:1d-3::0.112702 0.277778 +DEAL:1d-3::0.500000 0.444444 +DEAL:1d-3::0.887298 0.277778 DEAL:2d-1::0.333333 0.333333 1.00000 DEAL:2d-3::0.666667 0.166667 0.166667 DEAL:2d-3::0.166667 0.666667 0.166667 diff --git a/tests/simplex/step-12.cc b/tests/simplex/step-12.cc index 4c7fe7e2c1..29484a3aa2 100644 --- a/tests/simplex/step-12.cc +++ b/tests/simplex/step-12.cc @@ -391,11 +391,9 @@ namespace Step12 QGauss face_quad(degree + 1); #else - Simplex::QGauss quad(dim == 2 ? (degree == 1 ? 3 : 7) : - (degree == 1 ? 4 : 10)); + Simplex::QGauss quad(degree + 1); - Simplex::QGauss face_quad(dim == 2 ? (degree == 1 ? 2 : 3) : - (degree == 1 ? 3 : 7)); + Simplex::QGauss face_quad(degree + 1); #endif ScratchData scratch_data(mapping, fe, quad, face_quad); diff --git a/tests/simplex/step-18.cc b/tests/simplex/step-18.cc index 2bac10c4d1..3e2196005c 100644 --- a/tests/simplex/step-18.cc +++ b/tests/simplex/step-18.cc @@ -419,7 +419,7 @@ namespace Step18 : triangulation() , fe(Simplex::FE_P(degree), dim) , dof_handler(triangulation) - , quadrature_formula(Simplex::QGauss(fe.degree == 1 ? 4 : 10)) + , quadrature_formula(Simplex::QGauss(fe.degree + 1)) , mapping(Simplex::FE_P(1)) , present_time(0.0) , present_timestep(1.0) -- 2.39.5