dune-istl  2.2.1
Classes | Public Types | Public Member Functions | Friends | List of all members
Dune::BCRSMatrix< B, A > Class Template Reference

A sparse block matrix with compressed row storage. More...

#include <dune/istl/bcrsmatrix.hh>

Inheritance diagram for Dune::BCRSMatrix< B, A >:
Inheritance graph

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.
class  RealRowIterator
 Iterator access to matrix rows More...

Public Types

enum  { blocklevel = B::blocklevel+1 }
 increment block level counter More...
enum  BuildMode { row_wise, random, unknown }
 we support two modes More...
typedef B::field_type field_type
 export the type representing the field
typedef B block_type
 export the type representing the components
typedef A allocator_type
 export the allocator type
typedef
CompressedBlockVectorWindow< B,
A
row_type
 implement row_type with compressed vector
typedef A::size_type size_type
 The type for the index access and the size.
typedef RealRowIterator< row_typeiterator
 The iterator over the (mutable matrix rows.
typedef RealRowIterator< row_typeIterator
typedef Iterator RowIterator
 rename the iterators for easier access
typedef row_type::Iterator ColIterator
 Iterator for the entries of each row.
typedef RealRowIterator< const
row_type
const_iterator
 The const iterator over the matrix rows.
typedef RealRowIterator< const
row_type
ConstIterator
typedef ConstIterator ConstRowIterator
 rename the const row iterator for easier access
typedef row_type::ConstIterator ConstColIterator
 Const iterator to the entries of a row.

Public Member Functions

row_typeoperator[] (size_type i)
 random access to the rows
const row_typeoperator[] (size_type i) const
 same for read only access
Iterator begin ()
 Get iterator to first row.
Iterator end ()
 Get iterator to one beyond last row.
Iterator beforeEnd ()
Iterator beforeBegin ()
ConstIterator begin () const
 Get const iterator to first row.
ConstIterator end () const
 Get const iterator to one beyond last row.
ConstIterator beforeEnd () const
ConstIterator beforeBegin () const
 BCRSMatrix ()
 an empty matrix
 BCRSMatrix (size_type _n, size_type _m, size_type _nnz, BuildMode bm)
 matrix with known number of nonzeroes
 BCRSMatrix (size_type _n, size_type _m, BuildMode bm)
 matrix with unknown number of nonzeroes
 BCRSMatrix (const BCRSMatrix &Mat)
 copy constructor
 ~BCRSMatrix ()
 destructor
void setBuildMode (BuildMode bm)
 Sets the build mode of the matrix.
void setSize (size_type rows, size_type columns, size_type nnz=0)
 Set the size of the matrix.
BCRSMatrixoperator= (const BCRSMatrix &Mat)
 assignment
BCRSMatrixoperator= (const field_type &k)
 Assignment from a scalar.
CreateIterator createbegin ()
 get initial create iterator
CreateIterator createend ()
 get create iterator pointing to one after the last block
void setrowsize (size_type i, size_type s)
 set number of indices in row i to s
size_type getrowsize (size_type i) const
 get current number of indices in row i
void incrementrowsize (size_type i, size_type s=1)
 increment size of row i by s (1 by default)
void endrowsizes ()
 indicate that size of all rows is defined
void addindex (size_type row, size_type col)
 add index (row,col) to the matrix
void endindices ()
 indicate that all indices are defined, check consistency
BCRSMatrixoperator*= (const field_type &k)
 vector space multiplication with scalar
BCRSMatrixoperator/= (const field_type &k)
 vector space division by scalar
BCRSMatrixoperator+= (const BCRSMatrix &b)
 Add the entries of another matrix to this one.
BCRSMatrixoperator-= (const BCRSMatrix &b)
 Substract the entries of another matrix to this one.
BCRSMatrixaxpy (field_type alpha, const BCRSMatrix &b)
 Add the scaled entries of another matrix to this one.
template<class X , class Y >
void mv (const X &x, Y &y) const
 y = A x
template<class X , class Y >
void umv (const X &x, Y &y) const
 y += A x
template<class X , class Y >
void mmv (const X &x, Y &y) const
 y -= A x
template<class X , class Y >
void usmv (const field_type &alpha, const X &x, Y &y) const
 y += alpha A x
template<class X , class Y >
void mtv (const X &x, Y &y) const
 y = A^T x
template<class X , class Y >
void umtv (const X &x, Y &y) const
 y += A^T x
template<class X , class Y >
void mmtv (const X &x, Y &y) const
 y -= A^T x
template<class X , class Y >
void usmtv (const field_type &alpha, const X &x, Y &y) const
 y += alpha A^T x
template<class X , class Y >
void umhv (const X &x, Y &y) const
 y += A^H x
template<class X , class Y >
void mmhv (const X &x, Y &y) const
 y -= A^H x
template<class X , class Y >
void usmhv (const field_type &alpha, const X &x, Y &y) const
 y += alpha A^H x
double frobenius_norm2 () const
 square of frobenius norm, need for block recursion
double frobenius_norm () const
 frobenius norm: sqrt(sum over squared values of entries)
double infinity_norm () const
 infinity norm (row sum norm, how to generalize for blocks?)
double infinity_norm_real () const
 simplified infinity norm (uses Manhattan norm for complex values)
size_type N () const
 number of rows (counted in blocks)
size_type M () const
 number of columns (counted in blocks)
size_type nonzeroes () const
 number of blocks that are stored (the number of blocks that possibly are nonzero)
bool exists (size_type i, size_type j) const
 return true if (i,j) is in pattern

Friends

struct MatrixDimension< BCRSMatrix >
class CreateIterator
 allow CreateIterator to access internal data

Detailed Description

template<class B, class A = std::allocator<B>>
class Dune::BCRSMatrix< B, A >

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:

  1. Row-wise scheme
  2. Random scheme

Error checking: no error checking is provided normally. Setting the compile time switch DUNE_ISTL_WITH_CHECKING enables error checking.

Details:

  1. Row-wise scheme

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)).

