dune-istl
2.4
|
A sparse block matrix with compressed row storage. More...
#include <dune/istl/bcrsmatrix.hh>
Classes | |
class | CreateIterator |
Iterator class for sequential creation of blocks More... | |
class | Deallocator |
Class used by shared_ptr to deallocate memory using the proper allocator. More... | |
class | RealRowIterator |
Iterator access to matrix rows More... | |
Public Types | |
enum | BuildStage { notbuilt =0, notAllocated =0, building =1, rowSizesBuilt =2, built =3 } |
enum | { blocklevel = B::blocklevel+1 } |
increment block level counter More... | |
enum | BuildMode { row_wise, random, implicit, unknown } |
we support two modes More... | |
typedef B::field_type | field_type |
export the type representing the field More... | |
typedef B | block_type |
export the type representing the components More... | |
typedef A | allocator_type |
export the allocator type More... | |
typedef CompressedBlockVectorWindow< B, A > | row_type |
implement row_type with compressed vector More... | |
typedef A::size_type | size_type |
The type for the index access and the size. More... | |
typedef ::Dune::CompressionStatistics< size_type > | CompressionStatistics |
The type for the statistics object returned by compress() More... | |
typedef RealRowIterator< row_type > | iterator |
The iterator over the (mutable matrix rows. More... | |
typedef RealRowIterator< row_type > | Iterator |
typedef Iterator | RowIterator |
rename the iterators for easier access More... | |
typedef row_type::Iterator | ColIterator |
Iterator for the entries of each row. More... | |
typedef RealRowIterator< const row_type > | const_iterator |
The const iterator over the matrix rows. More... | |
typedef RealRowIterator< const row_type > | ConstIterator |
typedef ConstIterator | ConstRowIterator |
rename the const row iterator for easier access More... | |
typedef row_type::ConstIterator | ConstColIterator |
Const iterator to the entries of a row. More... | |
Public Member Functions | |
row_type & | operator[] (size_type i) |
random access to the rows More... | |
const row_type & | operator[] (size_type i) const |
same for read only access More... | |
Iterator | begin () |
Get iterator to first row. More... | |
Iterator | end () |
Get iterator to one beyond last row. More... | |
Iterator | beforeEnd () |
Iterator | beforeBegin () |
ConstIterator | begin () const |
Get const iterator to first row. More... | |
ConstIterator | end () const |
Get const iterator to one beyond last row. More... | |
ConstIterator | beforeEnd () const |
ConstIterator | beforeBegin () const |
BCRSMatrix () | |
an empty matrix More... | |
BCRSMatrix (size_type _n, size_type _m, size_type _nnz, BuildMode bm) | |
matrix with known number of nonzeroes More... | |
BCRSMatrix (size_type _n, size_type _m, BuildMode bm) | |
matrix with unknown number of nonzeroes More... | |
BCRSMatrix (size_type _n, size_type _m, size_type _avg, double _overflowsize, BuildMode bm) | |
construct matrix with a known average number of entries per row More... | |
BCRSMatrix (const BCRSMatrix &Mat) | |
copy constructor More... | |
~BCRSMatrix () | |
destructor More... | |
void | setBuildMode (BuildMode bm) |
Sets the build mode of the matrix. More... | |
void | setSize (size_type rows, size_type columns, size_type nnz=0) |
Set the size of the matrix. More... | |
void | setImplicitBuildModeParameters (size_type _avg, double _overflow) |
Set parameters needed for creation in implicit build mode. More... | |
BCRSMatrix & | operator= (const BCRSMatrix &Mat) |
assignment More... | |
BCRSMatrix & | operator= (const field_type &k) |
Assignment from a scalar. More... | |
CreateIterator | createbegin () |
get initial create iterator More... | |
CreateIterator | createend () |
get create iterator pointing to one after the last block More... | |
void | setrowsize (size_type i, size_type s) |
set number of indices in row i to s More... | |
size_type | getrowsize (size_type i) const |
get current number of indices in row i More... | |
void | incrementrowsize (size_type i, size_type s=1) |
increment size of row i by s (1 by default) More... | |
void | endrowsizes () |
indicate that size of all rows is defined More... | |
void | addindex (size_type row, size_type col) |
add index (row,col) to the matrix More... | |
template<typename It > | |
void | setIndices (size_type row, It begin, It end) |
Set all column indices for row from the given iterator range. More... | |
void | endindices () |
indicate that all indices are defined, check consistency More... | |
B & | entry (size_type row, size_type col) |
Returns reference to entry (row,col) of the matrix. More... | |
CompressionStatistics | compress () |
Finishes the buildstage in implicit mode. More... | |
BCRSMatrix & | operator*= (const field_type &k) |
vector space multiplication with scalar More... | |
BCRSMatrix & | operator/= (const field_type &k) |
vector space division by scalar More... | |
BCRSMatrix & | operator+= (const BCRSMatrix &b) |
Add the entries of another matrix to this one. More... | |
BCRSMatrix & | operator-= (const BCRSMatrix &b) |
Substract the entries of another matrix to this one. More... | |
BCRSMatrix & | axpy (field_type alpha, const BCRSMatrix &b) |
Add the scaled entries of another matrix to this one. More... | |
template<class X , class Y > | |
void | mv (const X &x, Y &y) const |
y = A x More... | |
template<class X , class Y > | |
void | umv (const X &x, Y &y) const |
y += A x More... | |
template<class X , class Y > | |
void | mmv (const X &x, Y &y) const |
y -= A x More... | |
template<class X , class Y > | |
void | usmv (const field_type &alpha, const X &x, Y &y) const |
y += alpha A x More... | |
template<class X , class Y > | |
void | mtv (const X &x, Y &y) const |
y = A^T x More... | |
template<class X , class Y > | |
void | umtv (const X &x, Y &y) const |
y += A^T x More... | |
template<class X , class Y > | |
void | mmtv (const X &x, Y &y) const |
y -= A^T x More... | |
template<class X , class Y > | |
void | usmtv (const field_type &alpha, const X &x, Y &y) const |
y += alpha A^T x More... | |
template<class X , class Y > | |
void | umhv (const X &x, Y &y) const |
y += A^H x More... | |
template<class X , class Y > | |
void | mmhv (const X &x, Y &y) const |
y -= A^H x More... | |
template<class X , class Y > | |
void | usmhv (const field_type &alpha, const X &x, Y &y) const |
y += alpha A^H x More... | |
FieldTraits< field_type >::real_type | frobenius_norm2 () const |
square of frobenius norm, need for block recursion More... | |
FieldTraits< field_type >::real_type | frobenius_norm () const |
frobenius norm: sqrt(sum over squared values of entries) More... | |
FieldTraits< field_type >::real_type | infinity_norm () const |
infinity norm (row sum norm, how to generalize for blocks?) More... | |
FieldTraits< field_type >::real_type | infinity_norm_real () const |
simplified infinity norm (uses Manhattan norm for complex values) More... | |
size_type | N () const |
number of rows (counted in blocks) More... | |
size_type | M () const |
number of columns (counted in blocks) More... | |
size_type | nonzeroes () const |
number of blocks that are stored (the number of blocks that possibly are nonzero) More... | |
BuildStage | buildStage () const |
The current build stage of the matrix. More... | |
BuildMode | buildMode () const |
The currently selected build mode of the matrix. More... | |
bool | exists (size_type i, size_type j) const |
return true if (i,j) is in pattern More... | |
Protected Types | |
typedef std::map< std::pair< size_type, size_type >, B > | OverflowType |
Protected Member Functions | |
void | setWindowPointers (ConstRowIterator row) |
void | setColumnPointers (ConstRowIterator row) |
Copy row sizes from iterator range starting at row and set column index pointers for all rows. More... | |
void | setDataPointers () |
Set data pointers for all rows. More... | |
void | copyWindowStructure (const BCRSMatrix &Mat) |
Copy the window structure from another matrix. More... | |
void | deallocate (bool deallocateRows=true) |
deallocate memory of the matrix. More... | |
void | allocate (size_type rows, size_type columns, size_type allocationSize_, bool allocateRows, bool allocate_data) |
Allocate memory for the matrix structure. More... | |
void | allocateData () |
void | implicit_allocate (size_type _n, size_type _m) |
organizes allocation implicit mode calculates correct array size to be allocated and sets the the window pointers to their correct positions for insertion. internally uses allocate() for the real allocation. More... | |
Protected Attributes | |
BuildMode | build_mode |
BuildStage | ready |
A::template rebind< B >::other | allocator_ |
A::template rebind< row_type >::other | rowAllocator_ |
A::template rebind< size_type >::other | sizeAllocator_ |
size_type | n |
size_type | m |
size_type | nnz |
size_type | allocationSize |
row_type * | r |
B * | a |
std::shared_ptr< size_type > | j |
size_type | avg |
double | overflowsize |
OverflowType | overflow |
A sparse block matrix with compressed row storage.
Implements a block compressed row storage scheme. The block type B can be any type implementing the matrix interface.
Different ways to build up a compressed row storage matrix are supported:
Error checking: no error checking is provided normally. Setting the compile time switch DUNE_ISTL_WITH_CHECKING enables error checking.
Details:
Rows are built up in sequential order. Size of the row and the column indices are defined. A row can be used as soon as it is initialized. With respect to memory there are two variants of this scheme: (a) number of non-zeroes known in advance (application finite difference schemes), (b) number of non-zeroes not known in advance (application: Sparse LU, ILU(n)).
For general finite element implementations the number of rows n is known, the number of non-zeroes might also be known (e.g. #edges + #nodes for P1) but the size of a row and the indices of a row can not be defined in sequential order.
implicit
mode tries to make use of that knowledge by allocating memory based on that average. Entries in rows with more non-zeroes than the average value are written to an overflow area during the initial assembly phase, up to a specified maximum number of overflow entries that must not be exceeded. After all indices are added a compression step optimizes the matrix and integrates any entries from the overflow area into the standard BCRS storage scheme.To use this mode use the following methods:
Construct the matrix via
Here, the parameter _avg
denotes the average number of matrix entries per row, while _overflowsize
reserves _n * _overflowsize * _avg
entries in the overflow area.
Start filling your matrix by calling entry(size_type row, size_type col), which returns the corresponding matrix entry, creating it on the fly if it did not exist yet. Please note that this method may be slightly slower than accessing entries via matrix[row][col]
after the initial assembly because of the additional overhead of searching the overflow area. The matrix pattern is created by implicitly by simply accessing nonzero entries during the initial matrix assembly.
After the entry-method has been called for each nonzero matrix entry at least once, you can call compress() to reorganize the data into the standard BCRS data layout, which sets the matrix state to built
. No matrix entries may be added after this step. compress() returns a value of type Dune::CompressionStatistics, which you can inspect to tune the construction parameters _avg
and _overflowsize
.
Use of copy constructor, assignment operator and matrix vector arithmetics are not supported until the matrix is fully built.
In the following sample code, an array with 28 entries will be reserved
typedef A Dune::BCRSMatrix< B, A >::allocator_type |
export the allocator type
typedef B Dune::BCRSMatrix< B, A >::block_type |
export the type representing the components
typedef row_type::Iterator Dune::BCRSMatrix< B, A >::ColIterator |
Iterator for the entries of each row.
typedef ::Dune::CompressionStatistics<size_type> Dune::BCRSMatrix< B, A >::CompressionStatistics |
The type for the statistics object returned by compress()
typedef RealRowIterator<const row_type> Dune::BCRSMatrix< B, A >::const_iterator |
The const iterator over the matrix rows.
typedef row_type::ConstIterator Dune::BCRSMatrix< B, A >::ConstColIterator |
Const iterator to the entries of a row.
typedef RealRowIterator<const row_type> Dune::BCRSMatrix< B, A >::ConstIterator |
typedef ConstIterator Dune::BCRSMatrix< B, A >::ConstRowIterator |
rename the const row iterator for easier access
typedef B::field_type Dune::BCRSMatrix< B, A >::field_type |
export the type representing the field
typedef RealRowIterator<row_type> Dune::BCRSMatrix< B, A >::iterator |
The iterator over the (mutable matrix rows.
typedef RealRowIterator<row_type> Dune::BCRSMatrix< B, A >::Iterator |
|
protected |
typedef CompressedBlockVectorWindow<B,A> Dune::BCRSMatrix< B, A >::row_type |
implement row_type with compressed vector
typedef Iterator Dune::BCRSMatrix< B, A >::RowIterator |
rename the iterators for easier access
typedef A::size_type Dune::BCRSMatrix< B, A >::size_type |
The type for the index access and the size.
|
inline |
an empty matrix
|
inline |
matrix with known number of nonzeroes
References Dune::BCRSMatrix< B, A >::allocate().
|
inline |
matrix with unknown number of nonzeroes
References Dune::BCRSMatrix< B, A >::allocate().
|
inline |
construct matrix with a known average number of entries per row
Constructs a matrix in implicit buildmode.
_n | number of rows of the matrix |
_m | number of columns of the matrix |
_avg | expected average number of entries per row |
_overflowsize | fraction of _n*_avg which is expected to be needed for elements that exceed _avg entries per row. |
References Dune::BCRSMatrix< B, A >::implicit, and Dune::BCRSMatrix< B, A >::implicit_allocate().
|
inline |
copy constructor
Does a deep copy as expected.
References Dune::BCRSMatrix< B, A >::allocate(), Dune::BCRSMatrix< B, A >::built, Dune::BCRSMatrix< B, A >::copyWindowStructure(), Dune::CompressedBlockVectorWindow< B, A >::getsize(), Dune::BCRSMatrix< B, A >::j, Dune::BCRSMatrix< B, A >::m, Dune::BCRSMatrix< B, A >::n, Dune::BCRSMatrix< B, A >::nnz, Dune::BCRSMatrix< B, A >::notAllocated, Dune::BCRSMatrix< B, A >::r, and Dune::BCRSMatrix< B, A >::ready.
|
inline |
destructor
References Dune::BCRSMatrix< B, A >::deallocate().
|
inline |
add index (row,col) to the matrix
This method can only be used when building the BCRSMatrix in random mode.
addindex adds a new column entry to the row. If this column entry already exists, nothing is done.
Don't call addindex after the setup phase is finished (after endindices is called).
References Dune::BCRSMatrix< B, A >::build_mode, Dune::BCRSMatrix< B, A >::building, Dune::BCRSMatrix< B, A >::built, col, Dune::BCRSMatrix< B, A >::end(), Dune::CompressedBlockVectorWindow< B, A >::getindexptr(), Dune::CompressedBlockVectorWindow< B, A >::getsize(), Dune::BCRSMatrix< B, A >::m, Dune::BCRSMatrix< B, A >::notAllocated, Dune::BCRSMatrix< B, A >::r, Dune::BCRSMatrix< B, A >::random, Dune::BCRSMatrix< B, A >::ready, and row.
Referenced by Dune::BTDMatrix< B, A >::BTDMatrix(), and test_IO().
|
inlineprotected |
Allocate memory for the matrix structure.
Sets the number of rows and columns of the matrix and allocates the memory needed for the storage of the matrix entries.
row | The number of rows the matrix should contain. |
columns | the number of columns the matrix should contain. |
allocationSize_ | The number of nonzero entries the matrix should hold (if omitted defaults to 0). |
allocateRow | Whether we have to allocate the row pointers, too. (Defaults to true) |
References Dune::BCRSMatrix< B, A >::allocateData(), and Dune::BCRSMatrix< B, A >::building.
Referenced by Dune::BCRSMatrix< B, A >::BCRSMatrix(), Dune::BCRSMatrix< B, A >::endrowsizes(), Dune::BCRSMatrix< B, A >::implicit_allocate(), Dune::BCRSMatrix< B, A >::operator=(), and Dune::BCRSMatrix< B, A >::setSize().
|
inlineprotected |
|
inline |
Add the scaled entries of another matrix to this one.
Matrix axpy operation: *this += alpha * b
alpha | Scaling factor. |
b | The matrix to add to this one. Its sparsity pattern has to be subset of the sparsity pattern of this matrix. |
References Dune::BCRSMatrix< B, A >::begin(), Dune::BCRSMatrix< B, A >::built, Dune::BCRSMatrix< B, A >::end(), Dune::BCRSMatrix< B, A >::j, Dune::BCRSMatrix< B, A >::M(), Dune::BCRSMatrix< B, A >::N(), and Dune::BCRSMatrix< B, A >::ready.
|
inline |
References Dune::BCRSMatrix< B, A >::r.
|
inline |
References Dune::BCRSMatrix< B, A >::r.
|
inline |
References Dune::BCRSMatrix< B, A >::n, and Dune::BCRSMatrix< B, A >::r.
|
inline |
References Dune::BCRSMatrix< B, A >::n, and Dune::BCRSMatrix< B, A >::r.
|
inline |
Get iterator to first row.
References Dune::BCRSMatrix< B, A >::r.
Referenced by Dune::SeqOverlappingSchwarzAssemblerHelper< S< BCRSMatrix< FieldMatrix< T, m, n >, A > >, true >::assembleLocalProblems(), Dune::BCRSMatrix< B, A >::axpy(), Dune::MatrixDimension< BCRSMatrix< B, TA > >::coldim(), Dune::BCRSMatrix< B, A >::compress(), Dune::BCRSMatrix< B, A >::copyWindowStructure(), Dune::BCRSMatrix< B, A >::endindices(), Dune::BCRSMatrix< B, A >::endrowsizes(), Dune::BCRSMatrix< B, A >::entry(), Dune::BCRSMatrix< B, A >::frobenius_norm2(), Dune::BCRSMatrix< B, A >::infinity_norm(), Dune::BCRSMatrix< B, A >::infinity_norm_real(), Dune::BCRSMatrix< B, A >::mmhv(), Dune::BCRSMatrix< B, A >::mmtv(), Dune::BCRSMatrix< B, A >::mmv(), Dune::BCRSMatrix< B, A >::mv(), Dune::BCRSMatrix< B, A >::operator*=(), Dune::BCRSMatrix< B, A >::operator+=(), Dune::BCRSMatrix< B, A >::operator-=(), Dune::BCRSMatrix< B, A >::operator/=(), test_IO(), test_Iter(), Dune::BCRSMatrix< B, A >::umhv(), Dune::BCRSMatrix< B, A >::umtv(), Dune::BCRSMatrix< B, A >::umv(), Dune::BCRSMatrix< B, A >::usmhv(), Dune::BCRSMatrix< B, A >::usmtv(), and Dune::BCRSMatrix< B, A >::usmv().
|
inline |
Get const iterator to first row.
References Dune::BCRSMatrix< B, A >::r.
|
inline |
The currently selected build mode of the matrix.
References Dune::BCRSMatrix< B, A >::build_mode.
|
inline |
The current build stage of the matrix.
References Dune::BCRSMatrix< B, A >::ready.
|
inline |
Finishes the buildstage in implicit mode.
Performs compression of index and data arrays with linear complexity in the number of nonzeroes.
After calling this method, the matrix is in the built state and no more entries can be added.
References Dune::BCRSMatrix< B, A >::a, Dune::BCRSMatrix< B, A >::allocationSize, Dune::CompressionStatistics< size_type >::avg, Dune::BCRSMatrix< B, A >::begin(), Dune::BCRSMatrix< B, A >::build_mode, Dune::BCRSMatrix< B, A >::building, Dune::BCRSMatrix< B, A >::built, Dune::CompressedBlockVectorWindow< B, A >::getindexptr(), Dune::CompressedBlockVectorWindow< B, A >::getsize(), Dune::BCRSMatrix< B, A >::implicit, Dune::BCRSMatrix< B, A >::j, Dune::CompressionStatistics< size_type >::maximum, Dune::CompressionStatistics< size_type >::mem_ratio, Dune::BCRSMatrix< B, A >::n, Dune::BCRSMatrix< B, A >::nnz, Dune::BCRSMatrix< B, A >::notAllocated, Dune::BCRSMatrix< B, A >::overflow, Dune::CompressionStatistics< size_type >::overflow_total, Dune::BCRSMatrix< B, A >::r, Dune::BCRSMatrix< B, A >::ready, Dune::CompressedBlockVectorWindow< B, A >::setindexptr(), Dune::CompressedBlockVectorWindow< B, A >::setptr(), and Dune::CompressedBlockVectorWindow< B, A >::setsize().
|
inlineprotected |
Copy the window structure from another matrix.
References Dune::BCRSMatrix< B, A >::begin(), Dune::BCRSMatrix< B, A >::built, Dune::BCRSMatrix< B, A >::n, Dune::BCRSMatrix< B, A >::r, Dune::BCRSMatrix< B, A >::row_wise, and Dune::BCRSMatrix< B, A >::setWindowPointers().
Referenced by Dune::BCRSMatrix< B, A >::BCRSMatrix(), and Dune::BCRSMatrix< B, A >::operator=().
|
inline |
get initial create iterator
References Dune::BCRSMatrix< B, A >::CreateIterator.
Referenced by test_IO(), and test_Iter().
|
inline |
get create iterator pointing to one after the last block
References Dune::BCRSMatrix< B, A >::CreateIterator, and Dune::BCRSMatrix< B, A >::n.
Referenced by test_IO(), and test_Iter().
|
inlineprotected |
deallocate memory of the matrix.
deallocateRows | Whether we have to deallocate the row pointers, too. If false they will not be touched. (Defaults to true). |
References col, Dune::BCRSMatrix< B, A >::n, Dune::BCRSMatrix< B, A >::notAllocated, and Dune::CompressedBlockVectorWindow< B, A >::set().
Referenced by Dune::BCRSMatrix< B, A >::operator=(), Dune::BCRSMatrix< B, A >::setSize(), and Dune::BCRSMatrix< B, A >::~BCRSMatrix().
|
inline |
Get iterator to one beyond last row.
References Dune::BCRSMatrix< B, A >::n, and Dune::BCRSMatrix< B, A >::r.
Referenced by Dune::BCRSMatrix< B, A >::addindex(), Dune::BCRSMatrix< B, A >::axpy(), Dune::MatrixDimension< BCRSMatrix< B, TA > >::coldim(), Dune::BCRSMatrix< B, A >::endindices(), Dune::BCRSMatrix< B, A >::entry(), Dune::BCRSMatrix< B, A >::exists(), Dune::BCRSMatrix< B, A >::frobenius_norm2(), Dune::BCRSMatrix< B, A >::infinity_norm(), Dune::BCRSMatrix< B, A >::infinity_norm_real(), Dune::BCRSMatrix< B, A >::mmhv(), Dune::BCRSMatrix< B, A >::mmtv(), Dune::BCRSMatrix< B, A >::mmv(), Dune::BCRSMatrix< B, A >::mv(), Dune::BCRSMatrix< B, A >::operator*=(), Dune::BCRSMatrix< B, A >::operator+=(), Dune::BCRSMatrix< B, A >::operator-=(), Dune::BCRSMatrix< B, A >::operator/=(), Dune::DynamicMatrixSubdomainSolver< BCRSMatrix< FieldMatrix< K, n, n >, Al >, X, Y >::setSubMatrix(), test_IO(), test_Iter(), Dune::BCRSMatrix< B, A >::umhv(), Dune::BCRSMatrix< B, A >::umtv(), Dune::BCRSMatrix< B, A >::umv(), Dune::BCRSMatrix< B, A >::usmhv(), Dune::BCRSMatrix< B, A >::usmtv(), and Dune::BCRSMatrix< B, A >::usmv().
|
inline |
Get const iterator to one beyond last row.
References Dune::BCRSMatrix< B, A >::n, and Dune::BCRSMatrix< B, A >::r.
|
inline |
indicate that all indices are defined, check consistency
References Dune::BCRSMatrix< B, A >::allocateData(), Dune::BCRSMatrix< B, A >::begin(), Dune::BCRSMatrix< B, A >::build_mode, Dune::BCRSMatrix< B, A >::building, Dune::BCRSMatrix< B, A >::built, Dune::BCRSMatrix< B, A >::end(), Dune::BCRSMatrix< B, A >::j, Dune::BCRSMatrix< B, A >::m, Dune::BCRSMatrix< B, A >::notAllocated, Dune::BCRSMatrix< B, A >::r, Dune::BCRSMatrix< B, A >::random, Dune::BCRSMatrix< B, A >::ready, Dune::BCRSMatrix< B, A >::setDataPointers(), and Dune::CompressedBlockVectorWindow< B, A >::setsize().
Referenced by Dune::BDMatrix< B, A >::BDMatrix(), Dune::BTDMatrix< B, A >::BTDMatrix(), and test_IO().
|
inline |
indicate that size of all rows is defined
References Dune::BCRSMatrix< B, A >::allocate(), Dune::BCRSMatrix< B, A >::begin(), Dune::BCRSMatrix< B, A >::build_mode, Dune::BCRSMatrix< B, A >::building, Dune::CompressedBlockVectorWindow< B, A >::getsize(), Dune::BCRSMatrix< B, A >::j, Dune::BCRSMatrix< B, A >::m, Dune::BCRSMatrix< B, A >::n, Dune::BCRSMatrix< B, A >::nnz, Dune::BCRSMatrix< B, A >::r, Dune::BCRSMatrix< B, A >::random, Dune::BCRSMatrix< B, A >::ready, Dune::BCRSMatrix< B, A >::rowSizesBuilt, and Dune::BCRSMatrix< B, A >::setColumnPointers().
Referenced by Dune::BDMatrix< B, A >::BDMatrix(), Dune::BTDMatrix< B, A >::BTDMatrix(), and test_IO().
|
inline |
Returns reference to entry (row,col) of the matrix.
This method can only be used when the matrix is in implicit building mode.
A reference to entry (row, col) of the matrix is returned. If entry (row, col) is accessed for the first time, it is created on the fly.
This method can only be used while building the matrix, after compression operator[] gives a much better performance.
References Dune::BCRSMatrix< B, A >::avg, Dune::BCRSMatrix< B, A >::begin(), Dune::BCRSMatrix< B, A >::build_mode, Dune::BCRSMatrix< B, A >::building, Dune::BCRSMatrix< B, A >::built, col, Dune::BCRSMatrix< B, A >::end(), Dune::CompressedBlockVectorWindow< B, A >::getindexptr(), Dune::CompressedBlockVectorWindow< B, A >::getptr(), Dune::CompressedBlockVectorWindow< B, A >::getsize(), Dune::BCRSMatrix< B, A >::implicit, Dune::BCRSMatrix< B, A >::m, Dune::BCRSMatrix< B, A >::n, Dune::BCRSMatrix< B, A >::notAllocated, Dune::BCRSMatrix< B, A >::overflow, Dune::BCRSMatrix< B, A >::r, Dune::BCRSMatrix< B, A >::ready, row, and Dune::CompressedBlockVectorWindow< B, A >::setsize().
|
inline |
return true if (i,j) is in pattern
References Dune::BCRSMatrix< B, A >::end(), Dune::BCRSMatrix< B, A >::m, Dune::BCRSMatrix< B, A >::n, and Dune::BCRSMatrix< B, A >::r.
|
inline |
frobenius norm: sqrt(sum over squared values of entries)
References Dune::BCRSMatrix< B, A >::frobenius_norm2().
|
inline |
square of frobenius norm, need for block recursion
References Dune::BCRSMatrix< B, A >::begin(), Dune::BCRSMatrix< B, A >::built, Dune::BCRSMatrix< B, A >::end(), Dune::BCRSMatrix< B, A >::j, and Dune::BCRSMatrix< B, A >::ready.
Referenced by Dune::BCRSMatrix< B, A >::frobenius_norm().
|
inline |
get current number of indices in row i
References Dune::CompressedBlockVectorWindow< B, A >::getsize(), Dune::BCRSMatrix< B, A >::n, and Dune::BCRSMatrix< B, A >::r.
|
inlineprotected |
organizes allocation implicit mode calculates correct array size to be allocated and sets the the window pointers to their correct positions for insertion. internally uses allocate() for the real allocation.
References Dune::BCRSMatrix< B, A >::allocate(), Dune::BCRSMatrix< B, A >::avg, Dune::BCRSMatrix< B, A >::building, Dune::BCRSMatrix< B, A >::implicit, Dune::BCRSMatrix< B, A >::n, Dune::BCRSMatrix< B, A >::notAllocated, and Dune::CompressedBlockVectorWindow< B, A >::set().
Referenced by Dune::BCRSMatrix< B, A >::BCRSMatrix(), and Dune::BCRSMatrix< B, A >::setSize().
|
inline |
increment size of row i by s (1 by default)
References Dune::BCRSMatrix< B, A >::build_mode, Dune::BCRSMatrix< B, A >::building, Dune::BCRSMatrix< B, A >::r, Dune::BCRSMatrix< B, A >::random, Dune::BCRSMatrix< B, A >::ready, and Dune::CompressedBlockVectorWindow< B, A >::setsize().
|
inline |
infinity norm (row sum norm, how to generalize for blocks?)
References Dune::BCRSMatrix< B, A >::begin(), Dune::BCRSMatrix< B, A >::built, Dune::BCRSMatrix< B, A >::end(), Dune::BCRSMatrix< B, A >::j, and Dune::BCRSMatrix< B, A >::ready.
|
inline |
simplified infinity norm (uses Manhattan norm for complex values)
References Dune::BCRSMatrix< B, A >::begin(), Dune::BCRSMatrix< B, A >::built, Dune::BCRSMatrix< B, A >::end(), Dune::BCRSMatrix< B, A >::j, and Dune::BCRSMatrix< B, A >::ready.
|
inline |
number of columns (counted in blocks)
References Dune::BCRSMatrix< B, A >::m.
Referenced by Dune::BCRSMatrix< B, A >::axpy(), Dune::MatrixDimension< BCRSMatrix< B, TA > >::coldim(), Dune::MatrixDimension< BCRSMatrix< FieldMatrix< B, n, m >, TA > >::coldim(), Dune::BCRSMatrix< B, A >::mmhv(), Dune::BCRSMatrix< B, A >::mmtv(), Dune::BCRSMatrix< B, A >::mmv(), Dune::BCRSMatrix< B, A >::mtv(), Dune::BCRSMatrix< B, A >::mv(), Dune::BCRSMatrix< B, A >::operator+=(), Dune::BCRSMatrix< B, A >::operator-=(), Dune::SuperLUMatrix< BCRSMatrix< FieldMatrix< B, n, m >, TA > >::setMatrix(), Dune::BCRSMatrix< B, A >::umhv(), Dune::BCRSMatrix< B, A >::umtv(), Dune::BCRSMatrix< B, A >::umv(), Dune::BCRSMatrix< B, A >::usmhv(), Dune::BCRSMatrix< B, A >::usmtv(), and Dune::BCRSMatrix< B, A >::usmv().
|
inline |
|
inline |
|
inline |
|
inline |
|
inline |
|
inline |
number of rows (counted in blocks)
References Dune::BCRSMatrix< B, A >::n.
Referenced by Dune::BCRSMatrix< B, A >::axpy(), Dune::MatrixDimension< BCRSMatrix< B, TA > >::coldim(), Dune::BDMatrix< B, A >::invert(), Dune::BCRSMatrix< B, A >::mmhv(), Dune::BCRSMatrix< B, A >::mmtv(), Dune::BCRSMatrix< B, A >::mmv(), Dune::BCRSMatrix< B, A >::mtv(), Dune::BCRSMatrix< B, A >::mv(), Dune::BCRSMatrix< B, A >::operator+=(), Dune::BCRSMatrix< B, A >::operator-=(), Dune::MatrixDimension< BCRSMatrix< B, TA > >::rowdim(), Dune::MatrixDimension< BCRSMatrix< FieldMatrix< B, n, m >, TA > >::rowdim(), Dune::SuperLUMatrix< BCRSMatrix< FieldMatrix< B, n, m >, TA > >::setMatrix(), Dune::BTDMatrix< B, A >::solve(), Dune::BCRSMatrix< B, A >::umhv(), Dune::BCRSMatrix< B, A >::umtv(), Dune::BCRSMatrix< B, A >::umv(), Dune::BCRSMatrix< B, A >::usmhv(), Dune::BCRSMatrix< B, A >::usmtv(), and Dune::BCRSMatrix< B, A >::usmv().
|
inline |
number of blocks that are stored (the number of blocks that possibly are nonzero)
References Dune::BCRSMatrix< B, A >::nnz.
|
inline |
vector space multiplication with scalar
References Dune::BCRSMatrix< B, A >::a, Dune::BCRSMatrix< B, A >::begin(), Dune::BCRSMatrix< B, A >::built, Dune::BCRSMatrix< B, A >::end(), Dune::BCRSMatrix< B, A >::j, Dune::BCRSMatrix< B, A >::nnz, and Dune::BCRSMatrix< B, A >::ready.
|
inline |
Add the entries of another matrix to this one.
b | The matrix to add to this one. Its sparsity pattern has to be subset of the sparsity pattern of this matrix. |
References Dune::BCRSMatrix< B, A >::begin(), Dune::BCRSMatrix< B, A >::built, Dune::BCRSMatrix< B, A >::end(), Dune::BCRSMatrix< B, A >::j, Dune::BCRSMatrix< B, A >::M(), Dune::BCRSMatrix< B, A >::N(), and Dune::BCRSMatrix< B, A >::ready.
|
inline |
Substract the entries of another matrix to this one.
b | The matrix to add to this one. Its sparsity pattern has to be subset of the sparsity pattern of this matrix. |
References Dune::BCRSMatrix< B, A >::begin(), Dune::BCRSMatrix< B, A >::built, Dune::BCRSMatrix< B, A >::end(), Dune::BCRSMatrix< B, A >::j, Dune::BCRSMatrix< B, A >::M(), Dune::BCRSMatrix< B, A >::N(), and Dune::BCRSMatrix< B, A >::ready.
|
inline |
vector space division by scalar
References Dune::BCRSMatrix< B, A >::a, Dune::BCRSMatrix< B, A >::begin(), Dune::BCRSMatrix< B, A >::built, Dune::BCRSMatrix< B, A >::end(), Dune::BCRSMatrix< B, A >::j, Dune::BCRSMatrix< B, A >::nnz, and Dune::BCRSMatrix< B, A >::ready.
|
inline |
assignment
Frees and reallocates space. Both sparsity pattern and values are set from Mat.
References Dune::BCRSMatrix< B, A >::allocate(), Dune::BCRSMatrix< B, A >::built, Dune::BCRSMatrix< B, A >::copyWindowStructure(), Dune::BCRSMatrix< B, A >::deallocate(), Dune::CompressedBlockVectorWindow< B, A >::getsize(), Dune::BCRSMatrix< B, A >::j, Dune::BCRSMatrix< B, A >::m, Dune::BCRSMatrix< B, A >::n, Dune::BCRSMatrix< B, A >::nnz, Dune::BCRSMatrix< B, A >::notAllocated, Dune::BCRSMatrix< B, A >::r, Dune::BCRSMatrix< B, A >::ready, and Dune::BCRSMatrix< B, A >::rowAllocator_.
Referenced by Dune::BDMatrix< B, A >::operator=(), and Dune::BTDMatrix< B, A >::operator=().
|
inline |
Assignment from a scalar.
References Dune::BCRSMatrix< B, A >::built, Dune::BCRSMatrix< B, A >::n, Dune::BCRSMatrix< B, A >::notAllocated, Dune::BCRSMatrix< B, A >::r, and Dune::BCRSMatrix< B, A >::ready.
|
inline |
random access to the rows
References Dune::BCRSMatrix< B, A >::build_mode, Dune::BCRSMatrix< B, A >::built, Dune::BCRSMatrix< B, A >::implicit, Dune::BCRSMatrix< B, A >::n, Dune::BCRSMatrix< B, A >::r, and Dune::BCRSMatrix< B, A >::ready.
|
inline |
same for read only access
References Dune::BCRSMatrix< B, A >::build_mode, Dune::BCRSMatrix< B, A >::built, Dune::BCRSMatrix< B, A >::implicit, Dune::BCRSMatrix< B, A >::n, Dune::BCRSMatrix< B, A >::r, and Dune::BCRSMatrix< B, A >::ready.
|
inline |
Sets the build mode of the matrix.
bm | The build mode to use. |
References Dune::BCRSMatrix< B, A >::build_mode, Dune::BCRSMatrix< B, A >::building, Dune::BCRSMatrix< B, A >::notAllocated, Dune::BCRSMatrix< B, A >::random, Dune::BCRSMatrix< B, A >::ready, Dune::BCRSMatrix< B, A >::row_wise, and Dune::BCRSMatrix< B, A >::unknown.
|
inlineprotected |
Copy row sizes from iterator range starting at row and set column index pointers for all rows.
This method does not modify the data pointers, as those are set only after building the pattern (to allow for a delayed allocation).
References Dune::BCRSMatrix< B, A >::n, row, Dune::CompressedBlockVectorWindow< B, A >::set(), Dune::CompressedBlockVectorWindow< B, A >::setindexptr(), and Dune::CompressedBlockVectorWindow< B, A >::setsize().
Referenced by Dune::BCRSMatrix< B, A >::endrowsizes().
|
inlineprotected |
Set data pointers for all rows.
This method assumes that column pointers and row sizes have been correctly set up by a prior call to setColumnPointers().
References Dune::BCRSMatrix< B, A >::a, Dune::CompressedBlockVectorWindow< B, A >::getsize(), Dune::BCRSMatrix< B, A >::n, Dune::CompressedBlockVectorWindow< B, A >::set(), and Dune::CompressedBlockVectorWindow< B, A >::setptr().
Referenced by Dune::BCRSMatrix< B, A >::endindices(), and Dune::BCRSMatrix< B, A >::CreateIterator::operator++().
|
inline |
Set parameters needed for creation in implicit build mode.
Use this method before setSize() to define storage behaviour of a matrix in implicit build mode
_avg | expected average number of entries per row |
_overflowsize | fraction of _n*_avg which is expected to be needed for elements that exceed _avg entries per row. |
References Dune::BCRSMatrix< B, A >::avg, Dune::BCRSMatrix< B, A >::notAllocated, Dune::BCRSMatrix< B, A >::overflowsize, and Dune::BCRSMatrix< B, A >::ready.
|
inline |
Set all column indices for row from the given iterator range.
The iterator range has to be of the same length as the previously set row size. The entries in the iterator range do not have to be in any particular order, but must not contain duplicate values.
Calling this method overwrites any previously set column indices!
References Dune::CompressedBlockVectorWindow< B, A >::getindexptr(), Dune::BCRSMatrix< B, A >::r, row, and Dune::compressed_base_array_unmanaged< B, A >::size().
|
inline |
set number of indices in row i to s
References Dune::BCRSMatrix< B, A >::build_mode, Dune::BCRSMatrix< B, A >::building, Dune::BCRSMatrix< B, A >::r, Dune::BCRSMatrix< B, A >::random, Dune::BCRSMatrix< B, A >::ready, and Dune::CompressedBlockVectorWindow< B, A >::setsize().
Referenced by Dune::BTDMatrix< B, A >::BTDMatrix(), and test_IO().
|
inline |
Set the size of the matrix.
Sets the number of rows and columns of the matrix and allocates the memory needed for the storage of the matrix entries.
rows | The number of rows the matrix should contain. |
columns | the number of columns the matrix should contain. |
nnz | The number of nonzero entries the matrix should hold (if omitted defaults to 0). Must be omitted in implicit mode. |
References Dune::BCRSMatrix< B, A >::allocate(), Dune::BCRSMatrix< B, A >::build_mode, Dune::BCRSMatrix< B, A >::deallocate(), Dune::BCRSMatrix< B, A >::implicit, Dune::BCRSMatrix< B, A >::implicit_allocate(), and Dune::BCRSMatrix< B, A >::nnz.
|
inlineprotected |
References Dune::BCRSMatrix< B, A >::n, row, and Dune::CompressedBlockVectorWindow< B, A >::set().
Referenced by Dune::BCRSMatrix< B, A >::copyWindowStructure().
|
inline |
|
inline |
|
inline |
|
inline |
|
inline |
|
inline |
|
protected |
Referenced by Dune::BCRSMatrix< B, A >::allocateData(), Dune::MatrixDimension< BCRSMatrix< B, TA > >::coldim(), Dune::BCRSMatrix< B, A >::compress(), Dune::BCRSMatrix< B, A >::operator*=(), Dune::BCRSMatrix< B, A >::CreateIterator::operator++(), Dune::BCRSMatrix< B, A >::operator/=(), and Dune::BCRSMatrix< B, A >::setDataPointers().
|
protected |
Referenced by Dune::BCRSMatrix< B, A >::compress().
|
protected |
Referenced by Dune::BCRSMatrix< B, A >::CreateIterator::operator++().
|
protected |
|
protected |
Referenced by Dune::BCRSMatrix< B, A >::addindex(), Dune::BCRSMatrix< B, A >::buildMode(), Dune::BCRSMatrix< B, A >::compress(), Dune::BCRSMatrix< B, A >::CreateIterator::CreateIterator(), Dune::BCRSMatrix< B, A >::endindices(), Dune::BCRSMatrix< B, A >::endrowsizes(), Dune::BCRSMatrix< B, A >::entry(), Dune::BCRSMatrix< B, A >::incrementrowsize(), Dune::BCRSMatrix< B, A >::operator[](), Dune::BCRSMatrix< B, A >::setBuildMode(), Dune::BCRSMatrix< B, A >::setrowsize(), and Dune::BCRSMatrix< B, A >::setSize().
|
protected |
Referenced by Dune::BCRSMatrix< B, A >::axpy(), Dune::BCRSMatrix< B, A >::BCRSMatrix(), Dune::MatrixDimension< BCRSMatrix< B, TA > >::coldim(), Dune::BCRSMatrix< B, A >::compress(), Dune::BCRSMatrix< B, A >::endindices(), Dune::BCRSMatrix< B, A >::endrowsizes(), Dune::BCRSMatrix< B, A >::frobenius_norm2(), Dune::BCRSMatrix< B, A >::infinity_norm(), Dune::BCRSMatrix< B, A >::infinity_norm_real(), Dune::BCRSMatrix< B, A >::mmhv(), Dune::BCRSMatrix< B, A >::mmtv(), Dune::BCRSMatrix< B, A >::mmv(), Dune::BCRSMatrix< B, A >::mv(), Dune::BCRSMatrix< B, A >::operator*=(), Dune::BCRSMatrix< B, A >::CreateIterator::operator++(), Dune::BCRSMatrix< B, A >::operator+=(), Dune::BCRSMatrix< B, A >::operator-=(), Dune::BCRSMatrix< B, A >::operator/=(), Dune::BCRSMatrix< B, A >::operator=(), Dune::BCRSMatrix< B, A >::umhv(), Dune::BCRSMatrix< B, A >::umtv(), Dune::BCRSMatrix< B, A >::umv(), Dune::BCRSMatrix< B, A >::usmhv(), Dune::BCRSMatrix< B, A >::usmtv(), and Dune::BCRSMatrix< B, A >::usmv().
|
protected |
Referenced by Dune::BCRSMatrix< B, A >::addindex(), Dune::BCRSMatrix< B, A >::BCRSMatrix(), Dune::BCRSMatrix< B, A >::endindices(), Dune::BCRSMatrix< B, A >::endrowsizes(), Dune::BCRSMatrix< B, A >::entry(), Dune::BCRSMatrix< B, A >::exists(), Dune::BCRSMatrix< B, A >::M(), and Dune::BCRSMatrix< B, A >::operator=().
|
protected |
Referenced by Dune::BCRSMatrix< B, A >::BCRSMatrix(), Dune::BCRSMatrix< B, A >::beforeEnd(), Dune::BCRSMatrix< B, A >::compress(), Dune::BCRSMatrix< B, A >::copyWindowStructure(), Dune::BCRSMatrix< B, A >::createend(), Dune::BCRSMatrix< B, A >::deallocate(), Dune::BCRSMatrix< B, A >::end(), Dune::BCRSMatrix< B, A >::endrowsizes(), Dune::BCRSMatrix< B, A >::entry(), Dune::BCRSMatrix< B, A >::exists(), Dune::BCRSMatrix< B, A >::getrowsize(), Dune::BCRSMatrix< B, A >::implicit_allocate(), Dune::BCRSMatrix< B, A >::N(), Dune::BCRSMatrix< B, A >::CreateIterator::operator++(), Dune::BCRSMatrix< B, A >::operator=(), Dune::BCRSMatrix< B, A >::operator[](), Dune::BCRSMatrix< B, A >::setColumnPointers(), Dune::BCRSMatrix< B, A >::setDataPointers(), and Dune::BCRSMatrix< B, A >::setWindowPointers().
|
protected |
Referenced by Dune::BCRSMatrix< B, A >::BCRSMatrix(), Dune::MatrixDimension< BCRSMatrix< B, TA > >::coldim(), Dune::BCRSMatrix< B, A >::compress(), Dune::BCRSMatrix< B, A >::endrowsizes(), Dune::BCRSMatrix< B, A >::nonzeroes(), Dune::BCRSMatrix< B, A >::operator*=(), Dune::BCRSMatrix< B, A >::CreateIterator::operator++(), Dune::BCRSMatrix< B, A >::operator/=(), Dune::BCRSMatrix< B, A >::operator=(), and Dune::BCRSMatrix< B, A >::setSize().
|
protected |
Referenced by Dune::BCRSMatrix< B, A >::compress(), and Dune::BCRSMatrix< B, A >::entry().
|
protected |
Referenced by Dune::BCRSMatrix< B, A >::setImplicitBuildModeParameters().
|
protected |
Referenced by Dune::BCRSMatrix< B, A >::addindex(), Dune::BCRSMatrix< B, A >::BCRSMatrix(), Dune::BCRSMatrix< B, A >::beforeBegin(), Dune::BCRSMatrix< B, A >::beforeEnd(), Dune::BCRSMatrix< B, A >::begin(), Dune::MatrixDimension< BCRSMatrix< B, TA > >::coldim(), Dune::BCRSMatrix< B, A >::compress(), Dune::BCRSMatrix< B, A >::copyWindowStructure(), Dune::BCRSMatrix< B, A >::end(), Dune::BCRSMatrix< B, A >::endindices(), Dune::BCRSMatrix< B, A >::endrowsizes(), Dune::BCRSMatrix< B, A >::entry(), Dune::BCRSMatrix< B, A >::exists(), Dune::BCRSMatrix< B, A >::getrowsize(), Dune::BCRSMatrix< B, A >::incrementrowsize(), Dune::BCRSMatrix< B, A >::CreateIterator::operator++(), Dune::BCRSMatrix< B, A >::operator=(), Dune::BCRSMatrix< B, A >::operator[](), Dune::MatrixDimension< BCRSMatrix< B, TA > >::rowdim(), Dune::BCRSMatrix< B, A >::setIndices(), and Dune::BCRSMatrix< B, A >::setrowsize().
|
protected |
Referenced by Dune::BCRSMatrix< B, A >::addindex(), Dune::BCRSMatrix< B, A >::axpy(), Dune::BCRSMatrix< B, A >::BCRSMatrix(), Dune::BCRSMatrix< B, A >::buildStage(), Dune::BCRSMatrix< B, A >::compress(), Dune::BCRSMatrix< B, A >::CreateIterator::CreateIterator(), Dune::BCRSMatrix< B, A >::endindices(), Dune::BCRSMatrix< B, A >::endrowsizes(), Dune::BCRSMatrix< B, A >::entry(), Dune::BCRSMatrix< B, A >::frobenius_norm2(), Dune::BCRSMatrix< B, A >::incrementrowsize(), Dune::BCRSMatrix< B, A >::infinity_norm(), Dune::BCRSMatrix< B, A >::infinity_norm_real(), Dune::BCRSMatrix< B, A >::mmhv(), Dune::BCRSMatrix< B, A >::mmv(), Dune::BCRSMatrix< B, A >::mtv(), Dune::BCRSMatrix< B, A >::mv(), Dune::BCRSMatrix< B, A >::operator*=(), Dune::BCRSMatrix< B, A >::CreateIterator::operator++(), Dune::BCRSMatrix< B, A >::operator+=(), Dune::BCRSMatrix< B, A >::operator-=(), Dune::BCRSMatrix< B, A >::operator/=(), Dune::BCRSMatrix< B, A >::operator=(), Dune::BCRSMatrix< B, A >::operator[](), Dune::BCRSMatrix< B, A >::setBuildMode(), Dune::BCRSMatrix< B, A >::setImplicitBuildModeParameters(), Dune::BCRSMatrix< B, A >::setrowsize(), Dune::BCRSMatrix< B, A >::umhv(), Dune::BCRSMatrix< B, A >::umtv(), Dune::BCRSMatrix< B, A >::umv(), Dune::BCRSMatrix< B, A >::usmhv(), Dune::BCRSMatrix< B, A >::usmtv(), and Dune::BCRSMatrix< B, A >::usmv().
|
protected |
Referenced by Dune::BCRSMatrix< B, A >::operator=().
|
protected |
Referenced by Dune::BCRSMatrix< B, A >::CreateIterator::operator++().