From: hartmann Date: Tue, 16 Mar 2004 16:38:30 +0000 (+0000) Subject: Make use of the new FE_Poly class. Move most of renumbering business to the TensorPro... X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=7de6f6d38f418cd31113905368508d655006253f;p=dealii-svn.git Make use of the new FE_Poly class. Move most of renumbering business to the TensorProductPolynomials class. git-svn-id: https://svn.dealii.org/trunk@8777 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/fe/fe_q.h b/deal.II/deal.II/include/fe/fe_q.h index f7c0d11b77..8491866a8d 100644 --- a/deal.II/deal.II/include/fe/fe_q.h +++ b/deal.II/deal.II/include/fe/fe_q.h @@ -2,7 +2,7 @@ // $Id$ // Version: $Name$ // -// Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003 by the deal.II authors +// Copyright (C) 1998, 1999, 2000, 2001, 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,12 +14,8 @@ #define __deal2__fe_q_h #include -#include #include -#include - -template class TensorProductPolynomials; -template class MappingQ; +#include /*!@addtogroup fe */ @@ -237,10 +233,10 @@ template class MappingQ; * Note the reverse ordering of degrees of freedom on the left and upper * line. * - * @author Wolfgang Bangerth, 1998, Ralf Hartmann, Guido Kanschat, 2001 + * @author Wolfgang Bangerth, 1998, Guido Kanschat, 2001, Ralf Hartmann, 2001, 2004 */ template -class FE_Q : public FiniteElement +class FE_Q : public FE_Poly,dim> { public: /** @@ -260,109 +256,6 @@ class FE_Q : public FiniteElement */ virtual std::string get_name () const; - /** - * Return the value 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. - */ - 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}. See the - * @ref{FiniteElementBase} base - * class for more information - * about the semantics of this - * function. - */ - 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. See the - * @ref{FiniteElementBase} base - * class for more information - * about the semantics of this - * function. - */ - 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, @@ -390,34 +283,6 @@ class FE_Q : public FiniteElement virtual void get_interpolation_matrix (const FiniteElementBase &source, FullMatrix &matrix) 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. @@ -489,59 +354,6 @@ class FE_Q : 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: @@ -680,159 +492,19 @@ class FE_Q : 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. + * Mapping from hierarchic to + * lexicographic numbering on + * first face. Hierarchic is the + * numbering of the shape + * functions. */ - 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; - }; + const std::vector face_index_map; /** * Allow access from other @@ -846,6 +518,8 @@ class FE_Q : public FiniteElement template friend class FE_Q; }; + + /*@}*/ /* -------------- declaration of explicit specializations ------------- */ diff --git a/deal.II/deal.II/source/fe/fe_q.cc b/deal.II/deal.II/source/fe/fe_q.cc index 53e7be3370..5b7fcfba56 100644 --- a/deal.II/deal.II/source/fe/fe_q.cc +++ b/deal.II/deal.II/source/fe/fe_q.cc @@ -2,7 +2,7 @@ // $Id$ // Version: $Name$ // -// Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003 by the deal.II authors +// Copyright (C) 1998, 1999, 2000, 2001, 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 @@ -12,8 +12,6 @@ //---------------------------------------------------------------- #include -#include -#include #include #include #include @@ -431,22 +429,21 @@ namespace FE_Q_Helper template FE_Q::FE_Q (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(FE_Q_Helper::invert_numbering (renumber)), - face_renumber(face_lexicographic_to_hierarchic_numbering (degree)), - polynomial_space(Polynomials::LagrangeEquidistant::generate_complete_basis(degree)) + FE_Poly, dim> ( + TensorProductPolynomials(Polynomials::LagrangeEquidistant::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_index_map(FE_Q_Helper::invert_numbering(face_lexicographic_to_hierarchic_numbering (degree))) { - - // copy constraint and embedding - // matrices if they are - // defined. otherwise leave them at - // invalid size + this->poly_space.set_numbering(FE_Q_Helper::invert_numbering( + lexicographic_to_hierarchic_numbering (*this, degree))); + + // compute constraint, embedding + // and restriction matrices initialize_constraints (); initialize_embedding (); initialize_restriction (); @@ -495,77 +492,6 @@ FE_Q::clone() const -template -double -FE_Q::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::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::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::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::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::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); -} - - - template void FE_Q:: @@ -592,7 +518,9 @@ get_interpolation_matrix (const FiniteElementBase &x_source_fe, Assert (interpolation_matrix.n() == source_fe.dofs_per_cell, ExcDimensionMismatch (interpolation_matrix.m(), source_fe.dofs_per_cell)); - + + const std::vector &index_map= + this->poly_space.get_numbering(); // compute the interpolation // matrices in much the same way as @@ -610,15 +538,13 @@ get_interpolation_matrix (const FiniteElementBase &x_source_fe, // cell and evaluate the // shape functions there const Point - p = FE_Q_Helper::generate_unit_point (j, this->dofs_per_cell, + p = FE_Q_Helper::generate_unit_point (index_map[j], this->dofs_per_cell, FE_Q_Helper::int2type()); for (unsigned int i=0; idofs_per_cell; ++i) - cell_interpolation(renumber[j],renumber[i]) - = polynomial_space.compute_value (i, p); + cell_interpolation(j,i) = this->poly_space.compute_value (i, p); for (unsigned int i=0; i::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(); const double step = 1./degree; Point p; @@ -682,7 +611,7 @@ void FE_Q::initialize_unit_support_points () if (dim>2) p(2) = iz * step; - this->unit_support_points[renumber[k++]] = p; + this->unit_support_points[index_map_inverse[k++]] = p; }; } @@ -709,6 +638,9 @@ void FE_Q::initialize_unit_face_support_points () n *= degree+1; this->unit_face_support_points.resize(n); + + const std::vector &face_index_map_inverse= + FE_Q_Helper::invert_numbering(face_index_map); const double step = 1./degree; Point p; @@ -724,7 +656,7 @@ void FE_Q::initialize_unit_face_support_points () if (codim>2) p(2) = iz * step; - this->unit_face_support_points[face_renumber[k++]] = p; + this->unit_face_support_points[face_index_map_inverse[k++]] = p; }; } @@ -1100,9 +1032,6 @@ FE_Q<2>::initialize_constraints () FullMatrix face_interpolation (n_small_functions, this->dofs_per_face); FullMatrix subface_interpolation (n_small_functions, n_small_functions); - - const std::vector - face_renumber_inverse (FE_Q_Helper::invert_numbering(face_renumber)); for (unsigned int i=0; i::initialize_constraints () // them in the order of // interpolation points. // - // face_renumber_inverse will + // face_index_map_inverse will // get us over this little // conversion for (unsigned int j=0; jdofs_per_face; ++j) { face_interpolation(i,j) - = face_polynomials.compute_value(face_renumber_inverse[j], p_face); + = face_polynomials.compute_value(face_index_map[j], p_face); // if the value is small up // to round-off, then // simply set it to zero to @@ -1361,6 +1290,9 @@ FE_Q::initialize_embedding () this->dofs_per_cell); FullMatrix subcell_interpolation (this->dofs_per_cell, this->dofs_per_cell); + const std::vector &index_map= + this->poly_space.get_numbering(); + for (unsigned int child=0; child::children_per_cell; ++child) this->prolongation[child].reinit (this->dofs_per_cell, this->dofs_per_cell); @@ -1373,7 +1305,7 @@ FE_Q::initialize_embedding () // evaluate the shape // functions there const Point p_subcell - = FE_Q_Helper::generate_unit_point (j, this->dofs_per_cell, + = FE_Q_Helper::generate_unit_point (index_map[j], this->dofs_per_cell, FE_Q_Helper::int2type()); const Point p_cell = GeometryInfo::child_to_cell_coordinates (p_subcell, child); @@ -1381,8 +1313,8 @@ FE_Q::initialize_embedding () for (unsigned int i=0; idofs_per_cell; ++i) { const double - cell_value = polynomial_space.compute_value (i, p_cell), - subcell_value = polynomial_space.compute_value (i, p_subcell); + cell_value = this->poly_space.compute_value (i, p_cell), + subcell_value = this->poly_space.compute_value (i, p_subcell); // cut off values that // are too small. note @@ -1417,16 +1349,16 @@ FE_Q::initialize_embedding () // degree*dim, times a // small constant. if (std::fabs(cell_value) < 2e-14*degree*dim) - cell_interpolation(renumber[j], renumber[i]) = 0.; + cell_interpolation(j, i) = 0.; else - cell_interpolation(renumber[j], renumber[i]) = cell_value; + cell_interpolation(j, i) = cell_value; if (std::fabs(subcell_value) < 2e-14*degree*dim) - subcell_interpolation(renumber[j], renumber[i]) = 0.; + subcell_interpolation(j, i) = 0.; else if (std::fabs(subcell_value-1) < 2e-14*degree*dim) - subcell_interpolation(renumber[j], renumber[i]) = 1.; - else + subcell_interpolation(j, i) = 1.; + else // we have put our // evaluation // points onto the @@ -1544,8 +1476,7 @@ FE_Q::initialize_restriction () for (; mother_dofdofs_per_cell; ++mother_dof) { const double val - = polynomial_space.compute_value(renumber_inverse[mother_dof], - p_cell); + = this->poly_space.compute_value(mother_dof, p_cell); if (std::fabs (val-1.) < 2e-14*degree*dim) // ok, this is the right // dof @@ -1560,8 +1491,7 @@ FE_Q::initialize_restriction () // check also the shape // functions after tat for (unsigned int j=mother_dof+1; jdofs_per_cell; ++j) - Assert (std::fabs (polynomial_space.compute_value(renumber_inverse[j], - p_cell)) + Assert (std::fabs (this->poly_space.compute_value(j, p_cell)) < 2e-14*degree*dim, ExcInternalError()); @@ -1593,8 +1523,7 @@ FE_Q::initialize_restriction () for (; child_dofdofs_per_cell; ++child_dof) { const double val - = polynomial_space.compute_value(renumber_inverse[child_dof], - p_subcell); + = this->poly_space.compute_value(child_dof, p_subcell); if (std::fabs (val-1.) < 2e-14*degree*dim) break; else @@ -1602,8 +1531,7 @@ FE_Q::initialize_restriction () ExcInternalError()); } for (unsigned int j=child_dof+1; jdofs_per_cell; ++j) - Assert (std::fabs (polynomial_space.compute_value(renumber_inverse[j], - p_subcell)) + Assert (std::fabs (this->poly_space.compute_value(j, p_subcell)) < 2e-14*degree*dim, ExcInternalError()); @@ -1618,294 +1546,10 @@ FE_Q::initialize_restriction () } - -template -UpdateFlags -FE_Q::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::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::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::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::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::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::n_base_elements () const -{ - return 1; -} - - - -template -const FiniteElement & -FE_Q::base_element (const unsigned int index) const -{ - Assert (index==0, ExcIndexRange(index, 0, 1)); - return *this; -} - - - -template -unsigned int -FE_Q::element_multiplicity (const unsigned int index) const -{ - Assert (index==0, ExcIndexRange(index, 0, 1)); - return 1; -} - - template bool