From: hartmann Date: Tue, 16 Mar 2004 16:39:39 +0000 (+0000) Subject: Make use of the new FE_Poly class. Move some of the renumbering business to the Tenso... X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=b38e3aecf78ac801f94bb008681bacdcabb3e2ad;p=dealii-svn.git Make use of the new FE_Poly class. Move some of the renumbering business to the TensorProductPolynomials class. git-svn-id: https://svn.dealii.org/trunk@8778 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/fe/fe_q_hierarchical.h b/deal.II/deal.II/include/fe/fe_q_hierarchical.h index 540fb2a363..2ddad0c47c 100644 --- a/deal.II/deal.II/include/fe/fe_q_hierarchical.h +++ b/deal.II/deal.II/include/fe/fe_q_hierarchical.h @@ -2,7 +2,7 @@ // $Id$ // Version: $Name$ // -// Copyright (C) 2002, 2003 by the deal.II authors +// Copyright (C) 2002, 2003, 2004 by the deal.II authors // // This file is subject to QPL and may not be distributed // without copyright and license information. Please refer @@ -14,10 +14,10 @@ #define __deal2__fe_q_hierarchical_h #include -#include #include +#include #include -#include +//#include template class TensorProductPolynomials; template class MappingQ; @@ -231,15 +231,15 @@ template class MappingQ; * Note the reverse ordering of degrees of freedom on the left and upper * line. * - * @author Brian Carnes, 2002 + * @author Brian Carnes, 2002, Ralf Hartmann 2004 */ template -class FE_Q_Hierarchical : public FiniteElement +class FE_Q_Hierarchical : public FE_Poly,dim> { public: /** * Constructor for tensor product - * polynomials of degree @p{p}. + * polynomials of degree @p p. */ FE_Q_Hierarchical (const unsigned int p); @@ -254,106 +254,6 @@ class FE_Q_Hierarchical : public FiniteElement */ virtual std::string get_name () const; - /** - * Return the value of the - * @p{i}th shape function at the - * point @p{p}. @p{p} is a point - * on the reference element. - */ - virtual double shape_value (const unsigned int i, - const Point &p) const; - - /** - * Return the value of the - * @p{component}th vector - * component of the @p{i}th shape - * function at the point - * @p{p}. See the - * @ref{FiniteElementBase} base - * class for more information - * about the semantics of this - * function. - * - * Since this element is scalar, - * the returned value is the same - * as if the function without the - * @p{_component} suffix were - * called, provided that the - * specified component is zero. - */ - virtual double shape_value_component (const unsigned int i, - const Point &p, - const unsigned int component) const; - - /** - * Return the gradient of the - * @p{i}th shape function at the - * point @p{p}. @p{p} is a point - * on the reference element, and - * likewise the gradient is the - * gradient on the unit cell with - * respect to unit cell - * coordinates. - */ - virtual Tensor<1,dim> shape_grad (const unsigned int i, - const Point &p) const; - - /** - * Return the gradient of the - * @p{component}th vector - * component of the @p{i}th shape - * function at the point - * @p{p}. See the - * @ref{FiniteElementBase} base - * class for more information - * about the semantics of this - * function. - * - * Since this element is scalar, - * the returned value is the same - * as if the function without the - * @p{_component} suffix were - * called, provided that the - * specified component is zero. - */ - virtual Tensor<1,dim> shape_grad_component (const unsigned int i, - const Point &p, - const unsigned int component) const; - - /** - * Return the tensor of second - * derivatives of the @p{i}th - * shape function at point @p{p} - * on the unit cell. The - * derivatives are derivatives on - * the unit cell with respect to - * unit cell coordinates. - */ - virtual Tensor<2,dim> shape_grad_grad (const unsigned int i, - const Point &p) const; - - /** - * Return the second derivative - * of the @p{component}th vector - * component of the @p{i}th shape - * function at the point - * @p{p}. See the - * @ref{FiniteElementBase} base - * class for more information - * about the semantics of this - * function. - * - * Since this element is scalar, - * the returned value is the same - * as if the function without the - * @p{_component} suffix were - * called, provided that the - * specified component is zero. - */ - virtual Tensor<2,dim> shape_grad_grad_component (const unsigned int i, - const Point &p, - const unsigned int component) const; - /** * Return the polynomial degree * of this finite element, @@ -362,34 +262,6 @@ class FE_Q_Hierarchical : public FiniteElement */ unsigned int get_degree () const; - /** - * Number of base elements in a - * mixed discretization. Since - * this is a scalar element, - * return one. - */ - virtual unsigned int n_base_elements () const; - - /** - * Access to base element - * objects. Since this element is - * scalar, @p{base_element(0)} is - * @p{this}, and all other - * indices throw an error. - */ - virtual const FiniteElement & - base_element (const unsigned int index) const; - - /** - * Multiplicity of base element - * @p{index}. Since this is a - * scalar element, - * @p{element_multiplicity(0)} - * returns one, and all other - * indices will throw an error. - */ - virtual unsigned int element_multiplicity (const unsigned int index) const; - /** * Check for non-zero values on a face. * @@ -438,59 +310,6 @@ class FE_Q_Hierarchical : public FiniteElement * constructors of @p{FESystem}. */ virtual FiniteElement * clone() const; - - /** - * Prepare internal data - * structures and fill in values - * independent of the cell. - */ - virtual - typename Mapping::InternalDataBase * - get_data (const UpdateFlags, - const Mapping& mapping, - const Quadrature& quadrature) const ; - - /** - * Implementation of the same - * function in - * @ref{FiniteElement}. - */ - virtual void - fill_fe_values (const Mapping &mapping, - const typename DoFHandler::cell_iterator &cell, - const Quadrature &quadrature, - typename Mapping::InternalDataBase &mapping_internal, - typename Mapping::InternalDataBase &fe_internal, - FEValuesData& data) const; - - /** - * Implementation of the same - * function in - * @ref{FiniteElement}. - */ - virtual void - fill_fe_face_values (const Mapping &mapping, - const typename DoFHandler::cell_iterator &cell, - const unsigned int face_no, - const Quadrature &quadrature, - typename Mapping::InternalDataBase &mapping_internal, - typename Mapping::InternalDataBase &fe_internal, - FEValuesData& data) const ; - - /** - * Implementation of the same - * function in - * @ref{FiniteElement}. - */ - virtual void - fill_fe_subface_values (const Mapping &mapping, - const typename DoFHandler::cell_iterator &cell, - const unsigned int face_no, - const unsigned int sub_no, - const Quadrature &quadrature, - typename Mapping::InternalDataBase &mapping_internal, - typename Mapping::InternalDataBase &fe_internal, - FEValuesData& data) const ; private: @@ -628,159 +447,16 @@ class FE_Q_Hierarchical : public FiniteElement */ void initialize_unit_face_support_points (); - /** - * Determine the values that need - * to be computed on the unit - * cell to be able to compute all - * values required by @p{flags}. - * - * For the purpuse of this - * function, refer to the - * documentation in - * @p{FiniteElement}. - * - * The effect in this element is - * as follows: if - * @p{update_values} is set in - * @p{flags}, copy it to the - * result. All other flags of the - * result are cleared, since - * everything else must be - * computed for each cell. - */ - virtual UpdateFlags update_once (const UpdateFlags flags) const; - - /** - * Determine the values that need - * to be computed on every - * cell to be able to compute all - * values required by @p{flags}. - * - * For the purpuse of this - * function, refer to the - * documentation in - * @p{FiniteElement}. - * - * The effect in this element is - * as follows: - * @begin{itemize} - * @item if @p{update_gradients} - * is set, the result will - * contain @p{update_gradients} - * and - * @p{update_covariant_transformation}. - * The latter is required to - * transform the gradient on the - * unit cell to the real - * cell. Remark, that the action - * required by - * @p{update_covariant_transformation} - * is actually performed by the - * @p{Mapping} object used in - * conjunction with this finite - * element. - * @item if - * @p{update_second_derivatives} - * is set, the result will - * contain - * @p{update_second_derivatives} - * and - * @p{update_covariant_transformation}. - * The rationale is the same as - * above and no higher - * derivatives of the - * transformation are required, - * since we use difference - * quotients for the actual - * computation. - * @end{itemize} - */ - virtual UpdateFlags update_each (const UpdateFlags flags) const; - /** * Degree of the polynomials. */ const unsigned int degree; - - /** - * Mapping from lexicographic to - * shape function numbering. - */ - const std::vector renumber; - - /** - * Inverse renumber - * vector. i.e. mapping from - * shape function numbering to - * lexicographic numbering. - */ - const std::vector renumber_inverse; /** * Mapping from lexicographic to * shape function numbering on first face. */ const std::vector face_renumber; - - /** - * Pointer to the tensor - * product polynomials. - */ - const TensorProductPolynomials polynomial_space; - - /** - * Fields of cell-independent data. - * - * For information about the - * general purpose of this class, - * see the documentation of the - * base class. - */ - class InternalData : public FiniteElementBase::InternalDataBase - { - public: - /** - * Array with shape function - * values in quadrature - * points. There is one - * row for each shape - * function, containing - * values for each quadrature - * point. - * - * In this array, we store - * the values of the shape - * function in the quadrature - * points on the unit - * cell. Since these values - * do not change under - * transformation to the real - * cell, we only need to copy - * them over when visiting a - * concrete cell. - */ - Table<2,double> shape_values; - - /** - * Array with shape function - * gradients in quadrature - * points. There is one - * row for each shape - * function, containing - * values for each quadrature - * point. - * - * We store the gradients in - * the quadrature points on - * the unit cell. We then - * only have to apply the - * transformation (which is a - * matrix-vector - * multiplication) when - * visiting an actual cell. - */ - Table<2,Tensor<1,dim> > shape_gradients; - }; /** * Allow access from other diff --git a/deal.II/deal.II/source/fe/fe_q_hierarchical.cc b/deal.II/deal.II/source/fe/fe_q_hierarchical.cc index edf384f409..7b9908a9d5 100644 --- a/deal.II/deal.II/source/fe/fe_q_hierarchical.cc +++ b/deal.II/deal.II/source/fe/fe_q_hierarchical.cc @@ -2,7 +2,7 @@ // $Id$ // Version: $Name$ // -// Copyright (C) 2002, 2003 by the deal.II authors +// Copyright (C) 2002, 2003, 2004 by the deal.II authors // // This file is subject to QPL and may not be distributed // without copyright and license information. Please refer @@ -47,17 +47,19 @@ namespace template FE_Q_Hierarchical::FE_Q_Hierarchical (const unsigned int degree) : - FiniteElement (FiniteElementData(get_dpo_vector(degree),1, degree), - std::vector (FiniteElementData(get_dpo_vector(degree),1, degree).dofs_per_cell, - false), - std::vector >(FiniteElementData(get_dpo_vector(degree),1, degree).dofs_per_cell, - std::vector(1,true))), - degree(degree), - renumber(lexicographic_to_hierarchic_numbering (*this, degree)), - renumber_inverse(invert_numbering(renumber)), - face_renumber(face_lexicographic_to_hierarchic_numbering (degree)), - polynomial_space(Polynomials::Hierarchical::generate_complete_basis(degree)) + FE_Poly, dim> ( + Polynomials::Hierarchical::generate_complete_basis(degree), + FiniteElementData(get_dpo_vector(degree),1, degree), + std::vector (FiniteElementData( + get_dpo_vector(degree),1, degree).dofs_per_cell, false), + std::vector >(FiniteElementData( + get_dpo_vector(degree),1, degree).dofs_per_cell, std::vector(1,true))), + degree(degree), + face_renumber(face_lexicographic_to_hierarchic_numbering (degree)) { + this->poly_space.set_numbering(invert_numbering( + lexicographic_to_hierarchic_numbering (*this, degree))); + // The matrix @p{dofs_cell} contains the // values of the linear functionals of // the master 1d cell applied to the @@ -133,78 +135,6 @@ FE_Q_Hierarchical::clone() const } - -template -double -FE_Q_Hierarchical::shape_value (const unsigned int i, - const Point &p) const -{ - Assert (idofs_per_cell, ExcIndexRange(i,0,this->dofs_per_cell)); - return polynomial_space.compute_value(renumber_inverse[i], p); -} - - - -template -double -FE_Q_Hierarchical::shape_value_component (const unsigned int i, - const Point &p, - const unsigned int component) const -{ - Assert (idofs_per_cell, ExcIndexRange(i,0,this->dofs_per_cell)); - Assert (component == 0, ExcIndexRange (component, 0, 1)); - return polynomial_space.compute_value(renumber_inverse[i], p); -} - - - -template -Tensor<1,dim> -FE_Q_Hierarchical::shape_grad (const unsigned int i, - const Point &p) const -{ - Assert (idofs_per_cell, ExcIndexRange(i,0,this->dofs_per_cell)); - return polynomial_space.compute_grad(renumber_inverse[i], p); -} - - - -template -Tensor<1,dim> -FE_Q_Hierarchical::shape_grad_component (const unsigned int i, - const Point &p, - const unsigned int component) const -{ - Assert (idofs_per_cell, ExcIndexRange(i,0,this->dofs_per_cell)); - Assert (component == 0, ExcIndexRange (component, 0, 1)); - return polynomial_space.compute_grad(renumber_inverse[i], p); -} - - - -template -Tensor<2,dim> -FE_Q_Hierarchical::shape_grad_grad (const unsigned int i, - const Point &p) const -{ - Assert (idofs_per_cell, ExcIndexRange(i,0,this->dofs_per_cell)); - return polynomial_space.compute_grad_grad(renumber_inverse[i], p); -} - - - -template -Tensor<2,dim> -FE_Q_Hierarchical::shape_grad_grad_component (const unsigned int i, - const Point &p, - const unsigned int component) const -{ - Assert (idofs_per_cell, ExcIndexRange(i,0,this->dofs_per_cell)); - Assert (component == 0, ExcIndexRange (component, 0, 1)); - return polynomial_space.compute_grad_grad(renumber_inverse[i], p); -} - - //---------------------------------------------------------------------- // Auxiliary functions //---------------------------------------------------------------------- @@ -414,6 +344,9 @@ initialize_embedding_and_restriction (const std::vector > &do const std::vector > &dofs_subcell) { const unsigned int dofs_1d = 2*this->dofs_per_vertex + this->dofs_per_line; + const std::vector &renumber= + this->poly_space.get_numbering(); + for (unsigned int c=0; c::children_per_cell; ++c) { @@ -455,16 +388,16 @@ initialize_embedding_and_restriction (const std::vector > &do unsigned int c1 = ((c==2) || (c==3)) ? 1 : 0; this->prolongation[c](j,i) = - dofs_subcell[c0](renumber_inverse[j] % dofs_1d, - renumber_inverse[i] % dofs_1d) * - dofs_subcell[c1]((renumber_inverse[j] - (renumber_inverse[j] % dofs_1d)) / dofs_1d, - (renumber_inverse[i] - (renumber_inverse[i] % dofs_1d)) / dofs_1d); + dofs_subcell[c0](renumber[j] % dofs_1d, + renumber[i] % dofs_1d) * + dofs_subcell[c1]((renumber[j] - (renumber[j] % dofs_1d)) / dofs_1d, + (renumber[i] - (renumber[i] % dofs_1d)) / dofs_1d); this->restriction[c](j,i) = - dofs_cell[c0](renumber_inverse[j] % dofs_1d, - renumber_inverse[i] % dofs_1d) * - dofs_cell[c1]((renumber_inverse[j] - (renumber_inverse[j] % dofs_1d)) / dofs_1d, - (renumber_inverse[i] - (renumber_inverse[i] % dofs_1d)) / dofs_1d); + dofs_cell[c0](renumber[j] % dofs_1d, + renumber[i] % dofs_1d) * + dofs_cell[c1]((renumber[j] - (renumber[j] % dofs_1d)) / dofs_1d, + (renumber[i] - (renumber[i] % dofs_1d)) / dofs_1d); } break; } @@ -478,20 +411,20 @@ initialize_embedding_and_restriction (const std::vector > &do unsigned int c2 = ((c==2) || (c==3) || (c==6) || (c==7)) ? 1 : 0; this->prolongation[c](j,i) = - dofs_subcell[c0](renumber_inverse[j] % dofs_1d, - renumber_inverse[i] % dofs_1d) * - dofs_subcell[c1](((renumber_inverse[j] - (renumber_inverse[j] % dofs_1d)) / dofs_1d) % dofs_1d, - ((renumber_inverse[i] - (renumber_inverse[i] % dofs_1d)) / dofs_1d) % dofs_1d) * - dofs_subcell[c2](((renumber_inverse[j] - (renumber_inverse[j] % dofs_1d)) / dofs_1d - (((renumber_inverse[j] - (renumber_inverse[j] % dofs_1d)) / dofs_1d ) % dofs_1d)) / dofs_1d, - ((renumber_inverse[i] - (renumber_inverse[i] % dofs_1d)) / dofs_1d - (((renumber_inverse[i] - (renumber_inverse[i] % dofs_1d)) / dofs_1d ) % dofs_1d)) / dofs_1d); + dofs_subcell[c0](renumber[j] % dofs_1d, + renumber[i] % dofs_1d) * + dofs_subcell[c1](((renumber[j] - (renumber[j] % dofs_1d)) / dofs_1d) % dofs_1d, + ((renumber[i] - (renumber[i] % dofs_1d)) / dofs_1d) % dofs_1d) * + dofs_subcell[c2](((renumber[j] - (renumber[j] % dofs_1d)) / dofs_1d - (((renumber[j] - (renumber[j] % dofs_1d)) / dofs_1d ) % dofs_1d)) / dofs_1d, + ((renumber[i] - (renumber[i] % dofs_1d)) / dofs_1d - (((renumber[i] - (renumber[i] % dofs_1d)) / dofs_1d ) % dofs_1d)) / dofs_1d); this->restriction[c](j,i) = - dofs_cell[c0](renumber_inverse[j] % dofs_1d, - renumber_inverse[i] % dofs_1d) * - dofs_cell[c1](((renumber_inverse[j] - (renumber_inverse[j] % dofs_1d)) / dofs_1d) % dofs_1d, - ((renumber_inverse[i] - (renumber_inverse[i] % dofs_1d)) / dofs_1d) % dofs_1d) * - dofs_cell[c2](((renumber_inverse[j] - (renumber_inverse[j] % dofs_1d)) / dofs_1d - (((renumber_inverse[j] - (renumber_inverse[j] % dofs_1d)) / dofs_1d ) % dofs_1d)) / dofs_1d, - ((renumber_inverse[i] - (renumber_inverse[i] % dofs_1d)) / dofs_1d - (((renumber_inverse[i] - (renumber_inverse[i] % dofs_1d)) / dofs_1d ) % dofs_1d)) / dofs_1d); + dofs_cell[c0](renumber[j] % dofs_1d, + renumber[i] % dofs_1d) * + dofs_cell[c1](((renumber[j] - (renumber[j] % dofs_1d)) / dofs_1d) % dofs_1d, + ((renumber[i] - (renumber[i] % dofs_1d)) / dofs_1d) % dofs_1d) * + dofs_cell[c2](((renumber[j] - (renumber[j] % dofs_1d)) / dofs_1d - (((renumber[j] - (renumber[j] % dofs_1d)) / dofs_1d ) % dofs_1d)) / dofs_1d, + ((renumber[i] - (renumber[i] % dofs_1d)) / dofs_1d - (((renumber[i] - (renumber[i] % dofs_1d)) / dofs_1d ) % dofs_1d)) / dofs_1d); } break; } @@ -512,6 +445,9 @@ void FE_Q_Hierarchical::initialize_unit_support_points () n *= degree+1; this->unit_support_points.resize(n); + + const std::vector &index_map_inverse= + this->poly_space.get_numbering_inverse(); Point p; // the method of numbering allows @@ -567,7 +503,7 @@ void FE_Q_Hierarchical::initialize_unit_support_points () else p(2) = .5; } - this->unit_support_points[renumber[k++]] = p; + this->unit_support_points[index_map_inverse[k++]] = p; }; } @@ -907,296 +843,6 @@ FE_Q_Hierarchical<1>::face_lexicographic_to_hierarchic_numbering (const unsigned #endif -template -UpdateFlags -FE_Q_Hierarchical::update_once (const UpdateFlags flags) const -{ - // for this kind of elements, only - // the values can be precomputed - // once and for all. set this flag - // if the values are requested at - // all - return (update_default | (flags & update_values)); -} - - - -template -UpdateFlags -FE_Q_Hierarchical::update_each (const UpdateFlags flags) const -{ - UpdateFlags out = update_default; - - if (flags & update_gradients) - out |= update_gradients | update_covariant_transformation; - if (flags & update_second_derivatives) - out |= update_second_derivatives | update_covariant_transformation; - - return out; -} - - - -//---------------------------------------------------------------------- -// Data field initialization -//---------------------------------------------------------------------- - -template -typename Mapping::InternalDataBase * -FE_Q_Hierarchical::get_data (const UpdateFlags update_flags, - const Mapping &mapping, - const Quadrature &quadrature) const -{ - // generate a new data object and - // initialize some fields - InternalData* data = new InternalData; - - // check what needs to be - // initialized only once and what - // on every cell/face/subface we - // visit - data->update_once = update_once(update_flags); - data->update_each = update_each(update_flags); - data->update_flags = data->update_once | data->update_each; - - const UpdateFlags flags(data->update_flags); - const unsigned int n_q_points = quadrature.n_quadrature_points; - - // some scratch arrays - std::vector values(0); - std::vector > grads(0); - std::vector > grad_grads(0); - - // initialize fields only if really - // necessary. otherwise, don't - // allocate memory - if (flags & update_values) - { - values.resize (this->dofs_per_cell); - data->shape_values.reinit (this->dofs_per_cell, - n_q_points); - } - - if (flags & update_gradients) - { - grads.resize (this->dofs_per_cell); - data->shape_gradients.reinit (this->dofs_per_cell, - n_q_points); - } - - // if second derivatives through - // finite differencing is required, - // then initialize some objects for - // that - if (flags & update_second_derivatives) - data->initialize_2nd (this, mapping, quadrature); - - // next already fill those fields - // of which we have information by - // now. note that the shape - // gradients are only those on the - // unit cell, and need to be - // transformed when visiting an - // actual cell - if (flags & (update_values | update_gradients)) - for (unsigned int i=0; idofs_per_cell; ++k) - data->shape_values[renumber[k]][i] = values[k]; - - if (flags & update_gradients) - for (unsigned int k=0; kdofs_per_cell; ++k) - data->shape_gradients[renumber[k]][i] = grads[k]; - } - return data; -} - - - - -//---------------------------------------------------------------------- -// Fill data of FEValues -//---------------------------------------------------------------------- - -template -void -FE_Q_Hierarchical::fill_fe_values ( - const Mapping &mapping, - const typename DoFHandler::cell_iterator &cell, - const Quadrature &quadrature, - typename Mapping::InternalDataBase &mapping_data, - typename Mapping::InternalDataBase &fedata, - FEValuesData &data) const -{ - // convert data object to internal - // data for this class. fails with - // an exception if that is not - // possible - InternalData &fe_data = dynamic_cast (fedata); - - const UpdateFlags flags(fe_data.current_update_flags()); - - for (unsigned int k=0; kdofs_per_cell; ++k) - { - if (flags & update_values) - for (unsigned int i=0; icompute_2nd (mapping, cell, - QProjector::DataSetDescriptor::cell(), - mapping_data, fe_data, data); -} - - - -template -void -FE_Q_Hierarchical::fill_fe_face_values ( - const Mapping &mapping, - const typename DoFHandler::cell_iterator &cell, - const unsigned int face, - const Quadrature &quadrature, - typename Mapping::InternalDataBase &mapping_data, - typename Mapping::InternalDataBase &fedata, - FEValuesData &data) const -{ - // convert data object to internal - // data for this class. fails with - // an exception if that is not - // possible - InternalData &fe_data = dynamic_cast (fedata); - - // offset determines which data set - // to take (all data sets for all - // faces are stored contiguously) - const typename QProjector::DataSetDescriptor offset - = (QProjector::DataSetDescriptor:: - face (face, cell->face_orientation(face), - quadrature.n_quadrature_points)); - - const UpdateFlags flags(fe_data.update_once | fe_data.update_each); - - for (unsigned int k=0; kdofs_per_cell; ++k) - { - if (flags & update_values) - for (unsigned int i=0; icompute_2nd (mapping, cell, offset, mapping_data, fe_data, data); -} - - - -template -void -FE_Q_Hierarchical::fill_fe_subface_values ( - const Mapping &mapping, - const typename DoFHandler::cell_iterator &cell, - const unsigned int face, - const unsigned int subface, - const Quadrature &quadrature, - typename Mapping::InternalDataBase &mapping_data, - typename Mapping::InternalDataBase &fedata, - FEValuesData &data) const -{ - // convert data object to internal - // data for this class. fails with - // an exception if that is not - // possible - InternalData &fe_data = dynamic_cast (fedata); - - // offset determines which data set - // to take (all data sets for all - // sub-faces are stored contiguously) - const typename QProjector::DataSetDescriptor offset - = (QProjector::DataSetDescriptor:: - sub_face (face, subface, cell->face_orientation(face), - quadrature.n_quadrature_points)); - - const UpdateFlags flags(fe_data.update_once | fe_data.update_each); - - for (unsigned int k=0; kdofs_per_cell; ++k) - { - if (flags & update_values) - for (unsigned int i=0; icompute_2nd (mapping, cell, offset, mapping_data, fe_data, data); -} - - - -template -unsigned int -FE_Q_Hierarchical::n_base_elements () const -{ - return 1; -} - - - -template -const FiniteElement & -FE_Q_Hierarchical::base_element (const unsigned int index) const -{ - Assert (index==0, ExcIndexRange(index, 0, 1)); - return *this; -} - - - -template -unsigned int -FE_Q_Hierarchical::element_multiplicity (const unsigned int index) const -{ - Assert (index==0, ExcIndexRange(index, 0, 1)); - return 1; -} - - template bool