KASKADE 7 development version
|
Classes to define finite element spaces and functions on them. More...
Modules | |
Working with finite element functions | |
Tools for representing and working with finite element functions. | |
Grid Management | |
Tools for managing grids and their refinement, and for iterating over grids. | |
Classes | |
class | Kaskade::BoundaryFace< Grid, Face, Displacement, dimw > |
A class for representing displaced/deformed boundary faces. More... | |
class | Kaskade::BoundaryLocator< GridView, Function, dimw > |
class | Kaskade::BoundaryMapper< Implementation_, Scalar_, GridView_ > |
A local to global mapper implementation for boundary spaces, with functions defined on the domain boundary (or parts of it) only. More... | |
struct | Kaskade::BoundaryMapper< Implementation_, Scalar_, GridView_ >::Element< m > |
Type of the FunctionSpaceElement, associated to the FEFunctionSpace. More... | |
class | Kaskade::CellLocator< GridView, dimw > |
class | Kaskade::DiagonalCombiner< Scalar > |
A base class for implementing combiners with diagonal structure. More... | |
struct | Kaskade::SpaceType< Spaces, Idx > |
Extracts the type of FE space with given index from set of spaces. More... | |
class | Kaskade::FunctionSpaceElement< FunctionSpace, m > |
A class for representing finite element functions. More... | |
class | Kaskade::FEFunctionSpace< LocalToGlobalMapper > |
A representation of a finite element function space defined over a domain covered by a grid. More... | |
class | Kaskade::HierarchicSimplexShapeFunction< ctype, dimension, Scalar > |
Scalar hierarchic shape functions on the unit simplex. More... | |
class | Kaskade::DiscontinuousHierarchicMapper< Scalar, GV > |
A degrees of freedom manager for FEFunctionSpace s of piecewise polynomials of order Order. More... | |
class | Kaskade::ContinuousHierarchicMapper< Scalar, GV > |
A degrees of freedom manager for globally continuous FEFunctionSpace s of piecewise polynomials. More... | |
class | Kaskade::DiscontinuousHierarchicExtensionMapper< Scalar, GV > |
A degrees of freedom manager for FEFunctionSpace s of piecewise polynomials of given order. More... | |
class | Kaskade::ContinuousHierarchicExtensionMapper< Scalar, GV > |
A degrees of freedom manager for globally continuous FEFunctionSpace s of piecewise polynomials. More... | |
class | Kaskade::EquidistantLagrange< dim, Real > |
Lagrange shape functions on an equidistant grid over the simplex. These are polynomial basis functions \( p_i \) over the unit simplex, together with interpolation nodes \( \xi_j \), such that \( p_i(\xi_j) = \delta_{ij} \) holds. More... | |
class | Kaskade::DiscontinuousLagrangeMapper< ScalarType, GV > |
A degrees of freedom manager for FEFunctionSpace s of piecewise polynomials of order Order. More... | |
class | Kaskade::ContinuousLagrangeMapper< ScalarType, GV > |
A degrees of freedom manager for globally continuous FEFunctionSpace s of piecewise polynomials. More... | |
class | Kaskade::ContinuousLagrangeMapperSubdomain< GV, SupportOracle, ScalarType > |
A local-to-global mapper for continuous finite elements on a subdomain. More... | |
class | Kaskade::MorleyMapper< ScalarType, GV > |
Degrees of freedom manager for Morley nonconforming elements. More... | |
class | Kaskade::VectorialConverterBase< GridView > |
A class mapping local vectorial shape function values and gradients to global shape function values and gradients. More... | |
class | Kaskade::HcurlConverter< GridView > |
A class mapping local vectorial shape function values and gradients to global shape function values and gradients in \( H(\text{curl}) \) conforming FE spaces. More... | |
class | Kaskade::HdivConverter< GridView > |
A class mapping local vectorial shape function values and gradients to global shape function values and gradients in \( H(\text{div}) \) conforming FE spaces. More... | |
class | Kaskade::HdivMapper< ScalarType, GV > |
Degrees of freedom manager for H(div) conforming elements. More... | |
class | Kaskade::NedelecMapper< ScalarType, GV > |
Degrees of freedom manager for Nedelec edge elements. More... | |
class | Kaskade::UniformPartitionedMapper< Implementation, Tagger, SFData > |
Base class for piecewise continuous mappers. More... | |
class | Kaskade::PiecewiseContinuousLagrangeMapper< Tagger, GV, ScalarType > |
A degrees of freedom manager for piecewise continuous FEFunctionSpace s of local polynomials. More... | |
class | Kaskade::ShapeFunction< ctype, dim, T, comp > |
Shape functions. More... | |
class | Kaskade::ShapeFunctionSet< ctype, dimension, T, comp > |
A set of shape functions. More... | |
struct | Kaskade::QuadratureTraits< QuadRule > |
A cache that stores suitable quadrature rules for quick retrieval. More... | |
struct | Kaskade::QuadratureTraits< Dune::QuadratureRule< ctype, dim > > |
A cache that stores suitable quadrature rules for quick retrieval. More... | |
class | Kaskade::UniformScalarMapper< Implementation, SFData > |
Base class for uniform scalar local to global mappers. More... | |
struct | Kaskade::VariationalArg< Scalar, dim, components > |
A class that stores values, gradients, and Hessians of evaluated shape functions. More... | |
class | Kaskade::ShapeFunctionCache< G, T, ComponentsEnd > |
This class caches values and derivatives of shape functions at quadrature points. More... | |
Typedefs | |
template<class Spaces > | |
using | Kaskade::Evaluators = typename boost::fusion::result_of::as_vector< typename boost::fusion::result_of::transform< Spaces, GetEvaluatorTypes >::type >::type |
the heterogeneous sequence type of Evaluators for the given spaces More... | |
template<int m> | |
using | Kaskade::DiscontinuousLagrangeMapperSubdomain< ScalarType, GV >::Element_t = FunctionSpaceElement< FEFunctionSpace< Self >, m > |
Type of the FunctionSpaceElement, associated to the FEFunctionSpace. More... | |
template<int m> | |
using | Kaskade::ContinuousLagrangeMapperSubdomain< GV, SupportOracle, ScalarType >::Element_t = FunctionSpaceElement< FEFunctionSpace< ContinuousLagrangeMapperSubdomain >, m > |
Type of the FunctionSpaceElement, associated to the FEFunctionSpace. More... | |
template<class Grid , class Scalar = double> | |
using | Kaskade::ConstantSpace = FEFunctionSpace< ConstantMapper< Scalar, typename Grid::LeafGridView > > |
A finite element space of globally constant functions on the leaf grid view. More... | |
template<class Grid , class Scalar = double> | |
using | Kaskade::H1Space = FEFunctionSpace< ContinuousLagrangeMapper< Scalar, typename Grid::LeafGridView > > |
An \( H^1 \) conforming finite element space on the leaf grid view with Lagrange basis. More... | |
template<class Grid , class Scalar = double> | |
using | Kaskade::H1HierarchicSpace = FEFunctionSpace< ContinuousHierarchicMapper< Scalar, typename Grid::LeafGridView > > |
An \( H^1 \) conforming finite element space on the leaf grid view with hierarchic basis. More... | |
template<class Grid , class Scalar = double> | |
using | Kaskade::H1HierarchicExtensionSpace = FEFunctionSpace< ContinuousHierarchicExtensionMapper< Scalar, typename Grid::LeafGridView > > |
An \( H^1 \) conforming finite element extension space on the leaf grid view with hierarchic basis. More... | |
template<class Grid , class Scalar = double> | |
using | Kaskade::L2Space = FEFunctionSpace< DiscontinuousLagrangeMapper< Scalar, typename Grid::LeafGridView > > |
An \( L^2 \) conforming finite element space on the leaf grid view with Lagrange basis. More... | |
template<class Grid , class Scalar = double> | |
using | Kaskade::L2HierarchicSpace = FEFunctionSpace< DiscontinuousHierarchicMapper< Scalar, typename Grid::LeafGridView > > |
An \( L^2 \) conforming finite element space on the leaf grid view with hierarchic basis. More... | |
template<class Grid , class Scalar = double> | |
using | Kaskade::L2HierarchicExtensionSpace = FEFunctionSpace< DiscontinuousHierarchicExtensionMapper< Scalar, typename Grid::LeafGridView > > |
An \( L^2 \) conforming finite element extension space on the leaf grid view with hierarchic basis. More... | |
Functions | |
template<size_t dim> | |
std::array< int, dim+1 > | Kaskade::barycentric (std::array< int, dim > const &x, int bsum) |
Converts integer coordinate index to barycentric indices. More... | |
template<class Scalar , int dim> | |
Vector< Scalar, dim > | Kaskade::Simplex::cellInterpolation (VertexPositions< Scalar, dim > const &vertices, EdgeDirections< Scalar, dim > const &boundaryEdges, FaceFlags< dim > const &isBoundaryFace, Barycentric< Scalar, dim > b) |
Computes a boundary displacement based on blended cubic rational Bezier curves. More... | |
Kaskade::BoundaryMapper< Implementation_, Scalar_, GridView_ >::BoundaryMapper (GridView const &gridView, int order) | |
BoundaryMapper: Create a boundary mapper over the entire boundary. More... | |
Kaskade::BoundaryMapper< Implementation_, Scalar_, GridView_ >::BoundaryMapper (GridView const &gridView, int order, std::vector< int > const &boundaryIds_, std::vector< int > const &usedIds_) | |
BoundaryMapper: Create a boundary mapper over parts of the boundary. More... | |
int | Kaskade::BoundaryMapper< Implementation_, Scalar_, GridView_ >::numFaces () |
numFaces: Counts number of faces on which this boundaryMapper is defined. This costs some time so use with care. More... | |
template<class Evaluator , class QuadratureRule > | |
void | Kaskade::useQuadratureRuleInEvaluators (Evaluator &evals, QuadratureRule const &qr, int subId) |
Tells all evaluators to use the given quadrature rule on the given subentity. More... | |
template<class Evaluators , class CoordType , int dim, int subDim> | |
void | Kaskade::moveEvaluatorsToIntegrationPoint (Evaluators &evals, Dune::FieldVector< CoordType, dim > const &x, Dune::QuadratureRule< CoordType, subDim > const &qr, int ip, int subId) |
Moves all provided evaluators to the given integration point, evaluating the shape functions there. More... | |
template<class Evaluators , class CoordType , int dim> | |
void | Kaskade::moveEvaluatorsToIntegrationPoint (Evaluators &evals, Dune::FieldVector< CoordType, dim > const &x) |
Moves all provided evaluators to the given integration point, evaluating the shape functions there. More... | |
template<class Evaluators > | |
int | Kaskade::maxOrder (Evaluators const &evals) |
Computes the maximum ansatz order used in a collection of evaluators. More... | |
template<class Map , class IndexSet , class LocalDofs > | |
size_t | Kaskade::computeDofStartIndices (int dim, Map &map, IndexSet const &is, LocalDofs &localDof) |
For each codimension (i.e., type of subentity) computes the number of global ansatz functions as well as an accumulated index into an array of all ansatz functions. More... | |
template<class FEFunction > | |
FEFunction::ValueType | Kaskade::integrateOverBoundary (FEFunction const &function) |
integrateOverBoundary computes the integral of an FE function over the whole boundary of the underlying grid. More... | |
template<template< class, class > class DomainMapper, typename Scalar , typename GridView , int m> | |
FunctionSpaceElement< FEFunctionSpace< BoundaryMapper< DomainMapper, Scalar, GridView > >, m >::ValueType | Kaskade::integrateOverBoundary (FunctionSpaceElement< FEFunctionSpace< BoundaryMapper< DomainMapper, Scalar, GridView > >, m > const &function) |
integrateOverBoundary computes the integral of an FE function which is restricted to (parts of) the boundary over the boundary of the underlying grid. More... | |
template<class Space , class Range > | |
std::vector< size_t > | Kaskade::algebraicPatch (Space const &space, Range const &cells) |
Computes the global indices of all ansatz functions whose support is a subset of the union of given cells. More... | |
template<class Space , class Range > | |
std::vector< size_t > | Kaskade::algebraicPatch (Space const &space, Range const &cells, size_t vertexIdx, std::vector< std::vector< size_t > > const &indicesOnFacets) |
Computes the global indices of all ansatz functions whose support is a subset of the union of given cells. More... | |
template<class GridView > | |
auto | Kaskade::computePatches (GridView const &gridview) |
Computes the patches around vertices of the grid view. More... | |
template<class Space , int m, class OutIter > | |
OutIter | vectorToSequence (FunctionSpaceElement< Space, m > const &v, OutIter iter) |
Writes the coefficients into a flat sequence. <Space,m> More... | |
template<class Spaces , class ShapeFunctionCache > | |
auto | getEvaluators (Spaces const &spaces, ShapeFunctionCache *cache, std::map< void const *, int > const &deriv=std::map< void const *, int >()) |
returns a heterogeneous sequence of evaluators for the given spaces More... | |
Classes to define finite element spaces and functions on them.
In this group the problem is addressed of creating and maintaining finite element spaces and functions that live there. The following tasks are fulfilled:
Objects of these classes are connected by a boost::signal in case that the grid is modified (adapted)
Standard FEM spaces are implemented:
Access to shape functions
using Kaskade::ConstantSpace = typedef FEFunctionSpace<ConstantMapper<Scalar,typename Grid::LeafGridView> > |
using Kaskade::DiscontinuousLagrangeMapperSubdomain< ScalarType, GV >::Element_t = FunctionSpaceElement<FEFunctionSpace<Self>,m> |
Type of the FunctionSpaceElement, associated to the FEFunctionSpace.
m | number of components |
Definition at line 364 of file lagrangespace.hh.
using Kaskade::ContinuousLagrangeMapperSubdomain< GV, SupportOracle, ScalarType >::Element_t = FunctionSpaceElement<FEFunctionSpace<ContinuousLagrangeMapperSubdomain>,m> |
Type of the FunctionSpaceElement, associated to the FEFunctionSpace.
m | number of components |
Definition at line 701 of file lagrangespace.hh.
using Kaskade::Evaluators = typedef typename boost::fusion::result_of::as_vector< typename boost::fusion::result_of::transform<Spaces,GetEvaluatorTypes>::type >::type |
the heterogeneous sequence type of Evaluators for the given spaces
Spaces | a heterogeneous sequence type of FEFunctionSpace s |
Definition at line 1856 of file functionspace.hh.
using Kaskade::H1HierarchicExtensionSpace = typedef FEFunctionSpace<ContinuousHierarchicExtensionMapper<Scalar,typename Grid::LeafGridView> > |
using Kaskade::H1HierarchicSpace = typedef FEFunctionSpace<ContinuousHierarchicMapper<Scalar,typename Grid::LeafGridView> > |
using Kaskade::H1Space = typedef FEFunctionSpace<ContinuousLagrangeMapper<Scalar,typename Grid::LeafGridView> > |
using Kaskade::L2HierarchicExtensionSpace = typedef FEFunctionSpace<DiscontinuousHierarchicExtensionMapper<Scalar,typename Grid::LeafGridView> > |
using Kaskade::L2HierarchicSpace = typedef FEFunctionSpace<DiscontinuousHierarchicMapper<Scalar,typename Grid::LeafGridView> > |
using Kaskade::L2Space = typedef FEFunctionSpace<DiscontinuousLagrangeMapper<Scalar,typename Grid::LeafGridView> > |
std::vector< size_t > Kaskade::algebraicPatch | ( | Space const & | space, |
Range const & | cells | ||
) |
Computes the global indices of all ansatz functions whose support is a subset of the union of given cells.
These ansatz functions span the space of FE functions on the patch with zero trace.
Space | a FEFunctionSpace instatiation |
Range | an STL range type with grid cells as value type |
space | the FE space considered |
cells | the cells forming the patch |
Definition at line 42 of file supports.hh.
Referenced by Kaskade::TaylorHoodPreconditioner< X >::TaylorHoodPreconditioner().
std::vector< size_t > Kaskade::algebraicPatch | ( | Space const & | space, |
Range const & | cells, | ||
size_t | vertexIdx, | ||
std::vector< std::vector< size_t > > const & | indicesOnFacets | ||
) |
Computes the global indices of all ansatz functions whose support is a subset of the union of given cells.
These ansatz functions span the space of FE functions on the patch with zero trace. This version also treats patches that touch the domain boundary correctly (the above version takes more indices on these patches).
Space | a FEFunctionSpace instatiation |
Range | an STL range type with grid cells as value type |
vertexIdx | the global index of patch middle vertex |
indicesOnFacets | contains for each facet (by global facet index) the indices of all ansatz function that do not vanish on that facet |
Definition at line 105 of file supports.hh.
std::array< int, dim+1 > Kaskade::barycentric | ( | std::array< int, dim > const & | x, |
int | bsum | ||
) |
Converts integer coordinate index to barycentric indices.
dim | the spatial dimension |
x | the coordinate index |
bsum | the barycentric sum |
The entries of x need to be between 0 and bsum, and their sum shall not exceed bsum. Then the returned index contains in the first dim entries just x, and in the last one bsum-sum.
Definition at line 148 of file barycentric.hh.
|
inline |
BoundaryMapper: Create a boundary mapper over the entire boundary.
gridView | The corresponding grid view. |
order | Order of the underlying FE mapper. |
Definition at line 193 of file boundaryspace.hh.
|
inline |
BoundaryMapper: Create a boundary mapper over parts of the boundary.
gridView | The corresponding grid view. |
order | Order of the underlying FE mapper. |
boundaryIds_ | Container which maps each boundary segment index, which occurs in the grid (consecutive from zero to #(boundary segments in level 0 grid view)-1) to an int property value. |
usedIds_ | States by int property value which boundary faces belong to the desired part of the boundary. |
Definition at line 207 of file boundaryspace.hh.
Vector< Scalar, dim > Kaskade::Simplex::cellInterpolation | ( | VertexPositions< Scalar, dim > const & | vertices, |
EdgeDirections< Scalar, dim > const & | boundaryEdges, | ||
FaceFlags< dim > const & | isBoundaryFace, | ||
Barycentric< Scalar, dim > | b | ||
) |
Computes a boundary displacement based on blended cubic rational Bezier curves.
This allows to define a G1 continuous boundary surface. On each vertex \( x_i \) of a single simplicial surface patch, i.e. a triangle in 3D and a line in 2D, a unit normal vector \( n_i \) is given.
vertices | the boundary simplex's vertices, as columns of this matrix |
normals | the unit normal vectors at the vertices, as columns of this matrix |
b | the barycentric coordinates of the point to interpolate, i.e. the weight of the corresponding vertices/normals. The entries in b shall sum up to 1. |
Referenced by Kaskade::BoundaryInterpolationDetail::getChildVertexDisplacements(), and Kaskade::BoundaryInterpolationDisplacement< GridView >::value().
size_t Kaskade::computeDofStartIndices | ( | int | dim, |
Map & | map, | ||
IndexSet const & | is, | ||
LocalDofs & | localDof | ||
) |
For each codimension (i.e., type of subentity) computes the number of global ansatz functions as well as an accumulated index into an array of all ansatz functions.
dim | the dimension of the grid |
map | an associative container with key type Dune::GeometryType and int values. |
is | the index set defining the entities |
localDof | a callable object returning the number of local dofs associated to an entity of given geometry type |
Definition at line 41 of file gridcombinatorics.hh.
auto Kaskade::computePatches | ( | GridView const & | gridview | ) |
Computes the patches around vertices of the grid view.
Definition at line 165 of file supports.hh.
Referenced by Kaskade::PatchDomainDecompositionPreconditioner< Space, m, StorageTag, SparseMatrixIndex >::PatchDomainDecompositionPreconditioner(), and Kaskade::TaylorHoodPreconditioner< X >::TaylorHoodPreconditioner().
|
related |
returns a heterogeneous sequence of evaluators for the given spaces
Spaces | a heterogeneous sequence of FEFunctionSpace s |
ShapeFunctionCache | a type for caching shape function values |
cache | pointer to a shape function cache - if null, no cache is used |
deriv | a map from space addresses to highest derivatives requested |
Definition at line 1838 of file functionspace.hh.
FEFunction::ValueType Kaskade::integrateOverBoundary | ( | FEFunction const & | function | ) |
integrateOverBoundary computes the integral of an FE function over the whole boundary of the underlying grid.
function | is the FE function to be integrated. |
FEFunction | is the type of the integrand. |
Definition at line 211 of file integration.hh.
Referenced by Kaskade::boundaryL2Norm().
FunctionSpaceElement< FEFunctionSpace< BoundaryMapper< DomainMapper, Scalar, GridView > >, m >::ValueType Kaskade::integrateOverBoundary | ( | FunctionSpaceElement< FEFunctionSpace< BoundaryMapper< DomainMapper, Scalar, GridView > >, m > const & | function | ) |
integrateOverBoundary computes the integral of an FE function which is restricted to (parts of) the boundary over the boundary of the underlying grid.
This specialized overload does integration only on those faces where the boundary FE function is nonzero.
function | is the FE function to be integrated. |
Definition at line 236 of file integration.hh.
int Kaskade::maxOrder | ( | Evaluators const & | evals | ) |
Computes the maximum ansatz order used in a collection of evaluators.
Definition at line 1954 of file functionspace.hh.
Referenced by Kaskade::edgeAveragingFull(), Kaskade::GridIterateDetail::gridIterateRange(), Kaskade::HierarchicExtensionShapeFunctionSetContainer< ctype, dimension, Scalar, restricted >::HierarchicExtensionShapeFunctionSetContainer(), and Kaskade::HierarchicShapeFunctionSetContainer< ctype, dimension, Scalar, restricted >::HierarchicShapeFunctionSetContainer().
void Kaskade::moveEvaluatorsToIntegrationPoint | ( | Evaluators & | evals, |
Dune::FieldVector< CoordType, dim > const & | x | ||
) |
Moves all provided evaluators to the given integration point, evaluating the shape functions there.
evals | the evaluators to be moved |
x | the evaluation point (local coordinates) |
Definition at line 1942 of file functionspace.hh.
void Kaskade::moveEvaluatorsToIntegrationPoint | ( | Evaluators & | evals, |
Dune::FieldVector< CoordType, dim > const & | x, | ||
Dune::QuadratureRule< CoordType, subDim > const & | qr, | ||
int | ip, | ||
int | subId | ||
) |
Moves all provided evaluators to the given integration point, evaluating the shape functions there.
evals | the evaluators to be moved |
x | the evaluation point (local coordinates) |
qr | the quadrature rule defining the integration point x |
ip | the number of the integration point in the quadrature rule qr |
subId | the sub-id of the entity on which the quadrature rule lives, w.r.t. the current cell |
Definition at line 1929 of file functionspace.hh.
Referenced by PointwiseCorrection< PolyEquation >::correct(), corrTangent(), Kaskade::edgeAveragingFull(), Kaskade::ErrorDistribution< Functional, ExtendedAnsatzVars >::BoundaryCache::evaluateAt(), Kaskade::ErrorDistribution< Functional, ExtendedAnsatzVars >::DomainCache::evaluateAt(), and Kaskade::GridIterateDetail::gridIterateRange().
|
inline |
numFaces: Counts number of faces on which this boundaryMapper is defined. This costs some time so use with care.
Definition at line 230 of file boundaryspace.hh.
void Kaskade::useQuadratureRuleInEvaluators | ( | Evaluator & | evals, |
QuadratureRule const & | qr, | ||
int | subId | ||
) |
Tells all evaluators to use the given quadrature rule on the given subentity.
Definition at line 1912 of file functionspace.hh.
Referenced by Kaskade::ErrorDistribution< Functional, ExtendedAnsatzVars >::BoundaryCache::evaluateAt(), Kaskade::ErrorDistribution< Functional, ExtendedAnsatzVars >::DomainCache::evaluateAt(), and Kaskade::GridIterateDetail::gridIterateRange().
|
related |
Writes the coefficients into a flat sequence. <Space,m>
Definition at line 1039 of file functionspace.hh.