]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Specialize NonMatching::FaceQuadratureGenerator in 1D 13597/head
authorSimon Sticko <simon@sticko.se>
Wed, 6 Apr 2022 15:52:58 +0000 (17:52 +0200)
committerSimon Sticko <simon@sticko.se>
Thu, 7 Apr 2022 17:04:54 +0000 (19:04 +0200)
This class must be specialized in 1D, because the general
FaceQuadratureGenerator<dim> class uses the QuadratureGenerator<dim-1>
class internally, which does not make sense when dim-1 = 0.

include/deal.II/non_matching/quadrature_generator.h
source/non_matching/quadrature_generator.cc
tests/non_matching/face_quadrature_generator.cc
tests/non_matching/face_quadrature_generator.output

index 2fa8214601f8cace4f272e2f5af776bd01143e17..97653b097a26bf4c68fee58e2cafc854466e3718 100644 (file)
@@ -362,6 +362,103 @@ namespace NonMatching
   };
 
 
+  /**
+   * Specialization of the FaceQuadratureGenerator class for the 1-dimensional
+   * case.
+   *
+   * In 1D, a face is only a point. Thus to generate the immersed
+   * quadrature rules we add a single 0-dimensional quadrature point to the
+   * inside or outside quadrature rule depending on if the level set function is
+   * positive or negative at the face. The added quadrature point will have
+   * weight equal to 1. The immersed surface quadrature over a face corresponds
+   * to integrating over a dim-1 dimensional curve. Thus, surface quadrature
+   * generated by this specialized class is always empty.
+   *
+   * This class must be specialized in 1D, because the general
+   * FaceQuadratureGenerator<dim> class uses the QuadratureGenerator<dim-1>
+   * class internally, which does not make sense when dim-1 = 0.
+   */
+  template <>
+  class FaceQuadratureGenerator<1>
+  {
+  public:
+    using AdditionalData = AdditionalQGeneratorData;
+
+    /**
+     * Constructor. The incoming hp::QCollection is not used. But this class
+     * must have the same signature as the non-specialized class.
+     */
+    FaceQuadratureGenerator(
+      const hp::QCollection<1> &quadratures1D,
+      const AdditionalData &    additional_data = AdditionalData());
+
+    /**
+     * Construct immersed quadratures rules for the incoming level set
+     * function on a given face of the BoundingBox.
+     *
+     * To get the constructed quadratures, use the functions
+     * get_inside_quadrature(),
+     * get_outside_quadrature(),
+     * get_surface_quadrature().
+     */
+    void
+    generate(const Function<1> &   level_set,
+             const BoundingBox<1> &box,
+             const unsigned int    face_index);
+
+    /**
+     * @copydoc FaceQuadratureGenerator<dim>::get_inside_quadrature()
+     */
+    const Quadrature<0> &
+    get_inside_quadrature() const;
+
+    /**
+     * @copydoc FaceQuadratureGenerator<dim>::get_outside_quadrature()
+     */
+    const Quadrature<0> &
+    get_outside_quadrature() const;
+
+    /**
+     * Return the quadrature rule for the region
+     * $\{x \in F : \psi(x) = 0 \}$
+     * where, $F$ is the face of the BoundingBox passed to generate().
+     *
+     * @note In 1D, this quadrature always contains 0 points.
+     */
+    const ImmersedSurfaceQuadrature<0, 1> &
+    get_surface_quadrature() const;
+
+    /**
+     * This function does nothing. It only exist to be compatible with
+     * FaceQuadratureGenerator<dim>.
+     */
+    void
+    set_1D_quadrature(const unsigned int q_index);
+
+  private:
+    /**
+     * Quadrature for the region
+     * $\{x \in F : \psi(x) < 0 \}$.
+     * Created in the last call to generate().
+     */
+    Quadrature<0> inside_quadrature;
+
+    /**
+     * Quadrature for the region
+     * $\{x \in F : \psi(x) > 0 \}$.
+     * Created in the last call to generate().
+     */
+    Quadrature<0> outside_quadrature;
+
+    /**
+     * Quadrature for the region
+     * $\{x \in F : \psi(x) = 0 \}$.
+     * This quadrature always contains zero points in 1D.
+     */
+    const ImmersedSurfaceQuadrature<0, 1> surface_quadrature;
+  };
+
+
   namespace internal
   {
     namespace QuadratureGeneratorImplementation
index 9f446a375db6e133655334ca8c67d3bd6b576165..da551e41a1212648eff450fcec5e86b7673de4fe 100644 (file)
@@ -1448,6 +1448,83 @@ namespace NonMatching
   {
     return surface_quadrature;
   }
+
+
+
+  FaceQuadratureGenerator<1>::FaceQuadratureGenerator(
+    const hp::QCollection<1> &quadratures1D,
+    const AdditionalData &    additional_data)
+  {
+    (void)quadratures1D;
+    (void)additional_data;
+  }
+
+
+
+  void
+  FaceQuadratureGenerator<1>::generate(const Function<1> &   level_set,
+                                       const BoundingBox<1> &box,
+                                       const unsigned int    face_index)
+  {
+    AssertIndexRange(face_index, GeometryInfo<1>::faces_per_cell);
+
+    // The only vertex the 1D-face has.
+    const Point<1> vertex =
+      box.vertex(GeometryInfo<1>::face_to_cell_vertices(face_index, 0));
+
+    const unsigned int          n_points = 1;
+    const double                weight   = 1;
+    const std::vector<Point<0>> points(n_points);
+    const std::vector<double>   weights(n_points, weight);
+
+    const double level_set_value = level_set.value(vertex);
+    if (level_set_value < 0)
+      {
+        inside_quadrature  = Quadrature<0>(points, weights);
+        outside_quadrature = Quadrature<0>();
+      }
+    else if (level_set_value > 0)
+      {
+        inside_quadrature  = Quadrature<0>();
+        outside_quadrature = Quadrature<0>(points, weights);
+      }
+    else
+      {
+        inside_quadrature  = Quadrature<0>();
+        outside_quadrature = Quadrature<0>();
+      }
+  }
+
+
+
+  void
+  FaceQuadratureGenerator<1>::set_1D_quadrature(const unsigned int q_index)
+  {
+    (void)q_index;
+  }
+
+
+
+  const Quadrature<0> &
+  FaceQuadratureGenerator<1>::get_inside_quadrature() const
+  {
+    return inside_quadrature;
+  }
+
+
+  const Quadrature<0> &
+  FaceQuadratureGenerator<1>::get_outside_quadrature() const
+  {
+    return outside_quadrature;
+  }
+
+
+
+  const ImmersedSurfaceQuadrature<0, 1> &
+  FaceQuadratureGenerator<1>::get_surface_quadrature() const
+  {
+    return surface_quadrature;
+  }
 } // namespace NonMatching
 #include "quadrature_generator.inst"
 DEAL_II_NAMESPACE_CLOSE