#include<dune/common/fmatrix.hh>
...
typedef FieldMatrix<double,2,2> M;
// third parameter is an optional upper bound for the number
// of nonzeros. If given the matrix will use one array for all values
// as opossed to one for each row.
BCRSMatrix<M> B(4,4,12,BCRSMatrix<M>::row_wise);
typedef BCRSMatrix<M>::CreateIterator Iter;
for(Iter row=B.createbegin(); row!=B.createend(); ++row){
// Add nonzeros for left neighbour, diagonal and right neighbour
if(row.index()>0)
row.insert(row.index()-1);
row.insert(row.index());
if(row.index()<B.N()-1)
row.insert(row.index()+1);
}
// Now the sparsity pattern is fully set up and we can add values
B[0][0]=2;
...
  1. Random scheme

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.

#include<dune/common/fmatrix.hh>
...
typedef FieldMatrix<double,2,2> M;
BCRSMatrix<M> B(4,4,BCRSMatrix<M>::random);
// initially set row size for each row
B.setrowsize(0,1);
B.setrowsize(3,4);
B.setrowsize(2,1);
B.setrowsize(1,1);
// increase row size for row 2
B.incrementrowsize(2)
// finalize row setup phase
B.endrowsizes();
// add column entries to rows
B.addindex(0,0);
B.addindex(3,1);
B.addindex(2,2);
B.addindex(1,1);
B.addindex(2,0);
B.addindex(3,2);
B.addindex(3,0);
B.addindex(3,3);
// finalize column setup phase
B.endindices();
// set entries using the random access operator
B[0][0] = 1;
B[1][1] = 2;
B[2][0] = 3;
B[2][2] = 4;
B[3][1] = 5;
B[3][2] = 6;
B[3][0] = 7;
B[3][3] = 8;

Member Typedef Documentation

template<class B, class A = std::allocator<B>>
typedef A Dune::BCRSMatrix< B, A >::allocator_type

export the allocator type

template<class B, class A = std::allocator<B>>
typedef B Dune::BCRSMatrix< B, A >::block_type

export the type representing the components

template<class B, class A = std::allocator<B>>
typedef row_type::Iterator Dune::BCRSMatrix< B, A >::ColIterator

Iterator for the entries of each row.

template<class B, class A = std::allocator<B>>
typedef RealRowIterator<const row_type> Dune::BCRSMatrix< B, A >::const_iterator

The const iterator over the matrix rows.

template<class B, class A = std::allocator<B>>
typedef row_type::ConstIterator Dune::BCRSMatrix< B, A >::ConstColIterator

Const iterator to the entries of a row.

template<class B, class A = std::allocator<B>>
typedef RealRowIterator<const row_type> Dune::BCRSMatrix< B, A >::ConstIterator
template<class B, class A = std::allocator<B>>
typedef ConstIterator Dune::BCRSMatrix< B, A >::ConstRowIterator

rename the const row iterator for easier access

template<class B, class A = std::allocator<B>>
typedef B::field_type Dune::BCRSMatrix< B, A >::field_type

export the type representing the field

