From a8eafbace9fe84b1c4f86aac02b11040183cfe73 Mon Sep 17 00:00:00 2001 From: hartmann Date: Tue, 16 Mar 2004 17:59:02 +0000 Subject: [PATCH] Make use of the new FE_Poly class. git-svn-id: https://svn.dealii.org/trunk@8782 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/include/fe/fe_dgq.h | 312 ++--------------------- deal.II/deal.II/source/fe/fe_dgq.cc | 372 +--------------------------- 2 files changed, 24 insertions(+), 660 deletions(-) diff --git a/deal.II/deal.II/include/fe/fe_dgq.h b/deal.II/deal.II/include/fe/fe_dgq.h index bb4af0c9fe..5bd213fe0b 100644 --- a/deal.II/deal.II/include/fe/fe_dgq.h +++ b/deal.II/deal.II/include/fe/fe_dgq.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,151 +14,47 @@ #define __deal2__fe_dgq_h #include -#include #include -#include +#include -template class TensorProductPolynomials; template class MappingQ; - /*!@addtogroup fe */ /*@{*/ /** - * Discontinuous tensor product elements based on equidistant support points. + * Discontinuous tensor product elements based on equidistant support + * points. * - * This is a discontinuous finite element using interpolating tensor - * product polynomials. The shape functions are Lagrangian + * This is a discontinuous finite element based on tensor products of + * Lagrangian polynomials. The shape functions are Lagrangian * interpolants of an equidistant grid of points on the unit cell. The - * points are numbered in lexicographical order, @p{x} running fastest. + * points are numbered in lexicographical order, with @p x running + * fastest. * - * @author Guido Kanschat, Ralf Hartmann, 2001 + * @author Ralf Hartmann, 2001, 2004, Guido Kanschat 2001 */ template -class FE_DGQ : public FiniteElement +class FE_DGQ : public FE_Poly,dim> { public: /** * Constructor for tensor product - * polynomials of degree @p{k}. + * polynomials of degree @p p. */ - FE_DGQ (const unsigned int k); + FE_DGQ (const unsigned int p); /** * Return a string that uniquely * identifies a finite * element. This class returns - * @p{FE_DGQ(degree)}, with - * @p{dim} and @p{degree} + * @p FE_DGQ(degree) , with + * @p dim and @p degree * replaced by appropriate * values. */ 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, @@ -186,34 +82,6 @@ class FE_DGQ : 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. @@ -302,59 +170,6 @@ class FE_DGQ : 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: /** @@ -368,45 +183,6 @@ class FE_DGQ : public FiniteElement * @p{FiniteElementData}. */ static std::vector get_dpo_vector(unsigned int degree); - - /** - * Given a set of flags indicating - * what quantities are requested - * from a @p{FEValues} object, - * return which of these can be - * precomputed once and for - * all. Often, the values of - * shape function at quadrature - * points can be precomputed, for - * example, in which case the - * return value of this function - * would be the logical and of - * the input @p{flags} and - * @p{update_values}. - * - * For the present kind of finite - * element, this is exactly the - * case. - */ - virtual UpdateFlags update_once (const UpdateFlags flags) const; - - /** - * This is the opposite to the - * above function: given a set of - * flags indicating what we want - * to know, return which of these - * need to be computed each time - * we visit a new cell. - * - * If for the computation of one - * quantity something else is - * also required (for example, we - * often need the covariant - * transformation when gradients - * need to be computed), include - * this in the result as well. - */ - virtual UpdateFlags update_each (const UpdateFlags flags) const; /** * Compute renumbering for rotation @@ -441,66 +217,6 @@ class FE_DGQ : public FiniteElement * Degree of the polynomials. */ const unsigned int degree; - - /** - * 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 dimensions. diff --git a/deal.II/deal.II/source/fe/fe_dgq.cc b/deal.II/deal.II/source/fe/fe_dgq.cc index b5cc39fb39..8d54bf1411 100644 --- a/deal.II/deal.II/source/fe/fe_dgq.cc +++ b/deal.II/deal.II/source/fe/fe_dgq.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 @@ -144,13 +144,13 @@ namespace template FE_DGQ::FE_DGQ (const unsigned int degree) : - FiniteElement (FiniteElementData(get_dpo_vector(degree), 1, degree), - std::vector(FiniteElementData(get_dpo_vector(degree),1, degree).dofs_per_cell, - true), - std::vector >(FiniteElementData(get_dpo_vector(degree),1, degree).dofs_per_cell, - std::vector(1,true))), - degree(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, true), + std::vector >(FiniteElementData( + get_dpo_vector(degree),1, degree).dofs_per_cell, std::vector(1,true))), + degree(degree) { // generate permutation/rotation // index sets to generate some @@ -361,77 +361,6 @@ FE_DGQ::clone() const } - -template -double -FE_DGQ::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(i, p); -} - - - -template -double -FE_DGQ::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(i, p); -} - - - -template -Tensor<1,dim> -FE_DGQ::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(i, p); -} - - -template -Tensor<1,dim> -FE_DGQ::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(i, p); -} - - - -template -Tensor<2,dim> -FE_DGQ::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(i, p); -} - - - -template -Tensor<2,dim> -FE_DGQ::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(i, p); -} - - //---------------------------------------------------------------------- // Auxiliary functions //---------------------------------------------------------------------- @@ -449,35 +378,6 @@ FE_DGQ::get_dpo_vector(unsigned int deg) } -template -UpdateFlags -FE_DGQ::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_DGQ::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; -} - - template void @@ -605,11 +505,11 @@ get_interpolation_matrix (const FiniteElementBase &x_source_fe, int2type()); for (unsigned int i=0; idofs_per_cell; ++i) cell_interpolation(j,i) - = polynomial_space.compute_value (i, p); + = this->poly_space.compute_value (i, p); for (unsigned int i=0; i &x_source_fe, - - -//---------------------------------------------------------------------- -// Data field initialization -//---------------------------------------------------------------------- - -template -typename Mapping::InternalDataBase * -FE_DGQ::get_data (const UpdateFlags update_flags, - const Mapping &mapping, - const Quadrature &quadrature) const -{ - // generate a new data object - 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; - - // have 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) - { - if (flags & update_values) - data->shape_values[k][i] = values[k]; - if (flags & update_gradients) - data->shape_gradients[k][i] = grads[k]; - } - } - return data; -} - - - //---------------------------------------------------------------------- // Fill data of FEValues //---------------------------------------------------------------------- -template -void -FE_DGQ::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_DGQ::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) - { - for (unsigned int i=0;icompute_2nd (mapping, cell, offset, mapping_data, fe_data, data); -} - - - -template -void -FE_DGQ::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) - { - for (unsigned int i=0;icompute_2nd (mapping, cell, offset, mapping_data, fe_data, data); -} - - - -template -unsigned int -FE_DGQ::n_base_elements () const -{ - return 1; -} - - - -template -const FiniteElement & -FE_DGQ::base_element (const unsigned int index) const -{ - Assert (index==0, ExcIndexRange(index, 0, 1)); - return *this; -} - - - -template -unsigned int -FE_DGQ::element_multiplicity (const unsigned int index) const -{ - Assert (index==0, ExcIndexRange(index, 0, 1)); - return 1; -} - - template bool -- 2.39.5