]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Implement domination information for at least some of the elements
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 1 Sep 2006 18:28:47 +0000 (18:28 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 1 Sep 2006 18:28:47 +0000 (18:28 +0000)
git-svn-id: https://svn.dealii.org/trunk@13789 0785d39b-7218-0410-832d-ea1e28bc413d

14 files changed:
deal.II/deal.II/include/fe/fe.h
deal.II/deal.II/include/fe/fe_base.h
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/include/fe/fe_dgq.h
deal.II/deal.II/include/fe/fe_q.h
deal.II/deal.II/include/fe/fe_system.h
deal.II/deal.II/source/fe/fe.cc
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
deal.II/deal.II/source/fe/fe_q.cc
deal.II/deal.II/source/fe/fe_system.cc

index 3be9e369f6970f2a81d089ba815c4e046340d2e9..286334e481331ab4f79b2473d5183e1dc57bcef1 100644 (file)
@@ -1059,6 +1059,21 @@ class FiniteElement : public Subscriptor,
     virtual
     std::vector<std::pair<unsigned int, unsigned int> >
     hp_quad_dof_identities (const FiniteElement<dim> &fe_other) const;
+
+                                    /**
+                                     * Return whether this element dominates
+                                     * the one given as argument, whether it
+                                     * is the other way around, whether
+                                     * neither dominates, or if either could
+                                     * dominate.
+                                     *
+                                     * For a definition of domination, see
+                                     * FiniteElementBase::Domination and in
+                                     * particular the hp paper.
+                                     */
+    virtual
+    typename FiniteElementData<dim>::Domination
+    compare_for_domination (const FiniteElement<dim> &fe_other) const;
     
                                     //@}
     
index 349df628d58d42a7b52df0784267f8ac734c5de2..066c71fd526df97bcf07bbc5b437df757c35a1ae 100644 (file)
@@ -155,6 +155,81 @@ class FiniteElementData
                                            */
          H2 = 0x0e
     };
+
+                                    /**
+                                     * An enum that describes the
+                                     * outcome of comparing two elements for
+                                     * mutual domination. If one element
+                                     * dominates another, then the
+                                     * restriction of the space described by
+                                     * the dominated element to a face of the
+                                     * cell is strictly larger than that of
+                                     * the dominating element. For example,
+                                     * in 2-d Q(2) elements dominate Q(4)
+                                     * elements, because the traces of Q(4)
+                                     * elements are quartic polynomials which
+                                     * is a space strictly larger than the
+                                     * quadratic polynomials (the restriction
+                                     * of the Q(2) element). In general, Q(k)
+                                     * dominates Q(k') if $k\le k'$.
+                                     *
+                                     * This enum is used in the
+                                     * FiniteElement::compare_fe_for_domination()
+                                     * function that is used in the context
+                                     * of hp finite element methods when
+                                     * determining what to do at faces where
+                                     * two different finite elements meet
+                                     * (see the hp paper for a more detailed
+                                     * description of the following). In that
+                                     * case, the degrees of freedom of one
+                                     * side need to be constrained to those
+                                     * on the other side. The determination
+                                     * which side is which is based on the
+                                     * outcome of a comparison for mutual
+                                     * domination: the dominated side is
+                                     * constrained to the dominating one.
+                                     *
+                                     * Note that there are situations where
+                                     * neither side dominates. The hp paper
+                                     * lists two case, with the simpler one
+                                     * being that a $Q_2\times Q_1$
+                                     * vector-valued element (i.e. a
+                                     * <code>FESystem(FE_Q(2),1,FE_Q(1),1)</code>)
+                                     * meets a $Q_1\times Q_2$ element: here,
+                                     * for each of the two vector-components,
+                                     * we can define a domination
+                                     * relationship, but it is different for
+                                     * the two components.
+                                     *
+                                     * It is clear that the concept of
+                                     * domination doesn't matter for
+                                     * discontinuous elements. However,
+                                     * discontinuous elements may be part of
+                                     * vector-valued elements and may
+                                     * therefore be compared against each
+                                     * other for domination. They should
+                                     * return
+                                     * <code>either_element_can_dominate</code>
+                                     * in that case. Likewise, when comparing
+                                     * two identical finite elements, they
+                                     * should return this code; the reason is
+                                     * that we can not decide which element
+                                     * will dominate at the time we look at
+                                     * the first component of, for example,
+                                     * two $Q_2\times Q_1$ and $Q_2\times
+                                     * Q_2$ elements, and have to keep our
+                                     * options open until we get to the
+                                     * second base element.
+                                     */
+    enum Domination
+    {
+         this_element_dominates,
+         other_element_dominates,
+         neither_element_dominates,
+         either_element_can_dominate
+    };
+    
+
     
                                     /**
                                      * Number of degrees of freedom on
@@ -472,6 +547,7 @@ class FiniteElementData
 };
 
 
+
 // --------- inline and template functions ---------------
 
 template <int dim>
@@ -483,6 +559,7 @@ FiniteElementData<dim>::n_dofs_per_vertex () const
 }
 
 
+
 template <int dim>
 inline
 unsigned int 
@@ -492,6 +569,7 @@ FiniteElementData<dim>::n_dofs_per_line () const
 }
 
 
+
 template <int dim>
 inline
 unsigned int 
@@ -501,6 +579,7 @@ FiniteElementData<dim>::n_dofs_per_quad () const
 }
 
 
+
 template <int dim>
 inline
 unsigned int 
@@ -510,6 +589,7 @@ FiniteElementData<dim>::n_dofs_per_hex () const
 }
 
 
+
 template <int dim>
 inline
 unsigned int 
@@ -519,6 +599,7 @@ FiniteElementData<dim>::n_dofs_per_face () const
 }
 
 
+
 template <int dim>
 inline
 unsigned int 
@@ -552,6 +633,7 @@ FiniteElementData<dim>::n_dofs_per_object () const
 }
 
 
+
 template <int dim>
 inline
 unsigned int 
index a1ffce4ac8343c7337ccbe8eb109c446bc2397cf..97a9dbd3401f8dca79e28f78037294af024ad32f 100644 (file)
@@ -115,6 +115,21 @@ class FE_DGP : public FE_Poly<PolynomialSpace<dim>,dim>
                                       */
     virtual bool hp_constraints_are_implemented () const;
 