index 4c0964ee7aba32539530a5095439e54850effe47..afb24dd8ce3eaa8788975e1f938da33935189cf2 100644 (file)
@@ -105,11 +105,40 @@ test_plane_cuts_through_center()
 
 
 
+/*
+ * Test the 1D-specialization of the FaceQuadratureGenerator class.
+ * Set up a 1D-box [0,1] and a level set function: psi(x) = x - 0.5,
+ * so that the level set function is negative at the left face and positive at
+ * the right. Generate quadrature rules at both faces and check that a single
+ * quadrature point at the inside/outside quadrature at the left/right face is
+ * generated.
+ */
+void
+test_1D()
+{
+  deallog << "test_1D" << std::endl;
+
+  const int                             dim = 1;
+  Point<dim>                            center(.5);
+  const Tensor<1, dim>                  normal = Point<dim>::unit_vector(0);
+  const Functions::LevelSet::Plane<dim> levelset(center, normal);
+
+  for (unsigned int f = 0; f < GeometryInfo<dim>::faces_per_cell; ++f)
+    {
+      deallog << "face = " << f << std::endl;
+      create_and_print_quadratures(levelset, f);
+      deallog << std::endl;
+    }
+}
+
+
+
 int
 main()
 {
   initlog();
 
+  test_1D();
   test_plane_cuts_through_center<2>();
   test_plane_cuts_through_center<3>();
 }
index c3a17d92597c7ab53e025fb1f7d9418aa36d544b..397a723acfa0dbd8a23d95bf2cba1c887a47a699 100644 (file)
@@ -1,4 +1,17 @@
 
+DEAL::test_1D
+DEAL::face = 0
+DEAL::Inside quadrature
+DEAL::1.00000
+DEAL::Outside quadrature
+DEAL::Surface quadrature
+DEAL::
+DEAL::face = 1
+DEAL::Inside quadrature
+DEAL::Outside quadrature
+DEAL::1.00000
+DEAL::Surface quadrature
+DEAL::
 DEAL::test_vertical_cuts_through_center
 DEAL::dim=2
 DEAL::plane direction = 0, face = 2

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.