dune-istl  2.6-git
Public Types | Public Member Functions | Protected Attributes | Friends | List of all members
Dune::Matrix< T, A > Class Template Reference

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. More...
 
typedef T block_type
 Export the type representing the components. More...
 
typedef A allocator_type
 Export the allocator. More...
 
typedef MatrixImp::DenseMatrixBase< T, A >::window_type row_type
 The type implementing a matrix row. More...
 
typedef A::size_type size_type
 Type for indices and sizes. More...
 
typedef MatrixImp::DenseMatrixBase< T, A >::Iterator RowIterator
 Iterator over the matrix rows. More...
 
typedef row_type::iterator ColIterator
 Iterator for the entries of each row. More...
 
typedef MatrixImp::DenseMatrixBase< T, A >::ConstIterator ConstRowIterator
 Const iterator over the matrix rows. More...
 
typedef row_type::const_iterator ConstColIterator
 Const iterator for the entries of each row. More...
 

Public Member Functions

 Matrix ()
 Create empty matrix. More...
 
 Matrix (size_type rows, size_type cols)
 Create uninitialized matrix of size rows x cols. More...
 
void setSize (size_type rows, size_type cols)
 Change the matrix size. More...
 
RowIterator begin ()
 Get iterator to first row. More...
 
RowIterator end ()
 Get iterator to one beyond last row. More...
 
RowIterator beforeEnd ()
 
RowIterator beforeBegin ()
 
ConstRowIterator begin () const
 Get const iterator to first row. More...
 
ConstRowIterator end () const
 Get const iterator to one beyond last row. More...
 
ConstRowIterator beforeEnd () const
 
ConstRowIterator beforeBegin () const
 
Matrixoperator= (const field_type &t)
 Assignment from scalar. More...
 
row_type operator[] (size_type row)
 The index operator. More...
 
const row_type operator[] (size_type row) const
 The const index operator. More...
 
size_type N () const
 Return the number of rows. More...
 
size_type M () const
 Return the number of columns. More...
 
Matrix< T > & operator*= (const field_type &scalar)
 Multiplication with a scalar. More...
 
Matrix< T > & operator/= (const field_type &scalar)
 Division by a scalar. More...
 
Matrixoperator+= (const Matrix &b)
 Add the entries of another matrix to this one. More...
 
Matrixoperator-= (const Matrix &b)
 Subtract the entries of another matrix from this one. More...
 
Matrix transpose () const
 Return the transpose of the matrix. More...
 
template<class X , class Y >
void mv (const X &x, Y &y) const
 y = 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 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 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_norm () const
 frobenius norm: sqrt(sum over squared values of entries) More...
 
FieldTraits< field_type >::real_type frobenius_norm2 () const
 square of frobenius norm, need for block recursion More...
 
template<typename ft = field_type, typename std::enable_if<!has_nan< ft >::value, int >::type = 0>
FieldTraits< ft >::real_type infinity_norm () const
 infinity norm (row sum norm, how to generalize for blocks?) More...
 
template<typename ft = field_type, typename std::enable_if<!has_nan< ft >::value, int >::type = 0>
FieldTraits< ft >::real_type infinity_norm_real () const
 simplified infinity norm (uses Manhattan norm for complex values) More...
 
template<typename ft = field_type, typename std::enable_if< has_nan< ft >::value, int >::type = 0>
FieldTraits< ft >::real_type infinity_norm () const
 infinity norm (row sum norm, how to generalize for blocks?) More...
 
template<typename ft = field_type, typename std::enable_if< has_nan< ft >::value, int >::type = 0>
FieldTraits< ft >::real_type infinity_norm_real () const
 simplified infinity norm (uses Manhattan norm for complex values) More...
 
bool exists (size_type i, size_type j) const
 return true if (i,j) is in pattern More...
 

Protected Attributes

MatrixImp::DenseMatrixBase< T, A > data_
 Abuse DenseMatrixBase as an engine for a 2d array ISTL-style. More...
 
