]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Implement ReferenceCell::face_measure().
authorDavid Wells <drwells@email.unc.edu>
Mon, 31 Mar 2025 13:14:02 +0000 (09:14 -0400)
committerDavid Wells <drwells@email.unc.edu>
Mon, 31 Mar 2025 22:52:02 +0000 (18:52 -0400)
We need these values in QProjector.

include/deal.II/grid/reference_cell.h
tests/grid/reference_cell_type_06.cc [new file with mode: 0644]
tests/grid/reference_cell_type_06.output [new file with mode: 0644]

index 046c926870740788d6e55d71f2d2f2be92987b81..17a1414aacb68ec7c136d4d7126198265d9bc262 100644 (file)
@@ -681,6 +681,24 @@ public:
   double
   volume() const;
 
+  /**
+   * Return the $d - 1$-dimensional measure of the face of a reference cell that
+   * corresponds to the current object, where $d$ is the dimension of the space
+   * it lives in.
+   *
+   * In this context the measure of a face depends on both the type of reference
+   * cell as well as the face number (i.e., `*this` and @p face_no). For
+   * example, the measures of the sides of a ReferenceCells::Triangle are $1$,
+   * $\sqrt{2}$, and $1$, respectively, whereas the measures of the sides of a
+   * ReferenceCells::Quadrilateral are all $1$.
+   *
+   * @note Like ReferenceCell::volume(), this function defines vertices to have
+   * a measure of $1$ to enable the definition of one-dimensional face
+   * quadratures.
+   */
+  double
+  face_measure(const unsigned int face_no) const;
+
   /**
    * Return the barycenter (i.e., the center of mass) of the reference
    * cell that corresponds to the current object. The function is not
@@ -3022,6 +3040,43 @@ ReferenceCell::volume() const
 
 
 
+inline double
+ReferenceCell::face_measure(const unsigned int face_no) const
+{
+  AssertIndexRange(face_no, n_faces());
+  constexpr double X = std::numeric_limits<double>::signaling_NaN();
+
+  static const std::array<double, 5> tri_faces{
+    {1.0, std::sqrt(2.0), 1.0, X, X}};
+  static const std::array<double, 5> tet_faces{
+    {0.5, 0.5, 0.5, std::sqrt(3.0) / 2.0, X}};
+  static const std::array<double, 5> pyramid_faces{
+    {4, std::sqrt(2.0), std::sqrt(2.0), std::sqrt(2.0), std::sqrt(2.0)}};
+  static const std::array<double, 5> wedge_faces{
+    {0.5, 0.5, 1.0, std::sqrt(2.0), 1.0}};
+
+  switch (this->kind)
+    {
+      case ReferenceCells::Vertex:
+      case ReferenceCells::Line:
+      case ReferenceCells::Quadrilateral:
+      case ReferenceCells::Hexahedron:
+        return 1.0;
+      case ReferenceCells::Triangle:
+        return tri_faces[face_no];
+      case ReferenceCells::Tetrahedron:
+        return tet_faces[face_no];
+      case ReferenceCells::Pyramid:
+        return pyramid_faces[face_no];
+      case ReferenceCells::Wedge:
+        return wedge_faces[face_no];
+    }
+  DEAL_II_ASSERT_UNREACHABLE();
+  return X;
+}
+
+
+
 template <int dim>
 inline Point<dim>
 ReferenceCell::barycenter() const
diff --git a/tests/grid/reference_cell_type_06.cc b/tests/grid/reference_cell_type_06.cc
new file mode 100644 (file)
index 0000000..6f7fd59
--- /dev/null
@@ -0,0 +1,108 @@
+// ------------------------------------------------------------------------
+//
+// SPDX-License-Identifier: LGPL-2.1-or-later
+// Copyright (C) 2022 - 2024 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.
+//
+// ------------------------------------------------------------------------
+
+
+// Test ReferenceCell::face_measure()
+
+#include <deal.II/base/quadrature_lib.h>
+
+#include <deal.II/fe/fe_nothing.h>
+#include <deal.II/fe/fe_values.h>
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/tria.h>
+
+#include <deal.II/hp/q_collection.h>
+
+#include "../tests.h"
+
+
+template <int dim>
+void
+test(const ReferenceCell &reference_cell)
+{
+  Triangulation<dim> triangulation;
+  GridGenerator::reference_cell(triangulation, reference_cell);
+
+  FE_Nothing<dim> fe(reference_cell);
+
+  hp::QCollection<dim - 1> quadratures;
+  // TODO: MappingQ asserts that face quadrature collections have exactly one
+  // entry. For now just special-case those.
+  if (reference_cell.is_hyper_cube())
+    quadratures.push_back(reference_cell.face_reference_cell(0)
+                            .template get_gauss_type_quadrature<dim - 1>(2));
+  else
+    for (unsigned int face_no = 0; face_no < reference_cell.n_faces();
+         ++face_no)
+      quadratures.push_back(reference_cell.face_reference_cell(face_no)
+                              .template get_gauss_type_quadrature<dim - 1>(2));
+
+  // Set up the objects to compute an integral on the reference cell
+  FEFaceValues<dim> fe_face_values(fe, quadratures, update_JxW_values);
+
+  for (unsigned int face_no = 0; face_no < reference_cell.n_faces(); ++face_no)
+    {
+      fe_face_values.reinit(triangulation.begin_active(), face_no);
+
+      double measure = 0;
+      for (unsigned int i = 0; i < fe_face_values.n_quadrature_points; ++i)
+        measure += fe_face_values.JxW(i);
+
+      deallog << "ReferenceCell: " << reference_cell.to_string() << std::endl;
+      deallog << "  computed face measure = " << measure << std::endl;
+      deallog << "  self-reported face measure = "
+              << reference_cell.face_measure(face_no) << std::endl;
+
+      // TriaAccessor<0, 1, 1>::measure() does not yet exist
+      if constexpr (dim > 1)
+        {
+          Assert(std::fabs(
+                   measure -
+                   triangulation.begin_active()->face(face_no)->measure()) <
+                   1e-12,
+                 ExcInternalError());
+        }
+      Assert(std::fabs(measure - reference_cell.face_measure(face_no)) < 1e-12,
+             ExcInternalError());
+    }
+}
+
+int
+main()
+{
+  initlog();
+
+  {
+    deallog.push("1D");
+    test<1>(ReferenceCells::Line);
+    deallog.pop();
+  }
+
+  {
+    deallog.push("2D");
+    test<2>(ReferenceCells::Quadrilateral);
+    test<2>(ReferenceCells::Triangle);
+    deallog.pop();
+  }
+
+  {
+    deallog.push("3D");
+    test<3>(ReferenceCells::Tetrahedron);
+    test<3>(ReferenceCells::Pyramid);
+    test<3>(ReferenceCells::Wedge);
+    test<3>(ReferenceCells::Hexahedron);
+    deallog.pop();
+  }
+}
diff --git a/tests/grid/reference_cell_type_06.output b/tests/grid/reference_cell_type_06.output
new file mode 100644 (file)
index 0000000..ff07dd1
--- /dev/null
@@ -0,0 +1,88 @@
+
+DEAL:1D::ReferenceCell: Line
+DEAL:1D::  computed face measure = 1.00000
+DEAL:1D::  self-reported face measure = 1.00000
+DEAL:1D::ReferenceCell: Line
+DEAL:1D::  computed face measure = 1.00000
+DEAL:1D::  self-reported face measure = 1.00000
+DEAL:2D::ReferenceCell: Quad
+DEAL:2D::  computed face measure = 1.00000
+DEAL:2D::  self-reported face measure = 1.00000
+DEAL:2D::ReferenceCell: Quad
+DEAL:2D::  computed face measure = 1.00000
+DEAL:2D::  self-reported face measure = 1.00000
+DEAL:2D::ReferenceCell: Quad
+DEAL:2D::  computed face measure = 1.00000
+DEAL:2D::  self-reported face measure = 1.00000
+DEAL:2D::ReferenceCell: Quad
+DEAL:2D::  computed face measure = 1.00000
+DEAL:2D::  self-reported face measure = 1.00000
+DEAL:2D::ReferenceCell: Tri
+DEAL:2D::  computed face measure = 1.00000
+DEAL:2D::  self-reported face measure = 1.00000
+DEAL:2D::ReferenceCell: Tri
+DEAL:2D::  computed face measure = 1.41421
+DEAL:2D::  self-reported face measure = 1.41421
+DEAL:2D::ReferenceCell: Tri
+DEAL:2D::  computed face measure = 1.00000
+DEAL:2D::  self-reported face measure = 1.00000
+DEAL:3D::ReferenceCell: Tet
+DEAL:3D::  computed face measure = 0.500000
+DEAL:3D::  self-reported face measure = 0.500000
+DEAL:3D::ReferenceCell: Tet
+DEAL:3D::  computed face measure = 0.500000
+DEAL:3D::  self-reported face measure = 0.500000
+DEAL:3D::ReferenceCell: Tet
+DEAL:3D::  computed face measure = 0.500000
+DEAL:3D::  self-reported face measure = 0.500000
+DEAL:3D::ReferenceCell: Tet
+DEAL:3D::  computed face measure = 0.866025
+DEAL:3D::  self-reported face measure = 0.866025
+DEAL:3D::ReferenceCell: Pyramid
+DEAL:3D::  computed face measure = 4.00000
+DEAL:3D::  self-reported face measure = 4.00000
+DEAL:3D::ReferenceCell: Pyramid
+DEAL:3D::  computed face measure = 1.41421
+DEAL:3D::  self-reported face measure = 1.41421
+DEAL:3D::ReferenceCell: Pyramid
+DEAL:3D::  computed face measure = 1.41421
+DEAL:3D::  self-reported face measure = 1.41421
+DEAL:3D::ReferenceCell: Pyramid
+DEAL:3D::  computed face measure = 1.41421
+DEAL:3D::  self-reported face measure = 1.41421
+DEAL:3D::ReferenceCell: Pyramid
+DEAL:3D::  computed face measure = 1.41421
+DEAL:3D::  self-reported face measure = 1.41421
+DEAL:3D::ReferenceCell: Wedge
+DEAL:3D::  computed face measure = 0.500000
+DEAL:3D::  self-reported face measure = 0.500000
+DEAL:3D::ReferenceCell: Wedge
+DEAL:3D::  computed face measure = 0.500000
+DEAL:3D::  self-reported face measure = 0.500000
+DEAL:3D::ReferenceCell: Wedge
+DEAL:3D::  computed face measure = 1.00000
+DEAL:3D::  self-reported face measure = 1.00000
+DEAL:3D::ReferenceCell: Wedge
+DEAL:3D::  computed face measure = 1.41421
+DEAL:3D::  self-reported face measure = 1.41421
+DEAL:3D::ReferenceCell: Wedge
+DEAL:3D::  computed face measure = 1.00000
+DEAL:3D::  self-reported face measure = 1.00000
+DEAL:3D::ReferenceCell: Hex
+DEAL:3D::  computed face measure = 1.00000
+DEAL:3D::  self-reported face measure = 1.00000
+DEAL:3D::ReferenceCell: Hex
+DEAL:3D::  computed face measure = 1.00000
+DEAL:3D::  self-reported face measure = 1.00000
+DEAL:3D::ReferenceCell: Hex
+DEAL:3D::  computed face measure = 1.00000
+DEAL:3D::  self-reported face measure = 1.00000
+DEAL:3D::ReferenceCell: Hex
+DEAL:3D::  computed face measure = 1.00000
+DEAL:3D::  self-reported face measure = 1.00000
+DEAL:3D::ReferenceCell: Hex
+DEAL:3D::  computed face measure = 1.00000
+DEAL:3D::  self-reported face measure = 1.00000
+DEAL:3D::ReferenceCell: Hex
+DEAL:3D::  computed face measure = 1.00000
+DEAL:3D::  self-reported face measure = 1.00000

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.