dune-pdelab  2.4-dev
Public Types | Public Member Functions | Protected Attributes | List of all members
Dune::PDELab::Laplace Class Reference

#include <dune/pdelab/localoperator/laplace.hh>

Inheritance diagram for Dune::PDELab::Laplace:
Inheritance graph

Public Types

enum  { doPatternVolume = true }
 
enum  { doAlphaVolume = true }
 
Flags for the sparsity pattern
enum  { doPatternVolume }
 Whether to assemble the pattern on the elements, i.e. whether or not pattern_volume() should be called. More...
 
enum  { doPatternVolumePostSkeleton }
 Whether to assemble the pattern on the elements after the skeleton has been handled, i.e. whether or not pattern_volume_post_skeleton() should be called. More...
 
enum  { doPatternSkeleton }
 Whether to assemble the pattern on the interior intersections, i.e. whether or not pattern_skeleton() should be called. More...
 
enum  { doPatternBoundary }
 Whether to assemble the pattern on the boundary intersections, i.e. whether or not pattern_boundary() should be called. More...
 
Flags for the non-constant part of the residual and the jacobian
enum  { doAlphaVolume }
 Whether to call the local operator's alpha_volume(), jacobian_apply_volume() and jacobian_volume(). More...
 
enum  { doAlphaVolumePostSkeleton }
 Whether to call the local operator's alpha_volume_post_skeleton(), jacobian_apply_volume_post_skeleton() and jacobian_volume_post_skeleton(). More...
 
enum  { doAlphaSkeleton }
 Whether to call the local operator's alpha_skeleton(), jacobian_apply_skeleton() and jacobian_skeleton(). More...
 
enum  { doAlphaBoundary }
 Whether to call the local operator's alpha_boundary(), jacobian_apply_boundary() and jacobian_boundary(). More...
 
Flags for the constant part of the residual
enum  { doLambdaVolume }
 Whether to call the local operator's lambda_volume(). More...
 
enum  { doLambdaVolumePostSkeleton }
 Whether to call the local operator's lambda_volume_post_skeleton(). More...
 
enum  { doLambdaSkeleton }
 Whether to call the local operator's lambda_skeleton(). More...
 
enum  { doLambdaBoundary }
 Whether to call the local operator's lambda_boundary(). More...
 
Special flags
enum  { doSkeletonTwoSided }
 Whether to visit the skeleton methods from both sides. More...
 

Public Member Functions

 Laplace (unsigned int quadOrder)
 Constructor. More...
 
template<typename EG , typename LFSU , typename X , typename LFSV , typename R >
void alpha_volume (const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, R &r) const
 Compute Laplace matrix times a given vector for one element. More...
 
template<typename EG , typename LFSU , typename X , typename LFSV , typename M >
void jacobian_volume (const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, M &matrix) const
 Compute the Laplace stiffness matrix for the element given in 'eg'. More...
 
template<typename LFSU , typename LFSV , typename LocalPattern >
void pattern_volume (const LFSU &lfsu, const LFSV &lfsv, LocalPattern &pattern) const
 

Protected Attributes

unsigned int quadOrder_
 

Detailed Description

a local operator for solving the Laplace equation

\begin{align*} - \Delta u &=& 0 \mbox{ in } \Omega, \\ -\nabla u \cdot \nu &=& 0 \mbox{ on } \partial\Omega_N \\ \end{align*}

with conforming finite elements on all types of grids in any dimension.

In other words, it only assembles the Laplace matrix.

Constructor & Destructor Documentation

Dune::PDELab::Laplace::Laplace ( unsigned int  quadOrder)
inline

Constructor.

Parameters
quadOrderOrder of the quadrature rule used for integrating over the element

Member Function Documentation

template<typename EG , typename LFSU , typename X , typename LFSV , typename R >
void Dune::PDELab::Laplace::alpha_volume ( const EG &  eg,
const LFSU &  lfsu,
const X &  x,
const LFSV &  lfsv,
R &  r 
) const
inline

Compute Laplace matrix times a given vector for one element.

This is used for matrix-free algorithms for the Laplace equation

Parameters
[in]egThe grid element we are assembling on
[in]lfsuLocal ansatz function space basis
[in]lfsvLocal test function space basis
[in]xInput vector
[out]rThe product of the Laplace matrix times x

References quadOrder_.

template<typename EG , typename LFSU , typename X , typename LFSV , typename M >
void Dune::PDELab::Laplace::jacobian_volume ( const EG &  eg,
const LFSU &  lfsu,
const X &  x,
const LFSV &  lfsv,
M &  matrix 
) const
inline

Compute the Laplace stiffness matrix for the element given in 'eg'.

Template Parameters
MType of the element stiffness matrix
Parameters
[in]egThe grid element we are assembling on
[in]lfsuLocal ansatz function space basis
[in]lfsvLocal test function space basis
[in]xCurrent configuration; gets ignored for linear problems like this one
[out]matrixElement stiffness matrix

References dim, and quadOrder_.

template<typename LFSU , typename LFSV , typename LocalPattern >
void Dune::PDELab::FullVolumePattern::pattern_volume ( const LFSU &  lfsu,
const LFSV &  lfsv,
LocalPattern &  pattern 
) const
inlineinherited

Member Data Documentation

unsigned int Dune::PDELab::Laplace::quadOrder_
protected

Referenced by alpha_volume(), and jacobian_volume().


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