Add a test with poor output.

This test passes on master, but shouldn't, since several of the printed
integrands are wrong. This will be corrected in the next commit.

diff --git a/tests/simplex/qgauss_01.cc b/tests/simplex/qgauss_01.cc
new file mode 100644 (file)
index 0000000..a517fbb
--- /dev/null
+++ b/tests/simplex/qgauss_01.cc
@@ -0,0 +1,82 @@
+// ---------------------------------------------------------------------
+//
+// 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.
+//
+// ---------------------------------------------------------------------
+
+#include <deal.II/base/function_lib.h>
+#include <deal.II/base/quadrature_lib.h>
+
+#include "../tests.h"
+
+template <int dim>
+void
+print(const unsigned int n_points_1D)
+{
+  deallog << "n_points_1D = " << n_points_1D << std::endl;
+  const QGaussSimplex<dim> quad(n_points_1D);
+
+  deallog << "quad size = " << quad.size() << std::endl;
+  for (unsigned int q = 0; q < quad.size(); ++q)
+    {
+      deallog << quad.point(q) << " ";
+      deallog << quad.weight(q) << " ";
+      deallog << std::endl;
+    }
+}
+
+template <int dim>
+void
+check_accuracy_1D(const unsigned int n_points_1D)
+{
+  const unsigned int accuracy = 2 * n_points_1D - 1;
+
+  Tensor<1, dim> monomial_powers;
+  unsigned int   sum = 0;
+  for (unsigned int d = 0; d < dim; ++d)
+    {
+      monomial_powers[d] += accuracy / dim;
+      sum += accuracy / dim;
+    }
+
+  // if we aren't at the correct degree then add the rest to the final
+  // component
+  monomial_powers[dim - 1] += accuracy - sum;
+
+  const Functions::Monomial<dim> func(monomial_powers);
+  const QGaussSimplex<dim>       quad(n_points_1D);
+
+  deallog << "Monomial powers = " << monomial_powers << std::endl;
+  double integrand = 0.0;
+  for (unsigned int q = 0; q < quad.size(); ++q)
+    integrand += quad.weight(q) * func.value(quad.point(q));
+  auto old_precision = deallog.precision(16);
+  deallog << "Integrand = " << integrand << std::endl;
+  deallog.precision(old_precision);
+}
+
+int
+main()
+{
+  initlog();
+
+  deallog << std::endl << std::endl;
+  check_accuracy_1D<2>(1);
+  check_accuracy_1D<2>(2);
+  check_accuracy_1D<2>(3);
+  check_accuracy_1D<2>(4);
+
+  check_accuracy_1D<3>(1);
+  check_accuracy_1D<3>(2);
+  check_accuracy_1D<3>(3);
+  check_accuracy_1D<3>(4);
+}

diff --git a/tests/simplex/qgauss_01.output b/tests/simplex/qgauss_01.output
new file mode 100644 (file)
index 0000000..d691178
--- /dev/null
+++ b/tests/simplex/qgauss_01.output
@@ -0,0 +1,19 @@
+
+DEAL::
+DEAL::
+DEAL::Monomial powers = 0.00000 1.00000
+DEAL::Integrand = 0.1666666666666667
+DEAL::Monomial powers = 1.00000 2.00000
+DEAL::Integrand = 0.01620370370370370
+DEAL::Monomial powers = 2.00000 3.00000
+DEAL::Integrand = 0.002380952380952380
+DEAL::Monomial powers = 3.00000 4.00000
+DEAL::Integrand = 0.0003968253968253969
+DEAL::Monomial powers = 0.00000 0.00000 1.00000
+DEAL::Integrand = 0.04166666666666666
+DEAL::Monomial powers = 1.00000 1.00000 1.00000
+DEAL::Integrand = 0.001507514161979123
+DEAL::Monomial powers = 1.00000 1.00000 3.00000
+DEAL::Integrand = 0.0001578591272523996
+DEAL::Monomial powers = 2.00000 2.00000 3.00000
+DEAL::Integrand = 6.613756613756609e-06