FE_Nedelec< dim, spacedim > Class Template Reference
[Finite element space descriptions]

Inheritance diagram for FE_Nedelec< dim, spacedim >:
Inheritance graph
[legend]

List of all members.

Classes

class  InternalData
struct  Matrices

Public Member Functions

 FE_Nedelec (const unsigned int p)
virtual std::string get_name () const
virtual double shape_value_component (const unsigned int i, const Point< dim > &p, const unsigned int component) const
virtual Tensor< 1, dim > shape_grad_component (const unsigned int i, const Point< dim > &p, const unsigned int component) const
virtual Tensor< 2, dim > shape_grad_grad_component (const unsigned int i, const Point< dim > &p, const unsigned int component) const
unsigned int get_degree () const
virtual unsigned int n_base_elements () const
virtual const FiniteElement
< dim, spacedim > & 
base_element (const unsigned int index) const
virtual unsigned int element_multiplicity (const unsigned int index) const
virtual bool has_support_on_face (const unsigned int shape_index, const unsigned int face_index) const
virtual unsigned int memory_consumption () const

Protected Member Functions

virtual FiniteElement< dim,
spacedim > * 
clone () const
virtual Mapping< dim, spacedim >
::InternalDataBase
get_data (const UpdateFlags, const Mapping< dim, spacedim > &mapping, const Quadrature< dim > &quadrature) const
virtual void fill_fe_values (const Mapping< dim, spacedim > &mapping, const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Quadrature< dim > &quadrature, typename Mapping< dim, spacedim >::InternalDataBase &mapping_internal, typename Mapping< dim, spacedim >::InternalDataBase &fe_internal, FEValuesData< dim, spacedim > &data, enum CellSimilarity::Similarity &cell_similarity) const
virtual void fill_fe_face_values (const Mapping< dim, spacedim > &mapping, const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const Quadrature< dim-1 > &quadrature, typename Mapping< dim, spacedim >::InternalDataBase &mapping_internal, typename Mapping< dim, spacedim >::InternalDataBase &fe_internal, FEValuesData< dim, spacedim > &data) const
virtual void fill_fe_subface_values (const Mapping< dim, spacedim > &mapping, const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const unsigned int sub_no, const Quadrature< dim-1 > &quadrature, typename Mapping< dim, spacedim >::InternalDataBase &mapping_internal, typename Mapping< dim, spacedim >::InternalDataBase &fe_internal, FEValuesData< dim, spacedim > &data) const

Private Member Functions

void initialize_constraints ()
void initialize_embedding ()
void initialize_restriction ()
void initialize_unit_support_points ()
void initialize_unit_face_support_points ()
virtual UpdateFlags update_once (const UpdateFlags flags) const
virtual UpdateFlags update_each (const UpdateFlags flags) const

Static Private Member Functions

static std::vector< unsigned intget_dpo_vector (const unsigned int degree)

Private Attributes

const unsigned int degree

Friends

class FE_Nedelec

Detailed Description

template<int dim, int spacedim = dim>
class FE_Nedelec< dim, spacedim >

Implementation of continuous Nedelec elements for the space H_curl. Note, however, that continuity only concerns the tangential component of the vector field.

The constructor of this class takes the degree p of this finite element. However, presently, only lowest order elements (i.e. p==1) are implemented. For a general overview of this element and its properties, see the report by Anna Schneebeli that is linked from the general documentation page of the library.

This class has not yet been implemented for the codimension one case (spacedim != dim).

Restriction on transformations

In some sense, the implementation of this element is not complete, but you will rarely notice. Here is the fact: since the element is vector-valued already on the unit cell, the Jacobian matrix (or its inverse) is needed already to generate the values of the shape functions on the cells in real space. This is in contrast to most other elements, where you only need the Jacobian for the gradients. Thus, to generate the gradients of Nedelec shape functions, one would need to have the derivatives of the inverse of the Jacobian matrix.

Basically, the Nedelec shape functions can be understood as the gradients of scalar shape functions on the real cell. They are thus the inverse Jacobian matrix times the gradients of scalar shape functions on the unit cell. The gradient of Nedelec shape functions is then, by the product rule, the sum of first the derivative (with respect to true coordinates) of the inverse Jacobian times the gradient (in unit coordinates) of the scalar shape function, plus second the inverse Jacobian times the derivative (in true coordinates) of the gradient (in unit coordinates) of the scalar shape functions. Note that each of the derivatives in true coordinates can be expressed as inverse Jacobian times gradient in unit coordinates.

