*/
FiniteElementData (const unsigned int dofs_per_vertex,
const unsigned int dofs_per_line,
- const unsigned int n_transform_functions) :
- dofs_per_vertex(dofs_per_vertex),
- dofs_per_line(dofs_per_line),
- dofs_per_face(dofs_per_vertex),
- total_dofs (2*dofs_per_vertex+dofs_per_line),
- n_transform_functions(n_transform_functions) {};
+ const unsigned int n_transform_functions);
+
+ /**
+ * Another frequently used useful
+ * constructor, copying the data from
+ * another object.
+ */
+ FiniteElementData (const FiniteElementData<1> &fe_data);
/**
* Comparison operator. It is not clear to
FiniteElementData (const unsigned int dofs_per_vertex,
const unsigned int dofs_per_line,
const unsigned int dofs_per_quad,
- const unsigned int n_transform_functions) :
- dofs_per_vertex(dofs_per_vertex),
- dofs_per_line(dofs_per_line),
- dofs_per_quad(dofs_per_quad),
- dofs_per_face(2*dofs_per_vertex+
- dofs_per_line),
- total_dofs (4*dofs_per_vertex+
- 4*dofs_per_line+
- dofs_per_quad),
- n_transform_functions (n_transform_functions) {};
+ const unsigned int n_transform_functions);
+
+ /**
+ * Another frequently used useful
+ * constructor, copying the data from
+ * another object.
+ */
+ FiniteElementData (const FiniteElementData<2> &fe_data);
/**
* Comparison operator. It is not clear to
* Construct an object of this type.
* You have to set the
* matrices explicitely after calling
- * this base class' constructor. For
- * #dim==1#, #dofs_per_quad# must be
- * zero.
+ * this base class' constructor.
*/
- FiniteElementBase (const unsigned int dofs_per_vertex,
- const unsigned int dofs_per_line,
- const unsigned int dofs_per_quad,
- const unsigned int n_transform_functions);
+ FiniteElementBase (const FiniteElementData<dim> &fe_data);
/**
* Return a readonly reference to the
* the assumption does not hold, then the distribution has to happen
* with all nodes on the small and the large cells. This is not
* implemented in the #DoFHandler# class as of now.
+ * \end{itemize}
+ *
+ * On a more general note, a concrete finite element class, derived from
+ * the #FiniteElement# class needs to fulfill at least the following criteria:
+ * \begin{itemize}
+ * \item Overload and define all pure virtual functions;
+ * \item If a finite element is to be used as part of a composed finite
+ * element, such as the #FESystem# of the #FEMixed# classes, it has
+ * to provide a \textit{static} function namend #get_fe_data#, which
+ * returns the data which also is stored in the #FiniteElementData# object
+ * at the bottom of the class hierarchy. This function needs to be static,
+ * since its output is needed at the time of construction of the composed
+ * finite element, at which time the subobjects denoting the parts of the
+ * composed element are not yet constructed. You need to make sure that
+ * each element has such a function; the composed elements have no way
+ * to make sure that the function is not inherited from a base class.
+ * You may only use an inherited such function, if the result would be
+ * the same.
+ *
+ * It is also a good idea to use the output of this function to construct
+ * the base class #FiniteElement# of a concrete element.
+ * \end{itemize}
*
* @author Wolfgang Bangerth, 1998
*/
/**
* Constructor
*/
- FiniteElement (const unsigned int dofs_per_vertex,
- const unsigned int dofs_per_line,
- const unsigned int dofs_per_quad,
- const unsigned int n_transform_functions) :
- FiniteElementBase<dim> (dofs_per_vertex,
- dofs_per_line,
- dofs_per_quad,
- n_transform_functions) {};
+ FiniteElement (const FiniteElementData<dim> &fe_data);
/**
* Destructor. Only declared to have a
*/
FECrissCross ();
+ /**
+ * Declare a static function which returns
+ * the number of degrees of freedom per
+ * vertex, line, face, etc. This function
+ * is used by the constructors, but is
+ * mainly needed for the composed finite
+ * elements, which assemble a FE object
+ * from several other FEs. See the
+ * documentation for the #FiniteElement#
+ * class for more information on this
+ * subject.
+ */
+ static const FiniteElementData<dim> get_fe_data ();
+
/**
* Return the value of the #i#th shape
* function at point #p# on the unit cell.
*/
FEDGConstant ();
- /**
+ /**
+ * Declare a static function which returns
+ * the number of degrees of freedom per
+ * vertex, line, face, etc. This function
+ * is used by the constructors, but is
+ * mainly needed for the composed finite
+ * elements, which assemble a FE object
+ * from several other FEs. See the
+ * documentation for the #FiniteElement#
+ * class for more information on this
+ * subject.
+ */
+ static const FiniteElementData<dim> get_fe_data ();
+
+ /**
* Return the value of the #i#th shape
* function at point #p# on the unit cell.
*/
* Constructor
*/
FEDGLinear();
+
+ /**
+ * Declare a static function which returns
+ * the number of degrees of freedom per
+ * vertex, line, face, etc. This function
+ * is used by the constructors, but is
+ * mainly needed for the composed finite
+ * elements, which assemble a FE object
+ * from several other FEs. See the
+ * documentation for the #FiniteElement#
+ * class for more information on this
+ * subject.
+ */
+ static const FiniteElementData<dim> get_fe_data ();
+
/**
* This function returns an error since the
* correct use of the restriction
* Constructor
*/
FEDGQuadraticSub();
+
+ /**
+ * Declare a static function which returns
+ * the number of degrees of freedom per
+ * vertex, line, face, etc. This function
+ * is used by the constructors, but is
+ * mainly needed for the composed finite
+ * elements, which assemble a FE object
+ * from several other FEs. See the
+ * documentation for the #FiniteElement#
+ * class for more information on this
+ * subject.
+ */
+ static const FiniteElementData<dim> get_fe_data ();
+
/**
* This function returns an error since the
* correct use of the restriction
* Constructor
*/
FEDGCubicSub();
+
+ /**
+ * Declare a static function which returns
+ * the number of degrees of freedom per
+ * vertex, line, face, etc. This function
+ * is used by the constructors, but is
+ * mainly needed for the composed finite
+ * elements, which assemble a FE object
+ * from several other FEs. See the
+ * documentation for the #FiniteElement#
+ * class for more information on this
+ * subject.
+ */
+ static const FiniteElementData<dim> get_fe_data ();
+
/**
* This function returns an error since the
* correct use of the restriction
* Constructor
*/
FEDGQuarticSub();
+
+ /**
+ * Declare a static function which returns
+ * the number of degrees of freedom per
+ * vertex, line, face, etc. This function
+ * is used by the constructors, but is
+ * mainly needed for the composed finite
+ * elements, which assemble a FE object
+ * from several other FEs. See the
+ * documentation for the #FiniteElement#
+ * class for more information on this
+ * subject.
+ */
+ static const FiniteElementData<dim> get_fe_data ();
+
/**
* This function returns an error since the
* correct use of the restriction
*/
FELinear (const int);
public:
+
+ /**
+ * Declare a static function which returns
+ * the number of degrees of freedom per
+ * vertex, line, face, etc. This function
+ * is used by the constructors, but is
+ * mainly needed for the composed finite
+ * elements, which assemble a FE object
+ * from several other FEs. See the
+ * documentation for the #FiniteElement#
+ * class for more information on this
+ * subject.
+ */
+ static const FiniteElementData<dim> get_fe_data ();
/**
* Return the value of the #i#th shape
* #FEDGQuadraticSub#.
*/
FEQuadraticSub (const int);
+
public:
+ /**
+ * Declare a static function which returns
+ * the number of degrees of freedom per
+ * vertex, line, face, etc. This function
+ * is used by the constructors, but is
+ * mainly needed for the composed finite
+ * elements, which assemble a FE object
+ * from several other FEs. See the
+ * documentation for the #FiniteElement#
+ * class for more information on this
+ * subject.
+ */
+ static const FiniteElementData<dim> get_fe_data ();
+
/**
* Return the value of the #i#th shape
* function at point #p# on the unit cell.
* Constructor
*/
FECubicSub ();
+
protected:
/**
* Constructor that is called by the
* #FEDGCubicSub#.
*/
FECubicSub (const int);
+
public:
+ /**
+ * Declare a static function which returns
+ * the number of degrees of freedom per
+ * vertex, line, face, etc. This function
+ * is used by the constructors, but is
+ * mainly needed for the composed finite
+ * elements, which assemble a FE object
+ * from several other FEs. See the
+ * documentation for the #FiniteElement#
+ * class for more information on this
+ * subject.
+ */
+ static const FiniteElementData<dim> get_fe_data ();
/**
* Return the value of the #i#th shape
* Constructor
*/
FEQuarticSub ();
+
protected:
/**
* Constructor that is called by the
* #FEDGQuarticSub#.
*/
FEQuarticSub (const int);
+
public:
+
+ /**
+ * Declare a static function which returns
+ * the number of degrees of freedom per
+ * vertex, line, face, etc. This function
+ * is used by the constructors, but is
+ * mainly needed for the composed finite
+ * elements, which assemble a FE object
+ * from several other FEs. See the
+ * documentation for the #FiniteElement#
+ * class for more information on this
+ * subject.
+ */
+ static const FiniteElementData<dim> get_fe_data ();
+
/**
* Return the value of the #i#th shape
* function at point #p# on the unit cell.
public:
/**
* Constructor. Simply passes through
- * its arguments to the base class.
+ * its arguments to the base class. For
+ * one space dimension, #dofs_per_quad#
+ * shall be zero.
*/
FELinearMapping (const unsigned int dofs_per_vertex,
const unsigned int dofs_per_line,
- const unsigned int dofs_per_quad=0) :
- FiniteElement<dim> (dofs_per_vertex,
- dofs_per_line,
- dofs_per_quad,
- GeometryInfo<dim>::vertices_per_cell) {};
+ const unsigned int dofs_per_quad=0);
/**
* Return the value of the #i#th shape
const dFMatrix &shape_values_transform,
const vector<vector<Tensor<1,dim> > > &shape_grad_transform,
const Boundary<dim> &boundary) const;
+
+ /**
+ * Exception
+ */
+ DeclException0 (ExcInvalidData);
};
+FiniteElementData<1>::FiniteElementData (const unsigned int dofs_per_vertex,
+ const unsigned int dofs_per_line,
+ const unsigned int n_transform_functions) :
+ dofs_per_vertex(dofs_per_vertex),
+ dofs_per_line(dofs_per_line),
+ dofs_per_face(dofs_per_vertex),
+ total_dofs (2*dofs_per_vertex+dofs_per_line),
+ n_transform_functions(n_transform_functions)
+{};
+
+
+
+FiniteElementData<1>::FiniteElementData (const FiniteElementData<1> &fe_data) :
+ dofs_per_vertex(fe_data.dofs_per_vertex),
+ dofs_per_line(fe_data.dofs_per_line),
+ dofs_per_face(fe_data.dofs_per_face),
+ total_dofs (fe_data.total_dofs),
+ n_transform_functions(fe_data.n_transform_functions)
+{};
+
+
+
bool FiniteElementData<1>::operator== (const FiniteElementData<1> &f) const {
return ((dofs_per_vertex == f.dofs_per_vertex) &&
(dofs_per_line == f.dofs_per_line) &&
+FiniteElementData<2>::FiniteElementData (const unsigned int dofs_per_vertex,
+ const unsigned int dofs_per_line,
+ const unsigned int dofs_per_quad,
+ const unsigned int n_transform_functions) :
+ dofs_per_vertex(dofs_per_vertex),
+ dofs_per_line(dofs_per_line),
+ dofs_per_quad(dofs_per_quad),
+ dofs_per_face(2*dofs_per_vertex+
+ dofs_per_line),
+ total_dofs (4*dofs_per_vertex+
+ 4*dofs_per_line+
+ dofs_per_quad),
+ n_transform_functions (n_transform_functions)
+{};
+
+
+
+FiniteElementData<2>::FiniteElementData (const FiniteElementData<2> &fe_data) :
+ dofs_per_vertex(fe_data.dofs_per_vertex),
+ dofs_per_line(fe_data.dofs_per_line),
+ dofs_per_quad(fe_data.dofs_per_quad),
+ dofs_per_face(fe_data.dofs_per_face),
+ total_dofs (fe_data.total_dofs),
+ n_transform_functions (fe_data.n_transform_functions)
+{};
+
+
+
bool FiniteElementData<2>::operator== (const FiniteElementData<2> &f) const {
return ((dofs_per_vertex == f.dofs_per_vertex) &&
(dofs_per_line == f.dofs_per_line) &&
#if deal_II_dimension == 1
template <>
-FiniteElementBase<1>::FiniteElementBase (const unsigned int dofs_per_vertex,
- const unsigned int dofs_per_line,
- const unsigned int dofs_per_quad,
- const unsigned int n_transform_funcs) :
- FiniteElementData<1> (dofs_per_vertex,
- dofs_per_line,
- n_transform_funcs)
+FiniteElementBase<1>::FiniteElementBase (const FiniteElementData<1> &fe_data) :
+ FiniteElementData<1> (fe_data)
{
- Assert (dofs_per_quad==0, ExcInternalError());
-
const unsigned int dim=1;
for (unsigned int i=0; i<GeometryInfo<dim>::children_per_cell; ++i)
{
#if deal_II_dimension == 2
template <>
-FiniteElementBase<2>::FiniteElementBase (const unsigned int dofs_per_vertex,
- const unsigned int dofs_per_line,
- const unsigned int dofs_per_quad,
- const unsigned int n_transform_funcs) :
- FiniteElementData<2> (dofs_per_vertex,
- dofs_per_line,
- dofs_per_quad,
- n_transform_funcs)
+FiniteElementBase<2>::FiniteElementBase (const FiniteElementData<2> &fe_data) :
+ FiniteElementData<2> (fe_data)
{
const unsigned int dim=2;
for (unsigned int i=0; i<GeometryInfo<dim>::children_per_cell; ++i)
/*------------------------------- FiniteElement ----------------------*/
+
+template <int dim>
+FiniteElement<dim>::FiniteElement (const FiniteElementData<dim> &fe_data) :
+ FiniteElementBase<dim> (fe_data) {};
+
+
+
#if deal_II_dimension == 1
// declare this function to be explicitely specialized before first use
#if deal_II_dimension == 1
+// declare explicit instantiation
+template <>
+const FiniteElementData<1> FECrissCross<1>::get_fe_data ();
+
+
+
template <>
FECrissCross<1>::FECrissCross () :
- FiniteElement<1> (0,0,0,0)
+ FiniteElement<1> (get_fe_data ())
{
Assert (false, ExcNotUseful());
};
+template <>
+const FiniteElementData<1>
+FECrissCross<1>::get_fe_data () {
+ // return more or less invalid data
+ return FiniteElementData<1> (0,0,0);
+};
+
+
+
template <>
double FECrissCross<1>::shape_value (const unsigned int, const Point<1> &) const {
Assert (false, ExcNotUseful());
#if deal_II_dimension == 2
+// declare explicit instantiation
+template <>
+const FiniteElementData<2> FECrissCross<2>::get_fe_data ();
+
+
template <>
FECrissCross<2>::FECrissCross () :
- FiniteElement<2> (1,0,1,5)
+ FiniteElement<2> (get_fe_data ())
{
interface_constraints(0,0) = 1./2.;
interface_constraints(0,1) = 1./2.;
+template <>
+const FiniteElementData<2>
+FECrissCross<2>::get_fe_data () {
+ return FiniteElementData<2> (1,0,1,5);
+};
+
+
+
+
template <>
inline
double FECrissCross<2>::shape_value (const unsigned int i,
+template <>
+const FiniteElementData<1>
+FECubicSub<1>::get_fe_data () {
+ return FiniteElementData<1> (1, 2, GeometryInfo<1>::vertices_per_cell);
+};
+
+
+
template <>
void FECubicSub<1>::initialize_matrices () {
+template <>
+const FiniteElementData<2>
+FECubicSub<2>::get_fe_data () {
+ return FiniteElementData<2> (1, 2, 4, GeometryInfo<2>::vertices_per_cell);
+};
+
+
+
template <>
void FECubicSub<2>::initialize_matrices () {
prolongation[0](0,0) = 1.0;
+#if deal_II_dimension == 1
+
+template <>
+const FiniteElementData<1>
+FEDGLinear<1>::get_fe_data () {
+ // no dofs at the vertices, all in the
+ // interior of the line
+ return FiniteElementData<1> (0,
+ FELinear<1>::get_fe_data().total_dofs,
+ FELinear<1>::get_fe_data().n_transform_functions);
+};
+
+
+template <>
+const FiniteElementData<1>
+FEDGQuadraticSub<1>::get_fe_data () {
+ // no dofs at the vertices, all in the
+ // interior of the line
+ return FiniteElementData<1> (0,
+ FEQuadraticSub<1>::get_fe_data().total_dofs,
+ FEQuadraticSub<1>::get_fe_data().n_transform_functions);
+};
+
+
+template <>
+const FiniteElementData<1>
+FEDGCubicSub<1>::get_fe_data () {
+ // no dofs at the vertices, all in the
+ // interior of the line
+ return FiniteElementData<1> (0,
+ FECubicSub<1>::get_fe_data().total_dofs,
+ FECubicSub<1>::get_fe_data().n_transform_functions);
+};
+
+
+template <>
+const FiniteElementData<1>
+FEDGQuarticSub<1>::get_fe_data () {
+ // no dofs at the vertices, all in the
+ // interior of the line
+ return FiniteElementData<1> (0,
+ FEQuarticSub<1>::get_fe_data().total_dofs,
+ FEQuarticSub<1>::get_fe_data().n_transform_functions);
+};
+
+#endif
+
+
+
+#if deal_II_dimension == 2
+
+template <>
+const FiniteElementData<2>
+FEDGLinear<2>::get_fe_data () {
+ // no dofs at the vertices or lines, all in the
+ // interior of the line
+ return FiniteElementData<2> (0, 0,
+ FELinear<2>::get_fe_data().total_dofs,
+ FELinear<2>::get_fe_data().n_transform_functions);
+};
+
+
+template <>
+const FiniteElementData<2>
+FEDGQuadraticSub<2>::get_fe_data () {
+ // no dofs at the vertices or lines, all in the
+ // interior of the line
+ return FiniteElementData<2> (0, 0,
+ FEQuadraticSub<2>::get_fe_data().total_dofs,
+ FEQuadraticSub<2>::get_fe_data().n_transform_functions);
+};
+
+
+template <>
+const FiniteElementData<2>
+FEDGCubicSub<2>::get_fe_data () {
+ // no dofs at the vertices or lines, all in the
+ // interior of the line
+ return FiniteElementData<2> (0, 0,
+ FECubicSub<2>::get_fe_data().total_dofs,
+ FECubicSub<2>::get_fe_data().n_transform_functions);
+};
+
+
+template <>
+const FiniteElementData<2>
+FEDGQuarticSub<2>::get_fe_data () {
+ // no dofs at the vertices or lines, all in the
+ // interior of the line
+ return FiniteElementData<2> (0, 0,
+ FEQuarticSub<2>::get_fe_data().total_dofs,
+ FEQuarticSub<2>::get_fe_data().n_transform_functions);
+};
+
+#endif
+
+
+
template <int dim>
const dFMatrix &
FEDGLinear<dim>::restrict (const unsigned int child) const {
};
+template <>
+const FiniteElementData<1>
+FEDGConstant<1>::get_fe_data () {
+ return FiniteElementData<1> (0, 1, GeometryInfo<1>::vertices_per_cell);
+};
+
+
+
template <>
void FEDGConstant<1>::get_face_support_points (const typename DoFHandler<1>::face_iterator &,
const Boundary<1> &,
prolongation[3](0,0) = 1.0;
};
-#endif
+template <>
+const FiniteElementData<2>
+FEDGConstant<2>::get_fe_data () {
+ return FiniteElementData<2> (0, 0, 1, GeometryInfo<2>::vertices_per_cell);
+};
+
+
+
+
+#endif
+
+template <>
+const FiniteElementData<1>
+FELinear<1>::get_fe_data () {
+ return FiniteElementData<1> (1, 0, GeometryInfo<1>::vertices_per_cell);
+};
+
+
+
template <>
void FELinear<1>::initialize_matrices () {
+template <>
+const FiniteElementData<2>
+FELinear<2>::get_fe_data () {
+ return FiniteElementData<2> (1, 0, 0, GeometryInfo<2>::vertices_per_cell);
+};
+
+
+
template <>
void FELinear<2>::initialize_matrices () {
restriction[0](0,0) = 1.0;
+template <>
+const FiniteElementData<1>
+FEQuadraticSub<1>::get_fe_data () {
+ return FiniteElementData<1> (1, 1, GeometryInfo<1>::vertices_per_cell);
+};
+
+
+
template <>
void FEQuadraticSub<1>::initialize_matrices () {
/*
+template <>
+const FiniteElementData<2>
+FEQuadraticSub<2>::get_fe_data () {
+ return FiniteElementData<2> (1, 1, 1, GeometryInfo<2>::vertices_per_cell);
+};
+
+
+
template <>
void FEQuadraticSub<2>::initialize_matrices () {
/*
};
+
+template <>
+const FiniteElementData<1>
+FEQuarticSub<1>::get_fe_data () {
+ return FiniteElementData<1> (1, 3, GeometryInfo<1>::vertices_per_cell);
+};
+
+
+
template <>
void FEQuarticSub<1>::initialize_matrices () {
prolongation[0](0,0) = 1.0;
+template <>
+const FiniteElementData<2>
+FEQuarticSub<2>::get_fe_data () {
+ return FiniteElementData<2> (1, 3, 9, GeometryInfo<2>::vertices_per_cell);
+};
+
+
+
template <>
void FEQuarticSub<2>::initialize_matrices () {
prolongation[0](0,0) = 1.0;
/*---------------------------- FELinearMapping ----------------------------------*/
+
+
#if deal_II_dimension == 1
+template <>
+FELinearMapping<1>::FELinearMapping (const unsigned int dofs_per_vertex,
+ const unsigned int dofs_per_line,
+ const unsigned int dofs_per_quad) :
+ FiniteElement<1> (FiniteElementData<1> (dofs_per_vertex,
+ dofs_per_line,
+ GeometryInfo<1>::vertices_per_cell))
+{
+ Assert (dofs_per_quad==0, ExcInvalidData());
+};
+
+
+
template <>
inline
double
#if deal_II_dimension == 2
+template <>
+FELinearMapping<2>::FELinearMapping (const unsigned int dofs_per_vertex,
+ const unsigned int dofs_per_line,
+ const unsigned int dofs_per_quad) :
+ FiniteElement<2> (FiniteElementData<2> (dofs_per_vertex,
+ dofs_per_line,
+ dofs_per_quad,
+ GeometryInfo<2>::vertices_per_cell))
+{};
+
+
+
template <>
inline
double