]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Implement hp identities for the DGP* element
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Sun, 3 Sep 2006 00:10:40 +0000 (00:10 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Sun, 3 Sep 2006 00:10:40 +0000 (00:10 +0000)
git-svn-id: https://svn.dealii.org/trunk@13807 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/fe/fe_dgp.h
deal.II/deal.II/include/fe/fe_dgp_monomial.h
deal.II/deal.II/include/fe/fe_dgp_nonparametric.h
deal.II/deal.II/source/fe/fe_dgp.cc
deal.II/deal.II/source/fe/fe_dgp_monomial.cc
deal.II/deal.II/source/fe/fe_dgp_nonparametric.cc

index d790e9689df0154dddca27ce01684bdca7b4bf13..1e9385d3ec7a20a45f33ef171b5b619fd4ad6f25 100644 (file)
@@ -101,6 +101,83 @@ class FE_DGP : public FE_Poly<PolynomialSpace<dim>,dim>
                                      */
     virtual std::string get_name () const;
     
+                                    /**
+                                     * @name Functions to support hp 
+                                     * @{
+                                     */
+    
+                                    /**
+                                     * If, on a vertex, several
+                                     * finite elements are active,
+                                     * the hp code first assigns the
+                                     * degrees of freedom of each of
+                                     * these FEs different global
+                                     * indices. It then calls this
+                                     * function to find out which of
+                                     * them should get identical
+                                     * values, and consequently can
+                                     * receive the same global DoF
+                                     * index. This function therefore
+                                     * returns a list of identities
+                                     * between DoFs of the present
+                                     * finite element object with the
+                                     * DoFs of @p fe_other, which is
+                                     * a reference to a finite
+                                     * element object representing
+                                     * one of the other finite
+                                     * elements active on this
+                                     * particular vertex. The
+                                     * function computes which of the
+                                     * degrees of freedom of the two
+                                     * finite element objects are
+                                     * equivalent, and returns a list
+                                     * of pairs of global dof indices
+                                     * in @p identities. The first
+                                     * index of each pair denotes one
+                                     * of the vertex dofs of the
+                                     * present element, whereas the
+                                     * second is the corresponding
+                                     * index of the other finite
+                                     * element.
+                                     *
+                                     * This being a discontinuous element,
+                                     * the set of such constraints is of
+                                     * course empty.
+                                     */
+    virtual
+    std::vector<std::pair<unsigned int, unsigned int> >
+    hp_vertex_dof_identities (const FiniteElement<dim> &fe_other) const;
+
+                                    /**
+                                     * Same as
+                                     * hp_vertex_dof_indices(),
+                                     * except that the function
+                                     * treats degrees of freedom on
+                                     * lines.
+                                     *
+                                     * This being a discontinuous element,
+                                     * the set of such constraints is of
+                                     * course empty.
+                                     */
+    virtual
+    std::vector<std::pair<unsigned int, unsigned int> >
+    hp_line_dof_identities (const FiniteElement<dim> &fe_other) const;
+
+                                    /**
+                                     * Same as
+                                     * hp_vertex_dof_indices(),
+                                     * except that the function
+                                     * treats degrees of freedom on
+                                     * quads.
+                                     *
+                                     * This being a discontinuous element,
+                                     * the set of such constraints is of
+                                     * course empty.
+                                     */
+    virtual
+    std::vector<std::pair<unsigned int, unsigned int> >
+    hp_quad_dof_identities (const FiniteElement<dim> &fe_other) const;
+
                                      /**
                                       * Return whether this element
                                       * implements its hanging node
@@ -130,6 +207,10 @@ class FE_DGP : public FE_Poly<PolynomialSpace<dim>,dim>
     FiniteElementDomination::Domination
     compare_for_domination (const FiniteElement<dim> &fe_other) const;
     
+                                    /**
+                                     * @}
+                                     */
+    
                                     /**
                                      * Return the matrix
                                      * interpolating from a face of
index e27bb6826c536712d1e08d58b250f16398e79d1a..cabf430e631fa7cf255e2bbd31f9b218ea8fcecd 100644 (file)
@@ -101,6 +101,83 @@ class FE_DGPMonomial : public FE_Poly<PolynomialsP<dim>,dim>
                                      */
     virtual std::string get_name () const;
     
+                                    /**
+                                     * @name Functions to support hp 
+                                     * @{
+                                     */
+    
+                                    /**
+                                     * If, on a vertex, several
+                                     * finite elements are active,
+                                     * the hp code first assigns the
+                                     * degrees of freedom of each of
+                                     * these FEs different global
+                                     * indices. It then calls this
+                                     * function to find out which of
+                                     * them should get identical
+                                     * values, and consequently can
+                                     * receive the same global DoF
+                                     * index. This function therefore
+                                     * returns a list of identities
+                                     * between DoFs of the present
+                                     * finite element object with the
+                                     * DoFs of @p fe_other, which is
+                                     * a reference to a finite
+                                     * element object representing
+                                     * one of the other finite
+                                     * elements active on this
+                                     * particular vertex. The
+                                     * function computes which of the
+                                     * degrees of freedom of the two
+                                     * finite element objects are
+                                     * equivalent, and returns a list
+                                     * of pairs of global dof indices
+                                     * in @p identities. The first
+                                     * index of each pair denotes one
+                                     * of the vertex dofs of the
+                                     * present element, whereas the
+                                     * second is the corresponding
+                                     * index of the other finite
+                                     * element.
+                                     *
+                                     * This being a discontinuous element,
+                                     * the set of such constraints is of
+                                     * course empty.
+                                     */
+    virtual
+    std::vector<std::pair<unsigned int, unsigned int> >
+    hp_vertex_dof_identities (const FiniteElement<dim> &fe_other) const;
+
+                                    /**
+                                     * Same as
+                                     * hp_vertex_dof_indices(),
+                                     * except that the function
+                                     * treats degrees of freedom on
+                                     * lines.
+                                     *
+                                     * This being a discontinuous element,
+                                     * the set of such constraints is of
+                                     * course empty.
+                                     */
+    virtual
+    std::vector<std::pair<unsigned int, unsigned int> >
+    hp_line_dof_identities (const FiniteElement<dim> &fe_other) const;
+
+                                    /**
+                                     * Same as
+                                     * hp_vertex_dof_indices(),
+                                     * except that the function
+                                     * treats degrees of freedom on
+                                     * quads.
+                                     *
+                                     * This being a discontinuous element,
+                                     * the set of such constraints is of
+                                     * course empty.
+                                     */
+    virtual
+    std::vector<std::pair<unsigned int, unsigned int> >
+    hp_quad_dof_identities (const FiniteElement<dim> &fe_other) const;
+
                                      /**
                                       * Return whether this element
                                       * implements its hanging node
@@ -131,6 +208,10 @@ class FE_DGPMonomial : public FE_Poly<PolynomialsP<dim>,dim>
     FiniteElementDomination::Domination
     compare_for_domination (const FiniteElement<dim> &fe_other) const;
 
+                                    /**
+                                     * @}
+                                     */
+    
                                     /**
                                      * Return the matrix
                                      * interpolating from the given
index 01a79fb1665ef3c96d79b2cd30b58675630eefa8..260ebc4bea55e25cb8bc05cc24b144f820c7227e 100644 (file)
@@ -256,6 +256,83 @@ class FE_DGPNonparametric : public FiniteElement<dim>
                                      const unsigned int        subface,
                                      FullMatrix<double>       &matrix) const;
 
+                                    /**
+                                     * @name Functions to support hp 
+                                     * @{
+                                     */
+    
+                                    /**
+                                     * If, on a vertex, several
+                                     * finite elements are active,
+                                     * the hp code first assigns the
+                                     * degrees of freedom of each of
+                                     * these FEs different global
+                                     * indices. It then calls this
+                                     * function to find out which of
+                                     * them should get identical
+                                     * values, and consequently can
+                                     * receive the same global DoF
+                                     * index. This function therefore
+                                     * returns a list of identities
+                                     * between DoFs of the present
+                                     * finite element object with the
+                                     * DoFs of @p fe_other, which is
+                                     * a reference to a finite
+                                     * element object representing
+                                     * one of the other finite
+                                     * elements active on this
+                                     * particular vertex. The
+                                     * function computes which of the
+                                     * degrees of freedom of the two
+                                     * finite element objects are
+                                     * equivalent, and returns a list
+                                     * of pairs of global dof indices
+                                     * in @p identities. The first
+                                     * index of each pair denotes one
+                                     * of the vertex dofs of the
+                                     * present element, whereas the
+                                     * second is the corresponding
+                                     * index of the other finite
+                                     * element.
+                                     *
+                                     * This being a discontinuous element,
+                                     * the set of such constraints is of
+                                     * course empty.
+                                     */
+    virtual
+    std::vector<std::pair<unsigned int, unsigned int> >
+    hp_vertex_dof_identities (const FiniteElement<dim> &fe_other) const;
+
+                                    /**
+                                     * Same as
+                                     * hp_vertex_dof_indices(),
+                                     * except that the function
+                                     * treats degrees of freedom on
+                                     * lines.
+                                     *
+                                     * This being a discontinuous element,
+                                     * the set of such constraints is of
+                                     * course empty.
+                                     */
+    virtual
+    std::vector<std::pair<unsigned int, unsigned int> >
+    hp_line_dof_identities (const FiniteElement<dim> &fe_other) const;
+
+                                    /**
+                                     * Same as
+                                     * hp_vertex_dof_indices(),
+                                     * except that the function
+                                     * treats degrees of freedom on
+                                     * quads.
+                                     *
+                                     * This being a discontinuous element,
+                                     * the set of such constraints is of
+                                     * course empty.
+                                     */
+    virtual
+    std::vector<std::pair<unsigned int, unsigned int> >
+    hp_quad_dof_identities (const FiniteElement<dim> &fe_other) const;
+
                                      /**
                                       * Return whether this element
                                       * implements its hanging node
@@ -286,6 +363,10 @@ class FE_DGPNonparametric : public FiniteElement<dim>
     FiniteElementDomination::Domination
     compare_for_domination (const FiniteElement<dim> &fe_other) const;
 
+                                    /**
+                                     * @}
+                                     */
+    
                                     /**
                                      * Check for non-zero values on a face.
                                      *
index 7d264e6734bd8af564d65d722248e064d8571995..bc3945796e35548c91211a631404b67341615f8a 100644 (file)
@@ -155,6 +155,63 @@ FE_DGP<dim>::hp_constraints_are_implemented () const
 
 
 
+template <int dim>
+std::vector<std::pair<unsigned int, unsigned int> >
+FE_DGP<dim>::
+hp_vertex_dof_identities (const FiniteElement<dim> &fe_other) const
+{
+                                  // there are no such constraints for DGP
+                                  // elements at all
+  if (dynamic_cast<const FE_DGP<dim>*>(&fe_other) != 0)
+    return
+      std::vector<std::pair<unsigned int, unsigned int> > ();
+  else
+    {
+      Assert (false, ExcNotImplemented());
+      return std::vector<std::pair<unsigned int, unsigned int> > ();
+    }
+}
+
+
+
+template <int dim>
+std::vector<std::pair<unsigned int, unsigned int> >
+FE_DGP<dim>::
+hp_line_dof_identities (const FiniteElement<dim> &fe_other) const
+{
+                                  // there are no such constraints for DGP
+                                  // elements at all
+  if (dynamic_cast<const FE_DGP<dim>*>(&fe_other) != 0)
+    return
+      std::vector<std::pair<unsigned int, unsigned int> > ();
+  else
+    {
+      Assert (false, ExcNotImplemented());
+      return std::vector<std::pair<unsigned int, unsigned int> > ();
+    }
+}
+
+
+
+template <int dim>
+std::vector<std::pair<unsigned int, unsigned int> >
+FE_DGP<dim>::
+hp_quad_dof_identities (const FiniteElement<dim>        &fe_other) const
+{
+                                  // there are no such constraints for DGP
+                                  // elements at all
+  if (dynamic_cast<const FE_DGP<dim>*>(&fe_other) != 0)
+    return
+      std::vector<std::pair<unsigned int, unsigned int> > ();
+  else
+    {
+      Assert (false, ExcNotImplemented());
+      return std::vector<std::pair<unsigned int, unsigned int> > ();
+    }
+}
+
+
+
 template <int dim>
 FiniteElementDomination::Domination
 FE_DGP<dim>::compare_for_domination (const FiniteElement<dim> &fe_other) const
index b0c3245d719692bb685cc1cc2a7616c2108b4330..c5950e0e163fe64f8f58f8cb3925c856107bf7d2 100644 (file)
@@ -316,6 +316,63 @@ FE_DGPMonomial<dim>::hp_constraints_are_implemented () const
 
 
 
+template <int dim>
+std::vector<std::pair<unsigned int, unsigned int> >
+FE_DGPMonomial<dim>::
+hp_vertex_dof_identities (const FiniteElement<dim> &fe_other) const
+{
+                                  // there are no such constraints for DGPMonomial
+                                  // elements at all
+  if (dynamic_cast<const FE_DGPMonomial<dim>*>(&fe_other) != 0)
+    return
+      std::vector<std::pair<unsigned int, unsigned int> > ();
+  else
+    {
+      Assert (false, ExcNotImplemented());
+      return std::vector<std::pair<unsigned int, unsigned int> > ();
+    }
+}
+
+
+
+template <int dim>
+std::vector<std::pair<unsigned int, unsigned int> >
+FE_DGPMonomial<dim>::
+hp_line_dof_identities (const FiniteElement<dim> &fe_other) const
+{
+                                  // there are no such constraints for DGPMonomial
+                                  // elements at all
+  if (dynamic_cast<const FE_DGPMonomial<dim>*>(&fe_other) != 0)
+    return
+      std::vector<std::pair<unsigned int, unsigned int> > ();
+  else
+    {
+      Assert (false, ExcNotImplemented());
+      return std::vector<std::pair<unsigned int, unsigned int> > ();
+    }
+}
+
+
+
+template <int dim>
+std::vector<std::pair<unsigned int, unsigned int> >
+FE_DGPMonomial<dim>::
+hp_quad_dof_identities (const FiniteElement<dim>        &fe_other) const
+{
+                                  // there are no such constraints for DGPMonomial
+                                  // elements at all
+  if (dynamic_cast<const FE_DGPMonomial<dim>*>(&fe_other) != 0)
+    return
+      std::vector<std::pair<unsigned int, unsigned int> > ();
+  else
+    {
+      Assert (false, ExcNotImplemented());
+      return std::vector<std::pair<unsigned int, unsigned int> > ();
+    }
+}
+
+
+
 template <int dim>
 FiniteElementDomination::Domination
 FE_DGPMonomial<dim>::
index 1d98ec84adced68f104be5ed4292985114ee27dd..4050ac8d2611aef2fcebab55616e4ce7954a290f 100644 (file)
@@ -485,6 +485,63 @@ FE_DGPNonparametric<dim>::hp_constraints_are_implemented () const
 
 
 
+template <int dim>
+std::vector<std::pair<unsigned int, unsigned int> >
+FE_DGPNonparametric<dim>::
+hp_vertex_dof_identities (const FiniteElement<dim> &fe_other) const
+{
+                                  // there are no such constraints for DGPNonparametric
+                                  // elements at all
+  if (dynamic_cast<const FE_DGPNonparametric<dim>*>(&fe_other) != 0)
+    return
+      std::vector<std::pair<unsigned int, unsigned int> > ();
+  else
+    {
+      Assert (false, ExcNotImplemented());
+      return std::vector<std::pair<unsigned int, unsigned int> > ();
+    }
+}
+
+
+
+template <int dim>
+std::vector<std::pair<unsigned int, unsigned int> >
+FE_DGPNonparametric<dim>::
+hp_line_dof_identities (const FiniteElement<dim> &fe_other) const
+{
+                                  // there are no such constraints for DGPNonparametric
+                                  // elements at all
+  if (dynamic_cast<const FE_DGPNonparametric<dim>*>(&fe_other) != 0)
+    return
+      std::vector<std::pair<unsigned int, unsigned int> > ();
+  else
+    {
+      Assert (false, ExcNotImplemented());
+      return std::vector<std::pair<unsigned int, unsigned int> > ();
+    }
+}
+
+
+
+template <int dim>
+std::vector<std::pair<unsigned int, unsigned int> >
+FE_DGPNonparametric<dim>::
+hp_quad_dof_identities (const FiniteElement<dim>        &fe_other) const
+{
+                                  // there are no such constraints for DGPNonparametric
+                                  // elements at all
+  if (dynamic_cast<const FE_DGPNonparametric<dim>*>(&fe_other) != 0)
+    return
+      std::vector<std::pair<unsigned int, unsigned int> > ();
+  else
+    {
+      Assert (false, ExcNotImplemented());
+      return std::vector<std::pair<unsigned int, unsigned int> > ();
+    }
+}
+
+
+
 template <int dim>
 FiniteElementDomination::Domination
 FE_DGPNonparametric<dim>::

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.