]> https://gitweb.dealii.org/ - dealii.git/commitdiff
QProjector: add another test to verify face quadrature rules. 18364/head
authorDavid Wells <drwells@email.unc.edu>
Mon, 21 Apr 2025 13:03:29 +0000 (09:03 -0400)
committerDavid Wells <drwells@email.unc.edu>
Mon, 21 Apr 2025 13:04:26 +0000 (09:04 -0400)
tests/base/qprojector_02.cc [new file with mode: 0644]
tests/base/qprojector_02.output [new file with mode: 0644]

diff --git a/tests/base/qprojector_02.cc b/tests/base/qprojector_02.cc
new file mode 100644 (file)
index 0000000..8a9d9d6
--- /dev/null
@@ -0,0 +1,194 @@
+// ------------------------------------------------------------------------
+//
+// SPDX-License-Identifier: LGPL-2.1-or-later
+// Copyright (C) 2025 - 2025 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// Part of the source code is dual licensed under Apache-2.0 WITH
+// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
+// governing the source code and code contributions can be found in
+// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
+//
+// ------------------------------------------------------------------------
+
+// Verify that the quadrature rules created by QProjector::project_to_face() and
+// QProjector::project_to_all_faces() are correct. I calculated the reference
+// values by-hand.
+
+#include <deal.II/base/qprojector.h>
+#include <deal.II/base/quadrature_lib.h>
+
+#include <deal.II/grid/reference_cell.h>
+
+#include "../tests.h"
+
+template <int dim, std::size_t N, typename F>
+void
+test(const ReferenceCell             reference_cell,
+     const hp::QCollection<dim - 1> &quadrature,
+     const std::array<double, N>    &integrals,
+     const F                        &integrand)
+{
+  Assert(reference_cell.n_faces() == N, ExcInternalError());
+  const auto all_quadratures =
+    QProjector<dim>::project_to_all_faces(reference_cell, quadrature);
+  for (const unsigned int face_no : reference_cell.face_indices())
+    for (types::geometric_orientation combined_orientation = 0;
+         combined_orientation < reference_cell.n_face_orientations(face_no);
+         ++combined_orientation)
+      {
+        const auto offset = QProjector<dim>::DataSetDescriptor::face(
+          reference_cell, face_no, combined_orientation, quadrature);
+        const Quadrature<dim> face_quadrature =
+          QProjector<dim>::project_to_face(
+            reference_cell,
+            quadrature[quadrature.size() == 1 ? 0 : face_no],
+            face_no,
+            combined_orientation);
+        double integral = 0;
+        for (unsigned int qp_n = 0; qp_n < face_quadrature.size(); ++qp_n)
+          {
+            integral += integrand(face_quadrature.point(qp_n)) *
+                        face_quadrature.weight(qp_n);
+            Assert(std::abs((all_quadratures.point(offset + qp_n) -
+                             face_quadrature.point(qp_n))
+                              .norm()) < 1e-14,
+                   ExcInternalError());
+            Assert(std::abs(all_quadratures.weight(offset + qp_n) -
+                            face_quadrature.weight(qp_n)) < 1e-14,
+                   ExcInternalError());
+          }
+        Assert(std::abs(integral - integrals[face_no]) <
+                 std::max(std::abs(integral), 1.0) * 1e-14,
+               ExcInternalError());
+        deallog << "(face_no, combined_orientation) = (" << face_no << ", "
+                << int(combined_orientation) << ")" << std::endl
+                << "value = " << integral << std::endl;
+      }
+}
+
+int
+main()
+{
+  initlog();
+
+  // test that we can correctly integrate (x + 1)^2 on each face
+  {
+    const QGauss<0> quadrature(1);
+    const auto      integrand = [&](const Point<1> &p) {
+      const double u = p[0] + 1;
+      return u * u;
+    };
+    deallog.push("Line");
+    {
+      const std::array<double, 2> integrals{{1.0, 4.0}};
+      test<1>(ReferenceCells::Line,
+              hp::QCollection<0>(quadrature),
+              integrals,
+              integrand);
+    }
+    deallog.pop();
+  }
+
+  // test that we can correctly integrate (x + 1)^2 * (y + 2)^2 on each face
+  {
+    const QGauss<1> quadrature(4);
+    const auto      integrand = [&](const Point<2> &p) {
+      const double u = p[0] + 1, v = p[1] + 2;
+      return u * u * v * v;
+    };
+    deallog.push("Triangle");
+    {
+      const std::array<double, 3> integrals{
+        {28.0 / 3.0, 203.0 / 15.0 * std::sqrt(2.0), 19.0 / 3.0}};
+      test<2>(ReferenceCells::Triangle,
+              hp::QCollection<1>(quadrature),
+              integrals,
+              integrand);
+    }
+    deallog.pop();
+
+    deallog.push("Quadrilateral");
+    {
+      const std::array<double, 4> integrals{
+        {19.0 / 3.0, 76.0 / 3.0, 28.0 / 3.0, 21.0}};
+      test<2>(ReferenceCells::Quadrilateral,
+              hp::QCollection<1>(quadrature),
+              integrals,
+              integrand);
+    }
+    deallog.pop();
+  }
+
+  // test that we can correctly integrate (x + 1)^2 (y + 2)^2 + (z + 3)^2 on
+  // each face
+  {
+    const QGauss<2>        quadrature_quadrilateral(4);
+    const QGaussSimplex<2> quadrature_triangle(4);
+    const auto             integrand = [&](const Point<3> &p) {
+      const double u = p[0] + 1.0, v = p[1] + 2.0, w = p[2] + 3.0;
+      return u * u * v * v + w * w;
+    };
+    deallog.push("Tetrahedron");
+    {
+      const std::array<double, 4> integrals{{421.0 / 45.0,
+                                             37.0 / 4.0,
+                                             25.0 / 3.0,
+                                             1879.0 / 180.0 * std::sqrt(3.0)}};
+      test<3>(ReferenceCells::Tetrahedron,
+              hp::QCollection<2>(quadrature_triangle),
+              integrals,
+              integrand);
+    }
+    deallog.pop();
+
+    deallog.push("Pyramid");
+    {
+      const std::array<double, 5> integrals{{532.0 / 9.0,
+                                             533.0 / 45.0 * std::sqrt(2.0),
+                                             1037.0 / 45.0 * std::sqrt(2.0),
+                                             596.0 / 45.0 * std::sqrt(2.0),
+                                             884.0 / 45.0 * std::sqrt(2.0)}};
+      test<3>(ReferenceCells::Pyramid,
+              hp::QCollection<2>(quadrature_quadrilateral,
+                                 quadrature_triangle,
+                                 quadrature_triangle,
+                                 quadrature_triangle,
+                                 quadrature_triangle),
+              integrals,
+              integrand);
+    }
+    deallog.pop();
+
+    const double X = 0.0;
+    deallog.push("Wedge");
+    {
+      const std::array<double, 5> integrals{{421.0 / 45.0,
+                                             1157.0 / 90.0,
+                                             65.0 / 3.0,
+                                             388.0 / 15.0 * std::sqrt(2.0),
+                                             56.0 / 3.0}};
+      test<3>(ReferenceCells::Wedge,
+              hp::QCollection<2>(quadrature_triangle,
+                                 quadrature_triangle,
+                                 quadrature_quadrilateral,
+                                 quadrature_quadrilateral,
+                                 quadrature_quadrilateral),
+              integrals,
+              integrand);
+    }
+    deallog.pop();
+
+    deallog.push("Hexahedron");
+    {
+      const std::array<double, 6> integrals{
+        {56 / 3.0, 113.0 / 3.0, 65 / 3.0, 100 / 3.0, 214.0 / 9.0, 277 / 9.0}};
+      test<3>(ReferenceCells::Hexahedron,
+              hp::QCollection<2>(quadrature_quadrilateral),
+              integrals,
+              integrand);
+    }
+    deallog.pop();
+  }
+}
diff --git a/tests/base/qprojector_02.output b/tests/base/qprojector_02.output
new file mode 100644 (file)
index 0000000..4deb8ee
--- /dev/null
@@ -0,0 +1,313 @@
+
+DEAL:Line::(face_no, combined_orientation) = (0, 0)
+DEAL:Line::value = 1.00000
+DEAL:Line::(face_no, combined_orientation) = (1, 0)
+DEAL:Line::value = 4.00000
+DEAL:Triangle::(face_no, combined_orientation) = (0, 0)
+DEAL:Triangle::value = 9.33333
+DEAL:Triangle::(face_no, combined_orientation) = (0, 1)
+DEAL:Triangle::value = 9.33333
+DEAL:Triangle::(face_no, combined_orientation) = (1, 0)
+DEAL:Triangle::value = 19.1390
+DEAL:Triangle::(face_no, combined_orientation) = (1, 1)
+DEAL:Triangle::value = 19.1390
+DEAL:Triangle::(face_no, combined_orientation) = (2, 0)
+DEAL:Triangle::value = 6.33333
+DEAL:Triangle::(face_no, combined_orientation) = (2, 1)
+DEAL:Triangle::value = 6.33333
+DEAL:Quadrilateral::(face_no, combined_orientation) = (0, 0)
+DEAL:Quadrilateral::value = 6.33333
+DEAL:Quadrilateral::(face_no, combined_orientation) = (0, 1)
+DEAL:Quadrilateral::value = 6.33333
+DEAL:Quadrilateral::(face_no, combined_orientation) = (1, 0)
+DEAL:Quadrilateral::value = 25.3333
+DEAL:Quadrilateral::(face_no, combined_orientation) = (1, 1)
+DEAL:Quadrilateral::value = 25.3333
+DEAL:Quadrilateral::(face_no, combined_orientation) = (2, 0)
+DEAL:Quadrilateral::value = 9.33333
+DEAL:Quadrilateral::(face_no, combined_orientation) = (2, 1)
+DEAL:Quadrilateral::value = 9.33333
+DEAL:Quadrilateral::(face_no, combined_orientation) = (3, 0)
+DEAL:Quadrilateral::value = 21.0000
+DEAL:Quadrilateral::(face_no, combined_orientation) = (3, 1)
+DEAL:Quadrilateral::value = 21.0000
+DEAL:Tetrahedron::(face_no, combined_orientation) = (0, 0)
+DEAL:Tetrahedron::value = 9.35556
+DEAL:Tetrahedron::(face_no, combined_orientation) = (0, 1)
+DEAL:Tetrahedron::value = 9.35556
+DEAL:Tetrahedron::(face_no, combined_orientation) = (0, 2)
+DEAL:Tetrahedron::value = 9.35556
+DEAL:Tetrahedron::(face_no, combined_orientation) = (0, 3)
+DEAL:Tetrahedron::value = 9.35556
+DEAL:Tetrahedron::(face_no, combined_orientation) = (0, 4)
+DEAL:Tetrahedron::value = 9.35556
+DEAL:Tetrahedron::(face_no, combined_orientation) = (0, 5)
+DEAL:Tetrahedron::value = 9.35556
+DEAL:Tetrahedron::(face_no, combined_orientation) = (1, 0)
+DEAL:Tetrahedron::value = 9.25000
+DEAL:Tetrahedron::(face_no, combined_orientation) = (1, 1)
+DEAL:Tetrahedron::value = 9.25000
+DEAL:Tetrahedron::(face_no, combined_orientation) = (1, 2)
+DEAL:Tetrahedron::value = 9.25000
+DEAL:Tetrahedron::(face_no, combined_orientation) = (1, 3)
+DEAL:Tetrahedron::value = 9.25000
+DEAL:Tetrahedron::(face_no, combined_orientation) = (1, 4)
+DEAL:Tetrahedron::value = 9.25000
+DEAL:Tetrahedron::(face_no, combined_orientation) = (1, 5)
+DEAL:Tetrahedron::value = 9.25000
+DEAL:Tetrahedron::(face_no, combined_orientation) = (2, 0)
+DEAL:Tetrahedron::value = 8.33333
+DEAL:Tetrahedron::(face_no, combined_orientation) = (2, 1)
+DEAL:Tetrahedron::value = 8.33333
+DEAL:Tetrahedron::(face_no, combined_orientation) = (2, 2)
+DEAL:Tetrahedron::value = 8.33333
+DEAL:Tetrahedron::(face_no, combined_orientation) = (2, 3)
+DEAL:Tetrahedron::value = 8.33333
+DEAL:Tetrahedron::(face_no, combined_orientation) = (2, 4)
+DEAL:Tetrahedron::value = 8.33333
+DEAL:Tetrahedron::(face_no, combined_orientation) = (2, 5)
+DEAL:Tetrahedron::value = 8.33333
+DEAL:Tetrahedron::(face_no, combined_orientation) = (3, 0)
+DEAL:Tetrahedron::value = 18.0807
+DEAL:Tetrahedron::(face_no, combined_orientation) = (3, 1)
+DEAL:Tetrahedron::value = 18.0807
+DEAL:Tetrahedron::(face_no, combined_orientation) = (3, 2)
+DEAL:Tetrahedron::value = 18.0807
+DEAL:Tetrahedron::(face_no, combined_orientation) = (3, 3)
+DEAL:Tetrahedron::value = 18.0807
+DEAL:Tetrahedron::(face_no, combined_orientation) = (3, 4)
+DEAL:Tetrahedron::value = 18.0807
+DEAL:Tetrahedron::(face_no, combined_orientation) = (3, 5)
+DEAL:Tetrahedron::value = 18.0807
+DEAL:Pyramid::(face_no, combined_orientation) = (0, 0)
+DEAL:Pyramid::value = 59.1111
+DEAL:Pyramid::(face_no, combined_orientation) = (0, 1)
+DEAL:Pyramid::value = 59.1111
+DEAL:Pyramid::(face_no, combined_orientation) = (0, 2)
+DEAL:Pyramid::value = 59.1111
+DEAL:Pyramid::(face_no, combined_orientation) = (0, 3)
+DEAL:Pyramid::value = 59.1111
+DEAL:Pyramid::(face_no, combined_orientation) = (0, 4)
+DEAL:Pyramid::value = 59.1111
+DEAL:Pyramid::(face_no, combined_orientation) = (0, 5)
+DEAL:Pyramid::value = 59.1111
+DEAL:Pyramid::(face_no, combined_orientation) = (0, 6)
+DEAL:Pyramid::value = 59.1111
+DEAL:Pyramid::(face_no, combined_orientation) = (0, 7)
+DEAL:Pyramid::value = 59.1111
+DEAL:Pyramid::(face_no, combined_orientation) = (1, 0)
+DEAL:Pyramid::value = 16.7506
+DEAL:Pyramid::(face_no, combined_orientation) = (1, 1)
+DEAL:Pyramid::value = 16.7506
+DEAL:Pyramid::(face_no, combined_orientation) = (1, 2)
+DEAL:Pyramid::value = 16.7506
+DEAL:Pyramid::(face_no, combined_orientation) = (1, 3)
+DEAL:Pyramid::value = 16.7506
+DEAL:Pyramid::(face_no, combined_orientation) = (1, 4)
+DEAL:Pyramid::value = 16.7506
+DEAL:Pyramid::(face_no, combined_orientation) = (1, 5)
+DEAL:Pyramid::value = 16.7506
+DEAL:Pyramid::(face_no, combined_orientation) = (2, 0)
+DEAL:Pyramid::value = 32.5898
+DEAL:Pyramid::(face_no, combined_orientation) = (2, 1)
+DEAL:Pyramid::value = 32.5898
+DEAL:Pyramid::(face_no, combined_orientation) = (2, 2)
+DEAL:Pyramid::value = 32.5898
+DEAL:Pyramid::(face_no, combined_orientation) = (2, 3)
+DEAL:Pyramid::value = 32.5898
+DEAL:Pyramid::(face_no, combined_orientation) = (2, 4)
+DEAL:Pyramid::value = 32.5898
+DEAL:Pyramid::(face_no, combined_orientation) = (2, 5)
+DEAL:Pyramid::value = 32.5898
+DEAL:Pyramid::(face_no, combined_orientation) = (3, 0)
+DEAL:Pyramid::value = 18.7305
+DEAL:Pyramid::(face_no, combined_orientation) = (3, 1)
+DEAL:Pyramid::value = 18.7305
+DEAL:Pyramid::(face_no, combined_orientation) = (3, 2)
+DEAL:Pyramid::value = 18.7305
+DEAL:Pyramid::(face_no, combined_orientation) = (3, 3)
+DEAL:Pyramid::value = 18.7305
+DEAL:Pyramid::(face_no, combined_orientation) = (3, 4)
+DEAL:Pyramid::value = 18.7305
+DEAL:Pyramid::(face_no, combined_orientation) = (3, 5)
+DEAL:Pyramid::value = 18.7305
+DEAL:Pyramid::(face_no, combined_orientation) = (4, 0)
+DEAL:Pyramid::value = 27.7814
+DEAL:Pyramid::(face_no, combined_orientation) = (4, 1)
+DEAL:Pyramid::value = 27.7814
+DEAL:Pyramid::(face_no, combined_orientation) = (4, 2)
+DEAL:Pyramid::value = 27.7814
+DEAL:Pyramid::(face_no, combined_orientation) = (4, 3)
+DEAL:Pyramid::value = 27.7814
+DEAL:Pyramid::(face_no, combined_orientation) = (4, 4)
+DEAL:Pyramid::value = 27.7814
+DEAL:Pyramid::(face_no, combined_orientation) = (4, 5)
+DEAL:Pyramid::value = 27.7814
+DEAL:Wedge::(face_no, combined_orientation) = (0, 0)
+DEAL:Wedge::value = 9.35556
+DEAL:Wedge::(face_no, combined_orientation) = (0, 1)
+DEAL:Wedge::value = 9.35556
+DEAL:Wedge::(face_no, combined_orientation) = (0, 2)
+DEAL:Wedge::value = 9.35556
+DEAL:Wedge::(face_no, combined_orientation) = (0, 3)
+DEAL:Wedge::value = 9.35556
+DEAL:Wedge::(face_no, combined_orientation) = (0, 4)
+DEAL:Wedge::value = 9.35556
+DEAL:Wedge::(face_no, combined_orientation) = (0, 5)
+DEAL:Wedge::value = 9.35556
+DEAL:Wedge::(face_no, combined_orientation) = (1, 0)
+DEAL:Wedge::value = 12.8556
+DEAL:Wedge::(face_no, combined_orientation) = (1, 1)
+DEAL:Wedge::value = 12.8556
+DEAL:Wedge::(face_no, combined_orientation) = (1, 2)
+DEAL:Wedge::value = 12.8556
+DEAL:Wedge::(face_no, combined_orientation) = (1, 3)
+DEAL:Wedge::value = 12.8556
+DEAL:Wedge::(face_no, combined_orientation) = (1, 4)
+DEAL:Wedge::value = 12.8556
+DEAL:Wedge::(face_no, combined_orientation) = (1, 5)
+DEAL:Wedge::value = 12.8556
+DEAL:Wedge::(face_no, combined_orientation) = (2, 0)
+DEAL:Wedge::value = 21.6667
+DEAL:Wedge::(face_no, combined_orientation) = (2, 1)
+DEAL:Wedge::value = 21.6667
+DEAL:Wedge::(face_no, combined_orientation) = (2, 2)
+DEAL:Wedge::value = 21.6667
+DEAL:Wedge::(face_no, combined_orientation) = (2, 3)
+DEAL:Wedge::value = 21.6667
+DEAL:Wedge::(face_no, combined_orientation) = (2, 4)
+DEAL:Wedge::value = 21.6667
+DEAL:Wedge::(face_no, combined_orientation) = (2, 5)
+DEAL:Wedge::value = 21.6667
+DEAL:Wedge::(face_no, combined_orientation) = (2, 6)
+DEAL:Wedge::value = 21.6667
+DEAL:Wedge::(face_no, combined_orientation) = (2, 7)
+DEAL:Wedge::value = 21.6667
+DEAL:Wedge::(face_no, combined_orientation) = (3, 0)
+DEAL:Wedge::value = 36.5810
+DEAL:Wedge::(face_no, combined_orientation) = (3, 1)
+DEAL:Wedge::value = 36.5810
+DEAL:Wedge::(face_no, combined_orientation) = (3, 2)
+DEAL:Wedge::value = 36.5810
+DEAL:Wedge::(face_no, combined_orientation) = (3, 3)
+DEAL:Wedge::value = 36.5810
+DEAL:Wedge::(face_no, combined_orientation) = (3, 4)
+DEAL:Wedge::value = 36.5810
+DEAL:Wedge::(face_no, combined_orientation) = (3, 5)
+DEAL:Wedge::value = 36.5810
+DEAL:Wedge::(face_no, combined_orientation) = (3, 6)
+DEAL:Wedge::value = 36.5810
+DEAL:Wedge::(face_no, combined_orientation) = (3, 7)
+DEAL:Wedge::value = 36.5810
+DEAL:Wedge::(face_no, combined_orientation) = (4, 0)
+DEAL:Wedge::value = 18.6667
+DEAL:Wedge::(face_no, combined_orientation) = (4, 1)
+DEAL:Wedge::value = 18.6667
+DEAL:Wedge::(face_no, combined_orientation) = (4, 2)
+DEAL:Wedge::value = 18.6667
+DEAL:Wedge::(face_no, combined_orientation) = (4, 3)
+DEAL:Wedge::value = 18.6667
+DEAL:Wedge::(face_no, combined_orientation) = (4, 4)
+DEAL:Wedge::value = 18.6667
+DEAL:Wedge::(face_no, combined_orientation) = (4, 5)
+DEAL:Wedge::value = 18.6667
+DEAL:Wedge::(face_no, combined_orientation) = (4, 6)
+DEAL:Wedge::value = 18.6667
+DEAL:Wedge::(face_no, combined_orientation) = (4, 7)
+DEAL:Wedge::value = 18.6667
+DEAL:Hexahedron::(face_no, combined_orientation) = (0, 0)
+DEAL:Hexahedron::value = 18.6667
+DEAL:Hexahedron::(face_no, combined_orientation) = (0, 1)
+DEAL:Hexahedron::value = 18.6667
+DEAL:Hexahedron::(face_no, combined_orientation) = (0, 2)
+DEAL:Hexahedron::value = 18.6667
+DEAL:Hexahedron::(face_no, combined_orientation) = (0, 3)
+DEAL:Hexahedron::value = 18.6667
+DEAL:Hexahedron::(face_no, combined_orientation) = (0, 4)
+DEAL:Hexahedron::value = 18.6667
+DEAL:Hexahedron::(face_no, combined_orientation) = (0, 5)
+DEAL:Hexahedron::value = 18.6667
+DEAL:Hexahedron::(face_no, combined_orientation) = (0, 6)
+DEAL:Hexahedron::value = 18.6667
+DEAL:Hexahedron::(face_no, combined_orientation) = (0, 7)
+DEAL:Hexahedron::value = 18.6667
+DEAL:Hexahedron::(face_no, combined_orientation) = (1, 0)
+DEAL:Hexahedron::value = 37.6667
+DEAL:Hexahedron::(face_no, combined_orientation) = (1, 1)
+DEAL:Hexahedron::value = 37.6667
+DEAL:Hexahedron::(face_no, combined_orientation) = (1, 2)
+DEAL:Hexahedron::value = 37.6667
+DEAL:Hexahedron::(face_no, combined_orientation) = (1, 3)
+DEAL:Hexahedron::value = 37.6667
+DEAL:Hexahedron::(face_no, combined_orientation) = (1, 4)
+DEAL:Hexahedron::value = 37.6667
+DEAL:Hexahedron::(face_no, combined_orientation) = (1, 5)
+DEAL:Hexahedron::value = 37.6667
+DEAL:Hexahedron::(face_no, combined_orientation) = (1, 6)
+DEAL:Hexahedron::value = 37.6667
+DEAL:Hexahedron::(face_no, combined_orientation) = (1, 7)
+DEAL:Hexahedron::value = 37.6667
+DEAL:Hexahedron::(face_no, combined_orientation) = (2, 0)
+DEAL:Hexahedron::value = 21.6667
+DEAL:Hexahedron::(face_no, combined_orientation) = (2, 1)
+DEAL:Hexahedron::value = 21.6667
+DEAL:Hexahedron::(face_no, combined_orientation) = (2, 2)
+DEAL:Hexahedron::value = 21.6667
+DEAL:Hexahedron::(face_no, combined_orientation) = (2, 3)
+DEAL:Hexahedron::value = 21.6667
+DEAL:Hexahedron::(face_no, combined_orientation) = (2, 4)
+DEAL:Hexahedron::value = 21.6667
+DEAL:Hexahedron::(face_no, combined_orientation) = (2, 5)
+DEAL:Hexahedron::value = 21.6667
+DEAL:Hexahedron::(face_no, combined_orientation) = (2, 6)
+DEAL:Hexahedron::value = 21.6667
+DEAL:Hexahedron::(face_no, combined_orientation) = (2, 7)
+DEAL:Hexahedron::value = 21.6667
+DEAL:Hexahedron::(face_no, combined_orientation) = (3, 0)
+DEAL:Hexahedron::value = 33.3333
+DEAL:Hexahedron::(face_no, combined_orientation) = (3, 1)
+DEAL:Hexahedron::value = 33.3333
+DEAL:Hexahedron::(face_no, combined_orientation) = (3, 2)
+DEAL:Hexahedron::value = 33.3333
+DEAL:Hexahedron::(face_no, combined_orientation) = (3, 3)
+DEAL:Hexahedron::value = 33.3333
+DEAL:Hexahedron::(face_no, combined_orientation) = (3, 4)
+DEAL:Hexahedron::value = 33.3333
+DEAL:Hexahedron::(face_no, combined_orientation) = (3, 5)
+DEAL:Hexahedron::value = 33.3333
+DEAL:Hexahedron::(face_no, combined_orientation) = (3, 6)
+DEAL:Hexahedron::value = 33.3333
+DEAL:Hexahedron::(face_no, combined_orientation) = (3, 7)
+DEAL:Hexahedron::value = 33.3333
+DEAL:Hexahedron::(face_no, combined_orientation) = (4, 0)
+DEAL:Hexahedron::value = 23.7778
+DEAL:Hexahedron::(face_no, combined_orientation) = (4, 1)
+DEAL:Hexahedron::value = 23.7778
+DEAL:Hexahedron::(face_no, combined_orientation) = (4, 2)
+DEAL:Hexahedron::value = 23.7778
+DEAL:Hexahedron::(face_no, combined_orientation) = (4, 3)
+DEAL:Hexahedron::value = 23.7778
+DEAL:Hexahedron::(face_no, combined_orientation) = (4, 4)
+DEAL:Hexahedron::value = 23.7778
+DEAL:Hexahedron::(face_no, combined_orientation) = (4, 5)
+DEAL:Hexahedron::value = 23.7778
+DEAL:Hexahedron::(face_no, combined_orientation) = (4, 6)
+DEAL:Hexahedron::value = 23.7778
+DEAL:Hexahedron::(face_no, combined_orientation) = (4, 7)
+DEAL:Hexahedron::value = 23.7778
+DEAL:Hexahedron::(face_no, combined_orientation) = (5, 0)
+DEAL:Hexahedron::value = 30.7778
+DEAL:Hexahedron::(face_no, combined_orientation) = (5, 1)
+DEAL:Hexahedron::value = 30.7778
+DEAL:Hexahedron::(face_no, combined_orientation) = (5, 2)
+DEAL:Hexahedron::value = 30.7778
+DEAL:Hexahedron::(face_no, combined_orientation) = (5, 3)
+DEAL:Hexahedron::value = 30.7778
+DEAL:Hexahedron::(face_no, combined_orientation) = (5, 4)
+DEAL:Hexahedron::value = 30.7778
+DEAL:Hexahedron::(face_no, combined_orientation) = (5, 5)
+DEAL:Hexahedron::value = 30.7778
+DEAL:Hexahedron::(face_no, combined_orientation) = (5, 6)
+DEAL:Hexahedron::value = 30.7778
+DEAL:Hexahedron::(face_no, combined_orientation) = (5, 7)
+DEAL:Hexahedron::value = 30.7778

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.