]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Convert FETools::compute_face_embedding_matrices() to use ArrayView.
authorWolfgang Bangerth <bangerth@colostate.edu>
Wed, 15 Jan 2025 21:35:04 +0000 (14:35 -0700)
committerWolfgang Bangerth <bangerth@colostate.edu>
Wed, 15 Jan 2025 21:35:04 +0000 (14:35 -0700)
include/deal.II/fe/fe_tools.h
include/deal.II/fe/fe_tools.templates.h
source/fe/fe_bdm.cc
source/fe/fe_raviart_thomas.cc
source/fe/fe_tools.inst.in

index ea447a2392ce1199f1e7bd46ee85873a74bcb6c4..69ed37fe1c4b5a6bbaaa458efea42f5c3aae2a39 100644 (file)
@@ -364,7 +364,7 @@ namespace FETools
    * @param fe The finite element for which to compute these matrices.
    *
    * @param matrices An array of <i>GeometryInfo<dim>::subfaces_per_face =
-   * 2<sup>dim-1</sup></i> FullMatrix objects,holding the embedding matrix for
+   * 2<sup>dim-1</sup></i> FullMatrix objects, holding the embedding matrix for
    * each subface.
    *
    * @param face_coarse The number of the face on the coarse side of the face
@@ -381,12 +381,11 @@ namespace FETools
    */
   template <int dim, typename number, int spacedim>
   void
-  compute_face_embedding_matrices(
-    const FiniteElement<dim, spacedim> &fe,
-    FullMatrix<number> (&matrices)[GeometryInfo<dim>::max_children_per_face],
-    const unsigned int face_coarse,
-    const unsigned int face_fine,
-    const double       threshold = 1.e-12);
+  compute_face_embedding_matrices(const FiniteElement<dim, spacedim>  &fe,
+                                  const ArrayView<FullMatrix<number>> &matrices,
+                                  const unsigned int face_coarse,
+                                  const unsigned int face_fine,
+                                  const double       threshold = 1.e-12);
 
   /**
    * For all possible (isotropic and anisotropic) refinement cases compute the
index e211523b699c8554e4713af440d512c5ce0e046f..b31cbe3020eeb791c244dc138dd5ebd55afb2f7c 100644 (file)
@@ -1811,13 +1811,13 @@ namespace FETools
 
   template <int dim, typename number, int spacedim>
   void
-  compute_face_embedding_matrices(
-    const FiniteElement<dim, spacedim> &fe,
-    FullMatrix<number> (&matrices)[GeometryInfo<dim>::max_children_per_face],
-    const unsigned int face_coarse,
-    const unsigned int face_fine,
-    const double       threshold)
+  compute_face_embedding_matrices(const FiniteElement<dim, spacedim>  &fe,
+                                  const ArrayView<FullMatrix<number>> &matrices,
+                                  const unsigned int face_coarse,
+                                  const unsigned int face_fine,
+                                  const double       threshold)
   {
+    AssertDimension(matrices.size(), GeometryInfo<dim>::max_children_per_face);
     const unsigned int face_no = face_coarse;
 
     Assert(face_coarse == 0, ExcNotImplemented());
index bbf49bbaa4a1b90a7997d203e4b97b7c77e3d210..569042773583b00cd4bc935307f81cf76982a20b 100644 (file)
@@ -90,7 +90,8 @@ FE_BDM<dim>::FE_BDM(const unsigned int deg)
   for (unsigned int i = 0; i < GeometryInfo<dim>::max_children_per_face; ++i)
     face_embeddings[i].reinit(this->n_dofs_per_face(face_no),
                               this->n_dofs_per_face(face_no));
-  FETools::compute_face_embedding_matrices(*this, face_embeddings, 0, 0, 1.);
+  FETools::compute_face_embedding_matrices(
+    *this, make_array_view(face_embeddings), 0, 0, 1.);
   this->interface_constraints.reinit((1 << (dim - 1)) *
                                        this->n_dofs_per_face(face_no),
                                      this->n_dofs_per_face(face_no));
index 163a3623fb88d97b22a504fed1b5416f27280ac7..c031f7ba567d541266668ff4c25732edec3cac57 100644 (file)
@@ -93,13 +93,13 @@ FE_RaviartThomas<dim>::FE_RaviartThomas(const unsigned int deg)
   // TODO[TL]: for anisotropic refinement we will probably need a table of
   // submatrices with an array for each refine case
   FullMatrix<double> face_embeddings[GeometryInfo<dim>::max_children_per_face];
-  for (unsigned int i = 0; i < GeometryInfo<dim>::max_children_per_face; ++i)
-    face_embeddings[i].reinit(this->n_dofs_per_face(face_no),
-                              this->n_dofs_per_face(face_no));
-  FETools::compute_face_embedding_matrices<dim, double>(*this,
-                                                        face_embeddings,
-                                                        0,
-                                                        0);
+  for (auto &face_embedding : face_embeddings)
+    face_embedding.reinit(this->n_dofs_per_face(face_no),
+                          this->n_dofs_per_face(face_no));
+  FETools::compute_face_embedding_matrices(*this,
+                                           make_array_view(face_embeddings),
+                                           0,
+                                           0);
   this->interface_constraints.reinit((1 << (dim - 1)) *
                                        this->n_dofs_per_face(face_no),
                                      this->n_dofs_per_face(face_no));
index 45b549d12512c484c779cad36e815c8e438745cb..7c9e101c04ee93aeaed06e244c6e1386179917ee 100644 (file)
@@ -228,10 +228,9 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS)
       template void
       compute_face_embedding_matrices<deal_II_dimension, double>(
         const FiniteElement<deal_II_dimension> &,
-        FullMatrix<double> (
-            &)[GeometryInfo<deal_II_dimension>::max_children_per_face],
-        unsigned int,
-        unsigned int,
+        const ArrayView<FullMatrix<double>> &,
+        const unsigned int,
+        const unsigned int,
         const double);
 
       template void

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.