Dune::EdgeS03DLocalBasis< D, R > Class Template Reference

Experimental lowest order edge elements for tetrahedrons. More...

#include <dune/localfunctions/whitney/edges03d/edges03dlocalbasis.hh>

List of all members.

Public Types

typedef LocalBasisTraits< D,
3, Dune::FieldVector< D, 3 >
, R, 3, Dune::FieldVector< R, 3 >
, Dune::FieldMatrix< R, 3, 3 > > 
Traits
 export type traits for function signature

Public Member Functions

 EdgeS03DLocalBasis ()
 contruct a local basis instance with default orientations
 EdgeS03DLocalBasis (unsigned int orientations)
unsigned int size () const
 number of shape functions
template<typename Geometry >
void evaluateFunctionGlobal (const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out, const Geometry &geometry) const
 Get global values of shape functions.
template<typename Geometry >
void evaluateJacobianGlobal (const typename Traits::DomainType &in, std::vector< typename Traits::JacobianType > &out, const Geometry &geometry) const
 Evaluate global Jacobian of all shape functions.
unsigned int order () const
 Polynomial order of the shape functions.

Detailed Description

template<class D, class R>
class Dune::EdgeS03DLocalBasis< D, R >

Experimental lowest order edge elements for tetrahedrons.

(S for simplex)

Template Parameters:
D Type to represent the field in the domain.
R Type to represent the field in the range.
Note:
This class does not implement the usual LocalBasisInterface since that does not make much sense for vector valued elements. An experimental interface providing global rather than local values is provided instead. Be warned that this interface is subject to change without notice, however.

Member Typedef Documentation

template<class D , class R >
typedef LocalBasisTraits< D,3,Dune::FieldVector<D,3>, R,3,Dune::FieldVector<R,3>, Dune::FieldMatrix<R,3,3> > Dune::EdgeS03DLocalBasis< D, R >::Traits

export type traits for function signature


Constructor & Destructor Documentation

template<class D , class R >
Dune::EdgeS03DLocalBasis< D, R >::EdgeS03DLocalBasis (  )  [inline]

contruct a local basis instance with default orientations

template<class D , class R >
Dune::EdgeS03DLocalBasis< D, R >::EdgeS03DLocalBasis ( unsigned int  orientations  )  [inline]

contruct a local basis instance with the given orientations

Parameters:
orientations Bit-map of orientations for each shape function; bit 0 = 0 means default orientation for the first shape function, bit 0 = 1 means inverted orientation for the first shape function.

Member Function Documentation

template<class D , class R >
template<typename Geometry >
void Dune::EdgeS03DLocalBasis< D, R >::evaluateFunctionGlobal ( const typename Traits::DomainType in,
std::vector< typename Traits::RangeType > &  out,
const Geometry &  geometry 
) const [inline]

Get global values of shape functions.

Experimental interface for getting global values of shape functions.

Vector valued shape functions often require a complicated transformation to get the global values from the local values:

\[ \psi_i(g(\mathbf{\hat x}))=L(\hat\psi_i(\mathbf x)) \]

where $g$ is transformation for coordinates, $L$ the tranformation for values and the hat $\hat{\phantom x}$ denotes local quantities. Worse, if we require for edge elements that

\[ \psi_i\cdot\mathbf t_j=\delta_{ij}\mbox{ on edge }j \]

and

\[ \hat\psi_i\cdot\mathbf{\hat t}_j=\delta_{ij}\mbox{ on edge }j \]

then we need a different $L$ for each shape function. This property is very desirable be eases debugging and the construction of constraints for grids with hanging nodes.

Template Parameters:
Geometry Type of geometry
Parameters:
[in] in Where to evaluate. NOTE: these are local coordinates
[out] out The global values of the shape functions.
[in] geometry The geometry used for the local to global mapping.
template<class D , class R >
template<typename Geometry >
void Dune::EdgeS03DLocalBasis< D, R >::evaluateJacobianGlobal ( const typename Traits::DomainType in,
std::vector< typename Traits::JacobianType > &  out,
const Geometry &  geometry 
) const [inline]

Evaluate global Jacobian of all shape functions.

Experimental interface for getting global values of the Jacobian of the shapefunctions.

This calculates the global Jacobian

\[ \mathrm J(\psi^i)= \left(\begin{array}{ccc} \partial\psi^i_0/\partial x_0 & \partial\psi^i_0/\partial x_1 & \cdots \\ \partial\psi^i_1/\partial x_0 & \partial\psi^i_1/\partial x_1 & \cdots \\ \vdots & \vdots & \ddots \end{array}\right) \]

Note that this are the derivatives of the global values by the global coordinates, evaluated at local coordinates.

Template Parameters:
Geometry Type of geometry
Parameters:
[in] in Where to evaluate. NOTE: these are local coordinates
[out] out The values of the global Jacobian.
[in] geometry The geometry used for the local to global mapping.
template<class D , class R >
unsigned int Dune::EdgeS03DLocalBasis< D, R >::order (  )  const [inline]

Polynomial order of the shape functions.

template<class D , class R >
unsigned int Dune::EdgeS03DLocalBasis< D, R >::size (  )  const [inline]

number of shape functions


The documentation for this class was generated from the following file:
Generated on Sat Apr 24 11:15:40 2010 for dune-localfunctions by  doxygen 1.6.3