]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Add functions to compute vertex dof equivalences
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 27 Jul 2006 00:39:44 +0000 (00:39 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 27 Jul 2006 00:39:44 +0000 (00:39 +0000)
git-svn-id: https://svn.dealii.org/trunk@13439 0785d39b-7218-0410-832d-ea1e28bc413d

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

index e5a792a7032e9ed820f4c07f7e493a9044c18734..7dd6811cdf7a9929e23eae1f27bf1fc02e8dea5a 100644 (file)
@@ -898,6 +898,73 @@ class FiniteElement : public Subscriptor,
     get_interpolation_matrix (const FiniteElement<dim> &source,
                              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.
+                                     */
+    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.
+                                     */
+    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.
+                                     */
+    virtual
+    std::vector<std::pair<unsigned int, unsigned int> >
+    hp_quad_dof_identities (const FiniteElement<dim> &fe_other) const;
+    
+                                    //@}
     
                                     /**
                                      * Comparison operator. We also
index fa44e97be7aff3aa4c54246f633649a2893af9fc..9ebe9dab1d8cd7694bb1907f168b2d7a39bd3b0e 100644 (file)
@@ -289,6 +289,73 @@ class FE_Q : public FE_Poly<TensorProductPolynomials<dim>,dim>
     virtual bool has_support_on_face (const unsigned int shape_index,
                                      const unsigned int face_index) 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.
+                                     */
+    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.
+                                     */
+    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.
+                                     */
+    virtual
+    std::vector<std::pair<unsigned int, unsigned int> >
+    hp_quad_dof_identities (const FiniteElement<dim> &fe_other) const;
+    
+                                    //@}
+
                                     /**
                                      * Determine an estimate for the
                                      * memory consumption (in bytes)
index 59081535efb3eabe68137c249e70f1a54230a658..13b1c8ea87e052ebaaa20a94248ebd4120aeb0b8 100644 (file)
@@ -440,7 +440,39 @@ get_interpolation_matrix (const FiniteElement<dim> &,
                ExcInterpolationNotImplemented());
 }
 
-                                   
+
+
+template <int dim>
+std::vector<std::pair<unsigned int, unsigned int> >
+FiniteElement<dim>::
+hp_vertex_dof_identities (const FiniteElement<dim> &) const
+{
+  Assert (false, ExcNotImplemented());
+  return std::vector<std::pair<unsigned int, unsigned int> > ();
+}
+
+
+
+template <int dim>
+std::vector<std::pair<unsigned int, unsigned int> >
+FiniteElement<dim>::
+hp_line_dof_identities (const FiniteElement<dim> &) const
+{
+  Assert (false, ExcNotImplemented());
+  return std::vector<std::pair<unsigned int, unsigned int> > ();
+}
+
+
+
+template <int dim>
+std::vector<std::pair<unsigned int, unsigned int> >
+FiniteElement<dim>::
+hp_quad_dof_identities (const FiniteElement<dim> &) const
+{
+  Assert (false, ExcNotImplemented());
+  return std::vector<std::pair<unsigned int, unsigned int> > ();
+}
+
 
 
 template <int dim>
index d21246db900dcd2d7aadc8a4c620be7e0c3321ac..3b2fbf8f0a6989944996de34d7e127132b9cb2ce 100644 (file)
@@ -324,6 +324,52 @@ get_interpolation_matrix (const FiniteElement<dim> &x_source_fe,
 }
 
 
+
+template <int dim>
+std::vector<std::pair<unsigned int, unsigned int> >
+FE_Q<dim>::
+hp_vertex_dof_identities (const FiniteElement<dim> &fe_other) const
+{
+                                  // we can presently only compute
+                                  // these identities if both FEs are
+                                  // FE_Qs. in that case, there
+                                  // should be exactly one single DoF
+                                  // of each FE at a vertex, and they
+                                  // should have identical value
+  const FE_Q<dim> *fe_q_other = dynamic_cast<const FE_Q<dim>*>(&fe_other);
+  if (fe_q_other != 0)
+    return
+      std::vector<std::pair<unsigned int, unsigned int> > (1, std::make_pair (0U, 0U));
+  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_Q<dim>::
+hp_line_dof_identities (const FiniteElement<dim> &/*fe_other*/) const
+{
+  Assert (false, ExcNotImplemented());
+  return std::vector<std::pair<unsigned int, unsigned int> > ();
+}
+
+
+
+template <int dim>
+std::vector<std::pair<unsigned int, unsigned int> >
+FE_Q<dim>::
+hp_quad_dof_identities (const FiniteElement<dim>        &/*fe_other*/) const
+{
+  Assert (false, ExcNotImplemented());
+  return std::vector<std::pair<unsigned int, unsigned int> > ();
+}
+
+
 //---------------------------------------------------------------------------
 // Auxiliary functions
 //---------------------------------------------------------------------------

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.