The problem is the derivative of the inverse Jacobian. This rank-3 tensor can actually be computed (and we did so in very early versions of the library), but is a large task and very time consuming, so we dropped it. Since it is not available, we simply drop this first term.

What this means for the present case: first the computation of gradients of Nedelec shape functions is wrong. Second, you will not notice this usually, for two reasons:

The first reason is that the gradient of the Jacobian vanishes if the cells are mapped by an affine mapping, to which the usual bilinear mapping reduces if the cell is a parallelogram. Then the gradient of the shape functions is computed exact, since the first term is zero.

Second, with the Nedelec elements, you will usually want to compute the curl, and extract and sum up the respective elements of the full gradient tensor. However, the curl of the Jacobian vanishes, so for the curl of shape functions the first term is irrelevant, and the curl will be computed correctly as well.

Interpolation to finer and coarser meshes

Each finite element class in deal.II provides matrices that are used to interpolate from coarser to finer meshes and the other way round. Interpolation from a mother cell to its children is usually trivial, since finite element spaces are normally nested and this kind of interpolation is therefore exact. On the other hand, when we interpolate from child cells to the mother cell, we usually have to throw away some information.

For continuous elements, this transfer usually happens by interpolating the values on the child cells at the support points of the shape functions of the mother cell. However, for discontinuous elements, we often use a projection from the child cells to the mother cell. The projection approach is only possible for discontinuous elements, since it cannot be guaranteed that the values of the projected functions on one cell and its neighbor match. In this case, only an interpolation can be used. (Internally, whether the values of a shape function are interpolated or projected, or better: whether the matrices the finite element provides are to be treated with the properties of a projection or of an interpolation, is controlled by the restriction_is_additive flag. See there for more information.)

Here, things are not so simple: since the element has some continuity requirements across faces, we can only resort to some kind of interpolation. On the other hand, for the lowest order elements, the values of generating functionals are the (constant) tangential values of the shape functions. We would therefore really like to take the mean value of the tangential values of the child faces, and make this the value of the mother face. Then, however, taking a mean value of two piecewise constant function is not an interpolation, but a restriction. Since this is not possible, we cannot use this.

To make a long story somewhat shorter, when interpolating from refined edges to a coarse one, we do not take the mean value, but pick only one (the one from the first child edge). While this is not optimal, it is certainly a valid choice (using an interpolation point that is not in the middle of the cell, but shifted to one side), and it also preserves the order of the interpolation.

Numbering of the degrees of freedom (DoFs)

Nedelec elements have their degrees of freedom on edges, with shape functions being vector valued and pointing in tangential direction. We use the standard enumeration and direction of edges in deal.II, yielding the following shape functions in 2d:

 *       3
 *    2-->--3
 *    |     |
 *   0^     ^1
 *    |     |
 *    0-->--1
 *       2
 *   

For the 3d case, the ordering follows the same scheme: the lines are numbered as described in the documentation of the Triangulation class, i.e.

 *       *---7---*        *---7---*
 *      /|       |       /       /|
 *     4 |       11     4       5 11
 *    /  10      |     /       /  |
 *   *   |       |    *---6---*   |
 *   |   *---3---*    |       |   *
 *   |  /       /     |       9  /
 *   8 0       1      8       | 1
 *   |/       /       |       |/
 *   *---2---*        *---2---*
 *   

and their directions are as follows:

 *         *--->---*        *--->---*
 *        /|       |       /       /|
 *       ^ |       ^      ^       ^ ^
 *      /  ^       |     /       /  |
 *     *   |       |    *--->---*   |
 *     |   *--->---*    |       |   *
 *     |  /       /     |       ^  /
 *     ^ ^       ^      ^       | ^
 *     |/       /       |       |/
 *     *--->---*        *--->---*
 *   

The element does not make much sense in 1d, so it is not implemented there.

Author:
Wolfgang Bangerth, Anna Schneebeli, 2002, 2003

Constructor & Destructor Documentation

template<int dim, int spacedim = dim>
FE_Nedelec< dim, spacedim >::FE_Nedelec ( const unsigned int  p  ) 

Constructor for the Nedelec element of degree p.


Member Function Documentation

