]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Implement hp identities for the DGQ element
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Sat, 2 Sep 2006 23:25:48 +0000 (23:25 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Sat, 2 Sep 2006 23:25:48 +0000 (23:25 +0000)
git-svn-id: https://svn.dealii.org/trunk@13801 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/fe/fe_dgq.h
deal.II/deal.II/source/fe/fe_dgq.cc

index 2eeca8b1ed116630f43a5fa41779005c09a645c7..112dc39c825baf64da8565f3a62d3c234bbfcfba 100644 (file)
@@ -173,6 +173,83 @@ class FE_DGQ : public FE_Poly<TensorProductPolynomials<dim>,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
@@ -202,6 +279,10 @@ class FE_DGQ : public FE_Poly<TensorProductPolynomials<dim>,dim>
     FiniteElementDomination::Domination
     compare_for_domination (const FiniteElement<dim> &fe_other) const;
 
+                                    /**
+                                     * @}
+                                     */
+    
                                     /**
                                      * Check for non-zero values on a face.
                                      *
index 8dd3b05dedb655f3ee4ff3bdd716b7f96ca71a31..42b433651c32a16d428bb9cd3143cf8b87ca0eb8 100644 (file)
@@ -492,6 +492,63 @@ FE_DGQ<dim>::hp_constraints_are_implemented () const
 
 
 
+template <int dim>
+std::vector<std::pair<unsigned int, unsigned int> >
+FE_DGQ<dim>::
+hp_vertex_dof_identities (const FiniteElement<dim> &fe_other) const
+{
+                                  // there are no such constraints for DGQ
+                                  // elements at all
+  if (dynamic_cast<const FE_DGQ<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_DGQ<dim>::
+hp_line_dof_identities (const FiniteElement<dim> &fe_other) const
+{
+                                  // there are no such constraints for DGQ
+                                  // elements at all
+  if (dynamic_cast<const FE_DGQ<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_DGQ<dim>::
+hp_quad_dof_identities (const FiniteElement<dim>        &fe_other) const
+{
+                                  // there are no such constraints for DGQ
+                                  // elements at all
+  if (dynamic_cast<const FE_DGQ<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_DGQ<dim>::compare_for_domination (const FiniteElement<dim> &fe_other) const

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.