dune-grid  2.2.0
Public Types | Public Member Functions | Static Public Attributes
Dune::SGeometry< mydim, cdim, GridImp > Class Template Reference

#include <dune/grid/sgrid.hh>

Inheritance diagram for Dune::SGeometry< mydim, cdim, GridImp >:
Dune::GeometryDefaultImplementation< mydim, cdim, GridImp, SGeometry >

List of all members.

Public Types

typedef GridImp::ctype ctype
 define type used for coordinates in grid module
typedef FieldVector< ctype, mydim > LocalCoordinate
typedef FieldVector< ctype, cdim > GlobalCoordinate
typedef FieldMatrix< ctype,
cdim, mydim > 
Jacobian
 type of jacobian (also of jacobian inverse transposed)
typedef FieldMatrix< ctype,
mydim, cdim > 
JacobianTransposed
 type of jacobian transposed

Public Member Functions

GeometryType type () const
 return the element type identifier
bool affine () const
 here we have always an affine geometry
int corners () const
 return the number of corners of this element. Corners are numbered 0...n-1
FieldVector< ctype, cdim > corner (const int i) const
 return i'th corner of the geometry
FieldVector< ctype, cdim > center () const
 return center of the geometry
FieldVector< ctype, cdim > global (const FieldVector< ctype, mydim > &local) const
 maps a local coordinate within reference element to global coordinate in element
FieldVector< ctype, mydim > local (const FieldVector< ctype, cdim > &global) const
 maps a global coordinate within the element to a local coordinate in its reference element
ctype integrationElement (const FieldVector< ctype, mydim > &local) const
ctype volume () const
 return volume of geometry
const FieldMatrix< ctype,
mydim, cdim > & 
jacobianTransposed (const FieldVector< ctype, mydim > &local) const
const FieldMatrix< ctype, cdim,
mydim > & 
jacobianInverseTransposed (const FieldVector< ctype, mydim > &local) const
void print (std::ostream &ss, int indent) const
 print internal data
void make (FieldMatrix< ctype, mydim+1, cdim > &__As)
 SGeometry ()
 constructor

Static Public Attributes

static const int mydimension
static const int coorddimension

Detailed Description

template<int mydim, int cdim, class GridImp>
class Dune::SGeometry< mydim, cdim, GridImp >

SGeometry realizes the concept of the geometric part of a mesh entity.

The geometric part of a mesh entity is a $d$-dimensional object in $\mathbf{R}^w$ where $d$ corresponds the template parameter dim and $w$ corresponds to the template parameter dimworld.

The $d$-dimensional object is a polyhedron given by a certain number of corners, which are vectors in $\mathbf{R}^w$.

The member function global provides a map from a topologically equivalent polyhedron ("reference element") in $\mathbf{R}^d$ to the given polyhedron. This map can be inverted by the member function local, where an appropriate projection is applied first, when $d\neq w$.

In the case of a structured mesh discretizing a generalized cube this map is linear and can be described as

\[ g(l) = s + \sum\limits_{i=0}^{d-1} l_ir^i\]

where $s\in\mathbf{R}^w$ is a given position vector, the $r^i\in\mathbf{R}^w$ are given direction vectors and $l\in\mathbf{R}^d$ is a local coordinate within the reference polyhedron. The direction vectors are assumed to be orthogonal with respect to the standard Eucliden inner product.

The $d$-dimensional reference polyhedron is given by the points $\{ (x_0,\ldots,x_{d-1}) \ | \ x_i\in\{0,1\}\ \}$.

In order to invert the map for a point $p$, we have to find a local coordinate $l$ such that $g(l)=p$. Of course this is only possible if $d=w$. In the general case $d\leq w$ we determine $l$ such that

\[(s,r^k) + \sum\limits_{i=0}^{d-1} l_i (r^i,r^k) = (p,r^k) \ \ \ \forall k=0,\ldots,d-1. \]

The resulting system is diagonal since the direction vectors are required to be orthogonal.


Member Typedef Documentation

template<int mydim, int cdim, class GridImp >
typedef GridImp::ctype Dune::SGeometry< mydim, cdim, GridImp >::ctype

define type used for coordinates in grid module

Reimplemented from Dune::GeometryDefaultImplementation< mydim, cdim, GridImp, SGeometry >.

typedef FieldVector< ctype, cdim > Dune::GeometryDefaultImplementation< mydim, cdim, GridImp, SGeometry >::GlobalCoordinate [inherited]
typedef FieldMatrix<ctype,cdim,mydim> Dune::GeometryDefaultImplementation< mydim, cdim, GridImp, SGeometry >::Jacobian [inherited]

type of jacobian (also of jacobian inverse transposed)

