dune-istl
2.2.1
|
A generic dynamic dense matrix. More...
#include <dune/istl/matrix.hh>
Public Types | |
enum | { blocklevel = T::blocklevel+1 } |
typedef T::field_type | field_type |
Export the type representing the underlying field. | |
typedef T | block_type |
Export the type representing the components. | |
typedef A | allocator_type |
Export the allocator. | |
typedef VariableBlockVector< T, A >::window_type | row_type |
The type implementing a matrix row. | |
typedef A::size_type | size_type |
Type for indices and sizes. | |
typedef VariableBlockVector< T, A >::Iterator | RowIterator |
Iterator over the matrix rows. | |
typedef row_type::iterator | ColIterator |
Iterator for the entries of each row. | |
typedef VariableBlockVector< T, A >::ConstIterator | ConstRowIterator |
Const iterator over the matrix rows. | |
typedef row_type::const_iterator | ConstColIterator |
Const iterator for the entries of each row. |
Public Member Functions | |
Matrix () | |
Create empty matrix. | |
Matrix (size_type rows, size_type cols) | |
Create uninitialized matrix of size rows x cols. | |
void | setSize (size_type rows, size_type cols) |
Change the matrix size. | |
RowIterator | begin () |
Get iterator to first row. | |
RowIterator | end () |
Get iterator to one beyond last row. | |
RowIterator | beforeEnd () |
RowIterator | beforeBegin () |
ConstRowIterator | begin () const |
Get const iterator to first row. | |
ConstRowIterator | end () const |
Get const iterator to one beyond last row. | |
ConstRowIterator | beforeEnd () const |
ConstRowIterator | beforeBegin () const |
Matrix & | operator= (const field_type &t) |
Assignment from scalar. | |
row_type & | operator[] (size_type row) |
The index operator. | |
const row_type & | operator[] (size_type row) const |
The const index operator. | |
size_type | N () const |
Return the number of rows. | |
size_type | M () const |
Return the number of columns. | |
size_type | rowdim () const |
The number of scalar rows. | |
size_type | coldim () const |
The number of scalar columns. | |
size_type | rowdim (size_type r) const |
The number of scalar rows. | |
size_type | coldim (size_type c) const |
The number of scalar columns. | |
Matrix< T > & | operator*= (const field_type &scalar) |
Multiplication with a scalar. | |
Matrix< T > & | operator/= (const field_type &scalar) |
Multiplication with a scalar. | |
Matrix & | operator+= (const Matrix &b) |
Add the entries of another matrix to this one. | |
Matrix & | operator-= (const Matrix &b) |
Subtract the entries of another matrix from this one. | |
Matrix | transpose () const |
Return the transpose of the matrix. | |
template<class X , class Y > | |
Y | transposedMult (const X &vec) |
Multiplication of the transposed matrix times a vector. | |
template<class X , class Y > | |
void | mv (const X &x, Y &y) const |
y = 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 | 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 |
![]() | |
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_norm () const |
frobenius norm: sqrt(sum over squared values of entries) | |
double | frobenius_norm2 () const |
square of frobenius norm, need for block recursion | |
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) | |
bool | exists (size_type i, size_type j) const |
return true if (i,j) is in pattern |
Protected Attributes | |
VariableBlockVector< T, A > | data_ |
Abuse VariableBlockVector as an engine for a 2d array ISTL-style. | |
size_type | cols_ |
Number of columns of the matrix. |
Friends | |
Matrix< T > | operator* (const Matrix< T > &m1, const Matrix< T > &m2) |
Generic matrix multiplication. | |
template<class X , class Y > | |
Y | operator* (const Matrix< T > &m, const X &vec) |
Generic matrix-vector multiplication. |
A generic dynamic dense matrix.
typedef A Dune::Matrix< T, A >::allocator_type |
Export the allocator.
typedef T Dune::Matrix< T, A >::block_type |
Export the type representing the components.
typedef row_type::iterator Dune::Matrix< T, A >::ColIterator |
Iterator for the entries of each row.
typedef row_type::const_iterator Dune::Matrix< T, A >::ConstColIterator |
Const iterator for the entries of each row.
typedef VariableBlockVector<T,A>::ConstIterator Dune::Matrix< T, A >::ConstRowIterator |
Const iterator over the matrix rows.
typedef T::field_type Dune::Matrix< T, A >::field_type |
Export the type representing the underlying field.
typedef VariableBlockVector<T,A>::window_type Dune::Matrix< T, A >::row_type |
The type implementing a matrix row.
typedef VariableBlockVector<T,A>::Iterator Dune::Matrix< T, A >::RowIterator |
Iterator over the matrix rows.
typedef A::size_type Dune::Matrix< T, A >::size_type |
Type for indices and sizes.
anonymous enum |
|
inline |
Create empty matrix.
|
inline |
Create uninitialized matrix of size rows x cols.
|
inline |
References Dune::VariableBlockVector< B, A >::beforeBegin(), and Dune::Matrix< T, A >::data_.
|
inline |
References Dune::VariableBlockVector< B, A >::beforeBegin(), and Dune::Matrix< T, A >::data_.
|
inline |
References Dune::VariableBlockVector< B, A >::beforeEnd(), and Dune::Matrix< T, A >::data_.
|
inline |
References Dune::VariableBlockVector< B, A >::beforeEnd(), and Dune::Matrix< T, A >::data_.
|
inline |
Get iterator to first row.
References Dune::VariableBlockVector< B, A >::begin(), and Dune::Matrix< T, A >::data_.
Referenced by Dune::SuperMatrixInitializer< BCRSMatrix< FieldMatrix< T, n, m >, A > >::addRowNnz(), Dune::CheckIfDiagonalPresent< Matrix, blocklevel, l >::check(), Dune::CheckIfDiagonalPresent< Matrix, 0, l >::check(), Dune::Matrix< T, A >::mmhv(), Dune::Matrix< T, A >::mmtv(), Dune::Matrix< T, A >::mmv(), Dune::Amg::AMG< M, X, S, PI, A >::pre(), Dune::printSparseMatrix(), Dune::Matrix< T, A >::umhv(), Dune::Matrix< T, A >::umtv(), Dune::Matrix< T, A >::usmhv(), and Dune::Matrix< T, A >::usmtv().
|
inline |
Get const iterator to first row.
References Dune::VariableBlockVector< B, A >::begin(), and Dune::Matrix< T, A >::data_.
|
inline |
The number of scalar columns.
References Dune::Matrix< T, A >::data_, Dune::Matrix< T, A >::N(), and Dune::base_array_unmanaged< B, A >::size().
|
inline |
The number of scalar columns.
References Dune::Matrix< T, A >::data_, Dune::Matrix< T, A >::M(), and Dune::Matrix< T, A >::N().
|
inline |
Get iterator to one beyond last row.
References Dune::Matrix< T, A >::data_, and Dune::VariableBlockVector< B, A >::end().
Referenced by Dune::CheckIfDiagonalPresent< Matrix, blocklevel, l >::check(), Dune::CheckIfDiagonalPresent< Matrix, 0, l >::check(), Dune::Matrix< T, A >::mmhv(), Dune::Matrix< T, A >::mmtv(), Dune::Matrix< T, A >::mmv(), Dune::Amg::AMG< M, X, S, PI, A >::pre(), Dune::printSparseMatrix(), Dune::Matrix< T, A >::umhv(), Dune::Matrix< T, A >::umtv(), Dune::Matrix< T, A >::usmhv(), and Dune::Matrix< T, A >::usmtv().
|
inline |
Get const iterator to one beyond last row.
References Dune::Matrix< T, A >::data_, and Dune::VariableBlockVector< B, A >::end().
|
inline |
return true if (i,j) is in pattern
References Dune::Matrix< T, A >::M(), and Dune::Matrix< T, A >::N().
|
inline |
frobenius norm: sqrt(sum over squared values of entries)
References Dune::Matrix< T, A >::frobenius_norm2().
|
inline |
square of frobenius norm, need for block recursion
References Dune::Matrix< T, A >::data_, Dune::Matrix< T, A >::M(), and Dune::Matrix< T, A >::N().
Referenced by Dune::Matrix< T, A >::frobenius_norm().
|
inline |
infinity norm (row sum norm, how to generalize for blocks?)
References Dune::Matrix< T, A >::data_, Dune::Matrix< T, A >::M(), and Dune::Matrix< T, A >::N().
|
inline |
simplified infinity norm (uses Manhattan norm for complex values)
References Dune::Matrix< T, A >::data_, Dune::Matrix< T, A >::M(), and Dune::Matrix< T, A >::N().
|
inline |
Return the number of columns.
References Dune::Matrix< T, A >::cols_.
Referenced by Dune::SuperLU< BCRSMatrix< FieldMatrix< T, n, m >, A > >::apply(), Dune::Matrix< T, A >::coldim(), Dune::MatrixDimension< Matrix< FieldMatrix< K, n, m >, TA > >::coldim(), Dune::Matrix< T, A >::exists(), Dune::Matrix< T, A >::frobenius_norm2(), Dune::Matrix< T, A >::infinity_norm(), Dune::Matrix< T, A >::infinity_norm_real(), Dune::Matrix< T, A >::mmhv(), Dune::Matrix< T, A >::mmtv(), Dune::Matrix< T, A >::mmv(), Dune::Matrix< T, A >::mtv(), Dune::Matrix< T, A >::mv(), Dune::Matrix< T, A >::operator+=(), Dune::Matrix< T, A >::operator-=(), Dune::printSparseMatrix(), Dune::Matrix< T, A >::rowdim(), Dune::SuperLU< BCRSMatrix< FieldMatrix< T, n, m >, A > >::setMatrix(), Dune::SuperLU< BCRSMatrix< FieldMatrix< T, n, m >, A > >::setSubMatrix(), Dune::Matrix< T, A >::transpose(), Dune::Matrix< T, A >::umhv(), Dune::Matrix< T, A >::umtv(), Dune::Matrix< T, A >::umv(), Dune::Matrix< T, A >::usmhv(), Dune::Matrix< T, A >::usmtv(), Dune::Matrix< T, A >::usmv(), and Dune::SuperLU< BCRSMatrix< FieldMatrix< T, n, m >, A > >::~SuperLU().
|
inline |
y -= A^H x
References Dune::Matrix< T, A >::begin(), Dune::Matrix< T, A >::end(), Dune::Matrix< T, A >::M(), and Dune::Matrix< T, A >::N().
|
inline |
y -= A^T x
References Dune::Matrix< T, A >::begin(), Dune::Matrix< T, A >::end(), Dune::Matrix< T, A >::M(), and Dune::Matrix< T, A >::N().
|
inline |
y -= A x
References Dune::Matrix< T, A >::begin(), Dune::Matrix< T, A >::end(), Dune::Matrix< T, A >::M(), and Dune::Matrix< T, A >::N().
|
inline |
y = A^T x
References Dune::Matrix< T, A >::M(), Dune::Matrix< T, A >::N(), and Dune::Matrix< T, A >::umtv().
|
inline |
|
inline |
Return the number of rows.
References Dune::Matrix< T, A >::data_, and Dune::VariableBlockVector< B, A >::N().
Referenced by Dune::SuperLU< BCRSMatrix< FieldMatrix< T, n, m >, A > >::apply(), Dune::Matrix< T, A >::coldim(), Dune::copyToSuperMatrix(), Dune::Matrix< T, A >::exists(), Dune::Matrix< T, A >::frobenius_norm2(), Dune::Matrix< T, A >::infinity_norm(), Dune::Matrix< T, A >::infinity_norm_real(), Dune::Matrix< T, A >::mmhv(), Dune::Matrix< T, A >::mmtv(), Dune::Matrix< T, A >::mmv(), Dune::Matrix< T, A >::mtv(), Dune::Matrix< T, A >::mv(), Dune::Matrix< T, A >::operator+=(), Dune::Matrix< T, A >::operator-=(), Dune::Matrix< T, A >::operator[](), Dune::printSparseMatrix(), Dune::Matrix< T, A >::rowdim(), Dune::MatrixDimension< Matrix< FieldMatrix< K, n, m >, TA > >::rowdim(), Dune::SuperLU< BCRSMatrix< FieldMatrix< T, n, m >, A > >::setMatrix(), Dune::SuperLU< BCRSMatrix< FieldMatrix< T, n, m >, A > >::setSubMatrix(), Dune::Matrix< T, A >::transpose(), Dune::Matrix< T, A >::umhv(), Dune::Matrix< T, A >::umtv(), Dune::Matrix< T, A >::umv(), Dune::Matrix< T, A >::usmhv(), Dune::Matrix< T, A >::usmtv(), Dune::Matrix< T, A >::usmv(), and Dune::SuperLU< BCRSMatrix< FieldMatrix< T, n, m >, A > >::~SuperLU().
|
inline |
Multiplication with a scalar.
References Dune::Matrix< T, A >::data_.
|
inline |
Add the entries of another matrix to this one.
b | The matrix to add to this one. Its size has to be the same as the size of this matrix. |
References Dune::Matrix< T, A >::data_, Dune::Matrix< T, A >::M(), and Dune::Matrix< T, A >::N().
|
inline |
Subtract the entries of another matrix from this one.
b | The matrix to add to this one. Its size has to be the same as the size of this matrix. |
References Dune::Matrix< T, A >::data_, Dune::Matrix< T, A >::M(), and Dune::Matrix< T, A >::N().
|
inline |
Multiplication with a scalar.
References Dune::Matrix< T, A >::data_.
|
inline |
Assignment from scalar.
References Dune::Matrix< T, A >::data_.
|
inline |
The index operator.
References Dune::Matrix< T, A >::data_, Dune::Matrix< T, A >::N(), and row.
|
inline |
The const index operator.
References Dune::Matrix< T, A >::data_, Dune::Matrix< T, A >::N(), and row.
|
inline |
The number of scalar rows.
References Dune::Matrix< T, A >::data_, Dune::Matrix< T, A >::M(), and Dune::VariableBlockVector< B, A >::N().
|
inline |
The number of scalar rows.
References Dune::Matrix< T, A >::data_, Dune::Matrix< T, A >::M(), and Dune::Matrix< T, A >::N().
|
inline |
Change the matrix size.
The way the data is handled is unpredictable.
References Dune::Matrix< T, A >::cols_, Dune::Matrix< T, A >::data_, and Dune::VariableBlockVector< B, A >::resize().
|
inline |
Return the transpose of the matrix.
References Dune::Matrix< T, A >::M(), and Dune::Matrix< T, A >::N().
|
inline |
y += A^H x
References Dune::Matrix< T, A >::begin(), Dune::Matrix< T, A >::end(), Dune::Matrix< T, A >::M(), and Dune::Matrix< T, A >::N().
|
inline |
y += A^T x
References Dune::Matrix< T, A >::begin(), Dune::Matrix< T, A >::end(), Dune::Matrix< T, A >::M(), and Dune::Matrix< T, A >::N().
Referenced by Dune::Matrix< T, A >::mtv().
|
inline |
|
inline |
y += alpha A^H x
References Dune::Matrix< T, A >::begin(), Dune::Matrix< T, A >::end(), Dune::Matrix< T, A >::M(), and Dune::Matrix< T, A >::N().
|
inline |
y += alpha A^T x
References Dune::Matrix< T, A >::begin(), Dune::Matrix< T, A >::end(), Dune::Matrix< T, A >::M(), and Dune::Matrix< T, A >::N().
|
inline |
|
friend |
Generic matrix multiplication.
|
friend |
Generic matrix-vector multiplication.
|
protected |
Number of columns of the matrix.
In general you can extract the same information from the data_ member. However if you want to be able to properly handle 0xn matrices then you need a separate member.
Referenced by Dune::Matrix< T, A >::M(), Dune::Matrix< T, A >::mv(), Dune::Matrix< T, A >::setSize(), Dune::Matrix< T, A >::umv(), and Dune::Matrix< T, A >::usmv().
|
protected |
Abuse VariableBlockVector as an engine for a 2d array ISTL-style.
This is almost as good as it can get. Further speedup may be possible by using the fact that all rows have the same length.
Referenced by Dune::Matrix< T, A >::beforeBegin(), Dune::Matrix< T, A >::beforeEnd(), Dune::Matrix< T, A >::begin(), Dune::Matrix< T, A >::coldim(), Dune::Matrix< T, A >::end(), Dune::Matrix< T, A >::frobenius_norm2(), Dune::Matrix< T, A >::infinity_norm(), Dune::Matrix< T, A >::infinity_norm_real(), Dune::Matrix< T, A >::mv(), Dune::Matrix< T, A >::N(), Dune::Matrix< T, A >::operator*=(), Dune::Matrix< T, A >::operator+=(), Dune::Matrix< T, A >::operator-=(), Dune::Matrix< T, A >::operator/=(), Dune::Matrix< T, A >::operator=(), Dune::Matrix< T, A >::operator[](), Dune::Matrix< T, A >::rowdim(), Dune::Matrix< T, A >::setSize(), Dune::Matrix< T, A >::umv(), and Dune::Matrix< T, A >::usmv().