template<class B, class A = std::allocator<B>>
typedef RealRowIterator<row_type> Dune::BCRSMatrix< B, A >::iterator

The iterator over the (mutable matrix rows.

template<class B, class A = std::allocator<B>>
typedef RealRowIterator<row_type> Dune::BCRSMatrix< B, A >::Iterator
template<class B, class A = std::allocator<B>>
typedef CompressedBlockVectorWindow<B,A> Dune::BCRSMatrix< B, A >::row_type

implement row_type with compressed vector

template<class B, class A = std::allocator<B>>
typedef Iterator Dune::BCRSMatrix< B, A >::RowIterator

rename the iterators for easier access

template<class B, class A = std::allocator<B>>
typedef A::size_type Dune::BCRSMatrix< B, A >::size_type

The type for the index access and the size.

Member Enumeration Documentation

template<class B, class A = std::allocator<B>>
anonymous enum

increment block level counter

Enumerator:
blocklevel 

The number of blocklevels the matrix contains.

template<class B, class A = std::allocator<B>>
enum Dune::BCRSMatrix::BuildMode

we support two modes

Enumerator:
row_wise 

Build in a row-wise manner.

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)).

random 

Build entries randomly.

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.

unknown 

Build mode not set!

Constructor & Destructor Documentation

template<class B, class A = std::allocator<B>>
Dune::BCRSMatrix< B, A >::BCRSMatrix ( )
inline

an empty matrix

template<class B, class A = std::allocator<B>>
Dune::BCRSMatrix< B, A >::BCRSMatrix ( size_type  _n,
size_type  _m,
size_type  _nnz,
BuildMode  bm 
)
inline

matrix with known number of nonzeroes

template<class B, class A = std::allocator<B>>
Dune::BCRSMatrix< B, A >::BCRSMatrix ( size_type  _n,
size_type  _m,
BuildMode  bm 
)
inline

matrix with unknown number of nonzeroes

template<class B, class A = std::allocator<B>>
Dune::BCRSMatrix< B, A >::BCRSMatrix ( const BCRSMatrix< B, A > &  Mat)
inline

copy constructor

Does a deep copy as expected.

References Dune::CompressedBlockVectorWindow< B, A >::getsize().

template<class B, class A = std::allocator<B>>
Dune::BCRSMatrix< B, A >::~BCRSMatrix ( )
inline

destructor

Member Function Documentation

template<class B, class A = std::allocator<B>>
void Dune::BCRSMatrix< B, A >::addindex ( size_type  row,
size_type  col 
)
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 col, Dune::BCRSMatrix< B, A >::end(), Dune::CompressedBlockVectorWindow< B, A >::getindexptr(), Dune::CompressedBlockVectorWindow< B, A >::getsize(), Dune::BCRSMatrix< B, A >::random, and row.

Referenced by test_IO().

template<class B, class A = std::allocator<B>>
BCRSMatrix& Dune::BCRSMatrix< B, A >::axpy ( field_type  alpha,
const BCRSMatrix< B, A > &  b 
)
inline

Add the scaled entries of another matrix to this one.

Matrix axpy operation: *this += alpha * b

Parameters
alphaScaling factor.
bThe 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 >::end(), Dune::BCRSMatrix< B, A >::M(), and Dune::BCRSMatrix< B, A >::N().

template<class B, class A = std::allocator<B>>
Iterator Dune::BCRSMatrix< B, A >::beforeBegin ( )
inline
Returns
an iterator that is positioned before the first row of the matrix.
template<class B, class A = std::allocator<B>>
ConstIterator Dune::BCRSMatrix< B, A >::beforeBegin ( ) const
inline
Returns
an iterator that is positioned before the first row of the matrix.
template<class B, class A = std::allocator<B>>
Iterator Dune::BCRSMatrix< B, A >::beforeEnd ( )
inline
Returns
an iterator that is positioned before the end iterator of the rows, i.e. at the last row.
template<class B, class A = std::allocator<B>>
ConstIterator Dune::BCRSMatrix< B, A >::beforeEnd ( ) const
inline
Returns
an iterator that is positioned before the end iterator of the rows. i.e. at the last row.
template<class B, class A = std::allocator<B>>
Iterator Dune::BCRSMatrix< B, A >::begin ( )
inline
template<class B, class A = std::allocator<B>>
ConstIterator Dune::BCRSMatrix< B, A >::begin ( ) const
inline

Get const iterator to first row.

template<class B, class A = std::allocator<B>>
CreateIterator Dune::BCRSMatrix< B, A >::createbegin ( )
inline

get initial create iterator

