1#ifndef DUNE_FEM_SPMATRIX_HH
2#define DUNE_FEM_SPMATRIX_HH
15#include <dune/fem/function/adaptivefunction/adaptivefunction.hh>
16#include <dune/fem/misc/functor.hh>
17#include <dune/fem/operator/common/localmatrix.hh>
18#include <dune/fem/operator/common/localmatrixwrapper.hh>
19#include <dune/fem/io/file/asciiparser.hh>
20#include <dune/fem/io/parameter.hh>
21#include <dune/fem/operator/common/operator.hh>
22#include <dune/fem/operator/matrix/columnobject.hh>
23#include <dune/fem/space/mapper/nonblockmapper.hh>
24#include <dune/fem/storage/objectstack.hh>
25#include <dune/fem/operator/common/stencil.hh>
26#include <dune/fem/operator/matrix/functor.hh>
27#include <dune/fem/solver/parameter.hh>
36 template <
class T,
class IndexT = std::size_t,
37 class ValuesVector = std::vector< T >,
38 class IndicesVector = std::vector< IndexT > >
53 static const int firstCol = 0;
59 values_(0), columns_(0), rows_(0), dim_({{0,0}}), maxNzPerRow_(0), compressed_(
false ), threading_( threading )
65 values_(0), columns_(0), rows_(0), dim_({{0,0}}), maxNzPerRow_(0), compressed_(
false ), threading_( threading )
79 template <
class Stencil>
85 for(
const auto& entry : sparsityPattern )
87 const auto& blockRow = entry.first;
88 const auto& blockColumnSet = entry.second;
91 const size_type nextRow = (blockRow + 1) * rowBlockSize;
92 for(
size_type row = blockRow * rowBlockSize; row < nextRow; ++row )
95 for(
const auto& blockColEntry : blockColumnSet )
97 size_type col = blockColEntry * colBlockSize;
98 for(
size_type c = 0; c<colBlockSize; ++c, ++col, ++column )
100 assert( column < endRow( row ) );
101 columns_[ column ] = col ;
123 assert((col>=0) && (col < dim_[1]));
124 assert((row>=0) && (row < dim_[0]));
127 assert( column != defaultCol && column != zeroCol );
129 values_ [ column ] = val;
130 columns_[ column ] = col;
136 assert((col>=0) && (col < dim_[1]));
137 assert((row>=0) && (row < dim_[0]));
140 assert( column != defaultCol && column != zeroCol );
142 values_ [ column ] += val;
143 columns_[ column ] = col;
147 template<
class ArgDFType,
class DestDFType>
148 void apply(
const ArgDFType& f, DestDFType& ret )
const
150 bool runParallel = threading_;
152 auto doApply = [
this, &f, &ret]()
154 constexpr auto blockSize = ArgDFType::DiscreteFunctionSpaceType::localBlockSize;
157 const auto slice = sliceBeginEnd( MPIManager::thread(), MPIManager::numThreads(), std::true_type() );
160 auto ret_it = ret.dofVector().find( slice.first );
162 const auto& fDofVec = f.dofVector();
163 for(
size_type row = slice.first; row<slice.second; ++row)
167 for(
size_type col = startRow( row ); col<endrow; ++col)
169 const auto realCol = columns_[ col ];
171 if( ! compressed_ && ((realCol == defaultCol) || (realCol == zeroCol)) )
174 const auto blockNr = realCol / blockSize ;
175 const auto dofNr = realCol % blockSize ;
176 (*ret_it) += values_[ col ] * fDofVec[ blockNr ][ dofNr ];
187 MPIManager :: run ( doApply );
205 assert((col>=0) && (col < dim_[1]));
206 assert((row>=0) && (row < dim_[0]));
209 for(
size_type i = startRow( row ); i<endrow; ++i )
211 if(columns_[ i ] == col)
220 std::fill( values_.begin(), values_.end(), 0 );
221 for (
auto &c : columns_) c = defaultCol;
227 assert((row>=0) && (row < dim_[0]));
230 for(
size_type idx = startRow( row ); idx<endrow; ++idx )
241 assert((row>=0) && (row <
rows()) );
242 assert((col>=0) && (col <
cols()) );
245 assert( column != defaultCol && column != zeroCol );
248 values_ [ column ] *= val;
262 return dim_[0] > 0 ? rows_[ dim_[0] ] : 0;
269 assert( row >= 0 && row < dim_[0] );
270 return endRow( row ) - startRow( row );
277 return std::pair<const field_type, size_type>(values_[index], columns_[index]);
281 void print(std::ostream& s=std::cout,
unsigned int offset=0)
const
286 for(
size_type pos = startRow( row ); pos<endrow; ++pos )
289 const auto column(rv.second);
290 const auto value(rv.first);
291 if((std::abs(value) > 1.e-15) && (column != defaultCol))
292 s << row+offset <<
" " << column+offset <<
" " << value << std::endl;
297 template <
class SizeT,
class NumericT >
298 void fillCSRStorage( std::vector< std::map<SizeT, NumericT> >& matrix )
const
300 matrix.resize( dim_[0] );
303 for(
size_type row = 0; row<dim_[ 0 ]; ++row )
305 auto& matRow = matrix[ row ];
307 for(
size_type col = startRow( row ); col<endrow; ++col)
309 const size_type realCol = columns_[ col ];
311 if( ! compressed_ && (realCol == defaultCol || realCol == zeroCol) )
314 matRow[ realCol ] = values_[ thisCol ];
321 if( ! compressed_ && (dim_[0] != 0) && (dim_[1] != 0))
325 for( newpos = startRow( 0 ); newpos < endRow( 0 ); ++newpos )
327 if( columns_[ newpos ] == defaultCol )
331 for(
size_type row = 1; row<dim_[0]; ++row )
336 rows_[ row ] = newpos;
337 for(; col<endrow; ++col, ++newpos )
339 if( columns_[ col ] == defaultCol )
342 assert( newpos <= col );
343 values_ [ newpos ] = values_ [ col ];
344 columns_[ newpos ] = columns_[ col ];
347 rows_[ dim_[0] ] = newpos ;
362 return rows_[ row+1 ];
365 std::tuple< ValuesVector&, IndicesVector&, IndicesVector& > exportCRS()
369 return std::tie(values_,columns_,rows_);
373 template<
class DiagType,
class ArgDFType,
class DestDFType,
class WType>
374 void forwardIterative(
const DiagType& diagInv,
const ArgDFType& b,
const DestDFType& xold, DestDFType& xnew,
const WType& w )
const
376 parallelIterative(diagInv, b, xold, xnew, w, std::true_type() );
380 template<
class DiagType,
class ArgDFType,
class DestDFType,
class WType>
381 void backwardIterative(
const DiagType& diagInv,
const ArgDFType& b,
const DestDFType& xold, DestDFType& xnew,
const WType& w )
const
383 parallelIterative(diagInv, b, xold, xnew, w, std::false_type() );
388 std::pair< size_type, size_type > sliceBeginEnd(
const size_type thread,
const size_type numThreads, std::true_type )
const
391 const size_type sliceBegin = (thread * sliceSize) ;
392 const size_type sliceEnd = (thread == numThreads-1 ) ? this->
rows(): (sliceBegin + sliceSize);
393 return std::make_pair( sliceBegin, sliceEnd );
397 std::pair< size_type, size_type > sliceBeginEnd(
const size_type thread,
const size_type numThreads, std::false_type )
const
399 std::pair< size_type, size_type > beginEnd = sliceBeginEnd( thread, numThreads, std::true_type() );
400 return std::make_pair( beginEnd.second-1, beginEnd.first-1 );
403 template<
class DiagType,
class ArgDFType,
class DestDFType,
class WType,
bool forward >
404 void parallelIterative(
const DiagType& diagInv,
const ArgDFType& b,
const DestDFType& xold, DestDFType& xnew,
405 const WType& w, std::integral_constant<bool, forward> direction )
const
407 bool runParallel = threading_ && (&xold != &xnew) ;
409 auto doIterate = [
this, &diagInv, &b, &xold, &xnew, &w, &direction] ()
412 const auto slice = sliceBeginEnd( MPIManager::thread(), MPIManager::numThreads(), direction );
415 b.dofVector().find( slice.first ),
417 xnew.dofVector().find( slice.first ),
428 MPIManager :: run ( doIterate );
430 catch (
const SingleThreadModeError& e )
444 template<
class DiagIt,
class ArgDFIt,
class DestDFType,
class DestDFIt,
445 class WType,
bool forward>
449 std::integral_constant<bool, forward> )
const
451 constexpr auto blockSize = DestDFType::DiscreteFunctionSpaceType::localBlockSize;
453 const auto nextRow = [&diag, &xit, &bit](
size_type &row)
455 if constexpr ( forward )
457 ++diag; ++xit; ++bit; ++row;
461 --diag; --xit; --bit; --row;
465 const auto& xOldVec = xold.dofVector();
466 for(; row != end; nextRow(row))
471 for(
size_type col = startRow( row ); col<endrow; ++col)
473 const auto realCol = columns_[ col ];
475 if( (realCol == row ) || (! compressed_ && ((realCol == defaultCol) || (realCol == zeroCol))) )
478 const auto blockNr = realCol / blockSize ;
479 const auto dofNr = realCol % blockSize ;
481 rhs -= values_[ col ] * xOldVec[ blockNr ][ dofNr ] ;
484 (*xit) = w * (rhs * (*diag));
490 constexpr auto colVal = defaultCol;
491 values_.resize(
rows*nz , 0 );
492 columns_.resize(
rows*nz , colVal );
493 rows_.resize(
rows+1 , 0 );
497 rows_[ i ] = rows_[ i-1 ] + nz ;
503 maxNzPerRow_ = nz+firstCol;
509 assert((col>=0) && (col < dim_[1]));
510 assert((row>=0) && (row < dim_[0]));
515 for( ; i < endR; ++i )
517 if( columns_[ i ] == defaultCol || columns_[ i ] == zeroCol || columns_[ i ] == col )
563 ValuesVector values_;
564 IndicesVector columns_;
567 std::array<size_type,2> dim_;
570 const bool threading_;
576 template<
class DomainSpace,
class RangeSpace,
581 template<
class MatrixObject >
584 template<
class MatrixObject >
588 typedef DomainSpace DomainSpaceType;
589 typedef RangeSpace RangeSpaceType;
590 typedef typename DomainSpaceType::EntityType DomainEntityType;
591 typedef typename RangeSpaceType::EntityType RangeEntityType;
592 typedef typename DomainSpaceType::EntityType ColumnEntityType;
593 typedef typename RangeSpaceType::EntityType RowEntityType;
595 typedef typename DomainSpaceType::BlockMapperType DomainBlockMapperType;
597 typedef typename RangeSpaceType::BlockMapperType RangeBlockMapperType;
605 static const size_type domainLocalBlockSize = DomainSpaceType::dimRange;
606 static const size_type rangeLocalBlockSize = RangeSpaceType::dimRange;
616 typedef LocalMatrixWrapper< LocalMatrixStackType > LocalMatrixType;
617 typedef ColumnObject< ThisType > LocalColumnObjectType;
622 const SolverParameter& param = SolverParameter() )
625 domainMapper_( domainSpace_.blockMapper() ),
626 rangeMapper_( rangeSpace_.blockMapper() ),
628 matrix_( param.threading() ),
629 localMatrixStack_( *this )
651 void finalizeAssembly()
const {
const_cast< ThisType&
> (*this).compress(); }
665 return new ObjectType( *
this, domainSpace_, rangeSpace_, domainMapper_, rangeMapper_ );
672 [[deprecated(
"Use TemporaryLocalMatrix,LocalContribution and {get,add,set}LocalMatrix")]]
673 LocalMatrixType
localMatrix(
const DomainEntityType &domainEntity,
const RangeEntityType &rangeEntity )
const
675 return LocalMatrixType( localMatrixStack_, domainEntity, rangeEntity );
682 [[deprecated(
"Use TemporaryLocalMatrix,LocalContribution and {get,add,set}LocalMatrix")]]
685 return LocalMatrixType( localMatrixStack_ );
689 LocalColumnObjectType
localColumn(
const DomainEntityType &domainEntity )
const
691 return LocalColumnObjectType( *
this, domainEntity );
694 void unitRow(
const size_type row )
696 for(
unsigned int i=0, r = row * domainLocalBlockSize; i<domainLocalBlockSize; ++i, ++r )
698 matrix_.clearRow( r );
699 matrix_.set( r, r, 1.0 );
703 template <
class LocalBlock>
704 void addBlock(
const size_type row,
const size_type col,
const LocalBlock& block )
706 std::array< size_type, rangeLocalBlockSize > rows;
707 std::array< size_type, domainLocalBlockSize > cols;
708 for(
unsigned int i=0, r = row * domainLocalBlockSize, c = col * domainLocalBlockSize; i<domainLocalBlockSize; ++i, ++r, ++c )
714 for(
unsigned int i=0; i<domainLocalBlockSize; ++i )
716 for(
unsigned int j=0; j<domainLocalBlockSize; ++j )
718 matrix_.add( rows[ i ], cols[ j ], block[ i ][ j ]);
723 template <
class LocalBlock>
724 void setBlock(
const size_type row,
const size_type col,
const LocalBlock& block )
726 std::array< size_type, rangeLocalBlockSize > rows;
727 std::array< size_type, domainLocalBlockSize > cols;
728 for(
unsigned int i=0, r = row * domainLocalBlockSize, c = col * domainLocalBlockSize; i<domainLocalBlockSize; ++i, ++r, ++c )
734 for(
unsigned int i=0; i<domainLocalBlockSize; ++i )
736 for(
unsigned int j=0; j<domainLocalBlockSize; ++j )
738 matrix_.set( rows[ i ], cols[ j ], block[ i ][ j ]);
743 template<
class LocalMatrix >
744 void addLocalMatrix (
const DomainEntityType &domainEntity,
const RangeEntityType &rangeEntity,
const LocalMatrix &localMat )
746 auto functor = [ &localMat, this ] ( std::pair< int, int > local,
const std::pair< size_type, size_type >& global )
748 matrix_.add( global.first, global.second, localMat.get( local.first, local.second ) );
751 rangeMapper_.mapEach( rangeEntity, makePairFunctor( domainMapper_, domainEntity, functor ) );
754 template<
class LocalMatrix,
class Scalar >
755 void addScaledLocalMatrix (
const DomainEntityType &domainEntity,
const RangeEntityType &rangeEntity,
const LocalMatrix &localMat,
const Scalar &s )
757 auto functor = [ &localMat, &s, this ] ( std::pair< int, int > local,
const std::pair< size_type, size_type >& global )
759 matrix_.add( global.first, global.second, s * localMat.get( local.first, local.second ) );
762 rangeMapper_.mapEach( rangeEntity, makePairFunctor( domainMapper_, domainEntity, functor ) );
765 template<
class LocalMatrix >
766 void setLocalMatrix (
const DomainEntityType &domainEntity,
const RangeEntityType &rangeEntity,
const LocalMatrix &localMat )
768 auto functor = [ &localMat, this ] ( std::pair< int, int > local,
const std::pair< size_type, size_type >& global )
770 matrix_.set( global.first, global.second, localMat.get( local.first, local.second ) );
773 rangeMapper_.mapEach( rangeEntity, makePairFunctor( domainMapper_, domainEntity, functor ) );
776 template<
class LocalMatrix >
777 void getLocalMatrix (
const DomainEntityType &domainEntity,
const RangeEntityType &rangeEntity, LocalMatrix &localMat )
const
779 auto functor = [ &localMat, this ] ( std::pair< int, int > local,
const std::pair< size_type, size_type >& global )
781 localMat.set( local.first, local.second, matrix_.get( global.first, global.second ) );
784 rangeMapper_.mapEach( rangeEntity, makePairFunctor( domainMapper_, domainEntity, functor ) );
797 void reserve (
const std::vector< Set >& sparsityPattern )
803 template <
class Stencil>
806 if( sequence_ != domainSpace_.sequence() )
809 if( (domainSpace_.begin() != domainSpace_.end()) && (rangeSpace_.begin() != rangeSpace_.end()) )
814 std::cout <<
"Reserve Matrix with (" << rangeSpace_.size() <<
"," << domainSpace_.size()<<
")" << std::endl;
815 std::cout <<
"Max number of base functions = (" << rangeMapper_.maxNumDofs() <<
","
816 << domainMapper_.maxNumDofs() <<
")" << std::endl;
820 const auto nonZeros = std::max(
static_cast<size_type
>(stencil.
maxNonZerosEstimate()*DomainSpaceType::localBlockSize),
821 matrix_.maxNzPerRow() );
822 matrix_.reserve( rangeSpace_.size(), domainSpace_.size(), nonZeros );
823 matrix_.fillPattern( stencil, RangeSpaceType::localBlockSize, DomainSpaceType::localBlockSize );
825 sequence_ = domainSpace_.sequence();
830 template<
class DomainFunction,
class RangeFunction >
831 void apply(
const DomainFunction &arg, RangeFunction &dest )
const
834 matrix_.apply( arg, dest );
841 template <
class DiscreteFunctionType >
844 assert( matrix_.rows() == matrix_.cols() );
845 const auto dofEnd = diag.
dend();
847 for(
auto dofIt = diag.
dbegin(); dofIt != dofEnd; ++dofIt, ++row )
849 assert( row < matrix_.rows() );
850 (*dofIt) = matrix_.get( row, row );
854 template <
class Container>
855 void setUnitRows(
const Container& unitRows,
const Container& auxRows )
857 for (
const auto& r : unitRows )
860 matrix_.set(r, r, 1.0);
863 for (
const auto& r : auxRows )
881 const DomainSpaceType &domainSpace_;
882 const RangeSpaceType &rangeSpace_;
883 DomainMapperType domainMapper_ ;
884 RangeMapperType rangeMapper_ ;
886 mutable MatrixType matrix_;
887 mutable LocalMatrixStackType localMatrixStack_;
893 template<
class DomainSpace,
class RangeSpace,
class Matrix >
894 template<
class MatrixObject >
897 typedef DomainSpace DomainSpaceType;
898 typedef RangeSpace RangeSpaceType;
904 typedef typename RangeSpaceType::RangeFieldType RangeFieldType;
905 typedef RangeFieldType LittleBlockType;
914 template<
class DomainSpace,
class RangeSpace,
class Matrix >
915 template<
class MatrixObject >
947 typedef std::vector< typename RangeMapperType::SizeType > RowIndicesType;
948 typedef std::vector< typename DomainMapperType::SizeType > ColumnIndicesType;
957 matrix_( matrixObject.
matrix() ),
958 domainMapper_( domainMapper ),
959 rangeMapper_( rangeMapper )
964 void init(
const DomainEntityType &domainEntity,
const RangeEntityType &rangeEntity )
967 BaseType::init( domainEntity, rangeEntity );
969 rowIndices_.resize( rangeMapper_.numDofs( rangeEntity ) );
970 rangeMapper_.mapEach( rangeEntity, AssignFunctor< RowIndicesType >( rowIndices_ ) );
972 columnIndices_.resize( domainMapper_.numDofs( domainEntity ) );
973 domainMapper_.mapEach( domainEntity, AssignFunctor< ColumnIndicesType >( columnIndices_ ) );
979 return rowIndices_.size();
985 return columnIndices_.size();
989 void add(size_type localRow, size_type localCol,
DofType value)
991 assert( value == value );
992 assert( (localRow >= 0) && (localRow < rows()) );
993 assert( (localCol >= 0) && (localCol < columns()) );
994 matrix_.add( rowIndices_[ localRow ], columnIndices_[ localCol ], value );
1000 assert( (localRow >= 0) && (localRow < rows()) );
1001 assert( (localCol >= 0) && (localCol < columns()) );
1002 return matrix_.get( rowIndices_[ localRow ], columnIndices_[ localCol ] );
1008 assert( (localRow >= 0) && (localRow < rows()) );
1009 assert( (localCol >= 0) && (localCol < columns()) );
1010 matrix_.set( rowIndices_[ localRow ], columnIndices_[ localCol ], value );
1016 assert( (localRow >= 0) && (localRow < rows()) );
1017 matrix_.unitRow( rowIndices_[ localRow ] );
1023 assert( (localRow >= 0) && (localRow < rows()) );
1024 matrix_.clearRow( rowIndices_[localRow]);
1030 assert( (localCol >= 0) && (localCol < columns()) );
1031 matrix_.clearCol( columnIndices_[localCol] );
1037 const size_type nrows = rows();
1038 for(size_type i=0; i < nrows; ++i )
1039 matrix_.clearRow( rowIndices_[ i ] );
1054 const size_type nrows = rows();
1055 const size_type ncols = columns();
1056 for(size_type i=0; i < nrows; ++i )
1058 for( size_type j=0; j < ncols; ++j )
1060 scale(i, j, value );
1069 assert( (localRow >= 0) && (localRow < rows()) );
1070 assert( (localCol >= 0) && (localCol < columns()) );
1071 matrix_.scale( rowIndices_[ localRow ], columnIndices_[ localCol ], value );
1078 RowIndicesType rowIndices_;
1079 ColumnIndicesType columnIndices_;
ConstDofIteratorType dbegin() const
Obtain the constant iterator pointing to the first dof.
Definition: discretefunction.hh:761
ConstDofIteratorType dend() const
Obtain the constant iterator pointing to the last dof.
Definition: discretefunction.hh:773
Default implementation for local matrix classes.
Definition: localmatrix.hh:287
Exception thrown when a code segment that is supposed to be only accessed in single thread mode is ac...
Definition: mpimanager.hh:43
LocalMatrix.
Definition: spmatrix.hh:918
void clearRow(size_type localRow)
set matrix row to zero
Definition: spmatrix.hh:1021
void add(size_type localRow, size_type localCol, DofType value)
add value to matrix entry
Definition: spmatrix.hh:989
MatrixObjectType::MatrixType MatrixType
type of matrix
Definition: spmatrix.hh:931
void clearCol(size_type localCol)
set matrix column to zero
Definition: spmatrix.hh:1028
void resort()
resort all global rows of matrix to have ascending numbering
Definition: spmatrix.hh:1043
size_type columns() const
return number of columns
Definition: spmatrix.hh:983
size_type rows() const
return number of rows
Definition: spmatrix.hh:977
Traits::RangeFieldType RangeFieldType
type of entries of little blocks
Definition: spmatrix.hh:934
DofType get(size_type localRow, size_type localCol) const
get matrix entry
Definition: spmatrix.hh:998
MatrixObject MatrixObjectType
type of matrix object
Definition: spmatrix.hh:921
Traits::DomainMapperType DomainMapperType
type of nonblocked domain mapper
Definition: spmatrix.hh:943
Traits::RangeMapperType RangeMapperType
type of nonblocked domain mapper
Definition: spmatrix.hh:945
LocalMatrixTraits< MatrixObjectType > Traits
type of the traits
Definition: spmatrix.hh:924
LocalMatrix(const MatrixObjectType &matrixObject, const DomainSpaceType &domainSpace, const RangeSpaceType &rangeSpace, const DomainMapperType &domainMapper, const RangeMapperType &rangeMapper)
constructor
Definition: spmatrix.hh:951
void scale(const DofType &value)
scale local matrix with a certain value
Definition: spmatrix.hh:1052
RangeFieldType DofType
type of the DoFs
Definition: spmatrix.hh:937
void scale(size_type localRow, size_type localCol, DofType value)
scale matrix entry with value
Definition: spmatrix.hh:1067
void clear()
clear all entries belonging to local matrix
Definition: spmatrix.hh:1035
Traits::LittleBlockType LittleBlockType
type of little blocks
Definition: spmatrix.hh:940
void set(size_type localRow, size_type localCol, DofType value)
set matrix entry to value
Definition: spmatrix.hh:1006
void unitRow(size_type localRow)
set matrix row to zero except diagonla entry
Definition: spmatrix.hh:1014
SparseRowMatrixObject.
Definition: spmatrix.hh:579
const DomainSpaceType & domainSpace() const
get domain space (i.e. space that builds the rows)
Definition: spmatrix.hh:633
LocalMatrixType localMatrix() const
Definition: spmatrix.hh:683
MatrixType & exportMatrix() const
get reference to storage object
Definition: spmatrix.hh:655
void apply(const DomainFunction &arg, RangeFunction &dest) const
apply matrix to discrete function
Definition: spmatrix.hh:831
SparseRowMatrixObject(const DomainSpaceType &domainSpace, const RangeSpaceType &rangeSpace, const SolverParameter ¶m=SolverParameter())
construct matrix object
Definition: spmatrix.hh:620
LocalColumnObjectType localColumn(const DomainEntityType &domainEntity) const
get local column
Definition: spmatrix.hh:689
void extractDiagonal(DiscreteFunctionType &diag) const
Definition: spmatrix.hh:842
LocalMatrixType localMatrix(const DomainEntityType &domainEntity, const RangeEntityType &rangeEntity) const
Definition: spmatrix.hh:673
ObjectType * newObject() const
interface method from LocalMatrixFactory
Definition: spmatrix.hh:663
void reserve(const Stencil &stencil, bool verbose=false)
reserve memory
Definition: spmatrix.hh:804
void compress()
compress matrix to a real CRS format
Definition: spmatrix.hh:794
void resort()
resort row numbering in matrix to have ascending numbering
Definition: spmatrix.hh:873
void clear()
clear matrix
Definition: spmatrix.hh:788
const RangeSpaceType & rangeSpace() const
get range space (i.e. space that builds the columns)
Definition: spmatrix.hh:639
MatrixType & matrix() const
get reference to storage object, for internal use
Definition: spmatrix.hh:646
SparseRowMatrix.
Definition: spmatrix.hh:40
void print(std::ostream &s=std::cout, unsigned int offset=0) const
print matrix
Definition: spmatrix.hh:281
size_type colIndex(size_type row, size_type col)
returns local col index for given global (row,col)
Definition: spmatrix.hh:507
std::pair< const field_type, size_type > realValue(size_type index) const
Definition: spmatrix.hh:275
SparseRowMatrix(const bool threading=true)
construct matrix of zero size
Definition: spmatrix.hh:58
void clearRow(const size_type row)
set all entries in row to zero
Definition: spmatrix.hh:225
IndexT size_type
matrix index type
Definition: spmatrix.hh:45
void add(const size_type row, const size_type col, const field_type val)
add value to row,col entry
Definition: spmatrix.hh:134
void apply(const ArgDFType &f, DestDFType &ret) const
ret = A*f
Definition: spmatrix.hh:148
size_type numNonZeros(size_type row) const
Definition: spmatrix.hh:267
void reserve(const size_type rows, const size_type cols, const size_type nz)
reserve memory for given rows, columns and number of non zeros
Definition: spmatrix.hh:71
void clear()
set all matrix entries to zero
Definition: spmatrix.hh:218
size_type numNonZeros() const
Definition: spmatrix.hh:260
void forwardIterative(const DiagType &diagInv, const ArgDFType &b, const DestDFType &xold, DestDFType &xnew, const WType &w) const
Apply Jacobi/SOR method.
Definition: spmatrix.hh:374
field_type get(const size_type row, const size_type col) const
return value of entry (row,col)
Definition: spmatrix.hh:203
void doParallelIterative(DiagIt diag, ArgDFIt bit, const DestDFType &xold, DestDFIt xit, const WType &w, size_type row, const size_type end, std::integral_constant< bool, forward >) const
Apply Jacobi/SOR method.
Definition: spmatrix.hh:446
void set(const size_type row, const size_type col, const field_type val)
set entry to value (also setting 0 will result in an entry)
Definition: spmatrix.hh:121
void resize(size_type rows, size_type cols, size_type nz)
resize matrix
Definition: spmatrix.hh:488
size_type rows() const
return number of rows
Definition: spmatrix.hh:109
size_type cols() const
return number of columns
Definition: spmatrix.hh:115
void scale(const size_type row, const size_type col, const field_type val)
scale all entries in row with a given value
Definition: spmatrix.hh:239
size_type maxNzPerRow() const
Definition: spmatrix.hh:253
void backwardIterative(const DiagType &diagInv, const ArgDFType &b, const DestDFType &xold, DestDFType &xnew, const WType &w) const
Apply Jacobi/SOR method.
Definition: spmatrix.hh:381
SparseRowMatrix(const size_type rows, const size_type cols, const size_type nz, const bool threading=true)
Definition: spmatrix.hh:64
ThisType MatrixBaseType
Definition: spmatrix.hh:49
void fillPattern(const Stencil &stencil, const size_type rowBlockSize, const size_type colBlockSize)
reserve memory for given rows, columns and number of non zeros
Definition: spmatrix.hh:80
T field_type
matrix field type
Definition: spmatrix.hh:43
a simple wrapper class for sparsity patterns provide as vector< set< size_t > >
Definition: stencil.hh:234
default implementation for a general operator stencil
Definition: stencil.hh:35
const GlobalStencilType & globalStencil() const
Return the full stencil.
Definition: stencil.hh:116
int maxNonZerosEstimate() const
Return an upper bound for the maximum number of non-zero entries in all rows.
Definition: stencil.hh:123
forward declaration
Definition: discretefunction.hh:51
A dense n x m matrix.
Definition: fmatrix.hh:117
Default exception if a function was called while the object is not in a valid state for that function...
Definition: exceptions.hh:281
A generic dynamic dense matrix.
Definition: matrix.hh:561
A::size_type size_type
Type for indices and sizes.
Definition: matrix.hh:577
typename Imp::BlockTraits< T >::field_type field_type
Export the type representing the underlying field.
Definition: matrix.hh:565
Default exception for dummy implementations.
Definition: exceptions.hh:263
#define DUNE_THROW(E, m)
Definition: exceptions.hh:218
constexpr auto max
Function object that returns the greater of the given values.
Definition: hybridutilities.hh:484
typename Overloads::ScalarType< std::decay_t< V > >::type Scalar
Element type of some SIMD type.
Definition: interface.hh:235
Dune namespace.
Definition: alignedallocator.hh:13
LocalMatrixTraits.
Definition: spmatrix.hh:896