dune-fem 2.8.0
Loading...
Searching...
No Matches
Public Types | Public Member Functions | Static Public Member Functions | Static Public Attributes | Protected Types | Protected Member Functions | Protected Attributes | List of all members
Dune::Fem::LumpingQuadrature< FieldImp, geometryId > Class Template Reference

#include <dune/fem/quadrature/lumpingquadrature.hh>

Inheritance diagram for Dune::Fem::LumpingQuadrature< FieldImp, geometryId >:
Inheritance graph

Public Types

typedef FieldImp FieldType
 
typedef BaseType::CoordinateType CoordinateType
 
enum  { codimension = 0 }
 to be revised, look at caching quad More...
 

Public Member Functions

 LumpingQuadrature (const GeometryType &gt, int order, int id)
 constructor filling the list of points and weights.
 
virtual GeometryType geometryType () const
 
virtual int order () const
 obtain order of the integration point list
 
const FieldTypeweight (size_t i) const
 obtain weight of i-th integration point
 
const CoordinateTypepoint (size_t i) const
 obtain coordinates of i-th integration point
 
size_t nop () const
 obtain the number of integration points
 
size_t id () const
 obtain the identifier of the integration point list
 
virtual std::vector< ElementCoordinateTypeinterpolationPoints (const int reqDim) const
 returns list of element interpolation points for a given face quadrature
 
virtual bool isFaceInterpolationQuadrature (const size_t numShapeFunctions) const
 return true if quadrature is also a set of interpolation points for a given number of shape functions
 

Static Public Member Functions

static std::size_t maxOrder ()
 maximal order of available quadratures
 

Static Public Attributes

static constexpr auto dimension = Dune::GeometryType(geometryId).dim()
 

Protected Types

typedef FieldVector< FieldType, dim+1 > ElementCoordinateType
 

Protected Member Functions

void addIntegrationPoint (const CoordinateType &point)
 Adds an integration point to the list.
 
void addQuadraturePoint (const CoordinateType &point, const FieldType weight)
 Adds a point-weight pair to the quadrature.
 
void setIntegrationPoints (std::vector< CoordinateType > &&points)
 Overwrites integration point list

 

Protected Attributes

std::vector< FieldTypeweights_
 
std::vector< CoordinateTypepoints_
 
const size_t id_
 

Detailed Description

template<class FieldImp, Dune::GeometryType::Id geometryId>
class Dune::Fem::LumpingQuadrature< FieldImp, geometryId >

Define a lumping quadrature for all geometries. Note, however, that this may not make sense for anything else than simplices or maybe hexagonal grids. For simplicial meshes the quadrature formula is exact on linear polynomials and hence the quadrature error is quadratic in the mesh-size. A mass-matrix assembled with the caching quadrature will be diagonal in the context of Lagrange spaces. Generally, it is a bad idea to use mass-lumping for anything else than linear (or maybe bilinear) finite elements.

There are probably much more efficient ways to perform mass-lumping. This "quadrature" approach is convenient, however, as it can simply be plugged into existing code by replacing the quadrature rules.

Member Typedef Documentation

◆ CoordinateType

template<class FieldImp , Dune::GeometryType::Id geometryId>
typedef BaseType::CoordinateType Dune::Fem::LumpingQuadrature< FieldImp, geometryId >::CoordinateType

◆ ElementCoordinateType

template<typename FieldImp , int dim>
typedef FieldVector< FieldType, dim+1 > Dune::Fem::IntegrationPointListImp< FieldImp, dim >::ElementCoordinateType
protectedinherited

◆ FieldType

template<class FieldImp , Dune::GeometryType::Id geometryId>
typedef FieldImp Dune::Fem::LumpingQuadrature< FieldImp, geometryId >::FieldType

Member Enumeration Documentation

◆ anonymous enum

template<typename FieldImp , int dim>
anonymous enum
inherited

to be revised, look at caching quad

Enumerator
codimension 

Constructor & Destructor Documentation

◆ LumpingQuadrature()

template<class FieldImp , Dune::GeometryType::Id geometryId>
Dune::Fem::LumpingQuadrature< FieldImp, geometryId >::LumpingQuadrature ( const GeometryType &  gt,
int  order,
int  id 
)
inline

constructor filling the list of points and weights.

Parameters
[in]gtgeometry type for which a quadrature is desired
[in]orderorder, ignored
[in]idunique identifier, ignored

Member Function Documentation

◆ addIntegrationPoint()

template<typename FieldImp , int dim>
void Dune::Fem::IntegrationPointListImp< FieldImp, dim >::addIntegrationPoint ( const CoordinateType point)
inlineprotectedinherited

Adds an integration point to the list.

This method allows derived classes to add integration points to the list. This mehtod should only be used within the constructor of the derived class.