+                                    /**
+                                     * Return whether this element dominates
+                                     * the one given as argument, whether it
+                                     * is the other way around, whether
+                                     * neither dominates, or if either could
+                                     * dominate.
+                                     *
+                                     * For a definition of domination, see
+                                     * FiniteElementBase::Domination and in
+                                     * particular the hp paper.
+                                     */
+    virtual
+    typename FiniteElementData<dim>::Domination
+    compare_for_domination (const FiniteElement<dim> &fe_other) const;
+    
                                     /**
                                      * Return the matrix
                                      * interpolating from a face of
index 9210b66e85b07cae76baf2627bda9af6d79e67df..7d2a27deab7fc37bfbcbcfa091453f4b0081c332 100644 (file)
@@ -116,6 +116,21 @@ class FE_DGPMonomial : public FE_Poly<PolynomialsP<dim>,dim>
                                       */
     virtual bool hp_constraints_are_implemented () const;
 
+                                    /**
+                                     * Return whether this element dominates
+                                     * the one given as argument, whether it
+                                     * is the other way around, whether
+                                     * neither dominates, or if either could
+                                     * dominate.
+                                     *
+                                     * For a definition of domination, see
+                                     * FiniteElementBase::Domination and in
+                                     * particular the hp paper.
+                                     */
+    virtual
+    typename FiniteElementData<dim>::Domination
+    compare_for_domination (const FiniteElement<dim> &fe_other) const;
+
                                     /**
                                      * Return the matrix
                                      * interpolating from the given
index 2c85dbd6d9f7c5489677efe54557e5ccf6c10050..9d6031478f03cdbe3d3edbfe467e4f9217421c56 100644 (file)
@@ -271,6 +271,21 @@ class FE_DGPNonparametric : public FiniteElement<dim>
                                       */
     virtual bool hp_constraints_are_implemented () const;
 
