From: wolf Date: Thu, 14 Jan 1999 09:53:52 +0000 (+0000) Subject: Modify the handling of the number of dofs in the finite element classes. This change... X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=aae65f5ee87ddfe890ac8144dee251bda2a8ad19;p=dealii-svn.git Modify the handling of the number of dofs in the finite element classes. This change of scheme is needed for the implementation of composed finite elements, such as systems or mixed elements. git-svn-id: https://svn.dealii.org/trunk@728 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/fe/fe.h b/deal.II/deal.II/include/fe/fe.h index 79726a139d..cf482fe16f 100644 --- a/deal.II/deal.II/include/fe/fe.h +++ b/deal.II/deal.II/include/fe/fe.h @@ -78,12 +78,14 @@ struct 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) {}; + 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 @@ -165,16 +167,14 @@ struct 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) {}; + 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 @@ -224,14 +224,9 @@ struct FiniteElementBase : * 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 &fe_data); /** * Return a readonly reference to the @@ -570,6 +565,28 @@ struct FiniteElementBase : * 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 */ @@ -579,14 +596,7 @@ class FiniteElement : public FiniteElementBase { /** * 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 (dofs_per_vertex, - dofs_per_line, - dofs_per_quad, - n_transform_functions) {}; + FiniteElement (const FiniteElementData &fe_data); /** * Destructor. Only declared to have a diff --git a/deal.II/deal.II/include/fe/fe_lib.criss_cross.h b/deal.II/deal.II/include/fe/fe_lib.criss_cross.h index 6e7d2c4931..1ffb9a07dc 100644 --- a/deal.II/deal.II/include/fe/fe_lib.criss_cross.h +++ b/deal.II/deal.II/include/fe/fe_lib.criss_cross.h @@ -170,6 +170,20 @@ class FECrissCross : public FiniteElement { */ 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 get_fe_data (); + /** * Return the value of the #i#th shape * function at point #p# on the unit cell. diff --git a/deal.II/deal.II/include/fe/fe_lib.dg.h b/deal.II/deal.II/include/fe/fe_lib.dg.h index db1ab4c585..c1e38e6e27 100644 --- a/deal.II/deal.II/include/fe/fe_lib.dg.h +++ b/deal.II/deal.II/include/fe/fe_lib.dg.h @@ -24,7 +24,21 @@ class FEDGConstant : public FELinearMapping { */ 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 get_fe_data (); + + /** * Return the value of the #i#th shape * function at point #p# on the unit cell. */ @@ -130,6 +144,21 @@ class FEDGLinear : public FELinear{ * 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 get_fe_data (); + /** * This function returns an error since the * correct use of the restriction @@ -169,6 +198,21 @@ class FEDGQuadraticSub : public FEQuadraticSub{ * 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 get_fe_data (); + /** * This function returns an error since the * correct use of the restriction @@ -209,6 +253,21 @@ class FEDGCubicSub : public FECubicSub{ * 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 get_fe_data (); + /** * This function returns an error since the * correct use of the restriction @@ -248,6 +307,21 @@ class FEDGQuarticSub : public FEQuarticSub{ * 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 get_fe_data (); + /** * This function returns an error since the * correct use of the restriction diff --git a/deal.II/deal.II/include/fe/fe_lib.lagrange.h b/deal.II/deal.II/include/fe/fe_lib.lagrange.h index fc3a28d624..60b5abac82 100644 --- a/deal.II/deal.II/include/fe/fe_lib.lagrange.h +++ b/deal.II/deal.II/include/fe/fe_lib.lagrange.h @@ -45,6 +45,20 @@ class FELinear : public FELinearMapping { */ 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 get_fe_data (); /** * Return the value of the #i#th shape @@ -144,7 +158,22 @@ class FEQuadraticSub : public FELinearMapping { * #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 get_fe_data (); + /** * Return the value of the #i#th shape * function at point #p# on the unit cell. @@ -246,6 +275,7 @@ class FECubicSub : public FELinearMapping { * Constructor */ FECubicSub (); + protected: /** * Constructor that is called by the @@ -258,8 +288,22 @@ class FECubicSub : public FELinearMapping { * #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 get_fe_data (); /** * Return the value of the #i#th shape @@ -363,6 +407,7 @@ class FEQuarticSub : public FELinearMapping { * Constructor */ FEQuarticSub (); + protected: /** * Constructor that is called by the @@ -375,7 +420,23 @@ class FEQuarticSub : public FELinearMapping { * #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 get_fe_data (); + /** * Return the value of the #i#th shape * function at point #p# on the unit cell. diff --git a/deal.II/deal.II/include/fe/fe_linear_mapping.h b/deal.II/deal.II/include/fe/fe_linear_mapping.h index 0b0af49d04..68250cd016 100644 --- a/deal.II/deal.II/include/fe/fe_linear_mapping.h +++ b/deal.II/deal.II/include/fe/fe_linear_mapping.h @@ -22,15 +22,13 @@ class FELinearMapping : public FiniteElement { 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 (dofs_per_vertex, - dofs_per_line, - dofs_per_quad, - GeometryInfo::vertices_per_cell) {}; + const unsigned int dofs_per_quad=0); /** * Return the value of the #i#th shape @@ -146,6 +144,11 @@ class FELinearMapping : public FiniteElement { const dFMatrix &shape_values_transform, const vector > > &shape_grad_transform, const Boundary &boundary) const; + + /** + * Exception + */ + DeclException0 (ExcInvalidData); }; diff --git a/deal.II/deal.II/source/fe/fe.cc b/deal.II/deal.II/source/fe/fe.cc index 7fe883f092..91a527b74d 100644 --- a/deal.II/deal.II/source/fe/fe.cc +++ b/deal.II/deal.II/source/fe/fe.cc @@ -27,6 +27,28 @@ FiniteElementData<1>::FiniteElementData () : +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) && @@ -52,6 +74,34 @@ FiniteElementData<2>::FiniteElementData () : +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) && @@ -69,16 +119,9 @@ bool FiniteElementData<2>::operator== (const FiniteElementData<2> &f) const { #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::children_per_cell; ++i) { @@ -95,14 +138,8 @@ FiniteElementBase<1>::FiniteElementBase (const unsigned int dofs_per_vertex, #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::children_per_cell; ++i) @@ -162,6 +199,13 @@ bool FiniteElementBase::operator == (const FiniteElementBase &f) const /*------------------------------- FiniteElement ----------------------*/ + +template +FiniteElement::FiniteElement (const FiniteElementData &fe_data) : + FiniteElementBase (fe_data) {}; + + + #if deal_II_dimension == 1 // declare this function to be explicitely specialized before first use diff --git a/deal.II/deal.II/source/fe/fe_lib.criss_cross.cc b/deal.II/deal.II/source/fe/fe_lib.criss_cross.cc index de29c860ed..adb3e7e29b 100644 --- a/deal.II/deal.II/source/fe/fe_lib.criss_cross.cc +++ b/deal.II/deal.II/source/fe/fe_lib.criss_cross.cc @@ -305,15 +305,30 @@ #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()); @@ -456,9 +471,14 @@ void FECrissCross<1>::fill_fe_values (const DoFHandler<1>::cell_iterator &, #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.; @@ -508,6 +528,15 @@ FECrissCross<2>::FECrissCross () : +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, diff --git a/deal.II/deal.II/source/fe/fe_lib.cubic.cc b/deal.II/deal.II/source/fe/fe_lib.cubic.cc index ad15a4653d..0384f217f6 100644 --- a/deal.II/deal.II/source/fe/fe_lib.cubic.cc +++ b/deal.II/deal.II/source/fe/fe_lib.cubic.cc @@ -36,6 +36,14 @@ FECubicSub<1>::FECubicSub (const int) : +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 () { @@ -274,6 +282,14 @@ FECubicSub<2>::FECubicSub (const int) : +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; diff --git a/deal.II/deal.II/source/fe/fe_lib.dg.cc b/deal.II/deal.II/source/fe/fe_lib.dg.cc index bf0d9454d1..0b8eb817c0 100644 --- a/deal.II/deal.II/source/fe/fe_lib.dg.cc +++ b/deal.II/deal.II/source/fe/fe_lib.dg.cc @@ -26,6 +26,104 @@ FEDGQuarticSub::FEDGQuarticSub(): +#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 const dFMatrix & FEDGLinear::restrict (const unsigned int child) const { diff --git a/deal.II/deal.II/source/fe/fe_lib.dg.constant.cc b/deal.II/deal.II/source/fe/fe_lib.dg.constant.cc index e42c6b7a87..aab9dcf5a3 100644 --- a/deal.II/deal.II/source/fe/fe_lib.dg.constant.cc +++ b/deal.II/deal.II/source/fe/fe_lib.dg.constant.cc @@ -42,6 +42,14 @@ FEDGConstant<1>::FEDGConstant () : }; +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> &, @@ -79,9 +87,19 @@ FEDGConstant<2>::FEDGConstant () : 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 + diff --git a/deal.II/deal.II/source/fe/fe_lib.linear.cc b/deal.II/deal.II/source/fe/fe_lib.linear.cc index c0ad349751..374efcce55 100644 --- a/deal.II/deal.II/source/fe/fe_lib.linear.cc +++ b/deal.II/deal.II/source/fe/fe_lib.linear.cc @@ -37,6 +37,14 @@ FELinear<1>::FELinear (const int) : +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 () { @@ -187,6 +195,14 @@ FELinear<2>::FELinear (const int) : +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; diff --git a/deal.II/deal.II/source/fe/fe_lib.quadratic.cc b/deal.II/deal.II/source/fe/fe_lib.quadratic.cc index 89a40fc343..3a982635d8 100644 --- a/deal.II/deal.II/source/fe/fe_lib.quadratic.cc +++ b/deal.II/deal.II/source/fe/fe_lib.quadratic.cc @@ -33,6 +33,14 @@ FEQuadraticSub<1>::FEQuadraticSub (const int) : +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 () { /* @@ -259,6 +267,14 @@ FEQuadraticSub<2>::FEQuadraticSub (const int) : +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 () { /* diff --git a/deal.II/deal.II/source/fe/fe_lib.quartic.cc b/deal.II/deal.II/source/fe/fe_lib.quartic.cc index 1656f834f1..d084d0f618 100644 --- a/deal.II/deal.II/source/fe/fe_lib.quartic.cc +++ b/deal.II/deal.II/source/fe/fe_lib.quartic.cc @@ -31,6 +31,15 @@ FEQuarticSub<1>::FEQuarticSub (const int) : }; + +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; @@ -259,6 +268,14 @@ FEQuarticSub<2>::FEQuarticSub (const int) : +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; diff --git a/deal.II/deal.II/source/fe/fe_linear_mapping.cc b/deal.II/deal.II/source/fe/fe_linear_mapping.cc index eb0dd4c6d5..19eb9d4314 100644 --- a/deal.II/deal.II/source/fe/fe_linear_mapping.cc +++ b/deal.II/deal.II/source/fe/fe_linear_mapping.cc @@ -14,8 +14,23 @@ /*---------------------------- 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 @@ -121,6 +136,18 @@ void FELinearMapping<1>::fill_fe_values (const DoFHandler<1>::cell_iterator &cel #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