Finite Element Module


Modules

 Equations

Classes

class  imaging::Assembler
 Assembles the stiffness matrix and force vector of a FE problem. More...
class  imaging::ElementIntegrator< N, N_NODES >
 Abstract base class of all integrator classes. More...
class  imaging::SquareIntegrator< N_NODES >
 Gaussian integrator on the square. More...
class  imaging::CubeIntegrator< N_NODES >
 Gaussian integrator on the square. More...
class  imaging::TriangleIntegrator< N_NODES >
 Gaussian integrator on the triangle. More...
class  imaging::TetrahedraIntegrator< N_NODES >
 Gaussian integrator on the tetrahedron. More...
class  imaging::IntervalIntegrator< N_NODES >
 Gaussian integrator on the symmetric unit interval. More...
class  imaging::PointIntegrator
 (Gaussian) integrator on a point. More...
class  imaging::FeFunctionInterface
 Abstract base class of FE approximations of functions. More...
class  imaging::fem_1d_types
 FE types for a 1-dimensional FE problem. More...
class  imaging::fem_2d_square_types
 FE types for a 2-dimensional FE problem with rectangular elements. More...
class  imaging::fem_2d_triangle_types
 FE types for a 2-dimensional FE problem with triangular elements. More...
class  imaging::fem_3d_cube_types
 FE types for a 3-dimensional FE problem with cube elements. More...
class  imaging::fem_3d_tetrahedra_types
 FE types for a 3-dimensional FE problem with tetrahedral elements. More...
class  imaging::FemKernel< fem_types >
 Computes the values of the shape functions and the element transformation in the integration nodes. More...
class  imaging::Grid< fem_types >
 Provides the data structure and accessor functions for a FE grid. More...
class  imaging::Image2Grid< fem_types >
 Constructs FE grids from multi-dimensional images and provides conversion functions between grid and image. More...
class  imaging::ShapeFunction< N_NODES, N_FACE_NODES, N >
 Abstract base class of all shape functions classes. More...
class  imaging::Bilinear2dShapeFunction
 Bilinear shape function on the square. More...
class  imaging::Linear2dShapeFunction
 Linear shape function on the triangle. More...
class  imaging::Linear3dShapeFunction
 Linear shape function on the tetrahedra. More...
class  imaging::Linear1dShapeFunction
 Linear shape function on the symmetric unit interval. More...
class  imaging::SimpleAssembler
 Assembles the stiffness matrix and force vector of a FE problem implementing SimpleEquationInterface. More...
class  imaging::Transform< N_VERTICES, N_FACES, N >
 Abstract base class of all transformations of the reference element to an element of the FE grid. More...
class  imaging::Square2dTransform
 Transformation of the square reference element. More...
class  imaging::Triangle2dTransform
 Transformation of the triangle reference element. More...
class  imaging::Tetrahedra3dTransform
 Transformation of the tetrahedra reference element. More...
class  imaging::Cube3dTransform
 Transformation of the symmetric unit interval reference element. More...
class  imaging::Interval1dTransform
 Transformation of the symmetric unit cube reference element. TODO:. More...

Functions

template<class fem_types, class function_t>
function_t::value_t imaging::integrate (const function_t &function, const Grid< fem_types > &grid)
template<class fem_types, class function_t>
function_t::value_t imaging::maximum (const function_t &function, const Grid< fem_types > &grid)
template<class fem_types, class function_t>
function_t::value_t imaging::minimum (const function_t &function, const Grid< fem_types > &grid)
template<class fem_types, class function_t>
function_t::value_t imaging::average (const function_t &function, const Grid< fem_types > &grid)
void imaging::uniform_grid (float_t lower_bound, float_t upper_bound, std::size_t n_elements, Grid< fem_1d_types > &grid, ublas::compressed_matrix< float_t > &stiffness_matrix_prototype, std::size_t system_size)
void imaging::circle_grid (float_t radius, std::size_t n_rings, Grid< fem_2d_triangle_types > &grid, ublas::compressed_matrix< float_t > &stiffness_matrix_prototype, std::size_t system_size)
void imaging::ellipse_grid (float_t a, float_t b, std::size_t n_rings, Grid< fem_2d_triangle_types > &grid, ublas::compressed_matrix< float_t > &stiffness_matrix_prototype, std::size_t system_size)
void imaging::triangulate_shape (const BoundaryDiscretizer< 2 > &shape_discretizer, float_t max_triangle_area, Grid< fem_2d_triangle_types > &grid, ublas::compressed_matrix< float_t > &stiffness_matrix_prototype, std::size_t system_size)

