1#ifndef DUNE_FEM_LOCALMASSMATRIX_HH
2#define DUNE_FEM_LOCALMASSMATRIX_HH
5#include <dune/common/dynmatrix.hh>
6#include <dune/common/fmatrix.hh>
9#include <dune/geometry/typeindex.hh>
32 template<
class DiscreteFunctionSpace,
class VolumeQuadrature >
39 typedef typename DiscreteFunctionSpaceType :: RangeFieldType
ctype;
41 typedef typename DiscreteFunctionSpaceType :: RangeType
RangeType;
43 enum {
dimRange = DiscreteFunctionSpaceType :: dimRange };
47 typedef Dune::FieldMatrix< ctype, dgNumDofs, dgNumDofs >
DGMatrixType;
50 typedef typename DiscreteFunctionSpaceType :: GridPartType
GridPartType;
52 typedef typename DiscreteFunctionSpaceType :: IndexSetType
IndexSetType;
53 typedef typename IndexSetType :: IndexType
IndexType;
57 typedef typename GridPartType :: GridType
GridType;
58 typedef typename DiscreteFunctionSpaceType :: EntityType
EntityType;
59 typedef typename EntityType :: Geometry
Geometry;
69 typedef typename GeometryInformationType :: DomainType
DomainType;
72 typedef Dune::DynamicMatrix< RangeFieldType >
MatrixType;
73 typedef Dune::DynamicVector< RangeFieldType >
VectorType;
76 std::shared_ptr< const DiscreteFunctionSpaceType >
spc_;
91 mutable std::vector< RangeType >
phi_;
94 typedef std::pair< std::unique_ptr< MatrixType >, std::unique_ptr< VectorType > >
MatrixPairType;
120 const int rows = matrix.rows();
121 for(
int r=0; r<rows; ++r )
123 for(
int c=0; c<r; ++c )
126 if(
std::abs(matrix[r][c]) > 1e-12 )
133 template<
class BasisFunctionSet >
138 const GeometryType geomType = geo.type();
139 typedef typename MassMatrixStorageType::iterator iterator;
142 auto it = massMap.find( numBasisFct );
143 if( it == massMap.end() )
145 std::pair< iterator, bool > insertPair = massMap.insert( std::make_pair( numBasisFct,
MatrixPairType(
nullptr,
nullptr) ) );
146 it = insertPair.first;
147 insertPair.first->second.first.reset(
new MatrixType( numBasisFct, numBasisFct, 0.0 ));
148 MatrixType& matrix = insertPair.first->second.first.operator *();
154 catch ( Dune::FMatrixError &e )
156 std::cerr <<
"Matrix is singular:" << std::endl << matrix << std::endl;
163 insertPair.first->second.second.reset(
new VectorType( matrix.rows() ) );
164 VectorType& diag = insertPair.first->second.second.operator *();
165 const int rows = matrix.rows();
166 for(
int row=0; row<rows; ++row )
168 diag[ row ] = matrix[ row ][ row ];
176 template<
class MassCaller,
class BasisFunctionSet >
180 const int numDofs = basisSet.
size();
185 if( numDofs !=
int(
matrix_.rows() ) )
186 matrix_.resize( numDofs, numDofs );
193 assert( numDofs ==
int(
matrix_.rows() ) );
264 template<
class BasisFunctionSet >
269 static const int basePointSetId = detail::SelectPointSetId< BasisFunctionSet >::value;
273 if constexpr ( quadPointSetId == basePointSetId )
275 const unsigned int numShapeFunctions = bfs.
size() /
dimRange;
277 .isInterpolationQuadrature(numShapeFunctions);
284 template<
class MassCaller,
class BasisFunctionSet,
class LocalFunction >
289 && !caller.hasMass() )
297 template<
class MassCaller,
class LocalFunction >
304 template<
class LocalFunction >
310 template<
class BasisFunctionSet,
class LocalFunction >
319 template<
class LocalFunction >
326 template<
class LocalMatrix >
329 const EntityType &entity = localMatrix.rangeEntity();
338 template<
class LocalMatrix >
341 const EntityType &entity = localMatrix.domainEntity();
359 template<
class MassCaller,
class BasisFunctionSet,
class LocalFunction >
367 const bool isAffine =
affine() || geo.affine();
369 assert(
affine() ? geo.affine() :
true );
372 if( isAffine && !caller.hasMass() )
378 lf[ l ] *= massVolInv;
396 template<
class LocalMatrix >
399 const EntityType &entity = localMatrix.rangeEntity();
401 assert(
dgNumDofs == localMatrix.columns() );
404 if(
affine() || geo.affine() )
414 const int rows = localMatrix.rows();
415 for(
int i = 0; i < rows; ++i )
418 dgRhs_[ j ] = localMatrix.get( i, j );
421 localMatrix.set( i, j,
dgX_[ j ] );
427 template<
class LocalMatrix >
430 const EntityType &entity = localMatrix.domainEntity();
432 assert(
dgNumDofs == localMatrix.columns() );
435 if(
affine() || geo.affine() )
445 const int rows = localMatrix.rows();
446 for(
int i = 0; i < rows; ++i )
449 dgRhs_[ j ] = localMatrix.get( i, j );
452 localMatrix.set( i, j,
dgX_[ j ] );
461 const int currentSequence =
space().sequence();
462 const unsigned int topologyId = entity.type().id();
487 template<
class MassCaller,
class BasisFunctionSet,
class LocalFunction >
496 const int numDofs = lf.
size();
497 rhs_.resize( numDofs );
498 for(
int l = 0; l < numDofs; ++l )
506 template<
class LocalMatrix >
513 const int cols = localMatrix.columns();
517 const int rows = localMatrix.rows();
518 for(
int i = 0; i < rows; ++i )
521 for(
int j = 0; j < cols; ++j )
522 rhs_[ j ] = localMatrix.get( i, j );
528 for(
int j = 0; j < cols; ++j )
529 localMatrix.set( i, j,
row_[ j ] );
534 template<
class LocalMatrix >
541 const int cols = localMatrix.columns();
545 const int rows = localMatrix.rows();
546 for(
int i = 0; i < rows; ++i )
549 for(
int j = 0; j < cols; ++j )
550 rhs_[ j ] = localMatrix.get( j, i );
556 for(
int j = 0; j < cols; ++j )
557 localMatrix.set( j, i,
row_[ j ] );
566 template<
class BasisFunctionSet,
class LocalFunction >
570 const int numDofs = lf.
size();
577 if( matrixPair.second )
579 const VectorType& diagonal = *matrixPair.second;
580 assert(
int(diagonal.size()) == numDofs );
584 const int nop = volQuad.nop();
586 for(
int l=0, qt = 0; qt < nop; ++qt )
588 const auto intel = geo.integrationElement( volQuad.point(qt) );
589 for (
int r = 0; r <
dimRange; ++r, ++l )
591 lf[ l ] *= diagonal[ l ] / intel;
600 rhs_.resize( numDofs );
601 for(
int l = 0; l < numDofs; ++l )
602 rhs_[ l ] = lf[ l ] * massVolInv;
604 const MatrixType& invMassMatrix = *matrixPair.first;
610 template <
class LocalMatrix>
614 LocalMatrix &localMatrix )
const
616 const int cols = localMatrix.columns();
618 assert(
int(refElemDiagonal.size()) == cols );
623 elementDiagonal.resize( cols );
625 const int nop = volQuad.nop();
627 for(
int l = 0, qt = 0; qt < nop; ++qt )
629 const auto intel = geo.integrationElement( volQuad.point(qt) );
630 for (
int r = 0; r <
dimRange; ++r,++l )
632 elementDiagonal[ l ] = refElemDiagonal[ l ] / intel;
635 return elementDiagonal;
638 template<
class LocalMatrix >
641 const int cols = localMatrix.columns();
647 if( matrixPair.second )
653 const int rows = localMatrix.rows();
654 for(
int i = 0; i < rows; ++i )
658 for(
int j = 0; j < cols; ++j )
659 row_[ j ] = elementDiagonal[ j ] * localMatrix.get( i, j );
662 for(
int j = 0; j < cols; ++j )
663 localMatrix.set( i, j,
row_[ j ] );
668 const MatrixType &invMassMatrix = *matrixPair.first;
675 const int rows = localMatrix.rows();
676 for(
int i = 0; i < rows; ++i )
680 for(
int j = 0; j < cols; ++j )
681 rhs_[ j ] = localMatrix.get( i, j ) * massVolInv;
687 for(
int j = 0; j < cols; ++j )
688 localMatrix.set( i, j,
row_[ j ] );
694 template<
class LocalMatrix >
697 const int cols = localMatrix.columns();
703 if( matrixPair.second )
709 const int rows = localMatrix.rows();
710 for(
int i = 0; i < rows; ++i )
714 for(
int j = 0; j < cols; ++j )
715 row_[ j ] = elementDiagonal[ j ] * localMatrix.get( j, i );
718 for(
int j = 0; j < cols; ++j )
719 localMatrix.set( j, i,
row_[ j ] );
724 const MatrixType &invMassMatrix = *matrixPair.first;
731 const int rows = localMatrix.rows();
732 for(
int i = 0; i < rows; ++i )
735 for(
int j = 0; j < cols; ++j )
736 rhs_[ j ] = localMatrix.get( j, i ) * massVolInv;
742 for(
int j = 0; j < cols; ++j )
743 localMatrix.set( j, i,
row_[ j ] );
760 if( (geomTypes.size() == 1) && geomTypes[0].isSimplex() )
770 template<
class MassCaller,
class Matrix >
773 std::size_t numDofs, Matrix &matrix )
const
775 assert( numDofs == set.size() );
783 if( caller.hasMass() )
790 template <
class Matrix>
798 const bool applyIntegrationElement =
true )
const
800 const int volNop = volQuad.nop();
801 for(
int qp=0; qp<volNop; ++qp)
804 const double intel = ( applyIntegrationElement ) ?
805 ( volQuad.weight(qp) * geo.integrationElement(volQuad.point(qp)) ) : volQuad.weight(qp) ;
808 set.evaluateAll(volQuad[qp],
phi_);
810 for(
int m=0; m<numDofs; ++m)
813 const ctype val = intel * (phi_m * phi_m);
816 for(
int k=m+1; k<numDofs; ++k)
818 const ctype val = intel * (phi_m *
phi_[k]);
827 template <
class MassCallerType,
class Matrix>
829 MassCallerType& caller,
835 Matrix& matrix)
const
837 typedef typename MassCallerType :: MassFactorType MassFactorType;
840 const int volNop = volQuad.nop();
841 for(
int qp=0; qp<volNop; ++qp)
844 const double intel = volQuad.weight(qp)
845 * geo.integrationElement(volQuad.point(qp));
848 set.evaluateAll( volQuad[qp],
phi_);
851 caller.mass( en, volQuad, qp, mass);
854 for(
int m=0; m<numDofs; ++m)
860 for(
int m=0; m<numDofs; ++m)
862 for(
int k=0; k<numDofs; ++k)
871 template <
class Matrix,
class Rhs,
class X>
873 const Matrix& matrix,
877 assert( (
int) matrix.rows() == size );
878 assert( (
int) matrix.cols() == size );
879 assert( (
int) rhs.size() == size );
881 for(
int row = 0; row < size; ++ row )
885 typedef typename Matrix :: const_row_reference MatRow;
886 MatRow matRow = matrix[ row ];
889 for(
int col = 0; col < size; ++ col )
891 sum += matRow[ col ] * rhs[ col ];
906 template<
class DiscreteFunctionSpace,
class VolumeQuadrature >
914 using BaseType :: LocalMassMatrixImplementation;
926 template<
class DiscreteFunctionSpace,
class VolumeQuadrature >
936 using BaseType :: LocalMassMatrixImplementation;
940 template <
class MassCallerType,
class BasisFunction,
class LocalFunctionType>
943 const BasisFunction &basisFunction,
944 LocalFunctionType& lf)
const
946 BaseType :: applyInverseDgOrthoNormalBasis( caller, en, basisFunction, lf );
948 template <
class MassCallerType,
class LocalFunctionType>
951 LocalFunctionType& lf)
const
953 BaseType :: applyInverseDgOrthoNormalBasis( caller, en, lf.basisFunctionSet(), lf );
957 template <
class LocalFunctionType>
959 LocalFunctionType& lf)
const
961 typename BaseType :: NoMassDummyCaller caller;
965 template <
class BasisFunction,
class LocalFunctionType>
967 const BasisFunction &basisFunction,
968 LocalFunctionType& lf)
const
970 typename BaseType :: NoMassDummyCaller caller;
975 template<
class LocalFunction >
982 template<
class LocalMatrix >
989 template<
class LocalMatrix >
Dune::Fem::Double abs(const Dune::Fem::Double &a)
Definition: double.hh:942
Definition: bindguard.hh:11
detail::SelectPointSetId< Quadrature, -Dune::QuadratureType::size > SelectQuadraturePointSetId
Definition: quadrature.hh:541
static std::shared_ptr< T > referenceToSharedPtr(T &t)
Definition: memory.hh:19
interface for local functions
Definition: localfunction.hh:77
const EntityType & entity() const
obtain the entity, this local function lives on
Definition: localfunction.hh:302
const BasisFunctionSetType & basisFunctionSet() const
obtain the basis function set for this local function
Definition: localfunction.hh:296
SizeType size() const
obtain the number of local DoFs
Definition: localfunction.hh:360
Helper class to check affinity of the grid's geometries.
Definition: checkgeomaffinity.hh:22
Local Mass Matrix inversion implementation, select the correct method in your implementation.
Definition: localmassmatrix.hh:34
DiscreteFunctionSpaceType::EntityType EntityType
Definition: localmassmatrix.hh:58
void applyInverse(const EntityType &entity, const BasisFunctionSet &basisFunctionSet, LocalFunction &lf) const
Definition: localmassmatrix.hh:311
void applyInverse(MassCaller &caller, const EntityType &entity, LocalFunction &lf) const
Definition: localmassmatrix.hh:298
std::vector< RangeType > phi_
Definition: localmassmatrix.hh:91
MatrixType & getLocalInverseMassMatrixDefault(MassCaller &caller, const EntityType &entity, const Geometry &geo, const BasisFunctionSet &basisSet) const
Definition: localmassmatrix.hh:177
AllGeomTypes< typename GridPartType ::IndexSetType, GridType > GeometryInformationType
Definition: localmassmatrix.hh:68
DiscreteFunctionSpaceType::IndexSetType IndexSetType
Definition: localmassmatrix.hh:52
bool checkInterpolationBFS(const BasisFunctionSet &bfs) const
Definition: localmassmatrix.hh:265
int volumeQuadratureOrder(const EntityType &entity) const
return appropriate quadrature order, default is 2 * order(entity)
Definition: localmassmatrix.hh:211
void applyInverseDefault(MassCaller &caller, const EntityType &entity, const Geometry &geo, const BasisFunctionSet &basisFunctionSet, LocalFunction &lf) const
Definition: localmassmatrix.hh:488
void applyInverseLocally(const EntityType &entity, const Geometry &geo, const BasisFunctionSet &basisFunctionSet, LocalFunction &lf) const
apply local mass matrix to local function lf
Definition: localmassmatrix.hh:567
const bool affine_
Definition: localmassmatrix.hh:81
const DiscreteFunctionSpaceType & space() const
Definition: localmassmatrix.hh:349
bool setup() const
setup and return affinity
Definition: localmassmatrix.hh:750
DiscreteFunctionSpaceType::RangeFieldType ctype
Definition: localmassmatrix.hh:39
GridPartType::GridType GridType
Definition: localmassmatrix.hh:57
Dune::FieldVector< ctype, dgNumDofs > DGVectorType
Definition: localmassmatrix.hh:48
const IndexSetType & indexSet_
Definition: localmassmatrix.hh:77
void applyInverse(const EntityType &entity, LocalFunction &lf) const
apply local dg mass matrix to local function lf without mass factor
Definition: localmassmatrix.hh:305
@ dimRange
Definition: localmassmatrix.hh:43
DGMatrixType dgMatrix_
Definition: localmassmatrix.hh:83
LocalMassMatrixImplementation(const DiscreteFunctionSpaceType &spc, int volQuadOrd)
constructor taking space and volume quadrature order
Definition: localmassmatrix.hh:217
Fem::GeometryAffinityCheck< VolumeQuadratureType > GeometryAffinityCheckType
Definition: localmassmatrix.hh:63
double getAffineMassFactor(const Geometry &geo) const
return mass factor for diagonal mass matrix
Definition: localmassmatrix.hh:259
@ StructuredGrid
Definition: localmassmatrix.hh:66
int maxNumDofs() const
Definition: localmassmatrix.hh:198
IndexSetType::IndexType IndexType
Definition: localmassmatrix.hh:53
VolumeQuadrature VolumeQuadratureType
Definition: localmassmatrix.hh:61
std::vector< MassMatrixStorageType > LocalInverseMassMatrixStorageType
Definition: localmassmatrix.hh:96
std::vector< RangeType > phiMass_
Definition: localmassmatrix.hh:92
LocalMassMatrixImplementation(const ThisType &other)
copy constructor
Definition: localmassmatrix.hh:239
const VectorType & setupInverseDiagonal(const EntityType &entity, const Geometry &geo, const VectorType &refElemDiagonal, LocalMatrix &localMatrix) const
Definition: localmassmatrix.hh:612
void applyInverse(LocalFunction &lf) const
apply local dg mass matrix to local function lf without mass factor
Definition: localmassmatrix.hh:320
Dune::DynamicVector< RangeFieldType > VectorType
Definition: localmassmatrix.hh:73
bool checkDiagonalMatrix(const MatrixType &matrix) const
Definition: localmassmatrix.hh:118
DGVectorType dgRhs_
Definition: localmassmatrix.hh:84
void applyInverse(MassCaller &caller, const EntityType &entity, const BasisFunctionSet &basisFunctionSet, LocalFunction &lf) const
Definition: localmassmatrix.hh:285
void leftMultiplyInverseDefault(const EntityType &entity, const Geometry &geo, LocalMatrix &localMatrix) const
compute M^-1 * localMatrix
Definition: localmassmatrix.hh:535
@ localBlockSize
Definition: localmassmatrix.hh:44
std::pair< std::unique_ptr< MatrixType >, std::unique_ptr< VectorType > > MatrixPairType
Definition: localmassmatrix.hh:94
MatrixPairType & getLocalInverseMassMatrix(const EntityType &entity, const Geometry &geo, const BasisFunctionSet &basisSet, int numBasisFct) const
Definition: localmassmatrix.hh:135
void multiply(const int size, const Matrix &matrix, const Rhs &rhs, X &x) const
Definition: localmassmatrix.hh:872
bool entityHasChanged(const EntityType &entity) const
returns true if the entity has been changed
Definition: localmassmatrix.hh:458
DiscreteFunctionSpaceType::RangeType RangeType
Definition: localmassmatrix.hh:41
DGVectorType dgX_
Definition: localmassmatrix.hh:84
LocalInverseMassMatrixStorageType localInverseMassMatrix_
Definition: localmassmatrix.hh:98
void buildMatrixWithMassFactor(MassCallerType &caller, const EntityType &en, const Geometry &geo, const BasisFunctionSetType &set, const VolumeQuadratureType &volQuad, const int numDofs, Matrix &matrix) const
build local mass matrix with mass factor
Definition: localmassmatrix.hh:828
DiscreteFunctionSpaceType::BasisFunctionSetType BasisFunctionSetType
Definition: localmassmatrix.hh:55
DiscreteFunctionSpace DiscreteFunctionSpaceType
Definition: localmassmatrix.hh:38
DiscreteFunctionSpaceType::RangeFieldType RangeFieldType
Definition: localmassmatrix.hh:40
MatrixType matrix_
Definition: localmassmatrix.hh:89
void rightMultiplyInverseDefault(const EntityType &entity, const Geometry &geo, LocalMatrix &localMatrix) const
compute localMatrix * M^-1
Definition: localmassmatrix.hh:507
void leftMultiplyInverseDgOrthoNormalBasis(LocalMatrix &localMatrix) const
compute M^-1 * localMatrix
Definition: localmassmatrix.hh:428
bool affine() const
returns true if geometry mapping is affine
Definition: localmassmatrix.hh:256
VectorType rhs_
Definition: localmassmatrix.hh:88
IndexType lastEntityIndex_
Definition: localmassmatrix.hh:101
@ dgNumDofs
Definition: localmassmatrix.hh:45
Dune::FieldMatrix< ctype, dgNumDofs, dgNumDofs > DGMatrixType
Definition: localmassmatrix.hh:47
void rightMultiplyInverseLocally(const EntityType &entity, const Geometry &geo, LocalMatrix &localMatrix) const
Definition: localmassmatrix.hh:639
GeometryInformationType::DomainType DomainType
Definition: localmassmatrix.hh:69
void applyInverseDgOrthoNormalBasis(MassCaller &caller, const EntityType &entity, const BasisFunctionSet &basisFunctionSet, LocalFunction &lf) const
Definition: localmassmatrix.hh:360
LocalMassMatrixImplementation(const DiscreteFunctionSpaceType &spc, std::function< int(const int)> volQuadOrderFct=[](const int order) { return Capabilities::DefaultQuadrature< DiscreteFunctionSpaceType >::volumeOrder(order);})
constructor taking space and volume quadrature order
Definition: localmassmatrix.hh:222
void rightMultiplyInverse(LocalMatrix &localMatrix) const
compute localMatrix * M^-1
Definition: localmassmatrix.hh:327
const std::function< int(const int)> volumeQuadratureOrder_
Definition: localmassmatrix.hh:80
unsigned int lastTopologyId_
Definition: localmassmatrix.hh:102
void leftMultiplyInverse(LocalMatrix &localMatrix) const
compute M^-1 * localMatrix
Definition: localmassmatrix.hh:339
void buildMatrix(MassCaller &caller, const EntityType &entity, const Geometry &geo, const BasisFunctionSetType &set, std::size_t numDofs, Matrix &matrix) const
build local mass matrix
Definition: localmassmatrix.hh:771
GeometryInformationType geoInfo_
Definition: localmassmatrix.hh:79
VectorType row_
Definition: localmassmatrix.hh:88
int maxVolumeQuadratureOrder() const
return appropriate quadrature order, default is 2 * order()
Definition: localmassmatrix.hh:204
void buildMatrixNoMassFactor(const EntityType &en, const Geometry &geo, const BasisFunctionSetType &set, const VolumeQuadratureType &volQuad, const int numDofs, Matrix &matrix, const bool applyIntegrationElement=true) const
build local mass matrix with mass factor
Definition: localmassmatrix.hh:791
std::shared_ptr< const DiscreteFunctionSpaceType > spc_
Definition: localmassmatrix.hh:76
EntityType::Geometry Geometry
Definition: localmassmatrix.hh:59
void leftMultiplyInverseLocally(const EntityType &entity, const Geometry &geo, LocalMatrix &localMatrix) const
compute M^-1 * localMatrix
Definition: localmassmatrix.hh:695
void rightMultiplyInverseDgOrthoNormalBasis(LocalMatrix &localMatrix) const
compute localMatrix * M^-1
Definition: localmassmatrix.hh:397
Dune::DynamicMatrix< RangeFieldType > MatrixType
Definition: localmassmatrix.hh:72
std::map< const int, MatrixPairType > MassMatrixStorageType
Definition: localmassmatrix.hh:95
int sequence_
Definition: localmassmatrix.hh:104
DiscreteFunctionSpaceType::GridPartType GridPartType
Definition: localmassmatrix.hh:50
Definition: localmassmatrix.hh:107
bool hasMass() const
Definition: localmassmatrix.hh:111
Dune::FieldMatrix< ctype, dimRange, dimRange > MassFactorType
Definition: localmassmatrix.hh:108
void mass(const EntityType &, const VolumeQuadratureType &, int, MassFactorType &) const
Definition: localmassmatrix.hh:113
Local Mass Matrix for arbitrary spaces.
Definition: localmassmatrix.hh:909
DG Local Mass Matrix for arbitrary spaces.
Definition: localmassmatrix.hh:929
BaseType::EntityType EntityType
Definition: localmassmatrix.hh:933
void applyInverse(LocalFunction &lf) const
apply local dg mass matrix to local function lf without mass factor
Definition: localmassmatrix.hh:976
void applyInverse(const EntityType &en, LocalFunctionType &lf) const
apply local dg mass matrix to local function lf without mass factor
Definition: localmassmatrix.hh:958
void applyInverse(MassCallerType &caller, const EntityType &en, LocalFunctionType &lf) const
Definition: localmassmatrix.hh:949
void applyInverse(MassCallerType &caller, const EntityType &en, const BasisFunction &basisFunction, LocalFunctionType &lf) const
Definition: localmassmatrix.hh:941
void applyInverse(const EntityType &en, const BasisFunction &basisFunction, LocalFunctionType &lf) const
Definition: localmassmatrix.hh:966
void rightMultiplyInverse(LocalMatrix &localMatrix) const
compute localMatrix * M^-1
Definition: localmassmatrix.hh:983
void leftMultiplyInverse(LocalMatrix &localMatrix) const
compute M^-1 * localMatrix
Definition: localmassmatrix.hh:990
Interface class for basis function sets.
Definition: basisfunctionset/basisfunctionset.hh:31
const EntityType & entity() const
return entity
std::size_t size() const
return size of basis function set
ctype referenceVolume(const GeometryType &type) const
return volume of reference element for geometry of type type
Definition: allgeomtypes.hh:72
default implementation uses method geomTypes of given index set. Used in DiscreteFunctionSpaces.
Definition: allgeomtypes.hh:99
const std ::vector< GeometryType > & geomTypes(unsigned int codim) const
returns vector with geometry tpyes this index set has indices for
Definition: allgeomtypes.hh:171
static int volumeOrder(const int k)
return quadrature order for volume quadratures for given polynomial order k
Definition: space/common/capabilities.hh:138