size_type cols_
 Number of columns of the matrix. More...
 

Friends

Matrix< T > operator* (const Matrix< T > &m1, const Matrix< T > &m2)
 Generic matrix multiplication. More...
 
template<class X , class Y >
operator* (const Matrix< T > &m, const X &vec)
 Generic matrix-vector multiplication. More...
 

Detailed Description

template<class T, class A = std::allocator<T>>
class Dune::Matrix< T, A >

A generic dynamic dense matrix.

Member Typedef Documentation

◆ allocator_type

template<class T, class A = std::allocator<T>>
typedef A Dune::Matrix< T, A >::allocator_type

Export the allocator.

◆ block_type

template<class T, class A = std::allocator<T>>
typedef T Dune::Matrix< T, A >::block_type

Export the type representing the components.

◆ ColIterator

template<class T, class A = std::allocator<T>>
typedef row_type::iterator Dune::Matrix< T, A >::ColIterator

Iterator for the entries of each row.

◆ ConstColIterator

template<class T, class A = std::allocator<T>>
typedef row_type::const_iterator Dune::Matrix< T, A >::ConstColIterator

Const iterator for the entries of each row.

◆ ConstRowIterator

template<class T, class A = std::allocator<T>>
typedef MatrixImp::DenseMatrixBase<T,A>::ConstIterator Dune::Matrix< T, A >::ConstRowIterator

Const iterator over the matrix rows.

◆ field_type

template<class T, class A = std::allocator<T>>
typedef T::field_type Dune::Matrix< T, A >::field_type

Export the type representing the underlying field.

◆ row_type

template<class T, class A = std::allocator<T>>
typedef MatrixImp::DenseMatrixBase<T,A>::window_type Dune::Matrix< T, A >::row_type

The type implementing a matrix row.

◆ RowIterator

template<class T, class A = std::allocator<T>>
typedef MatrixImp::DenseMatrixBase<T,A>::Iterator Dune::Matrix< T, A >::RowIterator

Iterator over the matrix rows.

◆ size_type

template<class T, class A = std::allocator<T>>
typedef A::size_type Dune::Matrix< T, A >::size_type

Type for indices and sizes.

Member Enumeration Documentation

◆ anonymous enum

template<class T, class A = std::allocator<T>>
anonymous enum
Enumerator
blocklevel 

The number of nesting levels the matrix contains.

Constructor & Destructor Documentation

◆ Matrix() [1/2]

template<class T, class A = std::allocator<T>>
Dune::Matrix< T, A >::Matrix ( )
inline

Create empty matrix.

◆ Matrix() [2/2]

template<class T, class A = std::allocator<T>>
Dune::Matrix< T, A >::Matrix ( size_type  rows,
size_type  cols 
)
inline

Create uninitialized matrix of size rows x cols.

Member Function Documentation

◆ beforeBegin() [1/2]

template<class T, class A = std::allocator<T>>
RowIterator Dune::Matrix< T, A >::beforeBegin ( )
inline
Returns
an iterator that is positioned before the first row of the matrix.

◆ beforeBegin() [2/2]

template<class T, class A = std::allocator<T>>
ConstRowIterator Dune::Matrix< T, A >::beforeBegin ( ) const
inline
Returns
an iterator that is positioned before the first row if the matrix.

◆ beforeEnd() [1/2]

template<class T, class A = std::allocator<T>>
RowIterator Dune::Matrix< T, A >::beforeEnd ( )
inline
Returns
an iterator that is positioned before the end iterator of the rows, i.e. at the last row.

◆ beforeEnd() [2/2]

template<class T, class A = std::allocator<T>>
ConstRowIterator Dune::Matrix< T, A >::beforeEnd ( ) const
inline
Returns
an iterator that is positioned before the end iterator of the rows. i.e. at the last row.

◆ begin() [1/2]

template<class T, class A = std::allocator<T>>
RowIterator Dune::Matrix< T, A >::begin ( )
inline

Get iterator to first row.

◆ begin() [2/2]

template<class T, class A = std::allocator<T>>
ConstRowIterator Dune::Matrix< T, A >::begin ( ) const
inline

Get const iterator to first row.

◆ end() [1/2]

template<class T, class A = std::allocator<T>>
RowIterator Dune::Matrix< T, A >::end ( )
inline

Get iterator to one beyond last row.

◆ end() [2/2]

template<class T, class A = std::allocator<T>>
ConstRowIterator Dune::Matrix< T, A >::end ( ) const
inline

Get const iterator to one beyond last row.

◆ exists()

template<class T, class A = std::allocator<T>>
bool Dune::Matrix< T, A >::exists ( size_type  i,
size_type  j 
) const
inline

return true if (i,j) is in pattern

◆ frobenius_norm()

template<class T, class A = std::allocator<T>>
FieldTraits<field_type>::real_type Dune::Matrix< T, A >::frobenius_norm ( ) const
inline

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

◆ frobenius_norm2()

template<class T, class A = std::allocator<T>>
FieldTraits<field_type>::real_type Dune::Matrix< T, A >::frobenius_norm2 ( ) const
inline

square of frobenius norm, need for block recursion

◆ infinity_norm() [1/2]

template<class T, class A = std::allocator<T>>
template<typename ft = field_type, typename std::enable_if<!has_nan< ft >::value, int >::type = 0>
FieldTraits<ft>::real_type Dune::Matrix< T, A >::infinity_norm ( ) const
inline

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

◆ infinity_norm() [2/2]

template<class T, class A = std::allocator<T>>
template<typename ft = field_type, typename std::enable_if< has_nan< ft >::value, int >::type = 0>
FieldTraits<ft>::real_type Dune::Matrix< T, A >::infinity_norm ( ) const
inline

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

◆ infinity_norm_real() [1/2]

template<class T, class A = std::allocator<T>>
template<typename ft = field_type, typename std::enable_if<!has_nan< ft >::value, int >::type = 0>
FieldTraits<ft>::real_type Dune::Matrix< T, A >::infinity_norm_real ( ) const
inline

simplified infinity norm (uses Manhattan norm for complex values)

◆ infinity_norm_real() [2/2]

template<class T, class A = std::allocator<T>>
template<typename ft = field_type, typename std::enable_if< has_nan< ft >::value, int >::type = 0>
FieldTraits<ft>::real_type Dune::Matrix< T, A >::infinity_norm_real ( ) const
inline

simplified infinity norm (uses Manhattan norm for complex values)

◆ M()

template<class T, class A = std::allocator<T>>
size_type Dune::Matrix< T, A >::M ( ) const
inline

Return the number of columns.

◆ mmhv()

template<class T, class A = std::allocator<T>>
template<class X , class Y >
void Dune::Matrix< T, A >::mmhv ( const X &  x,
Y &  y 
) const
inline

y -= A^H x

◆ mmtv()

template<class T, class A = std::allocator<T>>
template<class X , class Y >
void Dune::Matrix< T, A >::mmtv ( const X &  x,
Y &  y 
) const
inline

y -= A^T x

◆ mmv()

template<class T, class A = std::allocator<T>>
template<class X , class Y >
void Dune::Matrix< T, A >::mmv ( const X &  x,
Y &  y 
) const
inline

y -= A x

◆ mtv()

template<class T, class A = std::allocator<T>>
template<class X , class Y >
void Dune::Matrix< T, A >::mtv ( const X &  x,
Y &  y 
) const
inline

y = A^T x

◆ mv()

template<class T, class A = std::allocator<T>>
template<class X , class Y >
void Dune::Matrix< T, A >::mv ( const X &  x,
Y &  y 
) const
inline

y = A x

◆ N()

template<class T, class A = std::allocator<T>>
size_type Dune::Matrix< T, A >::N ( ) const
inline

Return the number of rows.

◆ operator*=()

template<class T, class A = std::allocator<T>>
Matrix<T>& Dune::Matrix< T, A >::operator*= ( const field_type scalar)
inline

Multiplication with a scalar.

◆ operator+=()

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

Add the entries of another matrix to this one.

Parameters
bThe matrix to add to this one. Its size has to be the same as the size of this matrix.

◆ operator-=()

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

Subtract the entries of another matrix from this one.

Parameters
bThe matrix to subtract from this one. Its size has to be the same as the size of this matrix.

◆ operator/=()

template<class T, class A = std::allocator<T>>
Matrix<T>& Dune::Matrix< T, A >::operator/= ( const field_type scalar)
inline

Division by a scalar.

◆ operator=()

template<class T, class A = std::allocator<T>>
Matrix& Dune::Matrix< T, A >::operator= ( const field_type t)
inline

Assignment from scalar.

◆ operator[]() [1/2]

template<class T, class A = std::allocator<T>>
row_type Dune::Matrix< T, A >::operator[] ( size_type  row)
inline

The index operator.

◆ operator[]() [2/2]

template<class T, class A = std::allocator<T>>
const row_type Dune::Matrix< T, A >::operator[] ( size_type  row) const
inline

The const index operator.

◆ setSize()

template<class T, class A = std::allocator<T>>
void Dune::Matrix< T, A >::setSize ( size_type  rows,
size_type  cols 
)
inline

Change the matrix size.

The way the data is handled is unpredictable.

◆ transpose()

template<class T, class A = std::allocator<T>>
Matrix Dune::Matrix< T, A >::transpose ( ) const
inline

Return the transpose of the matrix.

◆ umhv()

template<class T, class A = std::allocator<T>>
template<class X , class Y >
void Dune::Matrix< T, A >::umhv ( const X &  x,
Y &  y 
) const
inline

y += A^H x

◆ umtv()

template<class T, class A = std::allocator<T>>
template<class X , class Y >
void Dune::Matrix< T, A >::umtv ( const X &  x,
Y &  y 
) const
inline

y += A^T x

◆ umv()

template<class T, class A = std::allocator<T>>
template<class X , class Y >
void Dune::Matrix< T, A >::umv ( const X &  x,
Y &  y 
) const
inline

y += A x

◆ usmhv()

template<class T, class A = std::allocator<T>>
template<class X , class Y >
void Dune::Matrix< T, A >::usmhv ( const field_type alpha,
const X &  x,
Y &  y 
) const
inline

y += alpha A^H x

◆ usmtv()

template<class T, class A = std::allocator<T>>
template<class X , class Y >
void Dune::Matrix< T, A >::usmtv ( const field_type alpha,
const X &  x,
Y &  y 
) const
inline

y += alpha A^T x

◆ usmv()

template<class T, class A = std::allocator<T>>
template<class X , class Y >
void Dune::Matrix< T, A >::usmv ( const field_type alpha,
const X &  x,
Y &  y 
) const
inline

$ y += \alpha A x $

Friends And Related Function Documentation

◆ operator* [1/2]

template<class T, class A = std::allocator<T>>
Matrix<T> operator* ( const Matrix< T > &  m1,
const Matrix< T > &  m2 
)
friend

Generic matrix multiplication.

◆ operator* [2/2]

template<class T, class A = std::allocator<T>>
template<class X , class Y >
Y operator* ( const Matrix< T > &  m,
const X &  vec 
)
friend

Generic matrix-vector multiplication.

Member Data Documentation

◆ cols_

template<class T, class A = std::allocator<T>>
size_type Dune::Matrix< T, A >::cols_
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.

◆ data_

template<class T, class A = std::allocator<T>>
MatrixImp::DenseMatrixBase<T,A> Dune::Matrix< T, A >::data_
protected

Abuse DenseMatrixBase as an engine for a 2d array ISTL-style.


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