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;
//@}
*/
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
};
+
// --------- inline and template functions ---------------
template <int dim>
}
+
template <int dim>
inline
unsigned int
}
+
template <int dim>
inline
unsigned int
}
+
template <int dim>
inline
unsigned int
}
+
template <int dim>
inline
unsigned int
}
+
template <int dim>
inline
unsigned int
}
+
template <int dim>
inline
unsigned int
*/
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
*/
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
*/
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.
*
*/
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.
*
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;
//@}
/**
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;
//@}
/**
}
+
template <int dim>
bool
FiniteElement<dim>::prolongation_is_implemented () 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
+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,
+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 <>
+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,
}
+
+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
//---------------------------------------------------------------------------
+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,