dune-pdelab  2.4.1
Public Types | Public Member Functions | Protected Member Functions | Protected Attributes | Static Protected Attributes | List of all members
Dune::PDELab::LocalAssemblerBase< B, CU, CV > Class Template Reference

Base class for local assembler. More...

#include <dune/pdelab/gridoperator/common/assemblerutilities.hh>

Inheritance diagram for Dune::PDELab::LocalAssemblerBase< B, CU, CV >:
Inheritance graph

Public Types

typedef B::size_type SizeType
 

Public Member Functions

 LocalAssemblerBase ()
 construct AssemblerSpace More...
 
 LocalAssemblerBase (const CU &cu, const CV &cv)
 construct AssemblerSpace, with constraints More...
 
const CU & trialConstraints () const
 get the constraints on the trial grid function space More...
 
const CV & testConstraints () const
 get the constraints on the test grid function space More...
 
template<typename X >
enable_if< AlwaysTrue< X >::value &&!is_same< CV, EmptyTransformation >::value >::type forwardtransform (X &x, const bool postrestrict=false) const
 Transforms a vector $ \boldsymbol{x} $ from $ V$ to $ V'$. If postrestrict == true then $\boldsymbol{R}^T_{\boldsymbol{\tilde U}', \boldsymbol{U}'} \boldsymbol{S}_{\boldsymbol{\tilde V}}$ is applied instead of the full transformation. More...
 
template<typename X >
enable_if< AlwaysTrue< X >::value &&is_same< CV, EmptyTransformation >::value >::type forwardtransform (X &x, const bool postrestrict=false) const
 
template<typename X >
enable_if< AlwaysTrue< X >::value &&!is_same< CV, EmptyTransformation >::value >::type backtransform (X &x, const bool prerestrict=false) const
 Transforms a vector $ \boldsymbol{x} $ from $ V'$ to $ V$. If prerestrict == true then $\boldsymbol{S}^T_{\boldsymbol{\tilde U}}$ is applied instead of the full transformation. More...
 
template<typename X >
enable_if< AlwaysTrue< X >::value &&is_same< CV, EmptyTransformation >::value >::type backtransform (X &x, const bool prerestrict=false) const
 

Protected Member Functions

template<typename GCView , typename T >
void eread (const GCView &globalcontainer_view, LocalMatrix< T > &localcontainer) const
 read local stiffness matrix for entity More...
 
template<typename T , typename GCView >
void ewrite (const LocalMatrix< T > &localcontainer, GCView &globalcontainer_view) const
 write local stiffness matrix for entity More...
 
template<typename T , typename GCView >
void eadd (const LocalMatrix< T > &localcontainer, GCView &globalcontainer_view) const
 write local stiffness matrix for entity More...
 
template<typename M , typename GCView >
enable_if< AlwaysTrue< M >::value &&!is_same< CV, EmptyTransformation >::value >::type scatter_jacobian (M &local_container, GCView &global_container_view, bool symmetric_mode) const
 Scatter local jacobian to global container. More...
 
template<typename M , typename GCView >
enable_if< AlwaysTrue< M >::value &&is_same< CV, EmptyTransformation >::value >::type scatter_jacobian (M &local_container, GCView &global_container_view, bool symmetric_mode) const
 
template<typename M , typename GCView >
void etadd_symmetric (M &localcontainer, GCView &globalcontainer_view) const
 Add local matrix to global matrix, and apply Dirichlet constraints in a symmetric fashion. Apart from that, identical to etadd(). More...
 
template<typename M , typename GCView >
void etadd (const M &localcontainer, GCView &globalcontainer_view) const
 
template<typename Pattern , typename RI , typename CI >
enable_if< is_same< RI, CI >::value >::type add_diagonal_entry (Pattern &pattern, const RI &ri, const CI &ci) const
 
template<typename Pattern , typename RI , typename CI >
enable_if< !is_same< RI, CI >::value >::type add_diagonal_entry (Pattern &pattern, const RI &ri, const CI &ci) const
 
template<typename P , typename LFSVIndices , typename LFSUIndices , typename Index >
void add_entry (P &globalpattern, const LFSVIndices &lfsv_indices, Index i, const LFSUIndices &lfsu_indices, Index j) const
 Adding matrix entry to pattern with respect to the constraints contributions. This assembles the entries addressed by etadd(..). See the documentation there for more information about the matrix pattern. More...
 
template<typename GFSV , typename GC , typename C >
void set_trivial_rows (const GFSV &gfsv, GC &globalcontainer, const C &c) const
 insert dirichlet constraints for row and assemble T^T_U in constrained rows More...
 
template<typename GFSV , typename GC >
void set_trivial_rows (const GFSV &gfsv, GC &globalcontainer, const EmptyTransformation &c) const
 
template<typename GFSV , typename GC >
void handle_dirichlet_constraints (const GFSV &gfsv, GC &globalcontainer) const
 

Protected Attributes

const CU * pconstraintsu
 