+                                    /**
+                                     * Return whether this element dominates
+                                     * the one given as argument, whether it
+                                     * is the other way around, whether
+                                     * neither dominates, or if either could
+                                     * dominate.
+                                     *
+                                     * For a definition of domination, see
+                                     * FiniteElementBase::Domination and in
+                                     * particular the hp paper.
+                                     */
+    virtual
+    typename FiniteElementData<dim>::Domination
+    compare_for_domination (const FiniteElement<dim> &fe_other) const;
+
                                     /**
                                      * Check for non-zero values on a face.
                                      *
index 57f9e5b6eedd8a6d137ad1e863086b1e9fb3c70c..a9ecfb6dffd53c7abb72dea7a9126f3d99e7d0bf 100644 (file)
@@ -187,6 +187,21 @@ class FE_DGQ : public FE_Poly<TensorProductPolynomials<dim>,dim>
                                       */
     virtual bool hp_constraints_are_implemented () const;
 
+                                    /**
+                                     * Return whether this element dominates
+                                     * the one given as argument, whether it
+                                     * is the other way around, whether
+                                     * neither dominates, or if either could
+                                     * dominate.
+                                     *
+                                     * For a definition of domination, see
+                                     * FiniteElementBase::Domination and in
+                                     * particular the hp paper.
+                                     */
+    virtual
+    typename FiniteElementData<dim>::Domination
+    compare_for_domination (const FiniteElement<dim> &fe_other) const;
+
                                     /**
                                      * Check for non-zero values on a face.
                                      *
index 2b32672c7f0f62e6c77bc5a2505dbf58fd1e6117..c90ee992a3bbe15d9bd4db41165cfd79a4c42054 100644 (file)
@@ -421,6 +421,20 @@ class FE_Q : public FE_Poly<TensorProductPolynomials<dim>,dim>
     std::vector<std::pair<unsigned int, unsigned int> >
     hp_quad_dof_identities (const FiniteElement<dim> &fe_other) const;
     
+                                    /**
+                                     * Return whether this element dominates
+                                     * the one given as argument, whether it
+                                     * is the other way around, whether
+                                     * neither dominates, or if either could
+                                     * dominate.
+                                     *
+                                     * For a definition of domination, see
+                                     * FiniteElementBase::Domination and in
+                                     * particular the hp paper.
+                                     */
+    virtual
+    typename FiniteElementData<dim>::Domination
+    compare_for_domination (const FiniteElement<dim> &fe_other) const;
                                     //@}
 
                                     /**
index 7a6636a8fd37e0431a49c36766de640f477b6191..cb0eb81ac8d7cf91fcfb283c8aa742dbd1454606 100644 (file)
@@ -523,6 +523,20 @@ class FESystem : public FiniteElement<dim>
     std::vector<std::pair<unsigned int, unsigned int> >
     hp_quad_dof_identities (const FiniteElement<dim> &fe_other) const;
     
+                                    /**
+                                     * Return whether this element dominates
+                                     * the one given as argument, whether it
+                                     * is the other way around, whether
+                                     * neither dominates, or if either could
+                                     * dominate.
+                                     *
+                                     * For a definition of domination, see
+                                     * FiniteElementBase::Domination and in
+                                     * particular the hp paper.
+                                     */
+    virtual
+    typename FiniteElementData<dim>::Domination
+    compare_for_domination (const FiniteElement<dim> &fe_other) const;
                                     //@}
     
                                     /**
index 6cc38f3759e564dfe78d03ca7e5d48f4fe268cef..f5f620729aba849872bd4e61a5a9dbdcb3804396 100644 (file)
@@ -332,6 +332,7 @@ FiniteElement<dim>::component_to_block_index (const unsigned int index) const
 }
 
 
+
 template <int dim>
 bool
 FiniteElement<dim>::prolongation_is_implemented () const
@@ -517,6 +518,17 @@ hp_quad_dof_identities (const FiniteElement<dim> &) const
 
 
 
+template <int dim>
+typename FiniteElementData<dim>::Domination
+FiniteElement<dim>::
+compare_for_domination (const FiniteElement<dim> &) const
+{
+  Assert (false, ExcNotImplemented());
+  return FiniteElementData<dim>::neither_element_dominates;
+}
+
+
+
 template <int dim>
 bool
 FiniteElement<dim>::operator == (const FiniteElement<dim> &f) const
index 6736e464e0a0bb5c48bb39e5caee890da98b7c6e..3e99b0892d0e8b56931f57a41c5b21f15eeddc35 100644 (file)
@@ -155,6 +155,23 @@ FE_DGP<dim>::hp_constraints_are_implemented () const
 
 
 
+template <int dim>
+typename FiniteElementData<dim>::Domination
+FE_DGP<dim>::compare_for_domination (const FiniteElement<dim> &fe_other) const
+{
+                                  // check whether both are discontinuous
+                                  // elements and both could dominate, see
+                                  // the description of
+                                  // FiniteElementData<dim>::Domination
+  if (dynamic_cast<const FE_DGP<dim>*>(&fe_other) != 0)
+    return FiniteElementData<dim>::either_element_can_dominate;
+
+  Assert (false, ExcNotImplemented());
+  return FiniteElementData<dim>::neither_element_dominates;
+}
+
+
+
 template <int dim>
 bool
 FE_DGP<dim>::has_support_on_face (const unsigned int,
index 358542579151baa1b3bd1d2579402b926c281bb7..347ea12afc1bfbde14633bde154472902539a6fb 100644 (file)
@@ -316,6 +316,24 @@ FE_DGPMonomial<dim>::hp_constraints_are_implemented () const
 
 
 
+template <int dim>
+typename FiniteElementData<dim>::Domination
+FE_DGPMonomial<dim>::
+compare_for_domination (const FiniteElement<dim> &fe_other) const
+{
+                                  // check whether both are discontinuous
+                                  // elements and both could dominate, see
+                                  // the description of
+                                  // FiniteElementData<dim>::Domination
+  if (dynamic_cast<const FE_DGPMonomial<dim>*>(&fe_other) != 0)
+    return FiniteElementData<dim>::either_element_can_dominate;
+
+  Assert (false, ExcNotImplemented());
+  return FiniteElementData<dim>::neither_element_dominates;
+}
+
+
+
 #if deal_II_dimension == 1
 
 template <>
index 584121ae3547b5f91e9c9e9e0d57505415d685d5..73579bdfa24f9ca54359cecd66ecf9c0a35e0b9c 100644 (file)
@@ -485,6 +485,24 @@ FE_DGPNonparametric<dim>::hp_constraints_are_implemented () const
 
 
 
+template <int dim>
+typename FiniteElementData<dim>::Domination
+FE_DGPNonparametric<dim>::
+compare_for_domination (const FiniteElement<dim> &fe_other) const
+{
+                                  // check whether both are discontinuous
+                                  // elements and both could dominate, see
+                                  // the description of
+                                  // FiniteElementData<dim>::Domination
+  if (dynamic_cast<const FE_DGPNonparametric<dim>*>(&fe_other) != 0)
+    return FiniteElementData<dim>::either_element_can_dominate;
+
+  Assert (false, ExcNotImplemented());
+  return FiniteElementData<dim>::neither_element_dominates;
+}
+
+
+
 template <int dim>
 bool
 FE_DGPNonparametric<dim>::has_support_on_face (const unsigned int,
index 7fc62790ccf2b7c498ac36bc56d89a14affbecb6..d32ff982ac7e12754f8774581ad514428223d078 100644 (file)
@@ -677,6 +677,29 @@ hp_quad_dof_identities (const FiniteElement<dim>        &fe_other) const
 }
 
 
+
+template <int dim>
+typename FiniteElementData<dim>::Domination
+FE_Q<dim>::
+compare_for_domination (const FiniteElement<dim> &fe_other) const
+{
+  if (const FE_Q<dim> *fe_q_other
+      = dynamic_cast<const FE_Q<dim>*>(&fe_other))
+    {
+      if (this->degree < fe_q_other->degree)
+       return FiniteElementData<dim>::this_element_dominates;
+      else if (this->degree < fe_q_other->degree)
+       return FiniteElementData<dim>::either_element_can_dominate;
+      else
+       return FiniteElementData<dim>::other_element_dominates;
+    }
+  
+  Assert (false, ExcNotImplemented());
+  return FiniteElementData<dim>::neither_element_dominates;
+}
+
+
+
 //---------------------------------------------------------------------------
 // Auxiliary functions
 //---------------------------------------------------------------------------
index 96a41e1e31151309922e53bf40d66369a6be5ca5..d9c04abe0eb8f3b2390c1e927b1a6880bc834fa0 100644 (file)
@@ -2272,6 +2272,136 @@ FESystem<dim>::hp_quad_dof_identities (const FiniteElement<dim> &fe_other) const
 
 
 
+template <int dim>
+typename FiniteElementData<dim>::Domination
+FESystem<dim>::
+compare_for_domination (const FiniteElement<dim> &fe_other) const
+{
+                                  // at present all we can do is to compare
+                                  // with other FESystems that have the same
+                                  // number of components and bases
+  if (const FESystem<dim> *fe_sys_other
+      = dynamic_cast<const FESystem<dim>*>(&fe_other))
+    {
+      Assert (this->n_components() == fe_sys_other->n_components(),
+             ExcNotImplemented());
+      Assert (this->n_base_elements() == fe_sys_other->n_base_elements(),
+             ExcNotImplemented());
+
+      typename FiniteElementData<dim>::Domination
+       domination = FiniteElementData<dim>::either_element_can_dominate;
+
+                                      // loop over all base elements and do
+                                      // some sanity checks
+      for (unsigned int b=0; b<this->n_base_elements(); ++b)
+       {
+         Assert (this->base_element(b).n_components() ==
+                 fe_sys_other->base_element(b).n_components(),
+                 ExcNotImplemented());
+         Assert (this->element_multiplicity(b) ==
+                 fe_sys_other->element_multiplicity(b),
+                 ExcNotImplemented());
+
+                                          // for this pair of base elements,
+                                          // check who dominates
+         const typename FiniteElementData<dim>::Domination
+           base_domination
+           = (this->base_element(b)
+              .compare_for_domination (fe_sys_other->base_element(b)));
+
+                                          // now see what that means with
+                                          // regard to the previous state
+         switch (domination)
+           {
+             case FiniteElementData<dim>::either_element_can_dominate:
+             {
+                                                // we haven't made a decision
+                                                // yet, simply copy what this
+                                                // pair of bases have to say
+               domination = base_domination;
+               break;
+             }
+             
+             case FiniteElementData<dim>::this_element_dominates:
+             {
+                                                // the present element
+                                                // previously dominated. this
+                                                // will still be the case if
+                                                // either the present base
+                                                // dominates or if the two
+                                                // bases don't
+                                                // care. otherwise, there is
+                                                // a tie which we will not be
+                                                // able to escape from no
+                                                // matter what the other
+                                                // pairs of bases are going
+                                                // to say
+               switch (base_domination)
+                 {
+                   case FiniteElementData<dim>::either_element_can_dominate:
+                   case FiniteElementData<dim>::this_element_dominates:
+                   {
+                     break;
+                   }
+                   
+                   case FiniteElementData<dim>::other_element_dominates:
+                   case FiniteElementData<dim>::neither_element_dominates:
+                   {
+                     return FiniteElementData<dim>::neither_element_dominates;
+                   }
+
+                   default:
+                                                          // shouldn't get
+                                                          // here
+                         Assert (false, ExcInternalError());
+                 }
+               break;
+             }
+
+             case FiniteElementData<dim>::other_element_dominates:
+             {
+                                                // the  opposite case
+               switch (base_domination)
+                 {
+                   case FiniteElementData<dim>::either_element_can_dominate:
+                   case FiniteElementData<dim>::other_element_dominates:
+                   {
+                     break;
+                   }
+                   
+                   case FiniteElementData<dim>::this_element_dominates:
+                   case FiniteElementData<dim>::neither_element_dominates:
+                   {
+                     return FiniteElementData<dim>::neither_element_dominates;
+                   }
+
+                   default:
+                         Assert (false, ExcInternalError());
+                 }
+               break;
+             }
+             
+             default:
+                   Assert (false, ExcInternalError());
+           }
+       }
+
+                                      // if we've gotten here, then we've
+                                      // either found a winner or either
+                                      // element is fine being dominated
+      Assert (domination !=
+             FiniteElementData<dim>::neither_element_dominates,
+             ExcInternalError());
+
+      return domination;
+    }
+  
+  Assert (false, ExcNotImplemented());
+  return FiniteElementData<dim>::neither_element_dominates;
+}
+
+
+
 template <int dim>
 FiniteElementData<dim>
 FESystem<dim>::multiply_dof_numbers (const FiniteElementData<dim> &fe1,

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.