dune-fem
2.4.1-rc
|
#include </local/tomalk/somewhere/tmp/dune-fem/dune/fem/quadrature/lumpingquadrature.hh>
Public Types | |
enum | { dimension = TopologyType::dimension } |
typedef FieldImp | FieldType |
typedef Topology | TopologyType |
typedef BaseType::CoordinateType | CoordinateType |
enum | { codimension = 0 } |
to be revised, look at caching quad More... | |
Public Member Functions | |
LumpingQuadrature (const GeometryType &geometry, int ignored, int id) | |
constructor filling the list of points and weights. More... | |
virtual GeometryType | geometryType () const |
virtual int | order () const |
obtain order of the integration point list More... | |
const FieldType & | weight (size_t i) const |
obtain weight of i-th integration point More... | |
const CoordinateType & | point (size_t i) const |
obtain coordinates of i-th integration point More... | |
size_t | nop () const |
obtain the number of integration points More... | |
size_t | id () const |
obtain the identifier of the integration point list More... | |
Static Public Member Functions | |
static size_t | maxOrder () |
maximal order of available quadratures More... | |
Protected Types | |
typedef Dune::GenericGeometry::ReferenceDomain< TopologyType > | ReferenceDomain |
Protected Member Functions | |
void | addQuadraturePoint (const CoordinateType &point, const FieldType weight) |
Adds a point-weight pair to the quadrature. More... | |
Static Protected Attributes | |
static const unsigned int | topologyId = TopologyType::id |
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.
typedef BaseType::CoordinateType Dune::Fem::LumpingQuadrature< FieldImp, Topology >::CoordinateType |
typedef FieldImp Dune::Fem::LumpingQuadrature< FieldImp, Topology >::FieldType |
|
protected |
typedef Topology Dune::Fem::LumpingQuadrature< FieldImp, Topology >::TopologyType |
|
inherited |
|
inline |
constructor filling the list of points and weights.
[in] | gemoetry | geometry type for which a quadrature is desired |
[in] | ignored | order is ignored |
[in] | id | unique identifier, ignored |
References Dune::Fem::QuadratureImp< FieldImp, Topology::dimension >::addQuadraturePoint().
|
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.
Referenced by Dune::Fem::LumpingQuadrature< FieldImp, Topology >::LumpingQuadrature().
|
inlinevirtual |
|
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.
Referenced by Dune::Fem::CacheProvider< GridPart, 1 >::getMapper(), Dune::Fem::PointProvider< ct, dim, 1 >::getMappers(), Dune::Fem::PointProvider< ct, dim, 1 >::getPoints(), Dune::Fem::TwistProvider< ct, dim >::getTwistStorage(), and Dune::Fem::PointProvider< ct, dim, 0 >::registerQuadrature().
|
inlinestatic |
maximal order of available quadratures
|
inlineinherited |
obtain the number of integration points
Referenced by Dune::Fem::LocalDGPass< DiscreteModelImp, PreviousPassImp >::applyLocalNeighbor(), Dune::Fem::TwistMapperCreator< ct, dim >::createStorage(), Dune::Fem::PointProvider< ct, dim, 1 >::getMappers(), Dune::Fem::IntegrationPointListImp< ct, Topology::dimension >::point(), and Dune::Fem::PointProvider< ct, dim, 0 >::registerQuadrature().
|
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.
Implements Dune::Fem::IntegrationPointListImp< FieldImp, dim >.
|
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.
[in] | i | number of the integration point, 0 <= i < nop() |
Referenced by Dune::Fem::TwistMapperCreator< ct, dim >::createStorage(), Dune::Fem::PointProvider< ct, dim, 1 >::getMappers(), Dune::Fem::PointProvider< ct, dim, 0 >::registerQuadrature(), and Dune::Fem::QuadratureImp< FieldImp, 1 >::weight().
|
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.
[in] | i | number of the integration point, 0 <= i < nop() |
|
staticprotected |