]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Implement FE_WedgeDGP and FE_PyramidDGP 11405/head
authorPeter Munch <peterrmuench@gmail.com>
Sun, 20 Dec 2020 09:57:04 +0000 (10:57 +0100)
committerPeter Munch <peterrmuench@gmail.com>
Tue, 22 Dec 2020 18:37:19 +0000 (19:37 +0100)
include/deal.II/fe/fe_base.h
include/deal.II/simplex/fe_lib.h
source/fe/fe_data.cc
source/simplex/fe_lib.cc
source/simplex/fe_lib.inst.in
tests/simplex/fe_lib_01.cc
tests/simplex/fe_lib_01.output

index ac3deeac3e0e7084c8f349b8f47e8897a8d67daa..803776665e1cba01b5f2ab731913270b4fdc9d6e 100644 (file)
@@ -673,6 +673,18 @@ public:
   get_first_face_quad_index(const unsigned int face_no = 0) const;
 };
 
+namespace internal
+{
+  /**
+   * Utility function to convert "dofs per object" information
+   * of a @p dim dimensional reference cell @p cell_type.
+   */
+  internal::GenericDoFsPerObject
+  expand(const unsigned int               dim,
+         const std::vector<unsigned int> &dofs_per_object,
+         const ReferenceCell::Type        cell_type);
+} // namespace internal
+
 
 
 // --------- inline and template functions ---------------
index aa833316b9798bd3bc09e62b9eaf8f7caffef140..0d5fe5cc8e14e5c9e2046eb30e217c8eb65b0dae 100644 (file)
@@ -169,6 +169,24 @@ namespace Simplex
       const FiniteElement<dim, spacedim> &fe_other) const override;
   };
 
+  /**
+   * Base class of FE_WedgeP and FE_WedgeDGP.
+   *
+   * @note Only implemented for 3D.
+   *
+   * @ingroup simplex
+   */
+  template <int dim, int spacedim = dim>
+  class FE_Wedge : public dealii::FE_Poly<dim, spacedim>
+  {
+  public:
+    /**
+     * Constructor.
+     */
+    FE_Wedge(const unsigned int                                degree,
+             const internal::GenericDoFsPerObject &            dpos,
+             const typename FiniteElementData<dim>::Conformity conformity);
+  };
 
   /**
    * Implementation of a scalar Lagrange finite element on a wedge that yields
@@ -178,7 +196,7 @@ namespace Simplex
    * @ingroup simplex
    */
   template <int dim, int spacedim = dim>