◆ addQuadraturePoint()

void Dune::Fem::QuadratureImp< FieldImp, dim >::addQuadraturePoint ( const CoordinateType point,
const FieldType  weight 
)
inlineprotectedinherited

Adds a point-weight pair to the quadrature.

This method allows derived classes to add quadrature points (and their respective weights) to the list. This mehtod should only be used within the constructor of the derived class.

◆ geometryType()

template<class FieldImp , Dune::GeometryType::Id geometryId>
virtual GeometryType Dune::Fem::LumpingQuadrature< FieldImp, geometryId >::geometryType ( ) const
inlinevirtual

◆ id()

template<typename FieldImp , int dim>
size_t Dune::Fem::IntegrationPointListImp< FieldImp, dim >::id ( ) const
inlineinherited

obtain the identifier of the integration point list

The identifier of an integration point list must be globally unique. Even integration point lists for different dimensions must have different identifiers.

Note
Quadratures are considered distinct if they differ in one of the following points: geometry type, order, dimension or implementation.
Returns
globally unique identifier of the integration point list

◆ interpolationPoints()

template<typename FieldImp , int dim>
virtual std::vector< ElementCoordinateType > Dune::Fem::IntegrationPointListImp< FieldImp, dim >::interpolationPoints ( const int  reqDim) const
inlinevirtualinherited

returns list of element interpolation points for a given face quadrature

◆ isFaceInterpolationQuadrature()

template<typename FieldImp , int dim>
virtual bool Dune::Fem::IntegrationPointListImp< FieldImp, dim >::isFaceInterpolationQuadrature ( const size_t  numShapeFunctions) const
inlinevirtualinherited

return true if quadrature is also a set of interpolation points for a given number of shape functions

◆ maxOrder()

template<class FieldImp , Dune::GeometryType::Id geometryId>
static std::size_t Dune::Fem::LumpingQuadrature< FieldImp, geometryId >::maxOrder ( )
inlinestatic

maximal order of available quadratures

◆ nop()

template<typename FieldImp , int dim>
size_t Dune::Fem::IntegrationPointListImp< FieldImp, dim >::nop ( ) const
inlineinherited

obtain the number of integration points

Returns
number of integration points within this list

◆ order()

template<class FieldImp , Dune::GeometryType::Id geometryId>
virtual int Dune::Fem::LumpingQuadrature< FieldImp, geometryId >::order ( ) const
inlinevirtual

obtain order of the integration point list

The order of a quadrature is the maximal polynomial degree that is guaranteed to be integrated exactly by the quadrature.

In case of an integration point list, the definition of this value is left to the implementor.

Returns
the order of the integration point list

Implements Dune::Fem::IntegrationPointListImp< FieldImp, dim >.

◆ point()

template<typename FieldImp , int dim>
const CoordinateType & Dune::Fem::IntegrationPointListImp< FieldImp, dim >::point ( size_t  i) const
inlineinherited

obtain coordinates of i-th integration point

This method returns a reference to the coordinates of the i-th integration point for 0 <= i < nop(). The integration point is given in local coordinates, i.e., coordinates with respect to the reference element.

Parameters
[in]inumber of the integration point, 0 <= i < nop()
Returns
reference to i-th integration point

◆ setIntegrationPoints()

template<typename FieldImp , int dim>
void Dune::Fem::IntegrationPointListImp< FieldImp, dim >::setIntegrationPoints ( std::vector< CoordinateType > &&  points)
inlineprotectedinherited

Overwrites integration point list

◆ weight()

const FieldType & Dune::Fem::QuadratureImp< FieldImp, dim >::weight ( size_t  i) const
inlineinherited

obtain weight of i-th integration point

This method returns the weight of the i-th integration point for 0 <= i < nop() within the quadrature.

Note
The integration point can be obtained via the point() method.
The quadrature weights sum up to the volume of the reference element.
Parameters
[in]inumber of the integration point, 0 <= i < nop()
Returns
weight of the i-th integration point

Member Data Documentation

◆ dimension

template<class FieldImp , Dune::GeometryType::Id geometryId>
constexpr auto Dune::Fem::LumpingQuadrature< FieldImp, geometryId >::dimension = Dune::GeometryType(geometryId).dim()
staticconstexpr

◆ id_

template<typename FieldImp , int dim>
const size_t Dune::Fem::IntegrationPointListImp< FieldImp, dim >::id_
protectedinherited

◆ points_

template<typename FieldImp , int dim>
std::vector< CoordinateType > Dune::Fem::IntegrationPointListImp< FieldImp, dim >::points_
mutableprotectedinherited

◆ weights_

std::vector< FieldType > Dune::Fem::QuadratureImp< FieldImp, dim >::weights_
mutableprotectedinherited

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