// $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
#define __deal2__fe_q_h
#include <base/config.h>
-#include <base/polynomial.h>
#include <base/tensor_product_polynomials.h>
-#include <fe/fe.h>
-
-template <int dim> class TensorProductPolynomials;
-template <int dim> class MappingQ;
+#include <fe/fe_poly.h>
/*!@addtogroup fe */
* 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 <int dim>
-class FE_Q : public FiniteElement<dim>
+class FE_Q : public FE_Poly<TensorProductPolynomials<dim>,dim>
{
public:
/**
*/
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<dim> &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<dim> &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<dim> &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<dim> &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<dim> &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<dim> &p,
- const unsigned int component) const;
-
/**
* Return the polynomial degree
* of this finite element,
virtual void
get_interpolation_matrix (const FiniteElementBase<dim> &source,
FullMatrix<double> &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<dim> &
- 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.
* constructors of @p{FESystem}.
*/
virtual FiniteElement<dim> * clone() const;
-
- /**
- * Prepare internal data
- * structures and fill in values
- * independent of the cell.
- */
- virtual
- typename Mapping<dim>::InternalDataBase *
- get_data (const UpdateFlags,
- const Mapping<dim>& mapping,
- const Quadrature<dim>& quadrature) const ;
-
- /**
- * Implementation of the same
- * function in
- * @ref{FiniteElement}.
- */
- virtual void
- fill_fe_values (const Mapping<dim> &mapping,
- const typename DoFHandler<dim>::cell_iterator &cell,
- const Quadrature<dim> &quadrature,
- typename Mapping<dim>::InternalDataBase &mapping_internal,
- typename Mapping<dim>::InternalDataBase &fe_internal,
- FEValuesData<dim>& data) const;
-
- /**
- * Implementation of the same
- * function in
- * @ref{FiniteElement}.
- */
- virtual void
- fill_fe_face_values (const Mapping<dim> &mapping,
- const typename DoFHandler<dim>::cell_iterator &cell,
- const unsigned int face_no,
- const Quadrature<dim-1> &quadrature,
- typename Mapping<dim>::InternalDataBase &mapping_internal,
- typename Mapping<dim>::InternalDataBase &fe_internal,
- FEValuesData<dim>& data) const ;
-
- /**
- * Implementation of the same
- * function in
- * @ref{FiniteElement}.
- */
- virtual void
- fill_fe_subface_values (const Mapping<dim> &mapping,
- const typename DoFHandler<dim>::cell_iterator &cell,
- const unsigned int face_no,
- const unsigned int sub_no,
- const Quadrature<dim-1> &quadrature,
- typename Mapping<dim>::InternalDataBase &mapping_internal,
- typename Mapping<dim>::InternalDataBase &fe_internal,
- FEValuesData<dim>& data) const ;
private:
*/
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<unsigned int> renumber;
-
- /**
- * Inverse renumber
- * vector. i.e. mapping from
- * shape function numbering to
- * lexicographic numbering.
- */
- const std::vector<unsigned int> renumber_inverse;
/**
- * Mapping from lexicographic to
- * shape function numbering on first face.
- */
- const std::vector<unsigned int> face_renumber;
-
- /**
- * Pointer to the tensor
- * product polynomials.
- */
- const TensorProductPolynomials<dim> 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<dim>::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<unsigned int> face_index_map;
/**
* Allow access from other
template <int dim1> friend class FE_Q;
};
+
+
/*@}*/
/* -------------- declaration of explicit specializations ------------- */
// $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
//----------------------------------------------------------------
#include <base/quadrature.h>
-#include <base/polynomial.h>
-#include <base/tensor_product_polynomials.h>
#include <grid/tria.h>
#include <grid/tria_iterator.h>
#include <dofs/dof_accessor.h>
template <int dim>
FE_Q<dim>::FE_Q (const unsigned int degree)
:
- FiniteElement<dim> (FiniteElementData<dim>(get_dpo_vector(degree),1, degree),
- std::vector<bool> (FiniteElementData<dim>(get_dpo_vector(degree),1, degree).dofs_per_cell,
- false),
- std::vector<std::vector<bool> >(FiniteElementData<dim>(get_dpo_vector(degree),1, degree).dofs_per_cell,
- std::vector<bool>(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<TensorProductPolynomials<dim>, dim> (
+ TensorProductPolynomials<dim>(Polynomials::LagrangeEquidistant::generate_complete_basis(degree)),
+ FiniteElementData<dim>(get_dpo_vector(degree),1, degree),
+ std::vector<bool> (FiniteElementData<dim>(
+ get_dpo_vector(degree),1, degree).dofs_per_cell, false),
+ std::vector<std::vector<bool> >(FiniteElementData<dim>(
+ get_dpo_vector(degree),1, degree).dofs_per_cell, std::vector<bool>(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 ();
-template <int dim>
-double
-FE_Q<dim>::shape_value (const unsigned int i,
- const Point<dim> &p) const
-{
- Assert (i<this->dofs_per_cell, ExcIndexRange(i,0,this->dofs_per_cell));
- return polynomial_space.compute_value(renumber_inverse[i], p);
-}
-
-
-template <int dim>
-double
-FE_Q<dim>::shape_value_component (const unsigned int i,
- const Point<dim> &p,
- const unsigned int component) const
-{
- Assert (i<this->dofs_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 <int dim>
-Tensor<1,dim>
-FE_Q<dim>::shape_grad (const unsigned int i,
- const Point<dim> &p) const
-{
- Assert (i<this->dofs_per_cell, ExcIndexRange(i,0,this->dofs_per_cell));
- return polynomial_space.compute_grad(renumber_inverse[i], p);
-}
-
-
-
-template <int dim>
-Tensor<1,dim>
-FE_Q<dim>::shape_grad_component (const unsigned int i,
- const Point<dim> &p,
- const unsigned int component) const
-{
- Assert (i<this->dofs_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 <int dim>
-Tensor<2,dim>
-FE_Q<dim>::shape_grad_grad (const unsigned int i,
- const Point<dim> &p) const
-{
- Assert (i<this->dofs_per_cell, ExcIndexRange(i,0,this->dofs_per_cell));
- return polynomial_space.compute_grad_grad(renumber_inverse[i], p);
-}
-
-
-
-template <int dim>
-Tensor<2,dim>
-FE_Q<dim>::shape_grad_grad_component (const unsigned int i,
- const Point<dim> &p,
- const unsigned int component) const
-{
- Assert (i<this->dofs_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 <int dim>
void
FE_Q<dim>::
Assert (interpolation_matrix.n() == source_fe.dofs_per_cell,
ExcDimensionMismatch (interpolation_matrix.m(),
source_fe.dofs_per_cell));
-
+
+ const std::vector<unsigned int> &index_map=
+ this->poly_space.get_numbering();
// compute the interpolation
// matrices in much the same way as
// cell and evaluate the
// shape functions there
const Point<dim>
- 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<dim>());
for (unsigned int i=0; i<this->dofs_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<source_fe.dofs_per_cell; ++i)
- source_interpolation(renumber[j],source_fe.renumber[i])
- = source_fe.polynomial_space.compute_value (i, p);
+ source_interpolation(j,i) = source_fe.poly_space.compute_value (i, p);
}
// then compute the
n *= degree+1;
this->unit_support_points.resize(n);
+
+ const std::vector<unsigned int> &index_map_inverse=
+ this->poly_space.get_numbering_inverse();
const double step = 1./degree;
Point<dim> p;
if (dim>2)
p(2) = iz * step;
- this->unit_support_points[renumber[k++]] = p;
+ this->unit_support_points[index_map_inverse[k++]] = p;
};
}
n *= degree+1;
this->unit_face_support_points.resize(n);
+
+ const std::vector<unsigned int> &face_index_map_inverse=
+ FE_Q_Helper::invert_numbering(face_index_map);
const double step = 1./degree;
Point<codim> p;
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;
};
}
FullMatrix<double> face_interpolation (n_small_functions, this->dofs_per_face);
FullMatrix<double> subface_interpolation (n_small_functions, n_small_functions);
-
- const std::vector<unsigned int>
- face_renumber_inverse (FE_Q_Helper::invert_numbering(face_renumber));
for (unsigned int i=0; i<n_small_functions; ++i)
{
// 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; j<this->dofs_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
this->dofs_per_cell);
FullMatrix<double> subcell_interpolation (this->dofs_per_cell,
this->dofs_per_cell);
+ const std::vector<unsigned int> &index_map=
+ this->poly_space.get_numbering();
+
for (unsigned int child=0; child<GeometryInfo<dim>::children_per_cell; ++child)
this->prolongation[child].reinit (this->dofs_per_cell,
this->dofs_per_cell);
// evaluate the shape
// functions there
const Point<dim> 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<dim>());
const Point<dim> p_cell =
GeometryInfo<dim>::child_to_cell_coordinates (p_subcell, child);
for (unsigned int i=0; i<this->dofs_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
// 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
for (; mother_dof<this->dofs_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
// check also the shape
// functions after tat
for (unsigned int j=mother_dof+1; j<this->dofs_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());
for (; child_dof<this->dofs_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
ExcInternalError());
}
for (unsigned int j=child_dof+1; j<this->dofs_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());
}
-
-template <int dim>
-UpdateFlags
-FE_Q<dim>::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 <int dim>
-UpdateFlags
-FE_Q<dim>::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 <int dim>
-typename Mapping<dim>::InternalDataBase *
-FE_Q<dim>::get_data (const UpdateFlags update_flags,
- const Mapping<dim> &mapping,
- const Quadrature<dim> &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<double> values(0);
- std::vector<Tensor<1,dim> > grads(0);
- std::vector<Tensor<2,dim> > 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; i<n_q_points; ++i)
- {
- polynomial_space.compute(quadrature.point(i),
- values, grads, grad_grads);
-
- if (flags & update_values)
- for (unsigned int k=0; k<this->dofs_per_cell; ++k)
- data->shape_values[renumber[k]][i] = values[k];
-
- if (flags & update_gradients)
- for (unsigned int k=0; k<this->dofs_per_cell; ++k)
- data->shape_gradients[renumber[k]][i] = grads[k];
- }
- return data;
-}
-
-
-
-
-//----------------------------------------------------------------------
-// Fill data of FEValues
-//----------------------------------------------------------------------
-
-template <int dim>
-void
-FE_Q<dim>::fill_fe_values (const Mapping<dim> &mapping,
- const typename DoFHandler<dim>::cell_iterator &cell,
- const Quadrature<dim> &quadrature,
- typename Mapping<dim>::InternalDataBase &mapping_data,
- typename Mapping<dim>::InternalDataBase &fedata,
- FEValuesData<dim> &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<InternalData &> (fedata);
-
- const UpdateFlags flags(fe_data.current_update_flags());
-
- for (unsigned int k=0; k<this->dofs_per_cell; ++k)
- {
- if (flags & update_values)
- for (unsigned int i=0; i<quadrature.n_quadrature_points; ++i)
- data.shape_values(k,i) = fe_data.shape_values[k][i];
-
- if (flags & update_gradients)
- {
- Assert (data.shape_gradients[k].size() <=
- fe_data.shape_gradients[k].size(),
- ExcInternalError());
- mapping.transform_covariant(data.shape_gradients[k].begin(),
- data.shape_gradients[k].end(),
- fe_data.shape_gradients[k].begin(),
- mapping_data);
- };
- }
-
- if (flags & update_second_derivatives)
- this->compute_2nd (mapping, cell,
- QProjector<dim>::DataSetDescriptor::cell(),
- mapping_data, fe_data, data);
-}
-
-
-
-template <int dim>
-void
-FE_Q<dim>::fill_fe_face_values (const Mapping<dim> &mapping,
- const typename DoFHandler<dim>::cell_iterator &cell,
- const unsigned int face,
- const Quadrature<dim-1> &quadrature,
- typename Mapping<dim>::InternalDataBase &mapping_data,
- typename Mapping<dim>::InternalDataBase &fedata,
- FEValuesData<dim> &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<InternalData &> (fedata);
-
- // offset determines which data set
- // to take (all data sets for all
- // faces are stored contiguously)
- const typename QProjector<dim>::DataSetDescriptor offset
- = (QProjector<dim>::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; k<this->dofs_per_cell; ++k)
- {
- if (flags & update_values)
- for (unsigned int i=0; i<quadrature.n_quadrature_points; ++i)
- data.shape_values(k,i) = fe_data.shape_values[k][i+offset];
-
- if (flags & update_gradients)
- {
- Assert (data.shape_gradients[k].size() + offset <=
- fe_data.shape_gradients[k].size(),
- ExcInternalError());
- mapping.transform_covariant(data.shape_gradients[k].begin(),
- data.shape_gradients[k].end(),
- fe_data.shape_gradients[k].begin()+offset,
- mapping_data);
- };
- }
-
- if (flags & update_second_derivatives)
- this->compute_2nd (mapping, cell, offset, mapping_data, fe_data, data);
-}
-
-
-
-template <int dim>
-void
-FE_Q<dim>::fill_fe_subface_values (const Mapping<dim> &mapping,
- const typename DoFHandler<dim>::cell_iterator &cell,
- const unsigned int face,
- const unsigned int subface,
- const Quadrature<dim-1> &quadrature,
- typename Mapping<dim>::InternalDataBase &mapping_data,
- typename Mapping<dim>::InternalDataBase &fedata,
- FEValuesData<dim> &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<InternalData &> (fedata);
-
- // offset determines which data set
- // to take (all data sets for all
- // sub-faces are stored contiguously)
- const typename QProjector<dim>::DataSetDescriptor offset
- = (QProjector<dim>::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; k<this->dofs_per_cell; ++k)
- {
- if (flags & update_values)
- for (unsigned int i=0; i<quadrature.n_quadrature_points; ++i)
- data.shape_values(k,i) = fe_data.shape_values[k][i+offset];
-
- if (flags & update_gradients)
- {
- Assert (data.shape_gradients[k].size() + offset <=
- fe_data.shape_gradients[k].size(),
- ExcInternalError());
- mapping.transform_covariant(data.shape_gradients[k].begin(),
- data.shape_gradients[k].end(),
- fe_data.shape_gradients[k].begin()+offset,
- mapping_data);
- };
- }
-
- if (flags & update_second_derivatives)
- this->compute_2nd (mapping, cell, offset, mapping_data, fe_data, data);
-}
-
-
-
-template <int dim>
-unsigned int
-FE_Q<dim>::n_base_elements () const
-{
- return 1;
-}
-
-
-
-template <int dim>
-const FiniteElement<dim> &
-FE_Q<dim>::base_element (const unsigned int index) const
-{
- Assert (index==0, ExcIndexRange(index, 0, 1));
- return *this;
-}
-
-
-
-template <int dim>
-unsigned int
-FE_Q<dim>::element_multiplicity (const unsigned int index) const
-{
- Assert (index==0, ExcIndexRange(index, 0, 1));
- return 1;
-}
-
-
template <int dim>
bool