// $Id$
// Version: $Name$
//
-// Copyright (C) 2002, 2003, 2004, 2005 by the deal.II authors
+// Copyright (C) 2002, 2003, 2004, 2005, 2006 by the deal.II authors
//
// This file is subject to QPL and may not be distributed
// without copyright and license information. Please refer
*/
virtual std::string get_name () const;
+ /**
+ * Return whether this element
+ * implements its hanging node
+ * constraints in the new way,
+ * which has to be used to make
+ * elements "hp compatible".
+ *
+ * For the FE_DGP class the result is
+ * always true (independent of the degree
+ * of the element), as it has no hanging
+ * nodes (being a discontinuous element).
+ */
+ virtual bool hp_constraints_are_implemented () const;
+
+ /**
+ * Return the matrix
+ * interpolating from a face of
+ * of one element to the face of
+ * the neighboring element.
+ * The size of the matrix is
+ * then @p dofs_per_face times
+ * <tt>source.dofs_per_face</tt>.
+ *
+ * Derived elements will have to
+ * implement this function. They
+ * may only provide interpolation
+ * matrices for certain source
+ * finite elements, for example
+ * those from the same family. If
+ * they don't implement
+ * interpolation from a given
+ * element, then they must throw
+ * an exception of type
+ * FiniteElement<dim>::ExcInterpolationNotImplemented.
+ */
+ virtual void
+ get_face_interpolation_matrix (const FiniteElement<dim> &source,
+ FullMatrix<double> &matrix) const;
+
/**
+ * Return the matrix
+ * interpolating from a face of
+ * of one element to the face of
+ * the neighboring element.
+ * The size of the matrix is
+ * then @p dofs_per_face times
+ * <tt>source.dofs_per_face</tt>.
+ *
+ * Derived elements will have to
+ * implement this function. They
+ * may only provide interpolation
+ * matrices for certain source
+ * finite elements, for example
+ * those from the same family. If
+ * they don't implement
+ * interpolation from a given
+ * element, then they must throw
+ * an exception of type
+ * FiniteElement<dim>::ExcInterpolationNotImplemented.
+ */
+ virtual void
+ get_subface_interpolation_matrix (const FiniteElement<dim> &source,
+ const unsigned int subface,
+ FullMatrix<double> &matrix) const;
+
+ /**
* Check for non-zero values on a face.
*
* This function returns
// $Id$
// Version: $Name$
//
-// Copyright (C) 2004, 2005 by the deal.II authors
+// Copyright (C) 2004, 2005, 2006 by the deal.II authors
//
// This file is subject to QPL and may not be distributed
// without copyright and license information. Please refer
*/
virtual std::string get_name () const;
+ /**
+ * Return whether this element
+ * implements its hanging node
+ * constraints in the new way,
+ * which has to be used to make
+ * elements "hp compatible".
+ *
+ * For the FE_DGPMonomial class the
+ * result is always true (independent of
+ * the degree of the element), as it has
+ * no hanging nodes (being a
+ * discontinuous element).
+ */
+ virtual bool hp_constraints_are_implemented () const;
+
/**
* Return the matrix
* interpolating from the given
get_interpolation_matrix (const FiniteElement<dim> &source,
FullMatrix<double> &matrix) const;
+ /**
+ * Return the matrix
+ * interpolating from a face of
+ * of one element to the face of
+ * the neighboring element.
+ * The size of the matrix is
+ * then @p dofs_per_face times
+ * <tt>source.dofs_per_face</tt>.
+ *
+ * Derived elements will have to
+ * implement this function. They
+ * may only provide interpolation
+ * matrices for certain source
+ * finite elements, for example
+ * those from the same family. If
+ * they don't implement
+ * interpolation from a given
+ * element, then they must throw
+ * an exception of type
+ * FiniteElement<dim>::ExcInterpolationNotImplemented.
+ */
+ virtual void
+ get_face_interpolation_matrix (const FiniteElement<dim> &source,
+ FullMatrix<double> &matrix) const;
+
+ /**
+ * Return the matrix
+ * interpolating from a face of
+ * of one element to the face of
+ * the neighboring element.
+ * The size of the matrix is
+ * then @p dofs_per_face times
+ * <tt>source.dofs_per_face</tt>.
+ *
+ * Derived elements will have to
+ * implement this function. They
+ * may only provide interpolation
+ * matrices for certain source
+ * finite elements, for example
+ * those from the same family. If
+ * they don't implement
+ * interpolation from a given
+ * element, then they must throw
+ * an exception of type
+ * FiniteElement<dim>::ExcInterpolationNotImplemented.
+ */
+ virtual void
+ get_subface_interpolation_matrix (const FiniteElement<dim> &source,
+ const unsigned int subface,
+ FullMatrix<double> &matrix) const;
+
/**
* Check for non-zero values on a face.
*
// $Id$
// Version: $Name$
//
-// Copyright (C) 2002, 2003, 2004, 2005 by the deal.II authors
+// Copyright (C) 2002, 2003, 2004, 2005, 2006 by the deal.II authors
//
// This file is subject to QPL and may not be distributed
// without copyright and license information. Please refer
*/
virtual unsigned int element_multiplicity (const unsigned int index) const;
+ /**
+ * Return the matrix
+ * interpolating from a face of
+ * of one element to the face of
+ * the neighboring element.
+ * The size of the matrix is
+ * then @p dofs_per_face times
+ * <tt>source.dofs_per_face</tt>.
+ *
+ * Derived elements will have to
+ * implement this function. They
+ * may only provide interpolation
+ * matrices for certain source
+ * finite elements, for example
+ * those from the same family. If
+ * they don't implement
+ * interpolation from a given
+ * element, then they must throw
+ * an exception of type
+ * FiniteElement<dim>::ExcInterpolationNotImplemented.
+ */
+ virtual void
+ get_face_interpolation_matrix (const FiniteElement<dim> &source,
+ FullMatrix<double> &matrix) const;
+
+ /**
+ * Return the matrix
+ * interpolating from a face of
+ * of one element to the face of
+ * the neighboring element.
+ * The size of the matrix is
+ * then @p dofs_per_face times
+ * <tt>source.dofs_per_face</tt>.
+ *
+ * Derived elements will have to
+ * implement this function. They
+ * may only provide interpolation
+ * matrices for certain source
+ * finite elements, for example
+ * those from the same family. If
+ * they don't implement
+ * interpolation from a given
+ * element, then they must throw
+ * an exception of type
+ * FiniteElement<dim>::ExcInterpolationNotImplemented.
+ */
+ virtual void
+ get_subface_interpolation_matrix (const FiniteElement<dim> &source,
+ const unsigned int subface,
+ FullMatrix<double> &matrix) const;
+
+ /**
+ * Return whether this element
+ * implements its hanging node
+ * constraints in the new way,
+ * which has to be used to make
+ * elements "hp compatible".
+ *
+ * For the FE_DGPNonparametric class the
+ * result is always true (independent of
+ * the degree of the element), as it has
+ * no hanging nodes (being a
+ * discontinuous element).
+ */
+ virtual bool hp_constraints_are_implemented () const;
+
/**
* Check for non-zero values on a face.
*
virtual unsigned int memory_consumption () const;
+ private:
/**
* Declare a nested class which
* will hold static definitions of
// $Id$
// Version: $Name$
//
-// Copyright (C) 2001, 2002, 2003, 2004, 2005 by the deal.II authors
+// Copyright (C) 2001, 2002, 2003, 2004, 2005, 2006 by the deal.II authors
//
// This file is subject to QPL and may not be distributed
// without copyright and license information. Please refer
get_interpolation_matrix (const FiniteElement<dim> &source,
FullMatrix<double> &matrix) const;
+ /**
+ * Return the matrix
+ * interpolating from a face of
+ * of one element to the face of
+ * the neighboring element.
+ * The size of the matrix is
+ * then @p dofs_per_face times
+ * <tt>source.dofs_per_face</tt>.
+ *
+ * Derived elements will have to
+ * implement this function. They
+ * may only provide interpolation
+ * matrices for certain source
+ * finite elements, for example
+ * those from the same family. If
+ * they don't implement
+ * interpolation from a given
+ * element, then they must throw
+ * an exception of type
+ * FiniteElement<dim>::ExcInterpolationNotImplemented.
+ */
+ virtual void
+ get_face_interpolation_matrix (const FiniteElement<dim> &source,
+ FullMatrix<double> &matrix) const;
+
+ /**
+ * Return the matrix
+ * interpolating from a face of
+ * of one element to the face of
+ * the neighboring element.
+ * The size of the matrix is
+ * then @p dofs_per_face times
+ * <tt>source.dofs_per_face</tt>.
+ *
+ * Derived elements will have to
+ * implement this function. They
+ * may only provide interpolation
+ * matrices for certain source
+ * finite elements, for example
+ * those from the same family. If
+ * they don't implement
+ * interpolation from a given
+ * element, then they must throw
+ * an exception of type
+ * FiniteElement<dim>::ExcInterpolationNotImplemented.
+ */
+ virtual void
+ get_subface_interpolation_matrix (const FiniteElement<dim> &source,
+ const unsigned int subface,
+ FullMatrix<double> &matrix) const;
+
+ /**
+ * Return whether this element
+ * implements its hanging node
+ * constraints in the new way,
+ * which has to be used to make
+ * elements "hp compatible".
+ *
+ * For the FE_DGQ class the result is
+ * always true (independent of the degree
+ * of the element), as it has no hanging
+ * nodes (being a discontinuous element).
+ */
+ virtual bool hp_constraints_are_implemented () const;
+
/**
* Check for non-zero values on a face.
*
}
+
+template <int dim>
+void
+FE_DGP<dim>::
+get_face_interpolation_matrix (const FiniteElement<dim> &x_source_fe,
+ FullMatrix<double> &interpolation_matrix) const
+{
+ // this is only implemented, if the source
+ // FE is also a DGP element. in that case,
+ // both elements have no dofs on their
+ // faces and the face interpolation matrix
+ // is necessarily empty -- i.e. there isn't
+ // much we need to do here.
+ AssertThrow ((x_source_fe.get_name().find ("FE_DGP<") == 0)
+ ||
+ (dynamic_cast<const FE_DGP<dim>*>(&x_source_fe) != 0),
+ typename FiniteElement<dim>::
+ ExcInterpolationNotImplemented());
+
+ Assert (interpolation_matrix.m() == 0,
+ ExcDimensionMismatch (interpolation_matrix.m(),
+ 0));
+ Assert (interpolation_matrix.n() == 0,
+ ExcDimensionMismatch (interpolation_matrix.m(),
+ 0));
+}
+
+
+
+template <int dim>
+void
+FE_DGP<dim>::
+get_subface_interpolation_matrix (const FiniteElement<dim> &x_source_fe,
+ const unsigned int ,
+ FullMatrix<double> &interpolation_matrix) const
+{
+ // this is only implemented, if the source
+ // FE is also a DGP element. in that case,
+ // both elements have no dofs on their
+ // faces and the face interpolation matrix
+ // is necessarily empty -- i.e. there isn't
+ // much we need to do here.
+ AssertThrow ((x_source_fe.get_name().find ("FE_DGP<") == 0)
+ ||
+ (dynamic_cast<const FE_DGP<dim>*>(&x_source_fe) != 0),
+ typename FiniteElement<dim>::
+ ExcInterpolationNotImplemented());
+
+ Assert (interpolation_matrix.m() == 0,
+ ExcDimensionMismatch (interpolation_matrix.m(),
+ 0));
+ Assert (interpolation_matrix.n() == 0,
+ ExcDimensionMismatch (interpolation_matrix.m(),
+ 0));
+}
+
+
+
+template <int dim>
+bool
+FE_DGP<dim>::hp_constraints_are_implemented () const
+{
+ return true;
+}
+
+
+
template <int dim>
bool
FE_DGP<dim>::has_support_on_face (const unsigned int,
}
+template <int dim>
+void
+FE_DGPMonomial<dim>::
+get_face_interpolation_matrix (const FiniteElement<dim> &x_source_fe,
+ FullMatrix<double> &interpolation_matrix) const
+{
+ // this is only implemented, if the source
+ // FE is also a DGPMonomial element. in that case,
+ // both elements have no dofs on their
+ // faces and the face interpolation matrix
+ // is necessarily empty -- i.e. there isn't
+ // much we need to do here.
+ AssertThrow ((x_source_fe.get_name().find ("FE_DGPMonomial<") == 0)
+ ||
+ (dynamic_cast<const FE_DGPMonomial<dim>*>(&x_source_fe) != 0),
+ typename FiniteElement<dim>::
+ ExcInterpolationNotImplemented());
+
+ Assert (interpolation_matrix.m() == 0,
+ ExcDimensionMismatch (interpolation_matrix.m(),
+ 0));
+ Assert (interpolation_matrix.n() == 0,
+ ExcDimensionMismatch (interpolation_matrix.m(),
+ 0));
+}
+
+
+
+template <int dim>
+void
+FE_DGPMonomial<dim>::
+get_subface_interpolation_matrix (const FiniteElement<dim> &x_source_fe,
+ const unsigned int ,
+ FullMatrix<double> &interpolation_matrix) const
+{
+ // this is only implemented, if the source
+ // FE is also a DGPMonomial element. in that case,
+ // both elements have no dofs on their
+ // faces and the face interpolation matrix
+ // is necessarily empty -- i.e. there isn't
+ // much we need to do here.
+ AssertThrow ((x_source_fe.get_name().find ("FE_DGPMonomial<") == 0)
+ ||
+ (dynamic_cast<const FE_DGPMonomial<dim>*>(&x_source_fe) != 0),
+ typename FiniteElement<dim>::
+ ExcInterpolationNotImplemented());
+
+ Assert (interpolation_matrix.m() == 0,
+ ExcDimensionMismatch (interpolation_matrix.m(),
+ 0));
+ Assert (interpolation_matrix.n() == 0,
+ ExcDimensionMismatch (interpolation_matrix.m(),
+ 0));
+}
+
+
+
+template <int dim>
+bool
+FE_DGPMonomial<dim>::hp_constraints_are_implemented () const
+{
+ return true;
+}
+
+
+
#if deal_II_dimension == 1
template <>
+template <int dim>
+void
+FE_DGPNonparametric<dim>::
+get_face_interpolation_matrix (const FiniteElement<dim> &x_source_fe,
+ FullMatrix<double> &interpolation_matrix) const
+{
+ // this is only implemented, if the source
+ // FE is also a DGPNonparametric element. in that case,
+ // both elements have no dofs on their
+ // faces and the face interpolation matrix
+ // is necessarily empty -- i.e. there isn't
+ // much we need to do here.
+ AssertThrow ((x_source_fe.get_name().find ("FE_DGPNonparametric<") == 0)
+ ||
+ (dynamic_cast<const FE_DGPNonparametric<dim>*>(&x_source_fe) != 0),
+ typename FiniteElement<dim>::
+ ExcInterpolationNotImplemented());
+
+ Assert (interpolation_matrix.m() == 0,
+ ExcDimensionMismatch (interpolation_matrix.m(),
+ 0));
+ Assert (interpolation_matrix.n() == 0,
+ ExcDimensionMismatch (interpolation_matrix.m(),
+ 0));
+}
+
+
+
+template <int dim>
+void
+FE_DGPNonparametric<dim>::
+get_subface_interpolation_matrix (const FiniteElement<dim> &x_source_fe,
+ const unsigned int ,
+ FullMatrix<double> &interpolation_matrix) const
+{
+ // this is only implemented, if the source
+ // FE is also a DGPNonparametric element. in that case,
+ // both elements have no dofs on their
+ // faces and the face interpolation matrix
+ // is necessarily empty -- i.e. there isn't
+ // much we need to do here.
+ AssertThrow ((x_source_fe.get_name().find ("FE_DGPNonparametric<") == 0)
+ ||
+ (dynamic_cast<const FE_DGPNonparametric<dim>*>(&x_source_fe) != 0),
+ typename FiniteElement<dim>::
+ ExcInterpolationNotImplemented());
+
+ Assert (interpolation_matrix.m() == 0,
+ ExcDimensionMismatch (interpolation_matrix.m(),
+ 0));
+ Assert (interpolation_matrix.n() == 0,
+ ExcDimensionMismatch (interpolation_matrix.m(),
+ 0));
+}
+
+
+
+template <int dim>
+bool
+FE_DGPNonparametric<dim>::hp_constraints_are_implemented () const
+{
+ return true;
+}
+
+
+
template <int dim>
bool
FE_DGPNonparametric<dim>::has_support_on_face (const unsigned int,
-//---------------------------------------------------------------------------
-// Fill data of FEValues
-//---------------------------------------------------------------------------
+template <int dim>
+void
+FE_DGQ<dim>::
+get_face_interpolation_matrix (const FiniteElement<dim> &x_source_fe,
+ FullMatrix<double> &interpolation_matrix) const
+{
+ // this is only implemented, if the source
+ // FE is also a DGQ element. in that case,
+ // both elements have no dofs on their
+ // faces and the face interpolation matrix
+ // is necessarily empty -- i.e. there isn't
+ // much we need to do here.
+ AssertThrow ((x_source_fe.get_name().find ("FE_DGQ<") == 0)
+ ||
+ (dynamic_cast<const FE_DGQ<dim>*>(&x_source_fe) != 0),
+ typename FiniteElement<dim>::
+ ExcInterpolationNotImplemented());
+
+ Assert (interpolation_matrix.m() == 0,
+ ExcDimensionMismatch (interpolation_matrix.m(),
+ 0));
+ Assert (interpolation_matrix.n() == 0,
+ ExcDimensionMismatch (interpolation_matrix.m(),
+ 0));
+}
+
+
+
+template <int dim>
+void
+FE_DGQ<dim>::
+get_subface_interpolation_matrix (const FiniteElement<dim> &x_source_fe,
+ const unsigned int ,
+ FullMatrix<double> &interpolation_matrix) const
+{
+ // this is only implemented, if the source
+ // FE is also a DGQ element. in that case,
+ // both elements have no dofs on their
+ // faces and the face interpolation matrix
+ // is necessarily empty -- i.e. there isn't
+ // much we need to do here.
+ AssertThrow ((x_source_fe.get_name().find ("FE_DGQ<") == 0)
+ ||
+ (dynamic_cast<const FE_DGQ<dim>*>(&x_source_fe) != 0),
+ typename FiniteElement<dim>::
+ ExcInterpolationNotImplemented());
+
+ Assert (interpolation_matrix.m() == 0,
+ ExcDimensionMismatch (interpolation_matrix.m(),
+ 0));
+ Assert (interpolation_matrix.n() == 0,
+ ExcDimensionMismatch (interpolation_matrix.m(),
+ 0));
+}
+
+
+
+template <int dim>
+bool
+FE_DGQ<dim>::hp_constraints_are_implemented () const
+{
+ return true;
+}
+
+
template <int dim>
bool