const CV * pconstraintsv
 

Static Protected Attributes

static CU emptyconstraintsu
 
static CV emptyconstraintsv
 

Detailed Description

template<typename B, typename CU = EmptyTransformation, typename CV = EmptyTransformation>
class Dune::PDELab::LocalAssemblerBase< B, CU, CV >

Base class for local assembler.

This class provides some generic behavior required for local assembler implementations. This includes the access of the global vectors and matrices via local indices and local function spaces with regard to the constraint mappings.

Template Parameters
BThe matrix backend
CUConstraints maps for the individual dofs (trial space)
CVConstraints maps for the individual dofs (test space)

Member Typedef Documentation

template<typename B, typename CU = EmptyTransformation, typename CV = EmptyTransformation>
typedef B::size_type Dune::PDELab::LocalAssemblerBase< B, CU, CV >::SizeType

Constructor & Destructor Documentation

template<typename B, typename CU = EmptyTransformation, typename CV = EmptyTransformation>
Dune::PDELab::LocalAssemblerBase< B, CU, CV >::LocalAssemblerBase ( )
inline

construct AssemblerSpace

template<typename B, typename CU = EmptyTransformation, typename CV = EmptyTransformation>
Dune::PDELab::LocalAssemblerBase< B, CU, CV >::LocalAssemblerBase ( const CU &  cu,
const CV &  cv 
)
inline

construct AssemblerSpace, with constraints

Member Function Documentation

template<typename B, typename CU = EmptyTransformation, typename CV = EmptyTransformation>
template<typename Pattern , typename RI , typename CI >
enable_if< is_same<RI,CI>::value >::type Dune::PDELab::LocalAssemblerBase< B, CU, CV >::add_diagonal_entry ( Pattern &  pattern,
const RI &  ri,
const CI &  ci 
) const
inlineprotected
template<typename B, typename CU = EmptyTransformation, typename CV = EmptyTransformation>
template<typename Pattern , typename RI , typename CI >
enable_if< !is_same<RI,CI>::value >::type Dune::PDELab::LocalAssemblerBase< B, CU, CV >::add_diagonal_entry ( Pattern &  pattern,
const RI &  ri,
const CI &  ci 
) const
inlineprotected
template<typename B, typename CU = EmptyTransformation, typename CV = EmptyTransformation>
template<typename P , typename LFSVIndices , typename LFSUIndices , typename Index >
void Dune::PDELab::LocalAssemblerBase< B, CU, CV >::add_entry ( P &  globalpattern,
const LFSVIndices &  lfsv_indices,
Index  i,
const LFSUIndices &  lfsu_indices,
Index  j 
) const
inlineprotected

Adding matrix entry to pattern with respect to the constraints contributions. This assembles the entries addressed by etadd(..). See the documentation there for more information about the matrix pattern.

template<typename B, typename CU = EmptyTransformation, typename CV = EmptyTransformation>
template<typename X >
enable_if< AlwaysTrue<X>::value && !is_same< CV, EmptyTransformation >::value >::type Dune::PDELab::LocalAssemblerBase< B, CU, CV >::backtransform ( X &  x,
const bool  prerestrict = false 
) const
inline

Transforms a vector $ \boldsymbol{x} $ from $ V'$ to $ V$. If prerestrict == true then $\boldsymbol{S}^T_{\boldsymbol{\tilde U}}$ is applied instead of the full transformation.

template<typename B, typename CU = EmptyTransformation, typename CV = EmptyTransformation>
template<typename X >
enable_if< AlwaysTrue<X>::value && is_same< CV, EmptyTransformation >::value >::type Dune::PDELab::LocalAssemblerBase< B, CU, CV >::backtransform ( X &  x,
const bool  prerestrict = false 
) const
inline
template<typename B, typename CU = EmptyTransformation, typename CV = EmptyTransformation>
template<typename T , typename GCView >
void Dune::PDELab::LocalAssemblerBase< B, CU, CV >::eadd ( const LocalMatrix< T > &  localcontainer,
GCView &  globalcontainer_view 
) const
inlineprotected

write local stiffness matrix for entity

template<typename B, typename CU = EmptyTransformation, typename CV = EmptyTransformation>
template<typename GCView , typename T >
void Dune::PDELab::LocalAssemblerBase< B, CU, CV >::eread ( const GCView &  globalcontainer_view,
LocalMatrix< T > &  localcontainer 
) const
inlineprotected

read local stiffness matrix for entity

template<typename B, typename CU = EmptyTransformation, typename CV = EmptyTransformation>
template<typename M , typename GCView >
void Dune::PDELab::LocalAssemblerBase< B, CU, CV >::etadd ( const M &  localcontainer,
GCView &  globalcontainer_view 
) const
inlineprotected
template<typename B, typename CU = EmptyTransformation, typename CV = EmptyTransformation>
template<typename M , typename GCView >
void Dune::PDELab::LocalAssemblerBase< B, CU, CV >::etadd_symmetric ( M &  localcontainer,
GCView &  globalcontainer_view 
) const
inlineprotected

