1#ifndef DUNE_FEM_LOCALMASSMATRIX_HH
2#define DUNE_FEM_LOCALMASSMATRIX_HH
12#include <dune/fem/common/memory.hh>
13#include <dune/fem/common/utility.hh>
14#include <dune/fem/misc/checkgeomaffinity.hh>
15#include <dune/fem/quadrature/cachingquadrature.hh>
16#include <dune/fem/space/common/allgeomtypes.hh>
17#include <dune/fem/space/common/capabilities.hh>
18#include <dune/fem/version.hh>
32 template<
class DiscreteFunctionSpace,
class VolumeQuadrature,
bool refElemScaling = true >
39 typedef typename DiscreteFunctionSpaceType :: RangeFieldType ctype;
40 typedef typename DiscreteFunctionSpaceType :: RangeFieldType RangeFieldType;
41 typedef typename DiscreteFunctionSpaceType :: RangeType RangeType;
43 enum { dimRange = DiscreteFunctionSpaceType :: dimRange };
44 enum { localBlockSize = DiscreteFunctionSpaceType :: localBlockSize };
45 enum { dgNumDofs = localBlockSize };
50 typedef typename DiscreteFunctionSpaceType :: GridPartType GridPartType;
52 typedef typename DiscreteFunctionSpaceType :: IndexSetType IndexSetType;
55 typedef typename DiscreteFunctionSpaceType :: BasisFunctionSetType BasisFunctionSetType;
57 typedef typename GridPartType :: GridType GridType;
58 typedef typename DiscreteFunctionSpaceType :: EntityType EntityType;
59 typedef typename EntityType :: Geometry Geometry;
61 typedef VolumeQuadrature VolumeQuadratureType;
76 std::shared_ptr< const DiscreteFunctionSpaceType > spc_;
77 const IndexSetType& indexSet_;
80 const std::function<int(
const int)> volumeQuadratureOrder_;
91 mutable std::vector< RangeType > phi_;
92 mutable std::vector< RangeType > phiMass_;
94 typedef std::pair< std::unique_ptr< MatrixType >, std::unique_ptr< VectorType > > MatrixPairType;
95 typedef std::map< const int, MatrixPairType > MassMatrixStorageType;
96 typedef std::vector< MassMatrixStorageType > LocalInverseMassMatrixStorageType;
98 mutable LocalInverseMassMatrixStorageType localInverseMassMatrix_;
101 mutable IndexType lastEntityIndex_;
102 mutable unsigned int lastTopologyId_ ;
104 mutable int sequence_;
106 struct NoMassDummyCaller
111 bool hasMass ()
const {
return false; }
113 void mass (
const EntityType &,
const VolumeQuadratureType &,
int, MassFactorType & )
const
118 bool checkDiagonalMatrix(
const MatrixType& matrix )
const
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 >
135 getLocalInverseMassMatrix (
const EntityType &entity,
const Geometry &geo,
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 *();
156 std::cerr <<
"Matrix is singular:" << std::endl << matrix << std::endl;
159 const bool diagonal = checkDiagonalMatrix( matrix );
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 >
177 MatrixType &getLocalInverseMassMatrixDefault ( MassCaller &caller,
const EntityType &entity,
180 const int numDofs = basisSet.
size();
185 if( numDofs !=
int( matrix_.
rows() ) )
186 matrix_.
resize( numDofs, numDofs );
188 buildMatrix( caller, entity, geo, basisSet, numDofs, matrix_ );
193 assert( numDofs ==
int( matrix_.
rows() ) );
198 int maxNumDofs ()
const
200 return space().blockMapper().maxNumDofs() * localBlockSize;
206 return volumeQuadratureOrder_( space().order() );
213 return volumeQuadratureOrder_( space().order( entity ) );
224 : spc_( referenceToSharedPtr( spc ) )
225 , indexSet_( space().indexSet() )
226 , geoInfo_( indexSet_ )
227 , volumeQuadratureOrder_ ( volQuadOrderFct )
228 , affine_ (
setup() )
229 , rhs_(), row_(), matrix_()
230 , phi_( maxNumDofs() )
231 , phiMass_( maxNumDofs() )
233 , lastEntityIndex_( -1 )
234 , lastTopologyId_( ~0u )
241 indexSet_( space().indexSet() ),
242 geoInfo_( indexSet_ ),
243 volumeQuadratureOrder_( other.volumeQuadratureOrder_ ),
244 affine_( other.affine_ ),
245 rhs_( other.rhs_ ), row_( other.row_ ), matrix_( other.matrix_ ),
247 phiMass_( other.phiMass_ ),
249 lastEntityIndex_( other.lastEntityIndex_ ),
250 lastTopologyId_( other.lastTopologyId_ ),
251 sequence_( other.sequence_ )
263 if constexpr( refElemScaling )
270 return 1.0 / geo.volume();
274 template<
class BasisFunctionSet >
283 if constexpr ( quadPointSetId == basePointSetId )
285 const unsigned int numShapeFunctions = bfs.
size() / dimRange;
287 .isInterpolationQuadrature(numShapeFunctions);
294 template<
class MassCaller,
class BasisFunctionSet,
class LocalFunction >
297 Geometry geo = entity.geometry();
298 if( (
affine() || geo.affine() || checkInterpolationBFS(basisFunctionSet) )
299 && !caller.hasMass() )
307 template<
class MassCaller,
class LocalFunction >
314 template<
class LocalFunction >
317 NoMassDummyCaller caller;
320 template<
class BasisFunctionSet,
class LocalFunction >
323 NoMassDummyCaller caller;
329 template<
class LocalFunction >
336 template<
class LocalMatrix >
339 const EntityType &entity = localMatrix.rangeEntity();
340 Geometry geo = entity.geometry();
341 if( (
affine() || geo.affine() || checkInterpolationBFS(localMatrix.rangeBasisFunctionSet())) )
342 rightMultiplyInverseLocally( entity, geo, localMatrix );
348 template<
class LocalMatrix >
351 const EntityType &entity = localMatrix.domainEntity();
352 Geometry geo = entity.geometry();
353 if( (
affine() || geo.affine() || checkInterpolationBFS(localMatrix.rangeBasisFunctionSet())) )
359 const DiscreteFunctionSpaceType &space ()
const {
return *spc_; }
369 template<
class MassCaller,
class BasisFunctionSet,
class LocalFunction >
370 void applyInverseDgOrthoNormalBasis ( MassCaller &caller,
const EntityType &entity,
371 const BasisFunctionSet &basisFunctionSet, LocalFunction &lf )
const
374 assert( dgNumDofs == lf.size() );
382 if( isAffine && !caller.hasMass() )
387 for(
int l = 0; l < dgNumDofs; ++l )
388 lf[ l ] *= massVolInv;
393 for(
int l = 0; l < dgNumDofs; ++l )
394 dgRhs_[ l ] = lf[ l ];
396 buildMatrix( caller, entity, geo, basisFunctionSet, dgNumDofs, dgMatrix_ );
397 dgMatrix_.
solve( dgX_, dgRhs_ );
400 for(
int l = 0; l < dgNumDofs; ++l )
406 template<
class LocalMatrix >
409 const EntityType &entity = localMatrix.rangeEntity();
410 Geometry geo = entity.geometry();
411 assert( dgNumDofs == localMatrix.columns() );
414 if(
affine() || geo.affine() )
420 NoMassDummyCaller caller;
421 buildMatrix( caller, entity, geo, localMatrix.rangeBasisFunctionSet(), dgNumDofs, dgMatrix_ );
424 const int rows = localMatrix.rows();
425 for(
int i = 0; i < rows; ++i )
427 for(
int j = 0; j < dgNumDofs; ++j )
428 dgRhs_[ j ] = localMatrix.get( i, j );
429 dgMatrix_.
mtv( dgRhs_, dgX_ );
430 for(
int j = 0; j < dgNumDofs; ++j )
431 localMatrix.set( i, j, dgX_[ j ] );
437 template<
class LocalMatrix >
440 const EntityType &entity = localMatrix.domainEntity();
441 Geometry geo = entity.geometry();
442 assert( dgNumDofs == localMatrix.columns() );
445 if(
affine() || geo.affine() )
451 NoMassDummyCaller caller;
452 buildMatrix( caller, entity, geo, localMatrix.domainBasisFunctionSet(), dgNumDofs, dgMatrix_ );
455 const int rows = localMatrix.rows();
456 for(
int i = 0; i < rows; ++i )
458 for(
int j = 0; j < dgNumDofs; ++j )
459 dgRhs_[ j ] = localMatrix.get( i, j );
460 dgMatrix_.
mv( dgRhs_, dgX_ );
461 for(
int j = 0; j < dgNumDofs; ++j )
462 localMatrix.set( i, j, dgX_[ j ] );
471 const int currentSequence = space().sequence();
472 const unsigned int topologyId = entity.type().id();
473 const IndexType entityIndex = indexSet_.index( entity ) ;
476 if( sequence_ != currentSequence ||
477 lastEntityIndex_ != entityIndex ||
478 lastTopologyId_ != topologyId )
481 lastEntityIndex_ = entityIndex ;
482 sequence_ = currentSequence;
483 lastTopologyId_ = topologyId ;
497 template<
class MassCaller,
class BasisFunctionSet,
class LocalFunction >
503 = getLocalInverseMassMatrixDefault ( caller, entity, geo, basisFunctionSet );
506 const int numDofs = lf.
size();
507 rhs_.resize( numDofs );
508 for(
int l = 0; l < numDofs; ++l )
512 multiply( numDofs, invMassMatrix, rhs_, lf );
516 template<
class LocalMatrix >
519 NoMassDummyCaller caller;
521 = getLocalInverseMassMatrixDefault ( caller, entity, geo, localMatrix.rangeBasisFunctionSet() );
523 const int cols = localMatrix.columns();
527 const int rows = localMatrix.rows();
528 for(
int i = 0; i < rows; ++i )
531 for(
int j = 0; j < cols; ++j )
532 rhs_[ j ] = localMatrix.get( i, j );
535 invMassMatrix.
mtv( rhs_, row_ );
538 for(
int j = 0; j < cols; ++j )
539 localMatrix.set( i, j, row_[ j ] );
544 template<
class LocalMatrix >
547 NoMassDummyCaller caller;
549 = getLocalInverseMassMatrixDefault ( caller, entity, geo, localMatrix.rangeBasisFunctionSet() );
551 const int cols = localMatrix.columns();
555 const int rows = localMatrix.rows();
556 for(
int i = 0; i < rows; ++i )
559 for(
int j = 0; j < cols; ++j )
560 rhs_[ j ] = localMatrix.get( j, i );
563 invMassMatrix.
mv( rhs_, row_ );
566 for(
int j = 0; j < cols; ++j )
567 localMatrix.set( j, i, row_[ j ] );
576 template<
class BasisFunctionSet,
class LocalFunction >
580 const int numDofs = lf.
size();
583 MatrixPairType& matrixPair =
584 getLocalInverseMassMatrix( entity, geo, basisFunctionSet, numDofs );
587 if( matrixPair.second )
589 const VectorType& diagonal = *matrixPair.second;
590 assert(
int(diagonal.size()) == numDofs );
594 const int nop = volQuad.nop();
595 assert(nop*dimRange == numDofs);
596 for(
int l=0, qt = 0; qt < nop; ++qt )
598 const auto intel = geo.integrationElement( volQuad.point(qt) );
599 for (
int r = 0; r < dimRange; ++r, ++l )
601 lf[ l ] *= diagonal[ l ] / intel;
610 rhs_.resize( numDofs );
611 for(
int l = 0; l < numDofs; ++l )
612 rhs_[ l ] = lf[ l ] * massVolInv;
614 const MatrixType& invMassMatrix = *matrixPair.first;
616 multiply( numDofs, invMassMatrix, rhs_, lf );
620 template <
class LocalMatrix>
622 setupInverseDiagonal(
const EntityType &entity,
const Geometry &geo,
623 const VectorType& refElemDiagonal,
624 LocalMatrix &localMatrix )
const
626 const int cols = localMatrix.columns();
628 assert(
int(refElemDiagonal.size()) == cols );
632 VectorType& elementDiagonal = rhs_;
633 elementDiagonal.resize( cols );
635 const int nop = volQuad.nop();
636 assert(nop*dimRange == cols);
637 for(
int l = 0, qt = 0; qt < nop; ++qt )
640 for (
int r = 0; r < dimRange; ++r,++l )
642 elementDiagonal[ l ] = refElemDiagonal[ l ] / intel;
645 return elementDiagonal;
648 template<
class LocalMatrix >
649 void rightMultiplyInverseLocally (
const EntityType &entity,
const Geometry &geo, LocalMatrix &localMatrix )
const
651 const int cols = localMatrix.columns();
652 MatrixPairType& matrixPair =
653 getLocalInverseMassMatrix( entity, geo, localMatrix.rangeBasisFunctionSet(), cols );
657 if( matrixPair.second )
659 const VectorType& elementDiagonal =
660 setupInverseDiagonal( entity, geo, *matrixPair.second, localMatrix );
663 const int rows = localMatrix.rows();
664 for(
int i = 0; i < rows; ++i )
668 for(
int j = 0; j < cols; ++j )
669 row_[ j ] = elementDiagonal[ j ] * localMatrix.get( i, j );
672 for(
int j = 0; j < cols; ++j )
673 localMatrix.set( i, j, row_[ j ] );
678 const MatrixType &invMassMatrix = *matrixPair.first;
685 const int rows = localMatrix.rows();
686 for(
int i = 0; i < rows; ++i )
690 for(
int j = 0; j < cols; ++j )
691 rhs_[ j ] = localMatrix.get( i, j ) * massVolInv;
694 invMassMatrix.mtv( rhs_, row_ );
697 for(
int j = 0; j < cols; ++j )
698 localMatrix.set( i, j, row_[ j ] );
704 template<
class LocalMatrix >
707 const int cols = localMatrix.columns();
708 MatrixPairType& matrixPair =
709 getLocalInverseMassMatrix( entity, geo, localMatrix.rangeBasisFunctionSet(), cols );
713 if( matrixPair.second )
716 setupInverseDiagonal( entity, geo, *matrixPair.second, localMatrix );
719 const int rows = localMatrix.rows();
720 for(
int i = 0; i < rows; ++i )
724 for(
int j = 0; j < cols; ++j )
725 row_[ j ] = elementDiagonal[ j ] * localMatrix.get( j, i );
728 for(
int j = 0; j < cols; ++j )
729 localMatrix.set( j, i, row_[ j ] );
734 const MatrixType &invMassMatrix = *matrixPair.first;
741 const int rows = localMatrix.rows();
742 for(
int i = 0; i < rows; ++i )
745 for(
int j = 0; j < cols; ++j )
746 rhs_[ j ] = localMatrix.get( j, i ) * massVolInv;
749 invMassMatrix.
mv( rhs_, row_ );
752 for(
int j = 0; j < cols; ++j )
753 localMatrix.set( j, i, row_[ j ] );
767 const std::vector<Dune::GeometryType>& geomTypes = geoInfo_.
geomTypes(0);
770 if( (geomTypes.size() == 1) && geomTypes[0].isSimplex() )
780 template<
class MassCaller,
class Matrix >
782 const Geometry &geo,
const BasisFunctionSetType &set,
783 std::size_t numDofs,
Matrix &matrix )
const
785 assert( numDofs == set.size() );
793 if( caller.hasMass() )
800 template <
class Matrix>
802 const EntityType& en,
804 const BasisFunctionSetType& set,
805 const VolumeQuadratureType& volQuad,
808 const bool applyIntegrationElement =
true )
const
810 const int volNop = volQuad.nop();
811 for(
int qp=0; qp<volNop; ++qp)
814 const double intel = ( applyIntegrationElement ) ?
815 ( volQuad.weight(qp) * geo.integrationElement(volQuad.point(qp)) ) : volQuad.weight(qp) ;
818 set.evaluateAll(volQuad[qp], phi_);
820 for(
int m=0; m<numDofs; ++m)
822 const RangeType& phi_m = phi_[m];
823 const ctype val = intel * (phi_m * phi_m);
826 for(
int k=m+1; k<numDofs; ++k)
828 const ctype val = intel * (phi_m * phi_[k]);
837 template <
class MassCallerType,
class Matrix>
839 MassCallerType& caller,
840 const EntityType& en,
842 const BasisFunctionSetType& set,
843 const VolumeQuadratureType& volQuad,
847 typedef typename MassCallerType :: MassFactorType MassFactorType;
850 const int volNop = volQuad.nop();
851 for(
int qp=0; qp<volNop; ++qp)
854 const double intel = volQuad.weight(qp)
855 * geo.integrationElement(volQuad.point(qp));
858 set.evaluateAll( volQuad[qp], phi_);
861 caller.mass( en, volQuad, qp, mass);
864 for(
int m=0; m<numDofs; ++m)
866 mass.mv( phi_[m], phiMass_[m] );
870 for(
int m=0; m<numDofs; ++m)
872 for(
int k=0; k<numDofs; ++k)
874 matrix[m][k] += intel * (phiMass_[m] * phi_[k]);
881 template <
class Matrix,
class Rhs,
class X>
882 void multiply(
const int size,
887 assert( (
int) matrix.rows() ==
size );
888 assert( (
int) matrix.cols() ==
size );
889 assert( (
int) rhs.size() ==
size );
891 for(
int row = 0; row <
size; ++ row )
893 RangeFieldType sum = 0;
895 typedef typename Matrix :: const_row_reference MatRow;
896 MatRow matRow = matrix[ row ];
899 for(
int col = 0; col <
size; ++ col )
901 sum += matRow[ col ] * rhs[ col ];
916 template<
class DiscreteFunctionSpace,
class VolumeQuadrature >
936 template<
class DiscreteFunctionSpace,
class VolumeQuadrature,
bool refElemScaling = true >
943 typedef typename BaseType :: EntityType EntityType;
950 template <
class MassCallerType,
class BasisFunction,
class LocalFunctionType>
952 const EntityType& en,
953 const BasisFunction &basisFunction,
954 LocalFunctionType& lf)
const
956 BaseType :: applyInverseDgOrthoNormalBasis( caller, en, basisFunction, lf );
958 template <
class MassCallerType,
class LocalFunctionType>
960 const EntityType& en,
961 LocalFunctionType& lf)
const
963 BaseType :: applyInverseDgOrthoNormalBasis( caller, en, lf.basisFunctionSet(), lf );
967 template <
class LocalFunctionType>
969 LocalFunctionType& lf)
const
971 typename BaseType :: NoMassDummyCaller caller;
975 template <
class BasisFunction,
class LocalFunctionType>
977 const BasisFunction &basisFunction,
978 LocalFunctionType& lf)
const
980 typename BaseType :: NoMassDummyCaller caller;
985 template<
class LocalFunction >
992 template<
class LocalMatrix >
999 template<
class LocalMatrix >
void solve(V1 &x, const V2 &b, bool doPivoting=true) const
Solve system A x = b.
void mv(const X &x, Y &y) const
y = A x
Definition: densematrix.hh:368
void invert(bool doPivoting=true)
Compute inverse.
void mtv(const X &x, Y &y) const
y = A^T x
Definition: densematrix.hh:387
constexpr size_type rows() const
number of rows
Definition: densematrix.hh:709
void resize(size_type r, size_type c, value_type v=value_type())
resize matrix to r × c
Definition: dynmatrix.hh:106
Error thrown if operations of a FieldMatrix fail.
Definition: densematrix.hh:126
BaseType::IndexType IndexType
index type */
Definition: adaptiveleafindexset.hh:125
const std ::vector< GeometryType > & geomTypes(unsigned int codim) const
returns vector with geometry tpyes this index set has indices for
Definition: allgeomtypes.hh:171
Interface class for basis function sets.
Definition: basisfunctionset.hh:32
const EntityType & entity() const
return entity
std::size_t size() const
return size of basis function set
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
DG Local Mass Matrix for arbitrary spaces.
Definition: localmassmatrix.hh:939
void applyInverse(const EntityType &en, LocalFunctionType &lf) const
apply local dg mass matrix to local function lf without mass factor
Definition: localmassmatrix.hh:968
void rightMultiplyInverse(LocalMatrix &localMatrix) const
compute localMatrix * M^-1
Definition: localmassmatrix.hh:993
void applyInverse(MassCallerType &caller, const EntityType &en, const BasisFunction &basisFunction, LocalFunctionType &lf) const
Definition: localmassmatrix.hh:951
void applyInverse(LocalFunction &lf) const
apply local dg mass matrix to local function lf without mass factor
Definition: localmassmatrix.hh:986
void leftMultiplyInverse(LocalMatrix &localMatrix) const
compute M^-1 * localMatrix
Definition: localmassmatrix.hh:1000
Local Mass Matrix inversion implementation, select the correct method in your implementation.
Definition: localmassmatrix.hh:34
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:801
void applyInverseDefault(MassCaller &caller, const EntityType &entity, const Geometry &geo, const BasisFunctionSet &basisFunctionSet, LocalFunction &lf) const
Definition: localmassmatrix.hh:498
void leftMultiplyInverse(LocalMatrix &localMatrix) const
compute M^-1 * localMatrix
Definition: localmassmatrix.hh:349
void leftMultiplyInverseLocally(const EntityType &entity, const Geometry &geo, LocalMatrix &localMatrix) const
compute M^-1 * localMatrix
Definition: localmassmatrix.hh:705
bool entityHasChanged(const EntityType &entity) const
returns true if the entity has been changed
Definition: localmassmatrix.hh:468
LocalMassMatrixImplementation(const DiscreteFunctionSpaceType &spc, int volQuadOrd)
constructor taking space and volume quadrature order
Definition: localmassmatrix.hh:217
void applyInverse(MassCaller &caller, const EntityType &entity, const BasisFunctionSet &basisFunctionSet, LocalFunction &lf) const
Definition: localmassmatrix.hh:295
bool affine() const
returns true if geometry mapping is affine
Definition: localmassmatrix.hh:256
double getAffineMassFactor(const Geometry &geo) const
return mass factor for diagonal mass matrix
Definition: localmassmatrix.hh:259
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:577
void applyInverse(const EntityType &entity, LocalFunction &lf) const
apply local dg mass matrix to local function lf without mass factor
Definition: localmassmatrix.hh:315
void applyInverse(LocalFunction &lf) const
apply local dg mass matrix to local function lf without mass factor
Definition: localmassmatrix.hh:330
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:781
void rightMultiplyInverseDgOrthoNormalBasis(LocalMatrix &localMatrix) const
compute localMatrix * M^-1
Definition: localmassmatrix.hh:407
int maxVolumeQuadratureOrder() const
return appropriate quadrature order, default is 2 * order()
Definition: localmassmatrix.hh:204
LocalMassMatrixImplementation(const ThisType &other)
copy constructor
Definition: localmassmatrix.hh:239
int volumeQuadratureOrder(const EntityType &entity) const
return appropriate quadrature order, default is 2 * order(entity)
Definition: localmassmatrix.hh:211
void rightMultiplyInverseDefault(const EntityType &entity, const Geometry &geo, LocalMatrix &localMatrix) const
compute localMatrix * M^-1
Definition: localmassmatrix.hh:517
void applyInverse(MassCaller &caller, const EntityType &entity, LocalFunction &lf) const
Definition: localmassmatrix.hh:308
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:838
void rightMultiplyInverse(LocalMatrix &localMatrix) const
compute localMatrix * M^-1
Definition: localmassmatrix.hh:337
void leftMultiplyInverseDefault(const EntityType &entity, const Geometry &geo, LocalMatrix &localMatrix) const
compute M^-1 * localMatrix
Definition: localmassmatrix.hh:545
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
bool setup() const
setup and return affinity
Definition: localmassmatrix.hh:760
void leftMultiplyInverseDgOrthoNormalBasis(LocalMatrix &localMatrix) const
compute M^-1 * localMatrix
Definition: localmassmatrix.hh:438
Local Mass Matrix for arbitrary spaces.
Definition: localmassmatrix.hh:919
Unique label for each type of entities that can occur in DUNE grids.
Definition: type.hh:114
Wrapper class for geometries.
Definition: geometry.hh:71
Volume integrationElement(const LocalCoordinate &local) const
Return the factor appearing in the integral transformation formula.
Definition: geometry.hh:265
bool affine() const
Return true if the geometry mapping is affine and false otherwise.
Definition: geometry.hh:197
Compute indices for geometry types, taking the dimension into account.
Definition: typeindex.hh:90
static constexpr std::size_t index(const GeometryType >)
Compute the index for the given geometry type over all dimensions.
Definition: typeindex.hh:138
static constexpr std::size_t size(std::size_t maxdim)
Compute total number of geometry types up to and including the given dimension.
Definition: typeindex.hh:125
A generic dynamic dense matrix.
Definition: matrix.hh:561
This file implements a dense matrix with dynamic numbers of rows and columns.
Implements a matrix constructed from a given type representing a field and compile-time given number ...
Dune namespace.
Definition: alignedallocator.hh:13
constexpr std::integral_constant< std::size_t, sizeof...(II)> size(std::integer_sequence< T, II... >)
Return the size of the sequence.
Definition: integersequence.hh:75
Specialize with 'true' if the grid is a Cartesian grid. Cartesian grids satisfy the following propert...
Definition: capabilities.hh:48
static int volumeOrder(const int k)
return quadrature order for volume quadratures for given polynomial order k
Definition: capabilities.hh:138
Helper class to check affinity of the grid's geometries.
Definition: checkgeomaffinity.hh:22
selects Obj::pointSetId if available, otherwise defaultValue (default is -1)
Definition: utility.hh:174
Helper classes to provide indices for geometrytypes for use in a vector.