typedef FieldMatrix< ctype, mydim, cdim > Dune::GeometryDefaultImplementation< mydim, cdim, GridImp, SGeometry >::JacobianTransposed [inherited]

type of jacobian transposed

typedef FieldVector< ctype, mydim > Dune::GeometryDefaultImplementation< mydim, cdim, GridImp, SGeometry >::LocalCoordinate [inherited]

Constructor & Destructor Documentation

template<int mydim, int cdim, class GridImp >
Dune::SGeometry< mydim, cdim, GridImp >::SGeometry ( ) [inline]

constructor


Member Function Documentation

template<int mydim, int cdim, class GridImp >
bool Dune::SGeometry< mydim, cdim, GridImp >::affine ( ) const [inline]

here we have always an affine geometry

template<int mydim, int cdim, class GridImp >
FieldVector<ctype, cdim > Dune::SGeometry< mydim, cdim, GridImp >::center ( ) const [inline]

return center of the geometry

Reimplemented from Dune::GeometryDefaultImplementation< mydim, cdim, GridImp, SGeometry >.

template<int mydim, int cdim, class GridImp >
FieldVector< ctype, cdim > Dune::SGeometry< mydim, cdim, GridImp >::corner ( const int  i) const [inline]

return i'th corner of the geometry

Referenced by Dune::SGeometry< 0, cdim, GridImp >::global().

template<int mydim, int cdim, class GridImp >
int Dune::SGeometry< mydim, cdim, GridImp >::corners ( ) const [inline]

return the number of corners of this element. Corners are numbered 0...n-1

template<int mydim, int cdim, class GridImp >
FieldVector<ctype, cdim> Dune::SGeometry< mydim, cdim, GridImp >::global ( const FieldVector< ctype, mydim > &  local) const

maps a local coordinate within reference element to global coordinate in element

template<int mydim, int cdim, class GridImp >
ctype Dune::SGeometry< mydim, cdim, GridImp >::integrationElement ( const FieldVector< ctype, mydim > &  local) const [inline]

Integration over a general element is done by integrating over the reference element and using the transformation from the reference element to the global element as follows:

\[\int\limits_{\Omega_e} f(x) dx = \int\limits_{\Omega_{ref}} f(g(l)) A(l) dl \]

where $g$ is the local to global mapping and $A(l)$ is the integration element.

For a general map $g(l)$ involves partial derivatives of the map (surface element of the first kind if $d=2,w=3$, determinant of the Jacobian of the transformation for $d=w$, $\|dg/dl\|$ for $d=1$).

For linear elements, the derivatives of the map with respect to local coordinates do not depend on the local coordinates and are the same over the whole element.

For a structured mesh where all edges are parallel to the coordinate axes, the computation is the length, area or volume of the element is very simple to compute.

Each grid module implements the integration element with optimal efficieny. This will directly translate in substantial savings in the computation of finite element stiffness matrices.

References Dune::SGeometry< mydim, cdim, GridImp >::volume().

template<int mydim, int cdim, class GridImp >
const FieldMatrix<ctype,cdim,mydim>& Dune::SGeometry< mydim, cdim, GridImp >::jacobianInverseTransposed ( const FieldVector< ctype, mydim > &  local) const
template<int mydim, int cdim, class GridImp >
const FieldMatrix<ctype, mydim, cdim >& Dune::SGeometry< mydim, cdim, GridImp >::jacobianTransposed ( const FieldVector< ctype, mydim > &  local) const
template<int mydim, int cdim, class GridImp >
FieldVector<ctype, mydim> Dune::SGeometry< mydim, cdim, GridImp >::local ( const FieldVector< ctype, cdim > &  global) const

maps a global coordinate within the element to a local coordinate in its reference element

template<int mydim, int cdim, class GridImp >
void Dune::SGeometry< mydim, cdim, GridImp >::make ( FieldMatrix< ctype, mydim+1, cdim > &  __As)

The first dim columns of As contain the dim direction vectors. Column dim is the position vector. This format allows a consistent treatement of all dimensions, including 0 (the vertex).

template<int mydim, int cdim, class GridImp >
void Dune::SGeometry< mydim, cdim, GridImp >::print ( std::ostream &  ss,
int  indent 
) const

print internal data

template<int mydim, int cdim, class GridImp >
GeometryType Dune::SGeometry< mydim, cdim, GridImp >::type ( ) const [inline]

return the element type identifier

References Dune::cube.

template<int mydim, int cdim, class GridImp >
ctype Dune::SGeometry< mydim, cdim, GridImp >::volume ( ) const

Member Data Documentation

const int Dune::GeometryDefaultImplementation< mydim, cdim, GridImp, SGeometry >::coorddimension [static, inherited]
const int Dune::GeometryDefaultImplementation< mydim, cdim, GridImp, SGeometry >::mydimension [static, inherited]

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