template<int dim, int spacedim = dim>
virtual std::string FE_Nedelec< dim, spacedim >::get_name (  )  const [virtual]

Return a string that uniquely identifies a finite element. This class returns FE_Nedelec<dim>(degree), with dim and degree replaced by appropriate values.

Implements FiniteElement< dim, spacedim >.

template<int dim, int spacedim = dim>
virtual double FE_Nedelec< dim, spacedim >::shape_value_component ( const unsigned int  i,
const Point< dim > &  p,
const unsigned int  component 
) const [virtual]

Return the value of the componentth vector component of the ith shape function at the point p. See the FiniteElement base class for more information about the semantics of this function.

Reimplemented from FiniteElement< dim, spacedim >.

template<int dim, int spacedim = dim>
virtual Tensor<1,dim> FE_Nedelec< dim, spacedim >::shape_grad_component ( const unsigned int  i,
const Point< dim > &  p,
const unsigned int  component 
) const [virtual]

Return the gradient of the componentth vector component of the ith shape function at the point p. See the FiniteElement base class for more information about the semantics of this function.

Reimplemented from FiniteElement< dim, spacedim >.

template<int dim, int spacedim = dim>
virtual Tensor<2,dim> FE_Nedelec< dim, spacedim >::shape_grad_grad_component ( const unsigned int  i,
const Point< dim > &  p,
const unsigned int  component 
) const [virtual]

Return the second derivative of the componentth vector component of the ith shape function at the point p. See the FiniteElement base class for more information about the semantics of this function.

Reimplemented from FiniteElement< dim, spacedim >.

template<int dim, int spacedim = dim>
unsigned int FE_Nedelec< dim, spacedim >::get_degree (  )  const

Return the polynomial degree of this finite element, i.e. the value passed to the constructor.

template<int dim, int spacedim = dim>
virtual unsigned int FE_Nedelec< dim, spacedim >::n_base_elements (  )  const [virtual]

Number of base elements in a mixed discretization. Here, this is of course equal to one.

Implements FiniteElement< dim, spacedim >.

template<int dim, int spacedim = dim>
virtual const FiniteElement<dim,spacedim>& FE_Nedelec< dim, spacedim >::base_element ( const unsigned int  index  )  const [virtual]

Access to base element objects. Since this element is atomic, base_element(0) is this, and all other indices throw an error.

Implements FiniteElement< dim, spacedim >.

template<int dim, int spacedim = dim>
virtual unsigned int FE_Nedelec< dim, spacedim >::element_multiplicity ( const unsigned int  index  )  const [virtual]

Multiplicity of base element index. Since this is an atomic element, element_multiplicity(0) returns one, and all other indices will throw an error.

Implements FiniteElement< dim, spacedim >.

template<int dim, int spacedim = dim>
virtual bool FE_Nedelec< dim, spacedim >::has_support_on_face ( const unsigned int  shape_index,
const unsigned int  face_index 
) const [virtual]

This function returns true, if the shape function shape_index has non-zero values on the face face_index. For the lowest order Nedelec elements, this is actually the case for the one on which the shape function is defined and all neighboring ones.

Implementation of the interface in FiniteElement

Reimplemented from FiniteElement< dim, spacedim >.

template<int dim, int spacedim = dim>
virtual unsigned int FE_Nedelec< dim, spacedim >::memory_consumption (  )  const [virtual]

Determine an estimate for the memory consumption (in bytes) of this object.

This function is made virtual, since finite element objects are usually accessed through pointers to their base class, rather than the class itself.

Reimplemented from FiniteElement< dim, spacedim >.

template<int dim, int spacedim = dim>
virtual FiniteElement<dim,spacedim>* FE_Nedelec< dim, spacedim >::clone (  )  const [protected, virtual]

clone function instead of a copy constructor.

This function is needed by the constructors of FESystem.

Implements FiniteElement< dim, spacedim >.

template<int dim, int spacedim = dim>
virtual Mapping<dim,spacedim>::InternalDataBase* FE_Nedelec< dim, spacedim >::get_data ( const   UpdateFlags,
const Mapping< dim, spacedim > &  mapping,
const Quadrature< dim > &  quadrature 
) const [protected, virtual]

Prepare internal data structures and fill in values independent of the cell.

Implements FiniteElement< dim, spacedim >.