Detailed Description

Introduction

The Finite Element Module implements data structurs and methods to assembly the stiffness matrix and the force vector for a given equation and a given finite element grid. The classes of this module can be roughly divided into problem-centric and implementation-centric classes. In addition, there some utilities classes for frequent FE-related tasks.

The problem-centric classes provide data structure and functionality which are specific to a given FE problem and can only be partially provided by the imaging2 library. These cover

The implementation-centric classes provide the part of FE-code which does not depend on the actual geometry and equation. An implementation of implementation-centric classes is provided by the imaging2 library but it can be extended by the user. It consists of

The central idea of the FEM module is to be generic. I.e. the implementation is independent of the problem dimension and the actual type of elements and shape functions wherever possible. This is done by parametrizing all this information in classes called FEM traits. These trait classes encode all information regarding the dimension of the FE problem, the type of elements and shape functions used and the integrators for the numeric evaluation of the shape functions over the elements. An example of FEM traits for the linear approximation of a 1D FE problem is given below:

  class fem_1d_types
  {
  public:
    static const std::size_t data_dimension = 1;

    typedef Interval1dTransform transform_t;
    typedef Linear1dShapeFunction shape_function_t;
    typedef IntervalIntegrator<2> integrator_t;
    typedef PointIntegrator boundary_integrator_t;
  };

Classes which are designed to be independent of the actual choice of the above parameters are implemented as templates which will be instantiated for specific FEM traits. E.g. the code of the class Assembler does refer to the dimension of the FE problem as fem_types::data_dimension only. The compiler substitutes this expression, when Assembler is instantiated with specified FEM traits, e.g. Assembler<fem_1d_types>.

For a given FE problem the user first selects the right FEM traits class (e.g. fem_2d_square_types for a problem on a 2D image with square elements) and instantiates all components she needs for the FE computation with this trait class. This is shown in the following code example:

  img::Grid<img::fem_2d_square_types> grid;
  // construct grid...
  
  SomeEquation<img::fem_2d_square_types> equation;
  // provide data (i.e. boundary data, right hand side) to the equation...
  
  img::Assembler<img::fem_2d_square_types>
  img::ublas::compressed_matrix<img::float_t> stiffness_matrix;
  img::ublas::vector<img::float_t> force_vector;
  
  assembler.assemble(equation, grid, stiffness_matrix, force_vector);
  
  // solve the system of linear equations

This example also illustrates the differentiation into problem-centric and implementation-centric classes we introduced above. Although a complete data structure for the FE grid is provided, the actual construction of the grid depends on the geometry of the problem and has to be done by the user (although she can use helper functions provided by imaging2). Even more manual interaction is required for the equation. The user has to provide a class implementing the equation she wants to solve, and is in addition responsible for providing the equation with the data of the specific problem to solve. Both tasks are clearly problem-centric.

Instead, the class Assembler and its member assemble() is provided by the imaging2 library and can be used without modification for any choice of FEM traits. This class is implementation-centric.


Function Documentation

template<class fem_types, class function_t>
function_t::value_t imaging::average ( const function_t &  function,
const Grid< fem_types > &  grid 
) [inline]

void imaging::circle_grid ( float_t  radius,
std::size_t  n_rings,
Grid< fem_2d_triangle_types > &  grid,
ublas::compressed_matrix< float_t > &  stiffness_matrix_prototype,
std::size_t  system_size = 1 
)

#include <fem/utilities.hpp>