-  class FE_WedgeP : public dealii::FE_Poly<dim, spacedim>
+  class FE_WedgeP : public FE_Wedge<dim, spacedim>
   {
   public:
     /**
@@ -229,6 +247,55 @@ namespace Simplex
                            const unsigned int face_no = 0) const override;
   };
 
+  /**
+   * Implementation of a scalar Lagrange finite element on a wedge that yields
+   * the finite element space of discontinuous, piecewise polynomials of
+   * degree p.
+   *
+   * @ingroup simplex
+   */
+  template <int dim, int spacedim = dim>
+  class FE_WedgeDGP : public FE_Wedge<dim, spacedim>
+  {
+  public:
+    /**
+     * Constructor.
+     */
+    FE_WedgeDGP(const unsigned int degree);
+
+    /**
+     * @copydoc dealii::FiniteElement::clone()
+     */
+    std::unique_ptr<FiniteElement<dim, spacedim>>
+    clone() const override;
+
+    /**
+     * Return a string that uniquely identifies a finite element. This class
+     * returns <tt>Simplex::FE_WedgeDGP<dim>(degree)</tt>, with @p dim and @p degree
+     * replaced by appropriate values.
+     */
+    std::string
+    get_name() const override;
+  };
+
+  /**
+   * Base class of FE_PyramidP and FE_PyramidDGP.
+   *
+   * @note Only implemented for 3D.
+   *
+   * @ingroup simplex
+   */
+  template <int dim, int spacedim = dim>
+  class FE_Pyramid : public dealii::FE_Poly<dim, spacedim>
+  {
+  public:
+    /**
+     * Constructor.
+     */
+    FE_Pyramid(const unsigned int                                degree,
+               const internal::GenericDoFsPerObject &            dpos,
+               const typename FiniteElementData<dim>::Conformity conformity);
+  };
 
   /**
    * Implementation of a scalar Lagrange finite element on a pyramid that yields
@@ -238,7 +305,7 @@ namespace Simplex
    * @ingroup simplex
    */
   template <int dim, int spacedim = dim>
-  class FE_PyramidP : public dealii::FE_Poly<dim, spacedim>
+  class FE_PyramidP : public FE_Pyramid<dim, spacedim>
   {
   public:
     /**
@@ -289,6 +356,37 @@ namespace Simplex
                            const unsigned int face_no = 0) const override;
   };
 
+  /**
+   * Implementation of a scalar Lagrange finite element on a pyramid that yields
+   * the finite element space of discontinuous, piecewise polynomials of
+   * degree p.
+   *
+   * @ingroup simplex
+   */
+  template <int dim, int spacedim = dim>
+  class FE_PyramidDGP : public FE_Pyramid<dim, spacedim>
+  {
+  public:
+    /**
+     * Constructor.
+     */
+    FE_PyramidDGP(const unsigned int degree);
+
+    /**
+     * @copydoc dealii::FiniteElement::clone()
+     */
+    std::unique_ptr<FiniteElement<dim, spacedim>>
+    clone() const override;
+
+    /**
+     * Return a string that uniquely identifies a finite element. This class
+     * returns <tt>Simplex::FE_PyramidDGP<dim>(degree)</tt>, with @p dim and @p degree
+     * replaced by appropriate values.
+     */
+    std::string
+    get_name() const override;
+  };
+
 
 } // namespace Simplex
 
index 497f398444454b56103874880c5eb824fe00bb91..333baa7f4d5ea111ffad193ddb1c161f1cf40b37 100644 (file)
@@ -19,7 +19,7 @@
 
 DEAL_II_NAMESPACE_OPEN
 
-namespace
+namespace internal
 {
   internal::GenericDoFsPerObject
   expand(const unsigned int               dim,
@@ -118,7 +118,7 @@ namespace
 
     return result;
   }
-} // namespace
+} // namespace internal
 
 template <int dim>
 FiniteElementData<dim>::FiniteElementData(
@@ -147,7 +147,7 @@ FiniteElementData<dim>::FiniteElementData(
   const unsigned int               degree,
   const Conformity                 conformity,
   const BlockIndices &             block_indices)
-  : FiniteElementData(expand(dim, dofs_per_object, cell_type),
+  : FiniteElementData(internal::expand(dim, dofs_per_object, cell_type),
                       cell_type,
                       n_components,
                       degree,
index 0b1831d2bd2c94e5e256cf4a5bd193da27d0d59e..91dbf2e260061ac5e7685437064849bc5a7124c5 100644 (file)
@@ -85,11 +85,10 @@ namespace Simplex
      * Helper function to set up the dpo vector of FE_WedgeP for a given @p degree.
      */
     internal::GenericDoFsPerObject
-    get_dpo_vector_fe_wedge(const unsigned int degree)
+    get_dpo_vector_fe_wedge_p(const unsigned int degree)
     {
       internal::GenericDoFsPerObject dpo;
 
-      // all dofs are internal
       if (degree == 1)
         {
           dpo.dofs_per_object_exclusive  = {{1}, {0}, {0, 0, 0, 0, 0}, {0}};
@@ -117,14 +116,33 @@ namespace Simplex
     }
 
     /**
-     * Helper function to set up the dpo vector of FE_WedgeP for a given @p degree.
+     * Helper function to set up the dpo vector of FE_WedgeDGP for a given @p degree.
+     */
+    internal::GenericDoFsPerObject
+    get_dpo_vector_fe_wedge_dgp(const unsigned int degree)
+    {
+      unsigned int n_dofs = 0;
+
+      if (degree == 1)
+        n_dofs = 6;
+      else if (degree == 2)
+        n_dofs = 18;
+      else
+        Assert(false, ExcNotImplemented());
+
+      return internal::expand(3,
+                              {{0, 0, 0, n_dofs}},
+                              ReferenceCell::Type::Wedge);
+    }
+
+    /**
+     * Helper function to set up the dpo vector of FE_PyramidP for a given @p degree.
      */
     internal::GenericDoFsPerObject
-    get_dpo_vector_fe_pyramid(const unsigned int degree)
+    get_dpo_vector_fe_pyramid_p(const unsigned int degree)
     {
       internal::GenericDoFsPerObject dpo;
 
-      // all dofs are internal
       if (degree == 1)
         {
           dpo.dofs_per_object_exclusive  = {{1}, {0}, {0, 0, 0, 0, 0}, {0}};
@@ -141,6 +159,24 @@ namespace Simplex
 
       return dpo;
     }
+
+    /**
+     * Helper function to set up the dpo vector of FE_PyramidDGP for a given @p degree.
+     */
+    internal::GenericDoFsPerObject
+    get_dpo_vector_fe_pyramid_dgp(const unsigned int degree)
+    {
+      unsigned int n_dofs = 0;
+
+      if (degree == 1)
+        n_dofs = 5;
+      else
+        Assert(false, ExcNotImplemented());
+
+      return internal::expand(3,
+                              {{0, 0, 0, n_dofs}},
+                              ReferenceCell::Type::Pyramid);
+    }
   } // namespace
 
 
@@ -451,32 +487,28 @@ namespace Simplex
 
 
   template <int dim, int spacedim>
-  FE_WedgeP<dim, spacedim>::FE_WedgeP(const unsigned int degree)
+  FE_Wedge<dim, spacedim>::FE_Wedge(
+    const unsigned int                                degree,
+    const internal::GenericDoFsPerObject &            dpos,
+    const typename FiniteElementData<dim>::Conformity conformity)
     : dealii::FE_Poly<dim, spacedim>(
         Simplex::ScalarWedgePolynomial<dim>(degree),
-        FiniteElementData<dim>(get_dpo_vector_fe_wedge(degree),
+        FiniteElementData<dim>(dpos,
                                ReferenceCell::Type::Wedge,
                                1,
                                degree,
-                               FiniteElementData<dim>::H1),
-        std::vector<bool>(FiniteElementData<dim>(get_dpo_vector_fe_wedge(
-                                                   degree),
-                                                 ReferenceCell::Type::Wedge,
-                                                 1,
-                                                 degree)
-                            .dofs_per_cell,
-                          true),
+                               conformity),
+        std::vector<bool>(
+          FiniteElementData<dim>(dpos, ReferenceCell::Type::Wedge, 1, degree)
+            .dofs_per_cell,
+          true),
         std::vector<ComponentMask>(
-          FiniteElementData<dim>(get_dpo_vector_fe_wedge(degree),
-                                 ReferenceCell::Type::Wedge,
-                                 1,
-                                 degree)
+          FiniteElementData<dim>(dpos, ReferenceCell::Type::Wedge, 1, degree)
             .dofs_per_cell,
           std::vector<bool>(1, true)))
   {
     AssertDimension(dim, 3);
 
-
     if (degree == 1)
       {
         this->unit_support_points.emplace_back(1.0, 0.0, 0.0);
@@ -485,16 +517,20 @@ namespace Simplex
         this->unit_support_points.emplace_back(1.0, 0.0, 1.0);
         this->unit_support_points.emplace_back(0.0, 1.0, 1.0);
         this->unit_support_points.emplace_back(0.0, 0.0, 1.0);
-
-        // TODO
-        this->unit_face_support_points[0].emplace_back(1.0, 0.0);
-        this->unit_face_support_points[0].emplace_back(0.0, 1.0);
-        this->unit_face_support_points[0].emplace_back(0.0, 0.0);
       }
   }
 
 
 
+  template <int dim, int spacedim>
+  FE_WedgeP<dim, spacedim>::FE_WedgeP(const unsigned int degree)
+    : FE_Wedge<dim, spacedim>(degree,
+                              get_dpo_vector_fe_wedge_p(degree),
+                              FiniteElementData<dim>::H1)
+  {}
+
+
+
   template <int dim, int spacedim>
   std::unique_ptr<FiniteElement<dim, spacedim>>
   FE_WedgeP<dim, spacedim>::clone() const
@@ -579,7 +615,6 @@ namespace Simplex
   {
     (void)fe_other;
 
-
     AssertIndexRange(face_no, 5);
 
     if (face_no < 2)
@@ -604,26 +639,53 @@ namespace Simplex
 
 
   template <int dim, int spacedim>
-  FE_PyramidP<dim, spacedim>::FE_PyramidP(const unsigned int degree)
+  FE_WedgeDGP<dim, spacedim>::FE_WedgeDGP(const unsigned int degree)
+    : FE_Wedge<dim, spacedim>(degree,
+                              get_dpo_vector_fe_wedge_dgp(degree),
+                              FiniteElementData<dim>::L2)
+  {}
+
+
+
+  template <int dim, int spacedim>
+  std::unique_ptr<FiniteElement<dim, spacedim>>
+  FE_WedgeDGP<dim, spacedim>::clone() const
+  {
+    return std::make_unique<FE_WedgeDGP<dim, spacedim>>(*this);
+  }
+
+
+
+  template <int dim, int spacedim>
+  std::string
+  FE_WedgeDGP<dim, spacedim>::get_name() const
+  {
+    std::ostringstream namebuf;
+    namebuf << "FE_WedgeDGP<" << dim << ">(" << this->degree << ")";
+
+    return namebuf.str();
+  }
+
+
+
+  template <int dim, int spacedim>
+  FE_Pyramid<dim, spacedim>::FE_Pyramid(
+    const unsigned int                                degree,
+    const internal::GenericDoFsPerObject &            dpos,
+    const typename FiniteElementData<dim>::Conformity conformity)
     : dealii::FE_Poly<dim, spacedim>(
         Simplex::ScalarPyramidPolynomial<dim>(degree),
-        FiniteElementData<dim>(get_dpo_vector_fe_pyramid(degree),
+        FiniteElementData<dim>(dpos,
                                ReferenceCell::Type::Pyramid,
                                1,
                                degree,
-                               FiniteElementData<dim>::H1),
-        std::vector<bool>(FiniteElementData<dim>(get_dpo_vector_fe_pyramid(
-                                                   degree),
-                                                 ReferenceCell::Type::Pyramid,
-                                                 1,
-                                                 degree)
-                            .dofs_per_cell,
-                          true),
+                               conformity),
+        std::vector<bool>(
+          FiniteElementData<dim>(dpos, ReferenceCell::Type::Pyramid, 1, degree)
+            .dofs_per_cell,
+          true),
         std::vector<ComponentMask>(
-          FiniteElementData<dim>(get_dpo_vector_fe_pyramid(degree),
-                                 ReferenceCell::Type::Pyramid,
-                                 1,
-                                 degree)
+          FiniteElementData<dim>(dpos, ReferenceCell::Type::Pyramid, 1, degree)
             .dofs_per_cell,
           std::vector<bool>(1, true)))
   {
@@ -642,6 +704,15 @@ namespace Simplex
 
 
 
+  template <int dim, int spacedim>
+  FE_PyramidP<dim, spacedim>::FE_PyramidP(const unsigned int degree)
+    : FE_Pyramid<dim, spacedim>(degree,
+                                get_dpo_vector_fe_pyramid_p(degree),
+                                FiniteElementData<dim>::H1)
+  {}
+
+
+
   template <int dim, int spacedim>
   std::unique_ptr<FiniteElement<dim, spacedim>>
   FE_PyramidP<dim, spacedim>::clone() const
@@ -750,6 +821,36 @@ namespace Simplex
 
 
 
+  template <int dim, int spacedim>
+  FE_PyramidDGP<dim, spacedim>::FE_PyramidDGP(const unsigned int degree)
+    : FE_Pyramid<dim, spacedim>(degree,
+                                get_dpo_vector_fe_pyramid_dgp(degree),
+                                FiniteElementData<dim>::L2)
+  {}
+
+
+
+  template <int dim, int spacedim>
+  std::unique_ptr<FiniteElement<dim, spacedim>>
+  FE_PyramidDGP<dim, spacedim>::clone() const
+  {
+    return std::make_unique<FE_PyramidDGP<dim, spacedim>>(*this);
+  }
+
+
+
+  template <int dim, int spacedim>
+  std::string
+  FE_PyramidDGP<dim, spacedim>::get_name() const
+  {
+    std::ostringstream namebuf;
+    namebuf << "FE_PyramidDGP<" << dim << ">(" << this->degree << ")";
+
+    return namebuf.str();
+  }
+
+
+
 } // namespace Simplex
 
 // explicit instantiations
index 36b7279cd37aca0f5f98382620fd58fb203abcc7..be8135b313a764400684f95e78eb0db5868294b8 100644 (file)
@@ -21,9 +21,17 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : DIMENSIONS)
     template class Simplex::FE_Poly<deal_II_dimension, deal_II_space_dimension>;
     template class Simplex::FE_P<deal_II_dimension, deal_II_space_dimension>;
     template class Simplex::FE_DGP<deal_II_dimension, deal_II_space_dimension>;
+    template class Simplex::FE_Wedge<deal_II_dimension,
+                                     deal_II_space_dimension>;
     template class Simplex::FE_WedgeP<deal_II_dimension,
                                       deal_II_space_dimension>;
+    template class Simplex::FE_WedgeDGP<deal_II_dimension,
+                                        deal_II_space_dimension>;
+    template class Simplex::FE_Pyramid<deal_II_dimension,
+                                       deal_II_space_dimension>;
     template class Simplex::FE_PyramidP<deal_II_dimension,
                                         deal_II_space_dimension>;
+    template class Simplex::FE_PyramidDGP<deal_II_dimension,
+                                          deal_II_space_dimension>;
 #endif
   }
index e8f01a75bcf00ac738db96a3201f8d0f844b14c5..104d4d60dd7326e7224853ed877321f116a28ec4 100644 (file)
@@ -29,11 +29,24 @@ test(const FiniteElement<dim, spacedim> &fe)
 {
   deallog << fe.get_name() << ": " << std::endl;
 
+  const auto &reference_cell =
+    ReferenceCell::internal::Info::get_cell(fe.reference_cell_type());
+
   deallog << "  n_dofs_per_vertex(): " << fe.n_dofs_per_vertex() << std::endl;
   deallog << "  n_dofs_per_line():   " << fe.n_dofs_per_line() << std::endl;
-  deallog << "  n_dofs_per_quad():   " << fe.n_dofs_per_quad() << std::endl;
+
+  deallog << "  n_dofs_per_quad():   ";
+  for (unsigned int i = 0; i < (dim == 2 ? 1 : reference_cell.n_faces()); ++i)
+    deallog << fe.n_dofs_per_quad(i) << " ";
+  deallog << std::endl;
+
   deallog << "  n_dofs_per_hex():    " << fe.n_dofs_per_hex() << std::endl;
-  deallog << "  n_dofs_per_face():   " << fe.n_dofs_per_face() << std::endl;
+
+  deallog << "  n_dofs_per_face():   ";
+  for (unsigned int i = 0; i < reference_cell.n_faces(); ++i)
+    deallog << fe.n_dofs_per_face(i) << " ";
+  deallog << std::endl;
+
   deallog << "  n_dofs_per_cell():   " << fe.n_dofs_per_cell() << std::endl;
   deallog << "  tensor_degree():     " << fe.tensor_degree() << std::endl;
 
@@ -49,8 +62,19 @@ main()
   test(Simplex::FE_P<2>(2));
   test(Simplex::FE_P<3>(1));
   test(Simplex::FE_P<3>(2));
+
   test(Simplex::FE_DGP<2>(1));
   test(Simplex::FE_DGP<2>(2));
   test(Simplex::FE_DGP<3>(1));
   test(Simplex::FE_DGP<3>(2));
+
+  test(Simplex::FE_WedgeP<3>(1));
+  test(Simplex::FE_WedgeP<3>(2));
+
+  test(Simplex::FE_WedgeDGP<3>(1));
+  test(Simplex::FE_WedgeDGP<3>(2));
+
+  test(Simplex::FE_PyramidP<3>(1));
+
+  test(Simplex::FE_PyramidDGP<3>(1));
 }
index 4410d62f694c0321d122d920f64d33a09110f16d..ec3b0b86fdee1ff680b21fb3e59a94e7e0166658 100644 (file)
 DEAL::FE_P<2>(1): 
 DEAL::  n_dofs_per_vertex(): 1
 DEAL::  n_dofs_per_line():   0
-DEAL::  n_dofs_per_quad():   0
+DEAL::  n_dofs_per_quad():   0 
 DEAL::  n_dofs_per_hex():    0
-DEAL::  n_dofs_per_face():   2
+DEAL::  n_dofs_per_face():   2 2 2 
 DEAL::  n_dofs_per_cell():   3
 DEAL::  tensor_degree():     1
 DEAL::
 DEAL::FE_P<2>(2): 
 DEAL::  n_dofs_per_vertex(): 1
 DEAL::  n_dofs_per_line():   1
-DEAL::  n_dofs_per_quad():   0
+DEAL::  n_dofs_per_quad():   0 
 DEAL::  n_dofs_per_hex():    0
-DEAL::  n_dofs_per_face():   3
+DEAL::  n_dofs_per_face():   3 3 3 
 DEAL::  n_dofs_per_cell():   6
 DEAL::  tensor_degree():     2
 DEAL::
 DEAL::FE_P<3>(1): 
 DEAL::  n_dofs_per_vertex(): 1
 DEAL::  n_dofs_per_line():   0
-DEAL::  n_dofs_per_quad():   0
+DEAL::  n_dofs_per_quad():   0 0 0 0 
 DEAL::  n_dofs_per_hex():    0
-DEAL::  n_dofs_per_face():   3
+DEAL::  n_dofs_per_face():   3 3 3 3 
 DEAL::  n_dofs_per_cell():   4
 DEAL::  tensor_degree():     1
 DEAL::
 DEAL::FE_P<3>(2): 
 DEAL::  n_dofs_per_vertex(): 1
 DEAL::  n_dofs_per_line():   1
-DEAL::  n_dofs_per_quad():   0
+DEAL::  n_dofs_per_quad():   0 0 0 0 
 DEAL::  n_dofs_per_hex():    0
-DEAL::  n_dofs_per_face():   6
+DEAL::  n_dofs_per_face():   6 6 6 6 
 DEAL::  n_dofs_per_cell():   10
 DEAL::  tensor_degree():     2
 DEAL::
 DEAL::FE_DGP<2>(1): 
 DEAL::  n_dofs_per_vertex(): 0
 DEAL::  n_dofs_per_line():   0
-DEAL::  n_dofs_per_quad():   3
+DEAL::  n_dofs_per_quad():   3 
 DEAL::  n_dofs_per_hex():    0
-DEAL::  n_dofs_per_face():   0
+DEAL::  n_dofs_per_face():   0 0 0 
 DEAL::  n_dofs_per_cell():   3
 DEAL::  tensor_degree():     1
 DEAL::
 DEAL::FE_DGP<2>(2): 
 DEAL::  n_dofs_per_vertex(): 0
 DEAL::  n_dofs_per_line():   0
-DEAL::  n_dofs_per_quad():   6
+DEAL::  n_dofs_per_quad():   6 
 DEAL::  n_dofs_per_hex():    0
-DEAL::  n_dofs_per_face():   0
+DEAL::  n_dofs_per_face():   0 0 0 
 DEAL::  n_dofs_per_cell():   6
 DEAL::  tensor_degree():     2
 DEAL::
 DEAL::FE_DGP<3>(1): 
 DEAL::  n_dofs_per_vertex(): 0
 DEAL::  n_dofs_per_line():   0
-DEAL::  n_dofs_per_quad():   0
+DEAL::  n_dofs_per_quad():   0 0 0 0 
 DEAL::  n_dofs_per_hex():    4
-DEAL::  n_dofs_per_face():   0
+DEAL::  n_dofs_per_face():   0 0 0 0 
 DEAL::  n_dofs_per_cell():   4
 DEAL::  tensor_degree():     1
 DEAL::
 DEAL::FE_DGP<3>(2): 
 DEAL::  n_dofs_per_vertex(): 0
 DEAL::  n_dofs_per_line():   0
-DEAL::  n_dofs_per_quad():   0
+DEAL::  n_dofs_per_quad():   0 0 0 0 
 DEAL::  n_dofs_per_hex():    10
-DEAL::  n_dofs_per_face():   0
+DEAL::  n_dofs_per_face():   0 0 0 0 
 DEAL::  n_dofs_per_cell():   10
 DEAL::  tensor_degree():     2
 DEAL::
+DEAL::FE_WedgeP<3>(1): 
+DEAL::  n_dofs_per_vertex(): 1
+DEAL::  n_dofs_per_line():   0
+DEAL::  n_dofs_per_quad():   0 0 0 0 0 
+DEAL::  n_dofs_per_hex():    0
+DEAL::  n_dofs_per_face():   3 3 4 4 4 
+DEAL::  n_dofs_per_cell():   6
+DEAL::  tensor_degree():     1
+DEAL::
+DEAL::FE_WedgeP<3>(2): 
+DEAL::  n_dofs_per_vertex(): 1
+DEAL::  n_dofs_per_line():   1
+DEAL::  n_dofs_per_quad():   0 0 1 1 1 
+DEAL::  n_dofs_per_hex():    0
+DEAL::  n_dofs_per_face():   6 6 9 9 9 
+DEAL::  n_dofs_per_cell():   18
+DEAL::  tensor_degree():     2
+DEAL::
+DEAL::FE_WedgeDGP<3>(1): 
+DEAL::  n_dofs_per_vertex(): 0
+DEAL::  n_dofs_per_line():   0
+DEAL::  n_dofs_per_quad():   0 0 0 0 0 
+DEAL::  n_dofs_per_hex():    6
+DEAL::  n_dofs_per_face():   0 0 0 0 0 
+DEAL::  n_dofs_per_cell():   6
+DEAL::  tensor_degree():     1
+DEAL::
+DEAL::FE_WedgeDGP<3>(2): 
+DEAL::  n_dofs_per_vertex(): 0
+DEAL::  n_dofs_per_line():   0
+DEAL::  n_dofs_per_quad():   0 0 0 0 0 
+DEAL::  n_dofs_per_hex():    18
+DEAL::  n_dofs_per_face():   0 0 0 0 0 
+DEAL::  n_dofs_per_cell():   18
+DEAL::  tensor_degree():     2
+DEAL::
+DEAL::FE_PyramidP<3>(1): 
+DEAL::  n_dofs_per_vertex(): 1
+DEAL::  n_dofs_per_line():   0
+DEAL::  n_dofs_per_quad():   0 0 0 0 0 
+DEAL::  n_dofs_per_hex():    0
+DEAL::  n_dofs_per_face():   4 3 3 3 3 
+DEAL::  n_dofs_per_cell():   5
+DEAL::  tensor_degree():     1
+DEAL::
+DEAL::FE_PyramidDGP<3>(1): 
+DEAL::  n_dofs_per_vertex(): 0
+DEAL::  n_dofs_per_line():   0
+DEAL::  n_dofs_per_quad():   0 0 0 0 0 
+DEAL::  n_dofs_per_hex():    5
+DEAL::  n_dofs_per_face():   0 0 0 0 0 
+DEAL::  n_dofs_per_cell():   5
+DEAL::  tensor_degree():     1
+DEAL::

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.