template<int dim, int spacedim = dim>
virtual void FE_Nedelec< dim, spacedim >::fill_fe_values ( const Mapping< dim, spacedim > &  mapping,
const typename Triangulation< dim, spacedim >::cell_iterator &  cell,
const Quadrature< dim > &  quadrature,
typename Mapping< dim, spacedim >::InternalDataBase mapping_internal,
typename Mapping< dim, spacedim >::InternalDataBase fe_internal,
FEValuesData< dim, spacedim > &  data,
enum CellSimilarity::Similarity cell_similarity 
) const [protected, virtual]

Implementation of the same function in FiniteElement.

Implements FiniteElement< dim, spacedim >.

template<int dim, int spacedim = dim>
virtual void FE_Nedelec< dim, spacedim >::fill_fe_face_values ( const Mapping< dim, spacedim > &  mapping,
const typename Triangulation< dim, spacedim >::cell_iterator &  cell,
const unsigned int  face_no,
const Quadrature< dim-1 > &  quadrature,
typename Mapping< dim, spacedim >::InternalDataBase mapping_internal,
typename Mapping< dim, spacedim >::InternalDataBase fe_internal,
FEValuesData< dim, spacedim > &  data 
) const [protected, virtual]

Implementation of the same function in FiniteElement.

Implements FiniteElement< dim, spacedim >.

template<int dim, int spacedim = dim>
virtual void FE_Nedelec< dim, spacedim >::fill_fe_subface_values ( const Mapping< dim, spacedim > &  mapping,
const typename Triangulation< dim, spacedim >::cell_iterator &  cell,
const unsigned int  face_no,
const unsigned int  sub_no,
const Quadrature< dim-1 > &  quadrature,
typename Mapping< dim, spacedim >::InternalDataBase mapping_internal,
typename Mapping< dim, spacedim >::InternalDataBase fe_internal,
FEValuesData< dim, spacedim > &  data 
) const [protected, virtual]

Implementation of the same function in FiniteElement.

Implements FiniteElement< dim, spacedim >.

template<int dim, int spacedim = dim>
static std::vector<unsigned int> FE_Nedelec< dim, spacedim >::get_dpo_vector ( const unsigned int  degree  )  [static, private]

Only for internal use. Its full name is get_dofs_per_object_vector function and it creates the dofs_per_object vector that is needed within the constructor to be passed to the constructor of FiniteElementData.

template<int dim, int spacedim = dim>
void FE_Nedelec< dim, spacedim >::initialize_constraints (  )  [private]

Initialize the hanging node constraints matrices. Called from the constructor.

template<int dim, int spacedim = dim>
void FE_Nedelec< dim, spacedim >::initialize_embedding (  )  [private]

Initialize the embedding matrices. Called from the constructor.

template<int dim, int spacedim = dim>
void FE_Nedelec< dim, spacedim >::initialize_restriction (  )  [private]

Initialize the restriction matrices. Called from the constructor.

template<int dim, int spacedim = dim>
void FE_Nedelec< dim, spacedim >::initialize_unit_support_points (  )  [private]

Initialize the unit_support_points field of the FiniteElement class. Called from the constructor.

template<int dim, int spacedim = dim>
void FE_Nedelec< dim, spacedim >::initialize_unit_face_support_points (  )  [private]

Initialize the unit_face_support_points field of the FiniteElement class. Called from the constructor.

template<int dim, int spacedim = dim>
virtual UpdateFlags FE_Nedelec< dim, spacedim >::update_once ( const UpdateFlags  flags  )  const [private, virtual]

Given a set of flags indicating what quantities are requested from a 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 flags and update_values.

For the present kind of finite element, this is exactly the case.

Implements FiniteElement< dim, spacedim >.

template<int dim, int spacedim = dim>
virtual UpdateFlags FE_Nedelec< dim, spacedim >::update_each ( const UpdateFlags  flags  )  const [private, virtual]

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.

Implements FiniteElement< dim, spacedim >.


Friends And Related Function Documentation

template<int dim, int spacedim = dim>
friend class FE_Nedelec [friend]

Allow access from other dimensions.


Member Data Documentation

template<int dim, int spacedim = dim>
const unsigned int FE_Nedelec< dim, spacedim >::degree [private]

Degree of the polynomials.

Reimplemented from FiniteElementData< dim >.


The documentation for this class was generated from the following file:

deal.II documentation generated on Mon Nov 23 22:57:41 2009 by doxygen 1.6.1