References Dune::BCRSMatrix< B, A >::CreateIterator.

Referenced by test_IO(), and test_Iter().

template<class B, class A = std::allocator<B>>
CreateIterator Dune::BCRSMatrix< B, A >::createend ( )
inline

get create iterator pointing to one after the last block

References Dune::BCRSMatrix< B, A >::CreateIterator.

Referenced by test_IO(), and test_Iter().

template<class B, class A = std::allocator<B>>
Iterator Dune::BCRSMatrix< B, A >::end ( )
inline
template<class B, class A = std::allocator<B>>
ConstIterator Dune::BCRSMatrix< B, A >::end ( ) const
inline

Get const iterator to one beyond last row.

template<class B, class A = std::allocator<B>>
void Dune::BCRSMatrix< B, A >::endindices ( )
inline
template<class B, class A = std::allocator<B>>
void Dune::BCRSMatrix< B, A >::endrowsizes ( )
inline
template<class B, class A = std::allocator<B>>
bool Dune::BCRSMatrix< B, A >::exists ( size_type  i,
size_type  j 
) const
inline

return true if (i,j) is in pattern

References Dune::BCRSMatrix< B, A >::end().

template<class B, class A = std::allocator<B>>
double Dune::BCRSMatrix< B, A >::frobenius_norm ( ) const
inline

frobenius norm: sqrt(sum over squared values of entries)

References Dune::BCRSMatrix< B, A >::frobenius_norm2().

template<class B, class A = std::allocator<B>>
double Dune::BCRSMatrix< B, A >::frobenius_norm2 ( ) const
inline

square of frobenius norm, need for block recursion

References Dune::BCRSMatrix< B, A >::begin(), and Dune::BCRSMatrix< B, A >::end().

Referenced by Dune::BCRSMatrix< B, A >::frobenius_norm().

template<class B, class A = std::allocator<B>>
size_type Dune::BCRSMatrix< B, A >::getrowsize ( size_type  i) const
inline

get current number of indices in row i

References Dune::CompressedBlockVectorWindow< B, A >::getsize().

template<class B, class A = std::allocator<B>>
void Dune::BCRSMatrix< B, A >::incrementrowsize ( size_type  i,
size_type  s = 1 
)
inline

increment size of row i by s (1 by default)

References Dune::BCRSMatrix< B, A >::random, and Dune::CompressedBlockVectorWindow< B, A >::setsize().

template<class B, class A = std::allocator<B>>
double Dune::BCRSMatrix< B, A >::infinity_norm ( ) const
inline

infinity norm (row sum norm, how to generalize for blocks?)

References Dune::BCRSMatrix< B, A >::begin(), and Dune::BCRSMatrix< B, A >::end().

template<class B, class A = std::allocator<B>>
double Dune::BCRSMatrix< B, A >::infinity_norm_real ( ) const
inline

simplified infinity norm (uses Manhattan norm for complex values)

References Dune::BCRSMatrix< B, A >::begin(), and Dune::BCRSMatrix< B, A >::end().

template<class B, class A = std::allocator<B>>
size_type Dune::BCRSMatrix< B, A >::M ( ) const
inline
template<class B, class A = std::allocator<B>>
template<class X , class Y >
void Dune::BCRSMatrix< B, A >::mmhv ( const X &  x,
Y &  y 
) const
inline
template<class B, class A = std::allocator<B>>
template<class X , class Y >
void Dune::BCRSMatrix< B, A >::mmtv ( const X &  x,
Y &  y 
) const
inline
template<class B, class A = std::allocator<B>>
template<class X , class Y >
void Dune::BCRSMatrix< B, A >::mmv ( const X &  x,
Y &  y 
) const
inline
template<class B, class A = std::allocator<B>>
template<class X , class Y >
void Dune::BCRSMatrix< B, A >::mtv ( const X &  x,
Y &  y 
) const
inline
template<class B, class A = std::allocator<B>>
template<class X , class Y >
void Dune::BCRSMatrix< B, A >::mv ( const X &  x,
Y &  y 
) const
inline
template<class B, class A = std::allocator<B>>
size_type Dune::BCRSMatrix< B, A >::N ( ) const
inline
template<class B, class A = std::allocator<B>>
size_type Dune::BCRSMatrix< B, A >::nonzeroes ( ) const
inline

number of blocks that are stored (the number of blocks that possibly are nonzero)

template<class B, class A = std::allocator<B>>
BCRSMatrix& Dune::BCRSMatrix< B, A >::operator*= ( const field_type k)
inline

vector space multiplication with scalar

