Replace Simplex::QGauss(4) with a Witherden-Vincent rule.

diff --git a/include/deal.II/simplex/quadrature_lib.h b/include/deal.II/simplex/quadrature_lib.h
index 93f351f54a2d1960390b4618d394b8448950fa1c..ba554888716961a770f312c34d651546ab530cfa 100644 (file)
--- a/include/deal.II/simplex/quadrature_lib.h
+++ b/include/deal.II/simplex/quadrature_lib.h
@@ -31,10 +31,10 @@ namespace Simplex
    * 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
+   * a tensor product. The given value is translated for n_points_1D=1,2,3,4 to
    * following number of quadrature points for 2D and 3D:
-   *   - 2D: 1, 3, 7
-   *   - 3D: 1, 4, 10
+   *   - 2D: 1, 3, 7, 15
+   *   - 3D: 1, 4, 10, 35
    *
    * For 1D, the quadrature rule degenerates to a `QGauss<1>(n_points_1D)`.
    *

diff --git a/source/simplex/quadrature_lib.cc b/source/simplex/quadrature_lib.cc
index 624c3e8125e99a98ef5aa2e49387f5d988cb5ccf..cbbd1c9ad66f0e0edd6d8fdf44b885aaabbb0793 100644 (file)
--- a/source/simplex/quadrature_lib.cc
+++ b/source/simplex/quadrature_lib.cc
@@ -84,41 +84,7 @@ namespace Simplex
           }
         else if (n_points_1D == 4)
           {
-            const double q12 = 0.5;
-
-            // clang-format off
-            this->quadrature_points.emplace_back(0.06493051315916486, 0.06493051315916486);
-            this->quadrature_points.emplace_back(0.06493051315916486, 0.8701389736816703);
-            this->quadrature_points.emplace_back(0.8701389736816703, 0.06493051315916486);
-            this->quadrature_points.emplace_back(0.2845755842491703, 0.517039939069323);
-            this->quadrature_points.emplace_back(0.3135591843849315, 0.04386347179237249);
-            this->quadrature_points.emplace_back(0.19838447668150672, 0.2845755842491703);
-            this->quadrature_points.emplace_back(0.642577343822696, 0.3135591843849315);
-            this->quadrature_points.emplace_back(0.517039939069323, 0.19838447668150672);
-            this->quadrature_points.emplace_back(0.04386347179237249, 0.642577343822696);
-            this->quadrature_points.emplace_back(0.517039939069323, 0.2845755842491703);
-            this->quadrature_points.emplace_back(0.04386347179237249, 0.3135591843849315);
-            this->quadrature_points.emplace_back(0.19838447668150672, 0.517039939069323);
-            this->quadrature_points.emplace_back(0.642577343822696, 0.04386347179237249);
-            this->quadrature_points.emplace_back(0.2845755842491703, 0.19838447668150672);
-            this->quadrature_points.emplace_back(0.3135591843849315, 0.642577343822696);
-            // clang-format on
-
-            this->weights.emplace_back(q12 * 0.053077801790232484);
-            this->weights.emplace_back(q12 * 0.053077801790232484);
-            this->weights.emplace_back(q12 * 0.053077801790232484);
-            this->weights.emplace_back(q12 * 0.07085308369213365);
-            this->weights.emplace_back(q12 * 0.069274682079417);
-            this->weights.emplace_back(q12 * 0.07085308369213365);
-            this->weights.emplace_back(q12 * 0.069274682079417);
-            this->weights.emplace_back(q12 * 0.07085308369213365);
-            this->weights.emplace_back(q12 * 0.069274682079417);
-            this->weights.emplace_back(q12 * 0.07085308369213365);
-            this->weights.emplace_back(q12 * 0.069274682079417);
-            this->weights.emplace_back(q12 * 0.07085308369213365);
-            this->weights.emplace_back(q12 * 0.069274682079417);
-            this->weights.emplace_back(q12 * 0.07085308369213365);
-            this->weights.emplace_back(q12 * 0.069274682079417);
+            Quadrature<dim>::operator=(QWitherdenVincent<dim>(n_points_1D));
           }
       }
     else if (dim == 3)
@@ -176,81 +142,7 @@ namespace Simplex
           }
         else if (n_points_1D == 4)
           {
-            const double Q16 = 1.0 / 6.0;
-
-            // clang-format off
-            this->quadrature_points.emplace_back(0.25, 0.25, 0.25);
-            this->quadrature_points.emplace_back(0.3157011497782028, 0.3157011497782028, 0.3157011497782028);
-            this->quadrature_points.emplace_back(0.3157011497782028, 0.3157011497782028, 0.05289655066539156);
-            this->quadrature_points.emplace_back(0.3157011497782028, 0.05289655066539156, 0.3157011497782028);
-            this->quadrature_points.emplace_back(0.05289655066539156, 0.3157011497782028, 0.3157011497782028);
-            this->quadrature_points.emplace_back(0.05048982259839635, 0.05048982259839635, 0.44951017740160365);
-            this->quadrature_points.emplace_back(0.05048982259839635, 0.44951017740160365, 0.05048982259839635);
-            this->quadrature_points.emplace_back(0.44951017740160365, 0.05048982259839635, 0.05048982259839635);
-            this->quadrature_points.emplace_back(0.05048982259839635, 0.44951017740160365, 0.44951017740160365);
-            this->quadrature_points.emplace_back(0.44951017740160365, 0.05048982259839635, 0.44951017740160365);
-            this->quadrature_points.emplace_back(0.44951017740160365, 0.44951017740160365, 0.05048982259839635);
-            this->quadrature_points.emplace_back(0.188833831026001, 0.188833831026001, 0.575171637587);
-            this->quadrature_points.emplace_back(0.02126547254148314, 0.02126547254148314, 0.8108302410985486);
-            this->quadrature_points.emplace_back(0.188833831026001, 0.575171637587, 0.188833831026001);
-            this->quadrature_points.emplace_back(0.02126547254148314, 0.8108302410985486, 0.02126547254148314);
-            this->quadrature_points.emplace_back(0.575171637587, 0.188833831026001, 0.188833831026001);
-            this->quadrature_points.emplace_back(0.8108302410985486, 0.02126547254148314, 0.02126547254148314);
-            this->quadrature_points.emplace_back(0.188833831026001, 0.575171637587, 0.047160700360998065);
-            this->quadrature_points.emplace_back(0.02126547254148314, 0.8108302410985486, 0.1466388138184851);
-            this->quadrature_points.emplace_back(0.575171637587, 0.188833831026001, 0.047160700360998065);
-            this->quadrature_points.emplace_back(0.8108302410985486, 0.02126547254148314, 0.1466388138184851);
-            this->quadrature_points.emplace_back(0.575171637587, 0.047160700360998065, 0.188833831026001);
-            this->quadrature_points.emplace_back(0.8108302410985486, 0.1466388138184851, 0.02126547254148314);
-            this->quadrature_points.emplace_back(0.188833831026001, 0.188833831026001, 0.047160700360998065);
-            this->quadrature_points.emplace_back(0.02126547254148314, 0.02126547254148314, 0.1466388138184851);
-            this->quadrature_points.emplace_back(0.188833831026001, 0.047160700360998065, 0.188833831026001);
-            this->quadrature_points.emplace_back(0.02126547254148314, 0.1466388138184851, 0.02126547254148314);
-            this->quadrature_points.emplace_back(0.047160700360998065, 0.188833831026001, 0.188833831026001);
-            this->quadrature_points.emplace_back(0.1466388138184851, 0.02126547254148314, 0.02126547254148314);
-            this->quadrature_points.emplace_back(0.188833831026001, 0.047160700360998065, 0.575171637587);
-            this->quadrature_points.emplace_back(0.02126547254148314, 0.1466388138184851, 0.8108302410985486);
-            this->quadrature_points.emplace_back(0.047160700360998065, 0.188833831026001, 0.575171637587);
-            this->quadrature_points.emplace_back(0.1466388138184851, 0.02126547254148314, 0.8108302410985486);
-            this->quadrature_points.emplace_back(0.047160700360998065, 0.575171637587, 0.188833831026001);
-            this->quadrature_points.emplace_back(0.1466388138184851, 0.8108302410985486, 0.02126547254148314);
-            // clang-format on
-
-            this->weights.emplace_back(Q16 * 0.09548528946413085);
-            this->weights.emplace_back(Q16 * 0.04232958120996703);
-            this->weights.emplace_back(Q16 * 0.04232958120996703);
-            this->weights.emplace_back(Q16 * 0.04232958120996703);
-            this->weights.emplace_back(Q16 * 0.04232958120996703);
-            this->weights.emplace_back(Q16 * 0.03189692783285758);
-            this->weights.emplace_back(Q16 * 0.03189692783285758);
-            this->weights.emplace_back(Q16 * 0.03189692783285758);
-            this->weights.emplace_back(Q16 * 0.03189692783285758);
-            this->weights.emplace_back(Q16 * 0.03189692783285758);
-            this->weights.emplace_back(Q16 * 0.03189692783285758);
-            this->weights.emplace_back(Q16 * 0.03720713072833462);
-            this->weights.emplace_back(Q16 * 0.008110770829903342);
-            this->weights.emplace_back(Q16 * 0.03720713072833462);
-            this->weights.emplace_back(Q16 * 0.008110770829903342);
-            this->weights.emplace_back(Q16 * 0.03720713072833462);
-            this->weights.emplace_back(Q16 * 0.008110770829903342);
-            this->weights.emplace_back(Q16 * 0.03720713072833462);
-            this->weights.emplace_back(Q16 * 0.008110770829903342);
-            this->weights.emplace_back(Q16 * 0.03720713072833462);
-            this->weights.emplace_back(Q16 * 0.008110770829903342);
-            this->weights.emplace_back(Q16 * 0.03720713072833462);
-            this->weights.emplace_back(Q16 * 0.008110770829903342);
-            this->weights.emplace_back(Q16 * 0.03720713072833462);
-            this->weights.emplace_back(Q16 * 0.008110770829903342);
-            this->weights.emplace_back(Q16 * 0.03720713072833462);
-            this->weights.emplace_back(Q16 * 0.008110770829903342);
-            this->weights.emplace_back(Q16 * 0.03720713072833462);
-            this->weights.emplace_back(Q16 * 0.008110770829903342);
-            this->weights.emplace_back(Q16 * 0.03720713072833462);
-            this->weights.emplace_back(Q16 * 0.008110770829903342);
-            this->weights.emplace_back(Q16 * 0.03720713072833462);
-            this->weights.emplace_back(Q16 * 0.008110770829903342);
-            this->weights.emplace_back(Q16 * 0.03720713072833462);
-            this->weights.emplace_back(Q16 * 0.008110770829903342);
+            Quadrature<dim>::operator=(QWitherdenVincent<dim>(n_points_1D));
           }
       }