Add local matrix to global matrix, and apply Dirichlet constraints in a symmetric fashion. Apart from that, identical to etadd().

template<typename B, typename CU = EmptyTransformation, typename CV = EmptyTransformation>
template<typename T , typename GCView >
void Dune::PDELab::LocalAssemblerBase< B, CU, CV >::ewrite ( const LocalMatrix< T > &  localcontainer,
GCView &  globalcontainer_view 
) const
inlineprotected

write local stiffness matrix for entity

template<typename B, typename CU = EmptyTransformation, typename CV = EmptyTransformation>
template<typename X >
enable_if< AlwaysTrue<X>::value && !is_same< CV, EmptyTransformation >::value >::type Dune::PDELab::LocalAssemblerBase< B, CU, CV >::forwardtransform ( X &  x,
const bool  postrestrict = false 
) const
inline

Transforms a vector $ \boldsymbol{x} $ from $ V$ to $ V'$. If postrestrict == true then $\boldsymbol{R}^T_{\boldsymbol{\tilde U}', \boldsymbol{U}'} \boldsymbol{S}_{\boldsymbol{\tilde V}}$ is applied instead of the full transformation.

template<typename B, typename CU = EmptyTransformation, typename CV = EmptyTransformation>
template<typename X >
enable_if< AlwaysTrue<X>::value && is_same< CV, EmptyTransformation >::value >::type Dune::PDELab::LocalAssemblerBase< B, CU, CV >::forwardtransform ( X &  x,
const bool  postrestrict = false 
) const
inline
template<typename B, typename CU = EmptyTransformation, typename CV = EmptyTransformation>
template<typename GFSV , typename GC >
void Dune::PDELab::LocalAssemblerBase< B, CU, CV >::handle_dirichlet_constraints ( const GFSV &  gfsv,
GC &  globalcontainer 
) const
inlineprotected
template<typename B, typename CU = EmptyTransformation, typename CV = EmptyTransformation>
template<typename M , typename GCView >
enable_if< AlwaysTrue<M>::value && !is_same< CV, EmptyTransformation >::value >::type Dune::PDELab::LocalAssemblerBase< B, CU, CV >::scatter_jacobian ( M &  local_container,
GCView &  global_container_view,
bool  symmetric_mode 
) const
inlineprotected

Scatter local jacobian to global container.

template<typename B, typename CU = EmptyTransformation, typename CV = EmptyTransformation>
template<typename M , typename GCView >
enable_if< AlwaysTrue<M>::value && is_same< CV, EmptyTransformation >::value >::type Dune::PDELab::LocalAssemblerBase< B, CU, CV >::scatter_jacobian ( M &  local_container,
GCView &  global_container_view,
bool  symmetric_mode 
) const
inlineprotected
template<typename B, typename CU = EmptyTransformation, typename CV = EmptyTransformation>
template<typename GFSV , typename GC , typename C >
void Dune::PDELab::LocalAssemblerBase< B, CU, CV >::set_trivial_rows ( const GFSV &  gfsv,
GC &  globalcontainer,
const C &  c 
) const
inlineprotected

insert dirichlet constraints for row and assemble T^T_U in constrained rows

template<typename B, typename CU = EmptyTransformation, typename CV = EmptyTransformation>
template<typename GFSV , typename GC >
void Dune::PDELab::LocalAssemblerBase< B, CU, CV >::set_trivial_rows ( const GFSV &  gfsv,
GC &  globalcontainer,
const EmptyTransformation c 
) const
inlineprotected
template<typename B, typename CU = EmptyTransformation, typename CV = EmptyTransformation>
const CV& Dune::PDELab::LocalAssemblerBase< B, CU, CV >::testConstraints ( ) const
inline

get the constraints on the test grid function space

template<typename B, typename CU = EmptyTransformation, typename CV = EmptyTransformation>
const CU& Dune::PDELab::LocalAssemblerBase< B, CU, CV >::trialConstraints ( ) const
inline

get the constraints on the trial grid function space

Member Data Documentation

template<typename B, typename CU = EmptyTransformation, typename CV = EmptyTransformation>
CU Dune::PDELab::LocalAssemblerBase< B, CU, CV >::emptyconstraintsu
staticprotected
template<typename B, typename CU = EmptyTransformation, typename CV = EmptyTransformation>
CV Dune::PDELab::LocalAssemblerBase< B, CU, CV >::emptyconstraintsv
staticprotected
template<typename B, typename CU = EmptyTransformation, typename CV = EmptyTransformation>
const CU* Dune::PDELab::LocalAssemblerBase< B, CU, CV >::pconstraintsu
protected
template<typename B, typename CU = EmptyTransformation, typename CV = EmptyTransformation>
const CV* Dune::PDELab::LocalAssemblerBase< B, CU, CV >::pconstraintsv
protected

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