dune-fem 2.8-git
spmatrix.hh
Go to the documentation of this file.
1#ifndef DUNE_FEM_SPMATRIX_HH
2#define DUNE_FEM_SPMATRIX_HH
3
4// system includes
5#include <algorithm>
6#include <array>
7#include <cstdlib>
8#include <iostream>
9#include <limits>
10#include <string>
11#include <utility>
12#include <vector>
13
14// local includes
28
29namespace Dune
30{
31
32 namespace Fem
33 {
34
36 template <class T, class IndexT = std::size_t,
37 class ValuesVector = std::vector< T >,
38 class IndicesVector = std::vector< IndexT > >
40 {
41 public:
43 typedef T field_type;
45 typedef IndexT size_type;
50
53 static const int firstCol = 0;
54
55 SparseRowMatrix(const ThisType& ) = delete;
56
58 explicit SparseRowMatrix() :
59 values_(0), columns_(0), rows_(0), dim_({{0,0}}), maxNzPerRow_(0), compressed_( false )
60 {}
61
65 values_(0), columns_(0), rows_(0), dim_({{0,0}}), maxNzPerRow_(0), compressed_( false )
66 {
67 reserve(rows,cols,nz);
68 }
69
71 void reserve(const size_type rows, const size_type cols, const size_type nz)
72 {
73 // if( (rows != dim_[0]) || (cols != dim_[1]) || (nz != maxNzPerRow_))
74 resize(rows,cols,nz);
75 clear();
76 }
77
79 size_type rows () const
80 {
81 return dim_[0];
82 }
83
85 size_type cols () const
86 {
87 return dim_[1];
88 }
89
91 void set(const size_type row, const size_type col, const field_type val)
92 {
93 assert((col>=0) && (col < dim_[1]));
94 assert((row>=0) && (row < dim_[0]));
95
96 const size_type column = colIndex(row,col) ;
97 assert( column != defaultCol && column != zeroCol );
98
99 values_ [ column ] = val;
100 columns_[ column ] = col;
101 }
102
104 void add(const size_type row, const size_type col, const field_type val)
105 {
106 assert((col>=0) && (col < dim_[1]));
107 assert((row>=0) && (row < dim_[0]));
108
109 const size_type column = colIndex(row,col) ;
110 assert( column != defaultCol && column != zeroCol );
111
112 values_ [ column ] += val;
113 columns_[ column ] = col;
114 }
115
117 template<class ArgDFType, class DestDFType>
118 void apply(const ArgDFType& f, DestDFType& ret ) const
119 {
120 constexpr auto blockSize = ArgDFType::DiscreteFunctionSpaceType::localBlockSize;
121 auto ret_it = ret.dbegin();
122
123 for(size_type row = 0; row<dim_[0]; ++row)
124 {
125 const size_type endrow = endRow( row );
126 (*ret_it) = 0.0;
127 for(size_type col = startRow( row ); col<endrow; ++col)
128 {
129 const auto realCol = columns_[ col ];
130
131 if( ! compressed_ && ((realCol == defaultCol) || (realCol == zeroCol)) )
132 continue;
133
134 const auto blockNr = realCol / blockSize ;
135 const auto dofNr = realCol % blockSize ;
136 (*ret_it) += values_[ col ] * f.dofVector()[ blockNr ][ dofNr ];
137 }
138
139 ++ret_it;
140 }
141 }
142
144 field_type get(const size_type row, const size_type col) const
145 {
146 assert((col>=0) && (col < dim_[1]));
147 assert((row>=0) && (row < dim_[0]));
148
149 const size_type endrow = endRow( row );
150 for( size_type i = startRow( row ); i<endrow; ++i )
151 {
152 if(columns_[ i ] == col)
153 return values_[ i ];
154 }
155 return 0;
156 }
157
159 void clear()
160 {
161 std::fill( values_.begin(), values_.end(), 0 );
162 for (auto &c : columns_) c = defaultCol;
163 }
164
166 void clearRow(const size_type row)
167 {
168 assert((row>=0) && (row < dim_[0]));
169
170 const size_type endrow = endRow( row );
171 for(size_type idx = startRow( row ); idx<endrow; ++idx )
172 {
173 values_[ idx ] = 0;
174 // if ( !compressed_ )
175 // columns_[idx] = zeroCol;
176 }
177 }
178
180 void scale(const size_type row, const size_type col, const field_type val)
181 {
182 assert((row>=0) && (row < rows()) );
183 assert((col>=0) && (col < cols()) );
184
185 const size_type column = colIndex(row,col) ;
186 assert( column != defaultCol && column != zeroCol );
187
188 // scale value
189 values_ [ column ] *= val;
190 }
191
195 {
196 return maxNzPerRow_;
197 }
198
202 {
203 return dim_[0] > 0 ? rows_[ dim_[0] ] : 0;
204 }
205
209 {
210 assert( row >= 0 && row < dim_[0] );
211 return endRow( row ) - startRow( row );
212 }
213
216 std::pair<const field_type, size_type> realValue(size_type index) const
217 {
218 return std::pair<const field_type, size_type>(values_[index], columns_[index]);
219 }
220
222 void print(std::ostream& s=std::cout, unsigned int offset=0) const
223 {
224 for(size_type row=0; row<dim_[0]; ++row)
225 {
226 const size_type endrow = endRow( row );
227 for( size_type pos = startRow( row ); pos<endrow; ++pos )
228 {
229 const auto rv(realValue(pos));
230 const auto column(rv.second);
231 const auto value(rv.first);
232 if((std::abs(value) > 1.e-15) && (column != defaultCol))
233 s << row+offset << " " << column+offset << " " << value << std::endl;
234 }
235 }
236 }
237
238 template <class SizeT, class NumericT >
239 void fillCSRStorage( std::vector< std::map<SizeT, NumericT> >& matrix ) const
240 {
241 matrix.resize( dim_[0] );
242
243 size_type thisCol = 0;
244 for(size_type row = 0; row<dim_[ 0 ]; ++row )
245 {
246 auto& matRow = matrix[ row ];
247 const size_type endrow = endRow( row );
248 for(size_type col = startRow( row ); col<endrow; ++col)
249 {
250 const size_type realCol = columns_[ col ];
251
252 if( ! compressed_ && (realCol == defaultCol || realCol == zeroCol) )
253 continue;
254
255 matRow[ realCol ] = values_[ thisCol ];
256 }
257 }
258 }
259
260 void compress ()
261 {
262 if( ! compressed_ && (dim_[0] != 0) && (dim_[1] != 0))
263 {
264 // determine first row nonZeros
265 size_type newpos = 0 ;
266 for( newpos = startRow( 0 ); newpos < endRow( 0 ); ++newpos )
267 {
268 if( columns_[ newpos ] == defaultCol )
269 break;
270 }
271
272 for( size_type row = 1; row<dim_[0]; ++row )
273 {
274 const size_type endrow = endRow( row );
275 size_type col = startRow( row );
276 // update new row start position
277 rows_[ row ] = newpos;
278 for(; col<endrow; ++col, ++newpos )
279 {
280 if( columns_[ col ] == defaultCol )
281 break ;
282
283 assert( newpos <= col );
284 values_ [ newpos ] = values_ [ col ];
285 columns_[ newpos ] = columns_[ col ];
286 }
287 }
288 rows_[ dim_[0] ] = newpos ;
289
290 // values_.resize( newpos );
291 // columns_.resize( newpos );
292 compressed_ = true ;
293 }
294 }
295
296 size_type startRow ( const size_type row ) const
297 {
298 return rows_[ row ];
299 }
300
301 size_type endRow ( const size_type row ) const
302 {
303 return rows_[ row+1 ];
304 }
305
306 std::tuple< ValuesVector&, IndicesVector&, IndicesVector& > exportCRS()
307 {
308 // only return internal data in compressed status
309 compress();
310 return std::tie(values_,columns_,rows_);
311 }
312
313 protected:
316 {
317 constexpr auto colVal = defaultCol;
318 values_.resize( rows*nz , 0 );
319 columns_.resize( rows*nz , colVal );
320 rows_.resize( rows+1 , 0 );
321 rows_[ 0 ] = 0;
322 for( size_type i=1; i <= rows; ++i )
323 {
324 rows_[ i ] = rows_[ i-1 ] + nz ;
325 }
326 compressed_ = false;
327
328 dim_[0] = rows;
329 dim_[1] = cols;
331 }
332
335 {
336 assert((col>=0) && (col < dim_[1]));
337 assert((row>=0) && (row < dim_[0]));
338
339 const size_type endR = endRow( row );
340 size_type i = startRow( row );
341 // find local column or empty spot
342 for( ; i < endR; ++i )
343 {
344 if( columns_[ i ] == defaultCol || columns_[ i ] == zeroCol || columns_[ i ] == col )
345 {
346 return i;
347 }
348 }
349
350 assert(0);
351 DUNE_THROW( InvalidStateException, "Could not store entry in sparse matrix - no space available" );
352
353 // TODO: implement resize with 2*nz
354 std::abort();
355 return defaultCol;
356
357 /*
358 if(columns_[ i ] == col)
359 return i; // column already in matrix
360 else if( columns_[ i ] == defaultCol )
361 { // add this column at end of this row
362 ++nonZeros_[row];
363 return i;
364 }
365 else
366 {
367 std::abort();
368 // TODO re-implement
369 //
370 ++nonZeros_[row];
371 // must shift this row to add col at the position i
372 auto j = nz_-1; // last column
373 if (columns_[row*nz_+j] != defaultCol)
374 { // new space available - so resize
375 resize( rows(), cols(), (2 * nz_) );
376 j++;
377 }
378 for(;j>i;--j)
379 {
380 columns_[row*nz_+j] = columns_[row*nz_+j-1];
381 values_[row*nz_+j] = values_[row*nz_+j-1];
382 }
383 columns_[row*nz_+i] = col;
384 values_[row*nz_+i] = 0;
385 return i;
386 }
387 */
388 }
389
390 ValuesVector values_;
391 IndicesVector columns_;
392 IndicesVector rows_;
393
394 std::array<size_type,2> dim_;
397 };
398
399
400
402 template< class DomainSpace, class RangeSpace,
405 {
406 protected:
407 template< class MatrixObject >
408 struct LocalMatrixTraits;
409
410 template< class MatrixObject >
411 class LocalMatrix;
412
413 public:
414 typedef DomainSpace DomainSpaceType;
415 typedef RangeSpace RangeSpaceType;
416 typedef typename DomainSpaceType::EntityType DomainEntityType;
417 typedef typename RangeSpaceType::EntityType RangeEntityType;
418 typedef typename DomainSpaceType::EntityType ColumnEntityType;
419 typedef typename RangeSpaceType::EntityType RowEntityType;
420
421 typedef typename DomainSpaceType::BlockMapperType DomainBlockMapperType;
423 typedef typename RangeSpaceType::BlockMapperType RangeBlockMapperType;
425 typedef Matrix MatrixType;
426 typedef typename MatrixType::size_type size_type;
427 typedef typename MatrixType::field_type field_type;
428
430
431 static const size_type domainLocalBlockSize = DomainSpaceType::dimRange;
432 static const size_type rangeLocalBlockSize = RangeSpaceType::dimRange;
433
434 typedef Dune::FieldMatrix< field_type, rangeLocalBlockSize, domainLocalBlockSize > MatrixBlockType;
436
438
444
448 const SolverParameter& param = SolverParameter() )
451 domainMapper_( domainSpace_.blockMapper() ),
452 rangeMapper_( rangeSpace_.blockMapper() ),
453 sequence_( -1 ),
454 matrix_(),
456 param.preconditionMethod({SolverParameter::none,SolverParameter::jacobi})
458 localMatrixStack_( *this )
459 {}
460
463 {
464 return domainSpace_;
465 }
466
469 {
470 return rangeSpace_;
471 }
472
473 protected:
476 {
477 return matrix_;
478 }
479
480 void finalizeAssembly() const { const_cast< ThisType& > (*this).compress(); }
481
482 public:
485 {
487 return matrix_;
488 }
489
490
493 {
495 }
496
501 [[deprecated("Use TemporaryLocalMatrix,LocalContribution and {get,add,set}LocalMatrix")]]
502 LocalMatrixType localMatrix( const DomainEntityType &domainEntity, const RangeEntityType &rangeEntity ) const
503 {
504 return LocalMatrixType( localMatrixStack_, domainEntity, rangeEntity );
505 }
506
511 [[deprecated("Use TemporaryLocalMatrix,LocalContribution and {get,add,set}LocalMatrix")]]
513 {
515 }
516
519 {
520 return LocalColumnObjectType( *this, domainEntity );
521 }
522
523 void unitRow( const size_type row )
524 {
525 for( unsigned int i=0, r = row * domainLocalBlockSize; i<domainLocalBlockSize; ++i, ++r )
526 {
527 matrix_.clearRow( r );
528 matrix_.set( r, r, 1.0 );
529 }
530 }
531
532 template <class LocalBlock>
533 void addBlock( const size_type row, const size_type col, const LocalBlock& block )
534 {
535 std::array< size_type, rangeLocalBlockSize > rows;
536 std::array< size_type, domainLocalBlockSize > cols;
537 for( unsigned int i=0, r = row * domainLocalBlockSize, c = col * domainLocalBlockSize; i<domainLocalBlockSize; ++i, ++r, ++c )
538 {
539 rows[ i ] = r;
540 cols[ i ] = c;
541 }
542
543 for( unsigned int i=0; i<domainLocalBlockSize; ++i )
544 {
545 for( unsigned int j=0; j<domainLocalBlockSize; ++j )
546 {
547 matrix_.add( rows[ i ], cols[ j ], block[ i ][ j ]);
548 }
549 }
550 }
551
552 template <class LocalBlock>
553 void setBlock( const size_type row, const size_type col, const LocalBlock& block )
554 {
555 std::array< size_type, rangeLocalBlockSize > rows;
556 std::array< size_type, domainLocalBlockSize > cols;
557 for( unsigned int i=0, r = row * domainLocalBlockSize, c = col * domainLocalBlockSize; i<domainLocalBlockSize; ++i, ++r, ++c )
558 {
559 rows[ i ] = r;
560 cols[ i ] = c;
561 }
562
563 for( unsigned int i=0; i<domainLocalBlockSize; ++i )
564 {
565 for( unsigned int j=0; j<domainLocalBlockSize; ++j )
566 {
567 matrix_.set( rows[ i ], cols[ j ], block[ i ][ j ]);
568 }
569 }
570 }
571
572 template< class LocalMatrix >
573 void addLocalMatrix ( const DomainEntityType &domainEntity, const RangeEntityType &rangeEntity, const LocalMatrix &localMat )
574 {
575 auto functor = [ &localMat, this ] ( std::pair< int, int > local, const std::pair< size_type, size_type >& global )
576 {
577 matrix_.add( global.first, global.second, localMat.get( local.first, local.second ) );
578 };
579
580 rangeMapper_.mapEach( rangeEntity, makePairFunctor( domainMapper_, domainEntity, functor ) );
581 }
582
583 template< class LocalMatrix, class Scalar >
584 void addScaledLocalMatrix ( const DomainEntityType &domainEntity, const RangeEntityType &rangeEntity, const LocalMatrix &localMat, const Scalar &s )
585 {
586 auto functor = [ &localMat, &s, this ] ( std::pair< int, int > local, const std::pair< size_type, size_type >& global )
587 {
588 matrix_.add( global.first, global.second, s * localMat.get( local.first, local.second ) );
589 };
590
591 rangeMapper_.mapEach( rangeEntity, makePairFunctor( domainMapper_, domainEntity, functor ) );
592 }
593
594 template< class LocalMatrix >
595 void setLocalMatrix ( const DomainEntityType &domainEntity, const RangeEntityType &rangeEntity, const LocalMatrix &localMat )
596 {
597 auto functor = [ &localMat, this ] ( std::pair< int, int > local, const std::pair< size_type, size_type >& global )
598 {
599 matrix_.set( global.first, global.second, localMat.get( local.first, local.second ) );
600 };
601
602 rangeMapper_.mapEach( rangeEntity, makePairFunctor( domainMapper_, domainEntity, functor ) );
603 }
604
605 template< class LocalMatrix >
606 void getLocalMatrix ( const DomainEntityType &domainEntity, const RangeEntityType &rangeEntity, LocalMatrix &localMat ) const
607 {
608 auto functor = [ &localMat, this ] ( std::pair< int, int > local, const std::pair< size_type, size_type >& global )
609 {
610 localMat.set( local.first, local.second, matrix_.get( global.first, global.second ) );
611 };
612
613 rangeMapper_.mapEach( rangeEntity, makePairFunctor( domainMapper_, domainEntity, functor ) );
614 }
615
617 void clear()
618 {
619 matrix_.clear();
620 }
621
623 void compress() { matrix_.compress(); }
624
625 template <class Set>
626 void reserve (const std::vector< Set >& sparsityPattern )
627 {
629 }
630
632 template <class Stencil>
633 void reserve(const Stencil &stencil, bool verbose = false )
634 {
635 if( sequence_ != domainSpace_.sequence() )
636 {
637 // if empty grid do nothing (can appear in parallel runs)
638 if( (domainSpace_.begin() != domainSpace_.end()) && (rangeSpace_.begin() != rangeSpace_.end()) )
639 {
640 // output some info
641 if( verbose )
642 {
643 std::cout << "Reserve Matrix with (" << rangeSpace_.size() << "," << domainSpace_.size()<< ")" << std::endl;
644 std::cout << "Max number of base functions = (" << rangeMapper_.maxNumDofs() << ","
645 << domainMapper_.maxNumDofs() << ")" << std::endl;
646 }
647 // reserve matrix
648 const auto nonZeros = std::max( static_cast<size_type>(stencil.maxNonZerosEstimate()*DomainSpaceType::localBlockSize),
649 matrix_.maxNzPerRow() );
650 matrix_.reserve( rangeSpace_.size(), domainSpace_.size(), nonZeros );
651 }
652 sequence_ = domainSpace_.sequence();
653 }
654 }
655
657 template< class DomainFunction, class RangeFunction >
658 void apply( const DomainFunction &arg, RangeFunction &dest ) const
659 {
660 // do matrix vector multiplication
661 matrix_.apply( arg, dest );
662 // communicate data
663 dest.communicate();
664 }
665
668 template < class DiscreteFunctionType >
669 void extractDiagonal( DiscreteFunctionType& diag ) const
670 {
671 assert( matrix_.rows() == matrix_.cols() );
672 const auto dofEnd = diag.dend();
673 size_type row = 0;
674 for( auto dofIt = diag.dbegin(); dofIt != dofEnd; ++dofIt, ++row )
675 {
676 assert( row < matrix_.rows() );
677 (*dofIt) = matrix_.get( row, row );
678 }
679 }
680
681 template <class Vector>
682 void setUnitRows( const Vector &rows )
683 {
684 const auto &auxiliaryDofs = domainSpace().auxiliaryDofs();
685 for (auto r : rows)
686 {
687 matrix_.clearRow(r);
688 matrix_.set(r,r,auxiliaryDofs.contains( r )? 0.0 : 1.0);
689 }
690 }
691
693 void resort()
694 {
695 DUNE_THROW(NotImplemented,"SpMatrixObject::resort is not implemented");
696 // this method does not even exist on SpMatrix!!!
697 // matrix_.resort();
698 }
699
700 protected:
709 };
710
711
712
714 template< class DomainSpace, class RangeSpace, class Matrix >
715 template< class MatrixObject >
716 struct SparseRowMatrixObject< DomainSpace, RangeSpace, Matrix >::LocalMatrixTraits
717 {
718 typedef DomainSpace DomainSpaceType;
719 typedef RangeSpace RangeSpaceType;
720
722
723 typedef typename SparseRowMatrixObjectType::template LocalMatrix< MatrixObject > LocalMatrixType;
724
725 typedef typename RangeSpaceType::RangeFieldType RangeFieldType;
727
730 };
731
732
733
735 template< class DomainSpace, class RangeSpace, class Matrix >
736 template< class MatrixObject >
737 class SparseRowMatrixObject< DomainSpace, RangeSpace, Matrix >::LocalMatrix
738 : public LocalMatrixDefault< LocalMatrixTraits< MatrixObject > >
739 {
740 public:
742 typedef MatrixObject MatrixObjectType;
743
746
747 private:
749
750 public:
752 typedef typename MatrixObjectType::MatrixType MatrixType;
753
756
759
762
767
768 typedef std::vector< typename RangeMapperType::SizeType > RowIndicesType;
769 typedef std::vector< typename DomainMapperType::SizeType > ColumnIndicesType;
770
772 LocalMatrix( const MatrixObjectType &matrixObject,
775 const DomainMapperType& domainMapper,
776 const RangeMapperType& rangeMapper )
778 matrix_( matrixObject.matrix() ),
779 domainMapper_( domainMapper ),
780 rangeMapper_( rangeMapper )
781 {}
782
783 LocalMatrix( const LocalMatrix & ) = delete;
784
785 void init( const DomainEntityType &domainEntity, const RangeEntityType &rangeEntity )
786 {
787 // initialize base functions sets
788 BaseType::init( domainEntity, rangeEntity );
789 // rows are determined by the range space
790 rowIndices_.resize( rangeMapper_.numDofs( rangeEntity ) );
791 rangeMapper_.mapEach( rangeEntity, AssignFunctor< RowIndicesType >( rowIndices_ ) );
792 // columns are determined by the domain space
793 columnIndices_.resize( domainMapper_.numDofs( domainEntity ) );
794 domainMapper_.mapEach( domainEntity, AssignFunctor< ColumnIndicesType >( columnIndices_ ) );
795 }
796
799 {
800 return rowIndices_.size();
801 }
802
805 {
806 return columnIndices_.size();
807 }
808
810 void add(size_type localRow, size_type localCol, DofType value)
811 {
812 assert( value == value );
813 assert( (localRow >= 0) && (localRow < rows()) );
814 assert( (localCol >= 0) && (localCol < columns()) );
815 matrix_.add( rowIndices_[ localRow ], columnIndices_[ localCol ], value );
816 }
817
819 DofType get(size_type localRow, size_type localCol) const
820 {
821 assert( (localRow >= 0) && (localRow < rows()) );
822 assert( (localCol >= 0) && (localCol < columns()) );
823 return matrix_.get( rowIndices_[ localRow ], columnIndices_[ localCol ] );
824 }
825
827 void set(size_type localRow, size_type localCol, DofType value)
828 {
829 assert( (localRow >= 0) && (localRow < rows()) );
830 assert( (localCol >= 0) && (localCol < columns()) );
831 matrix_.set( rowIndices_[ localRow ], columnIndices_[ localCol ], value );
832 }
833
835 void unitRow(size_type localRow)
836 {
837 assert( (localRow >= 0) && (localRow < rows()) );
838 matrix_.unitRow( rowIndices_[ localRow ] );
839 }
840
842 void clearRow( size_type localRow )
843 {
844 assert( (localRow >= 0) && (localRow < rows()) );
845 matrix_.clearRow( rowIndices_[localRow]);
846 }
847
849 void clearCol( size_type localCol )
850 {
851 assert( (localCol >= 0) && (localCol < columns()) );
852 matrix_.clearCol( columnIndices_[localCol] );
853 }
854
856 void clear()
857 {
858 const size_type nrows = rows();
859 for(size_type i=0; i < nrows; ++i )
860 matrix_.clearRow( rowIndices_[ i ] );
861 }
862
864 void resort()
865 {
866 DUNE_THROW(NotImplemented,"SpMatrixObject::LocalMatrix::resort is not implemented");
867 //const size_type nrows = rows();
868 //for(size_type i=0; i < nrows; ++i )
869 //matrix_.resortRow( rowIndices_[ i ] );
870 }
871
873 void scale( const DofType& value )
874 {
875 const size_type nrows = rows();
876 const size_type ncols = columns();
877 for(size_type i=0; i < nrows; ++i )
878 {
879 for( size_type j=0; j < ncols; ++j )
880 {
881 scale(i, j, value );
882 }
883 }
884 }
885
886 protected:
888 void scale(size_type localRow, size_type localCol, DofType value)
889 {
890 assert( (localRow >= 0) && (localRow < rows()) );
891 assert( (localCol >= 0) && (localCol < columns()) );
892 matrix_.scale( rowIndices_[ localRow ], columnIndices_[ localCol ], value );
893 }
894
895 protected:
901 };
902
903 } // namespace Fem
904
905} // namespace Dune
906
907#endif // #ifndef DUNE_FEM_SPMATRIX_HH
Dune::Fem::Double abs(const Dune::Fem::Double &a)
Definition: double.hh:942
double max(const Dune::Fem::Double &v, const double p)
Definition: double.hh:965
Definition: bindguard.hh:11
PairFunctor< Mapper, Entity, Functor > makePairFunctor(const Mapper &mapper, const Entity &entity, Functor functor)
Definition: operator/matrix/functor.hh:65
static constexpr T max(T a)
Definition: utility.hh:77
Definition: misc/functor.hh:31
RangeSpaceType::EntityType RangeEntityType
Definition: localmatrix.hh:66
DomainSpaceType::EntityType DomainEntityType
Definition: localmatrix.hh:65
Default implementation for local matrix classes.
Definition: localmatrix.hh:285
Definition: localmatrixwrapper.hh:48
default implementation for a general operator stencil
Definition: stencil.hh:35
int maxNonZerosEstimate() const
Return an upper bound for the maximum number of non-zero entries in all row.
Definition: stencil.hh:104
a simple wrapper class for sparsity patterns provide as vector< set< size_t > >
Definition: stencil.hh:191
Definition: columnobject.hh:12
SparseRowMatrix.
Definition: spmatrix.hh:40
void print(std::ostream &s=std::cout, unsigned int offset=0) const
print matrix
Definition: spmatrix.hh:222
size_type colIndex(size_type row, size_type col)
returns local col index for given global (row,col)
Definition: spmatrix.hh:334
std::array< size_type, 2 > dim_
Definition: spmatrix.hh:394
std::pair< const field_type, size_type > realValue(size_type index) const
Definition: spmatrix.hh:216
void fillCSRStorage(std::vector< std::map< SizeT, NumericT > > &matrix) const
Definition: spmatrix.hh:239
size_type endRow(const size_type row) const
Definition: spmatrix.hh:301
void clearRow(const size_type row)
set all entries in row to zero
Definition: spmatrix.hh:166
IndexT size_type
matrix index type
Definition: spmatrix.hh:45
size_type maxNzPerRow_
Definition: spmatrix.hh:395
void add(const size_type row, const size_type col, const field_type val)
add value to row,col entry
Definition: spmatrix.hh:104
std::tuple< ValuesVector &, IndicesVector &, IndicesVector & > exportCRS()
Definition: spmatrix.hh:306
void apply(const ArgDFType &f, DestDFType &ret) const
ret = A*f
Definition: spmatrix.hh:118
size_type numNonZeros(size_type row) const
Definition: spmatrix.hh:208
IndicesVector columns_
Definition: spmatrix.hh:391
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:159
static const size_type zeroCol
Definition: spmatrix.hh:52
SparseRowMatrix< field_type, size_type, ValuesVector, IndicesVector > ThisType
Definition: spmatrix.hh:46
static const size_type defaultCol
Definition: spmatrix.hh:51
size_type numNonZeros() const
Definition: spmatrix.hh:201
field_type get(const size_type row, const size_type col) const
return value of entry (row,col)
Definition: spmatrix.hh:144
static const int firstCol
Definition: spmatrix.hh:53
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:91
size_type startRow(const size_type row) const
Definition: spmatrix.hh:296
void resize(size_type rows, size_type cols, size_type nz)
resize matrix
Definition: spmatrix.hh:315
size_type rows() const
return number of rows
Definition: spmatrix.hh:79
IndicesVector rows_
Definition: spmatrix.hh:392
size_type cols() const
return number of columns
Definition: spmatrix.hh:85
SparseRowMatrix(const size_type rows, const size_type cols, const size_type nz)
Definition: spmatrix.hh:64
ValuesVector values_
Definition: spmatrix.hh:390
bool compressed_
Definition: spmatrix.hh:396
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:180
void compress()
Definition: spmatrix.hh:260
SparseRowMatrix(const ThisType &)=delete
size_type maxNzPerRow() const
Definition: spmatrix.hh:194
SparseRowMatrix()
construct matrix of zero size
Definition: spmatrix.hh:58
ThisType MatrixBaseType
Definition: spmatrix.hh:49
T field_type
matrix field type
Definition: spmatrix.hh:43
SparseRowMatrixObject.
Definition: spmatrix.hh:405
NonBlockMapper< DomainBlockMapperType, DomainSpaceType::localBlockSize > DomainMapperType
Definition: spmatrix.hh:422
RangeSpaceType::EntityType RangeEntityType
Definition: spmatrix.hh:417
const DomainSpaceType & domainSpace() const
get domain space (i.e. space that builds the rows)
Definition: spmatrix.hh:462
const DomainSpaceType & domainSpace_
Definition: spmatrix.hh:701
DomainMapperType domainMapper_
Definition: spmatrix.hh:703
LocalMatrixType localMatrix() const
Definition: spmatrix.hh:512
DomainSpaceType::EntityType DomainEntityType
Definition: spmatrix.hh:416
RangeSpaceType::BlockMapperType RangeBlockMapperType
Definition: spmatrix.hh:423
Matrix MatrixType
Definition: spmatrix.hh:425
MatrixType::size_type size_type
Definition: spmatrix.hh:426
LocalMatrixStackType localMatrixStack_
Definition: spmatrix.hh:708
MatrixType & exportMatrix() const
get reference to storage object
Definition: spmatrix.hh:484
ThisType LocalMatrixFactoryType
Definition: spmatrix.hh:440
LocalMatrixWrapper< LocalMatrixStackType > LocalMatrixType
Definition: spmatrix.hh:442
void setUnitRows(const Vector &rows)
Definition: spmatrix.hh:682
void addScaledLocalMatrix(const DomainEntityType &domainEntity, const RangeEntityType &rangeEntity, const LocalMatrix &localMat, const Scalar &s)
Definition: spmatrix.hh:584
SparseRowMatrixObject< DomainSpaceType, RangeSpaceType, MatrixType > ThisType
Definition: spmatrix.hh:429
RangeSpace RangeSpaceType
Definition: spmatrix.hh:415
void setBlock(const size_type row, const size_type col, const LocalBlock &block)
Definition: spmatrix.hh:553
RangeSpaceType::EntityType RowEntityType
Definition: spmatrix.hh:419
static const size_type rangeLocalBlockSize
Definition: spmatrix.hh:432
void reserve(const std::vector< Set > &sparsityPattern)
Definition: spmatrix.hh:626
LocalMatrix< ThisType > ObjectType
Definition: spmatrix.hh:439
ColumnObject< ThisType > LocalColumnObjectType
Definition: spmatrix.hh:443
void apply(const DomainFunction &arg, RangeFunction &dest) const
apply matrix to discrete function
Definition: spmatrix.hh:658
Fem::ObjectStack< LocalMatrixFactoryType > LocalMatrixStackType
Definition: spmatrix.hh:441
int sequence_
Definition: spmatrix.hh:705
SparseRowMatrixObject(const DomainSpaceType &domainSpace, const RangeSpaceType &rangeSpace, const SolverParameter &param=SolverParameter())
construct matrix object
Definition: spmatrix.hh:446
static const size_type domainLocalBlockSize
Definition: spmatrix.hh:431
LocalColumnObjectType localColumn(const DomainEntityType &domainEntity) const
get local column
Definition: spmatrix.hh:518
void extractDiagonal(DiscreteFunctionType &diag) const
Definition: spmatrix.hh:669
void setLocalMatrix(const DomainEntityType &domainEntity, const RangeEntityType &rangeEntity, const LocalMatrix &localMat)
Definition: spmatrix.hh:595
MatrixType::field_type field_type
Definition: spmatrix.hh:427
LocalMatrixType localMatrix(const DomainEntityType &domainEntity, const RangeEntityType &rangeEntity) const
Definition: spmatrix.hh:502
MatrixType PreconditionMatrixType
Definition: spmatrix.hh:437
DomainSpace DomainSpaceType
Definition: spmatrix.hh:411
DomainSpaceType::BlockMapperType DomainBlockMapperType
Definition: spmatrix.hh:421
DomainSpaceType::EntityType ColumnEntityType
Definition: spmatrix.hh:418
ObjectType * newObject() const
interface method from LocalMatrixFactory
Definition: spmatrix.hh:492
NonBlockMapper< RangeBlockMapperType, RangeSpaceType::localBlockSize > RangeMapperType
Definition: spmatrix.hh:424
Dune::FieldMatrix< field_type, rangeLocalBlockSize, domainLocalBlockSize > MatrixBlockType
Definition: spmatrix.hh:434
void reserve(const Stencil &stencil, bool verbose=false)
reserve memory
Definition: spmatrix.hh:633
void addLocalMatrix(const DomainEntityType &domainEntity, const RangeEntityType &rangeEntity, const LocalMatrix &localMat)
Definition: spmatrix.hh:573
MatrixType matrix_
Definition: spmatrix.hh:706
void compress()
compress matrix to a real CRS format
Definition: spmatrix.hh:623
void resort()
resort row numbering in matrix to have ascending numbering
Definition: spmatrix.hh:693
void getLocalMatrix(const DomainEntityType &domainEntity, const RangeEntityType &rangeEntity, LocalMatrix &localMat) const
Definition: spmatrix.hh:606
void unitRow(const size_type row)
Definition: spmatrix.hh:523
void addBlock(const size_type row, const size_type col, const LocalBlock &block)
Definition: spmatrix.hh:533
RangeMapperType rangeMapper_
Definition: spmatrix.hh:704
MatrixBlockType block_type
Definition: spmatrix.hh:435
void clear()
clear matrix
Definition: spmatrix.hh:617
bool preconditioning_
Definition: spmatrix.hh:707
const RangeSpaceType & rangeSpace_
Definition: spmatrix.hh:702
const RangeSpaceType & rangeSpace() const
get range space (i.e. space that builds the columns)
Definition: spmatrix.hh:468
MatrixType & matrix() const
get reference to storage object, for internal use
Definition: spmatrix.hh:475
void finalizeAssembly() const
Definition: spmatrix.hh:480
LocalMatrixTraits.
Definition: spmatrix.hh:717
SparseRowMatrixObject< DomainSpaceType, RangeSpaceType, Matrix > SparseRowMatrixObjectType
Definition: spmatrix.hh:721
RangeFieldType LittleBlockType
Definition: spmatrix.hh:726
RangeSpaceType::RangeFieldType RangeFieldType
Definition: spmatrix.hh:725
DomainSpace DomainSpaceType
Definition: spmatrix.hh:718
RangeSpace RangeSpaceType
Definition: spmatrix.hh:719
SparseRowMatrixObjectType::template LocalMatrix< MatrixObject > LocalMatrixType
Definition: spmatrix.hh:723
SparseRowMatrixObjectType::RangeMapperType RangeMapperType
Definition: spmatrix.hh:729
SparseRowMatrixObjectType::DomainMapperType DomainMapperType
Definition: spmatrix.hh:728
LocalMatrix.
Definition: spmatrix.hh:739
MatrixType & matrix_
Definition: spmatrix.hh:896
LocalMatrix(const LocalMatrix &)=delete
void clearRow(size_type localRow)
set matrix row to zero
Definition: spmatrix.hh:842
const DomainMapperType & domainMapper_
Definition: spmatrix.hh:897
void add(size_type localRow, size_type localCol, DofType value)
add value to matrix entry
Definition: spmatrix.hh:810
void init(const DomainEntityType &domainEntity, const RangeEntityType &rangeEntity)
Definition: spmatrix.hh:785
MatrixObjectType::MatrixType MatrixType
type of matrix
Definition: spmatrix.hh:752
void clearCol(size_type localCol)
set matrix column to zero
Definition: spmatrix.hh:849
void resort()
resort all global rows of matrix to have ascending numbering
Definition: spmatrix.hh:864
size_type columns() const
return number of columns
Definition: spmatrix.hh:804
ColumnIndicesType columnIndices_
Definition: spmatrix.hh:900
size_type rows() const
return number of rows
Definition: spmatrix.hh:798
RowIndicesType rowIndices_
Definition: spmatrix.hh:899
std::vector< typename DomainMapperType::SizeType > ColumnIndicesType
Definition: spmatrix.hh:769
Traits::RangeFieldType RangeFieldType
type of entries of little blocks
Definition: spmatrix.hh:755
DofType get(size_type localRow, size_type localCol) const
get matrix entry
Definition: spmatrix.hh:819
MatrixObject MatrixObjectType
type of matrix object
Definition: spmatrix.hh:742
std::vector< typename RangeMapperType::SizeType > RowIndicesType
Definition: spmatrix.hh:768
Traits::DomainMapperType DomainMapperType
type of nonblocked domain mapper
Definition: spmatrix.hh:764
Traits::RangeMapperType RangeMapperType
type of nonblocked domain mapper
Definition: spmatrix.hh:766
LocalMatrixTraits< MatrixObjectType > Traits
type of the traits
Definition: spmatrix.hh:745
LocalMatrix(const MatrixObjectType &matrixObject, const DomainSpaceType &domainSpace, const RangeSpaceType &rangeSpace, const DomainMapperType &domainMapper, const RangeMapperType &rangeMapper)
constructor
Definition: spmatrix.hh:772
const RangeMapperType & rangeMapper_
Definition: spmatrix.hh:898
void scale(const DofType &value)
scale local matrix with a certain value
Definition: spmatrix.hh:873
RangeFieldType DofType
type of the DoFs
Definition: spmatrix.hh:758
void scale(size_type localRow, size_type localCol, DofType value)
scale matrix entry with value
Definition: spmatrix.hh:888
void clear()
clear all entries belonging to local matrix
Definition: spmatrix.hh:856
Traits::LittleBlockType LittleBlockType
type of little blocks
Definition: spmatrix.hh:761
void set(size_type localRow, size_type localCol, DofType value)
set matrix entry to value
Definition: spmatrix.hh:827
void unitRow(size_type localRow)
set matrix row to zero except diagonla entry
Definition: spmatrix.hh:835
Definition: solver/parameter.hh:15
static const int none
Definition: solver/parameter.hh:40
static const int jacobi
Definition: solver/parameter.hh:45