Constructs a circular grid of radius. The triangular elements are laid out in n_rings rings. The i-th ring contains 4 * i elements. In addition, stiffness_matrix_prototype is resized to the correct size for grid. In case the PDE to be solved is not scalar but a system of equation, the user has to pass the number of equations (system_size) to ensure that stiffness_matrix_prototype is sized correctly.

References imaging::ellipse_grid().

void imaging::ellipse_grid ( float_t  a,
float_t  b,
std::size_t  n_rings,
Grid< fem_2d_triangle_types > &  grid,
ublas::compressed_matrix< float_t > &  stiffness_matrix_prototype,
std::size_t  system_size = 1 
)

#include <fem/utilities.hpp>

Constructs a elliptic grid with axis of length a and b. The triangular elements are laid out in n_rings rings. The i-th ring contains 4 * i elements. In addition, stiffness_matrix_prototype is resized to the correct size for grid. In case the PDE to be solved is not scalar but a system of equation, the user has to pass the number of equations (system_size) to ensure that stiffness_matrix_prototype is sized correctly.

References imaging::max(), imaging::PI, imaging::Grid< fem_types >::set_boundary_element(), imaging::Grid< fem_types >::set_dimensions(), imaging::Grid< fem_types >::set_global_node_index(), imaging::Grid< fem_types >::set_global_vertex_index(), and imaging::Grid< fem_types >::set_vertex().

Referenced by imaging::circle_grid().

template<class fem_types, class function_t>
function_t::value_t imaging::integrate ( const function_t &  function,
const Grid< fem_types > &  grid 
) [inline]

template<class fem_types, class function_t>
function_t::value_t imaging::maximum ( const function_t &  function,
const Grid< fem_types > &  grid 
) [inline]

template<class fem_types, class function_t>
function_t::value_t imaging::minimum ( const function_t &  function,
const Grid< fem_types > &  grid 
) [inline]

void imaging::triangulate_shape ( const BoundaryDiscretizer< 2 > &  shape_discretizer,
float_t  max_triangle_area,
Grid< fem_2d_triangle_types > &  grid,
ublas::compressed_matrix< float_t > &  stiffness_matrix_prototype,
std::size_t  system_size = 1 
)

#include <fem/utilities.hpp>

Triangulates an arbitrary, planar shape given by shape_discretizer and computes a grid from the triangulation. The maximum area of each triangle will not exceed max_triangle_area. In addition, stiffness_matrix_prototype is resized to the correct size for grid. In case the PDE to be solved is not scalar but a system of equation, the user has to pass the number of equations (system_size) to ensure that stiffness_matrix_prototype is sized correctly.

References imaging::Grid< fem_types >::n_elements(), imaging::BoundaryDiscretizer< N >::n_points(), imaging::Grid< fem_types >::n_vertices(), imaging::Grid< fem_types >::set_boundary_element(), imaging::Grid< fem_types >::set_dimensions(), imaging::Grid< fem_types >::set_global_node_index(), imaging::Grid< fem_types >::set_global_vertex_index(), and imaging::Grid< fem_types >::set_vertex().

void imaging::uniform_grid ( float_t  lower_bound,
float_t  upper_bound,
std::size_t  n_elements,
Grid< fem_1d_types > &  grid,
ublas::compressed_matrix< float_t > &  stiffness_matrix_prototype,
std::size_t  system_size = 1 
)

#include <fem/utilities.hpp>

Constructs a regular grid (i.e. with elements of constant size) on the interval defined by lower_bound and upper_bound for a system of equations. In addition, stiffness_matrix_prototype is resized to the correct size for grid. In case the PDE to be solved is not scalar but a system of equation, the user has to pass the number of equations (system_size) to ensure that stiffness_matrix_prototype is sized correctly.

References imaging::Grid< fem_types >::set_boundary_element(), imaging::Grid< fem_types >::set_dimensions(), imaging::Grid< fem_types >::set_global_node_index(), imaging::Grid< fem_types >::set_global_vertex_index(), and imaging::Grid< fem_types >::set_vertex().


Generated on Tue Feb 10 10:01:31 2009 for imaging2 by  doxygen 1.5.5