References Dune::BCRSMatrix< B, A >::begin(), and Dune::BCRSMatrix< B, A >::end().

template<class B, class A = std::allocator<B>>
BCRSMatrix& Dune::BCRSMatrix< B, A >::operator+= ( const BCRSMatrix< B, A > &  b)
inline

Add the entries of another matrix to this one.

Parameters
bThe 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 >::end(), Dune::BCRSMatrix< B, A >::M(), and Dune::BCRSMatrix< B, A >::N().

template<class B, class A = std::allocator<B>>
BCRSMatrix& Dune::BCRSMatrix< B, A >::operator-= ( const BCRSMatrix< B, A > &  b)
inline

Substract the entries of another matrix to this one.

Parameters
bThe 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 >::end(), Dune::BCRSMatrix< B, A >::M(), and Dune::BCRSMatrix< B, A >::N().

template<class B, class A = std::allocator<B>>
BCRSMatrix& Dune::BCRSMatrix< B, A >::operator/= ( const field_type k)
inline

vector space division by scalar

References Dune::BCRSMatrix< B, A >::begin(), and Dune::BCRSMatrix< B, A >::end().

template<class B, class A = std::allocator<B>>
BCRSMatrix& Dune::BCRSMatrix< B, A >::operator= ( const BCRSMatrix< B, A > &  Mat)
inline

assignment

Frees and reallocates space. Both sparsity pattern and values are set from Mat.

References Dune::CompressedBlockVectorWindow< B, A >::getsize().

template<class B, class A = std::allocator<B>>
BCRSMatrix& Dune::BCRSMatrix< B, A >::operator= ( const field_type k)
inline

Assignment from a scalar.

Reimplemented in Dune::BTDMatrix< B, A >, and Dune::BDMatrix< B, A >.

template<class B, class A = std::allocator<B>>
row_type& Dune::BCRSMatrix< B, A >::operator[] ( size_type  i)
inline

random access to the rows

template<class B, class A = std::allocator<B>>
const row_type& Dune::BCRSMatrix< B, A >::operator[] ( size_type  i) const
inline

same for read only access

template<class B, class A = std::allocator<B>>
void Dune::BCRSMatrix< B, A >::setBuildMode ( BuildMode  bm)
inline

Sets the build mode of the matrix.

Parameters
bmThe build mode to use.
template<class B, class A = std::allocator<B>>
void Dune::BCRSMatrix< B, A >::setrowsize ( size_type  i,
size_type  s 
)
inline

set number of indices in row i to s

References Dune::BCRSMatrix< B, A >::random, and Dune::CompressedBlockVectorWindow< B, A >::setsize().

Referenced by test_IO().

template<class B, class A = std::allocator<B>>
void Dune::BCRSMatrix< B, A >::setSize ( size_type  rows,
size_type  columns,
size_type  nnz = 0 
)
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.

Warning
After calling this methods on an already allocated (and probably setup matrix) results in all the structure and data being deleted. I.~e. one has to setup the matrix again.
Parameters
rowsThe number of rows the matrix should contain.
columnsthe number of columns the matrix should contain.
nnzThe number of nonzero entries the matrix should hold (if omitted defaults to 0).
template<class B, class A = std::allocator<B>>
template<class X , class Y >
void Dune::BCRSMatrix< B, A >::umhv ( const X &  x,
Y &  y 
) const
inline
template<class B, class A = std::allocator<B>>
template<class X , class Y >
void Dune::BCRSMatrix< B, A >::umtv ( const X &  x,
Y &  y 
) const
inline
template<class B, class A = std::allocator<B>>
template<class X , class Y >
void Dune::BCRSMatrix< B, A >::umv ( const X &  x,
Y &  y 
) const
inline
template<class B, class A = std::allocator<B>>
template<class X , class Y >
void Dune::BCRSMatrix< B, A >::usmhv ( const field_type alpha,
const X &  x,
Y &  y 
) const
inline
template<class B, class A = std::allocator<B>>
template<class X , class Y >
void Dune::BCRSMatrix< B, A >::usmtv ( const field_type alpha,
const X &  x,
Y &  y 
) const
inline
template<class B, class A = std::allocator<B>>
template<class X , class Y >
void Dune::BCRSMatrix< B, A >::usmv ( const field_type alpha,
const X &  x,
Y &  y 
) const
inline

Friends And Related Function Documentation

template<class B, class A = std::allocator<B>>
friend class CreateIterator
friend
template<class B, class A = std::allocator<B>>
friend struct MatrixDimension< BCRSMatrix >
friend

The documentation for this class was generated from the following file: