Cantera 2.6.0
Public Member Functions | Protected Attributes | List of all members
BandMatrix Class Reference

A class for banded matrices, involving matrix inversion processes. More...

#include <BandMatrix.h>

Inheritance diagram for BandMatrix:
[legend]
Collaboration diagram for BandMatrix:
[legend]

Public Member Functions

 BandMatrix ()
 Base Constructor. More...
 
 BandMatrix (size_t n, size_t kl, size_t ku, doublereal v=0.0)
 Creates a banded matrix and sets all elements to zero. More...
 
 BandMatrix (const BandMatrix &y)
 
BandMatrixoperator= (const BandMatrix &y)
 
void resize (size_t n, size_t kl, size_t ku, doublereal v=0.0)
 Resize the matrix problem. More...
 
void bfill (doublereal v=0.0)
 Fill or zero the matrix. More...
 
doublereal & operator() (size_t i, size_t j)
 Index into the (i,j) element. More...
 
doublereal operator() (size_t i, size_t j) const
 Constant Index into the (i,j) element. More...
 
doublereal & value (size_t i, size_t j)
 Return a changeable reference to element (i,j). More...
 
doublereal value (size_t i, size_t j) const
 Return the value of element (i,j). More...
 
size_t index (size_t i, size_t j) const
 Returns the location in the internal 1D array corresponding to the (i,j) element in the banded array. More...
 
doublereal _value (size_t i, size_t j) const
 Return the value of the (i,j) element for (i,j) within the bandwidth. More...
 
virtual size_t nRows () const
 Return the number of rows in the matrix. More...
 
size_t nColumns () const
 Number of columns. More...
 
size_t nSubDiagonals () const
 Number of subdiagonals. More...
 
size_t nSuperDiagonals () const
 Number of superdiagonals. More...
 
size_t ldim () const
 Return the number of rows of storage needed for the band storage. More...
 
virtual void mult (const doublereal *b, doublereal *prod) const
 Multiply A*b and write result to prod. More...
 
virtual void leftMult (const doublereal *const b, doublereal *const prod) const
 Multiply b*A and write result to prod. More...
 
int factor ()
 Perform an LU decomposition, the LAPACK routine DGBTRF is used. More...
 
int solve (const doublereal *const b, doublereal *const x)
 Solve the matrix problem Ax = b. More...
 
int solve (doublereal *b, size_t nrhs=1, size_t ldb=0)
 Solve the matrix problem Ax = b. More...
 
virtual vector_fp::iterator begin ()
 Returns an iterator for the start of the band storage data. More...
 
vector_fp::iterator end ()
 Returns an iterator for the end of the band storage data. More...
 
vector_fp::const_iterator begin () const
 Returns a const iterator for the start of the band storage data. More...
 
vector_fp::const_iterator end () const
 Returns a const iterator for the end of the band storage data. More...
 
virtual void zero ()
 Zero the matrix elements. More...
 
virtual doublereal rcond (doublereal a1norm)
 Returns an estimate of the inverse of the condition number for the matrix. More...
 
virtual int factorAlgorithm () const
 Returns the factor algorithm used. More...
 
virtual doublereal oneNorm () const
 Returns the one norm of the matrix. More...
 
virtual doublereal * ptrColumn (size_t j)
 Return a pointer to the top of column j. More...
 
virtual doublereal *const * colPts ()
 Return a vector of const pointers to the columns. More...
 
virtual size_t checkRows (doublereal &valueSmall) const
 Check to see if we have any zero rows in the Jacobian. More...
 
virtual size_t checkColumns (doublereal &valueSmall) const
 Check to see if we have any zero columns in the Jacobian. More...
 
int info () const
 LAPACK "info" flag after last factor/solve operation. More...
 
- Public Member Functions inherited from GeneralMatrix
 GeneralMatrix ()
 Base Constructor. More...
 
virtual int factorQR ()
 Factors the A matrix using the QR algorithm, overwriting A. More...
 
virtual doublereal rcondQR ()
 Returns an estimate of the inverse of the condition number for the matrix. More...
 
virtual void useFactorAlgorithm (int fAlgorithm)
 Change the way the matrix is factored. More...
 
virtual void clearFactorFlag ()
 clear the factored flag More...
 
virtual bool factored () const
 true if the current factorization is up to date with the matrix More...
 

Protected Attributes

vector_fp data
 Matrix data. More...
 
vector_fp ludata
 Factorized data. More...
 
size_t m_n
 Number of rows and columns of the matrix. More...
 
size_t m_kl
 Number of subdiagonals of the matrix. More...
 
size_t m_ku
 Number of super diagonals of the matrix. More...
 
doublereal m_zero
 value of zero More...
 
std::unique_ptr< PivData > m_ipiv
 Pivot vector. More...
 
std::vector< doublereal * > m_colPtrs
 Vector of column pointers. More...
 
std::vector< double * > m_lu_col_ptrs
 
vector_int iwork_
 Extra work array needed - size = n. More...
 
vector_fp work_
 Extra dp work array needed - size = 3n. More...
 
int m_info
 
- Protected Attributes inherited from GeneralMatrix
int m_factored
 Indicates whether the matrix is factored. More...
 

Detailed Description

A class for banded matrices, involving matrix inversion processes.

The class is based upon the LAPACK banded storage matrix format.

An important issue with this class is that it stores both the original data and the LU factorization of the data. This means that the banded matrix typically will take up twice the room that it is expected to take.

QR factorizations of banded matrices are not included in the original LAPACK work. Add-ons are available. However, they are not included here. Instead we just use the stock LU decompositions.

This class is a derived class of the base class GeneralMatrix. However, within the oneD directory, the class is used as is, without reference to the GeneralMatrix base type.

Definition at line 34 of file BandMatrix.h.

Constructor & Destructor Documentation

◆ BandMatrix() [1/3]

Base Constructor.

Create an 0 by 0 matrix, and initialize all elements to 0.

Definition at line 43 of file BandMatrix.cpp.

◆ ~BandMatrix()

~BandMatrix ( )

Definition at line 53 of file BandMatrix.cpp.

◆ BandMatrix() [2/3]

BandMatrix ( size_t  n,
size_t  kl,
size_t  ku,
doublereal  v = 0.0 
)

Creates a banded matrix and sets all elements to zero.

Create an n by n banded matrix, and initialize all elements to v.

Parameters
nsize of the square matrix
klband size on the lower portion of the matrix
kuband size on the upper portion of the matrix
vinitial value of all matrix components.

Definition at line 58 of file BandMatrix.cpp.

References BandMatrix::data, BandMatrix::ludata, BandMatrix::m_colPtrs, BandMatrix::m_ipiv, and BandMatrix::m_n.

◆ BandMatrix() [3/3]

BandMatrix ( const BandMatrix y)

Definition at line 80 of file BandMatrix.cpp.

Member Function Documentation

◆ operator=()

BandMatrix & operator= ( const BandMatrix y)

Definition at line 104 of file BandMatrix.cpp.

◆ resize()

void resize ( size_t  n,
size_t  kl,
size_t  ku,
doublereal  v = 0.0 
)

Resize the matrix problem.

All data is lost

Parameters
nsize of the square matrix
klband size on the lower portion of the matrix
kuband size on the upper portion of the matrix
vinitial value of all matrix components.

Definition at line 127 of file BandMatrix.cpp.

References BandMatrix::data, BandMatrix::ludata, BandMatrix::m_colPtrs, GeneralMatrix::m_factored, BandMatrix::m_ipiv, BandMatrix::m_kl, BandMatrix::m_ku, and BandMatrix::m_n.

◆ bfill()

void bfill ( doublereal  v = 0.0)

Fill or zero the matrix.

Parameters
vFill value, defaults to zero.

Definition at line 146 of file BandMatrix.cpp.

References BandMatrix::data, and GeneralMatrix::m_factored.

◆ operator()() [1/2]

doublereal & operator() ( size_t  i,
size_t  j 
)
virtual

Index into the (i,j) element.

Parameters
irow
jcolumn
Returns
a changeable reference to the matrix entry

Implements GeneralMatrix.

Definition at line 158 of file BandMatrix.cpp.

References BandMatrix::value().

◆ operator()() [2/2]

doublereal operator() ( size_t  i,
size_t  j 
) const
virtual

Constant Index into the (i,j) element.

Parameters
irow
jcolumn
Returns
an unchangeable reference to the matrix entry

Implements GeneralMatrix.

Definition at line 163 of file BandMatrix.cpp.

References BandMatrix::value().

◆ value() [1/2]

doublereal & value ( size_t  i,
size_t  j 
)

Return a changeable reference to element (i,j).

Since this method may alter the element value, it may need to be refactored, so the flag m_factored is set to false.

Parameters
irow
jcolumn
Returns
a reference to the value of the matrix entry

Definition at line 168 of file BandMatrix.cpp.

References BandMatrix::data, BandMatrix::index(), GeneralMatrix::m_factored, BandMatrix::m_kl, and BandMatrix::m_zero.

Referenced by BandMatrix::oneNorm(), and BandMatrix::operator()().

◆ value() [2/2]

doublereal value ( size_t  i,
size_t  j 
) const

Return the value of element (i,j).

This method does not alter the array.

Parameters
irow
jcolumn
Returns
the value of the matrix entry

Definition at line 177 of file BandMatrix.cpp.

References BandMatrix::data, BandMatrix::index(), and BandMatrix::m_kl.

◆ index()

size_t index ( size_t  i,
size_t  j 
) const

Returns the location in the internal 1D array corresponding to the (i,j) element in the banded array.

Parameters
irow
jcolumn
Returns
the index of the matrix entry

Definition at line 185 of file BandMatrix.cpp.

References BandMatrix::m_kl, and BandMatrix::m_ku.

Referenced by BandMatrix::_value(), and BandMatrix::value().

◆ _value()

doublereal _value ( size_t  i,
size_t  j 
) const

Return the value of the (i,j) element for (i,j) within the bandwidth.

For efficiency, this method does not check that (i,j) are within the bandwidth; it is up to the calling program to insure that this is true.

Parameters
irow
jcolumn
Returns
the value of the matrix entry

Definition at line 190 of file BandMatrix.cpp.

References BandMatrix::data, and BandMatrix::index().

Referenced by BandMatrix::checkColumns(), BandMatrix::checkRows(), BandMatrix::leftMult(), BandMatrix::mult(), and BandMatrix::oneNorm().

◆ nRows()

size_t nRows ( ) const
virtual

Return the number of rows in the matrix.

Implements GeneralMatrix.

Definition at line 195 of file BandMatrix.cpp.

References BandMatrix::m_n.

Referenced by BandMatrix::factor(), and Cantera::operator<<().

◆ nColumns()

size_t nColumns ( ) const

Number of columns.

Definition at line 200 of file BandMatrix.cpp.

References BandMatrix::m_n.

Referenced by BandMatrix::factor(), and BandMatrix::solve().

◆ nSubDiagonals()

size_t nSubDiagonals ( ) const

Number of subdiagonals.

Definition at line 205 of file BandMatrix.cpp.

References BandMatrix::m_kl.

Referenced by BandMatrix::factor().

◆ nSuperDiagonals()

size_t nSuperDiagonals ( ) const

Number of superdiagonals.

Definition at line 210 of file BandMatrix.cpp.

References BandMatrix::m_ku.

Referenced by BandMatrix::factor().

◆ ldim()

size_t ldim ( ) const

Return the number of rows of storage needed for the band storage.

Definition at line 215 of file BandMatrix.cpp.

References BandMatrix::m_kl, and BandMatrix::m_ku.

Referenced by BandMatrix::factor().

◆ mult()

void mult ( const doublereal *  b,
doublereal *  prod 
) const
virtual

Multiply A*b and write result to prod.

Implements GeneralMatrix.

Definition at line 220 of file BandMatrix.cpp.

References BandMatrix::_value(), BandMatrix::m_kl, BandMatrix::m_ku, and BandMatrix::m_n.

◆ leftMult()

void leftMult ( const doublereal *const  b,
doublereal *const  prod 
) const
virtual

Multiply b*A and write result to prod.

Parameters
bVector to do the lh multiplication
prodOUTPUT vector to receive the result

Implements GeneralMatrix.

Definition at line 233 of file BandMatrix.cpp.

References BandMatrix::_value(), BandMatrix::m_kl, BandMatrix::m_ku, and BandMatrix::m_n.

◆ factor()

int factor ( )
virtual

Perform an LU decomposition, the LAPACK routine DGBTRF is used.

The factorization is saved in ludata.

Returns
a success flag. 0 indicates a success; ~0 indicates some error occurred, see the LAPACK documentation

Implements GeneralMatrix.

Definition at line 246 of file BandMatrix.cpp.

References BandMatrix::data, BandMatrix::ldim(), BandMatrix::ludata, GeneralMatrix::m_factored, BandMatrix::m_ipiv, BandMatrix::nColumns(), BandMatrix::nRows(), BandMatrix::nSubDiagonals(), and BandMatrix::nSuperDiagonals().

Referenced by BandMatrix::solve().

◆ solve() [1/2]

int solve ( const doublereal *const  b,
doublereal *const  x 
)

Solve the matrix problem Ax = b.

Parameters
bINPUT RHS of the problem
xOUTPUT solution to the problem
Returns
a success flag. 0 indicates a success; ~0 indicates some error occurred, see the LAPACK documentation

Definition at line 267 of file BandMatrix.cpp.

References BandMatrix::m_n, and BandMatrix::solve().

Referenced by BandMatrix::solve(), and MultiNewton::step().

◆ solve() [2/2]

int solve ( doublereal *  b,
size_t  nrhs = 1,
size_t  ldb = 0 
)
virtual

Solve the matrix problem Ax = b.

Parameters
bINPUT RHS of the problem OUTPUT solution to the problem
nrhsNumber of right hand sides to solve
ldbLeading dimension of b. Default is nColumns()
Returns
a success flag. 0 indicates a success; ~0 indicates some error occurred, see the LAPACK documentation

Implements GeneralMatrix.

Definition at line 273 of file BandMatrix.cpp.

References BandMatrix::factor(), GeneralMatrix::m_factored, and BandMatrix::nColumns().

◆ begin() [1/2]

vector_fp::iterator begin ( )
virtual

Returns an iterator for the start of the band storage data.

Iterator points to the beginning of the data, and it is changeable.

Implements GeneralMatrix.

Definition at line 301 of file BandMatrix.cpp.

References BandMatrix::data, and GeneralMatrix::m_factored.

◆ end() [1/2]

vector_fp::iterator end ( )

Returns an iterator for the end of the band storage data.

Iterator points to the end of the data, and it is changeable.

Definition at line 307 of file BandMatrix.cpp.

References BandMatrix::data, and GeneralMatrix::m_factored.

◆ begin() [2/2]

vector_fp::const_iterator begin ( ) const
virtual

Returns a const iterator for the start of the band storage data.

Iterator points to the beginning of the data, and it is not changeable.

Implements GeneralMatrix.

Definition at line 313 of file BandMatrix.cpp.

References BandMatrix::data.

◆ end() [2/2]

vector_fp::const_iterator end ( ) const

Returns a const iterator for the end of the band storage data.

Iterator points to the end of the data, and it is not changeable.

Definition at line 318 of file BandMatrix.cpp.

References BandMatrix::data.

◆ zero()

void zero ( )
virtual

Zero the matrix elements.

Implements GeneralMatrix.

Definition at line 152 of file BandMatrix.cpp.

References BandMatrix::data, and GeneralMatrix::m_factored.

◆ rcond()

doublereal rcond ( doublereal  a1norm)
virtual

Returns an estimate of the inverse of the condition number for the matrix.

The matrix must have been previously factored using the LU algorithm

Parameters
a1normNorm of the matrix
Returns
the inverse of the condition number

Implements GeneralMatrix.

Definition at line 335 of file BandMatrix.cpp.

References BandMatrix::iwork_, BandMatrix::ludata, GeneralMatrix::m_factored, BandMatrix::m_ipiv, BandMatrix::m_kl, BandMatrix::m_ku, BandMatrix::m_n, BandMatrix::rcond(), and BandMatrix::work_.

Referenced by BandMatrix::rcond().

◆ factorAlgorithm()

int factorAlgorithm ( ) const
virtual

Returns the factor algorithm used.

This method will always return 0 (LU) for band matrices.

Implements GeneralMatrix.

Definition at line 358 of file BandMatrix.cpp.

◆ oneNorm()

doublereal oneNorm ( ) const
virtual

Returns the one norm of the matrix.

Implements GeneralMatrix.

Definition at line 363 of file BandMatrix.cpp.

References BandMatrix::_value(), BandMatrix::m_kl, BandMatrix::m_ku, BandMatrix::m_n, and BandMatrix::value().

◆ ptrColumn()

doublereal * ptrColumn ( size_t  j)
virtual

Return a pointer to the top of column j.

Column values are assumed to be contiguous in memory (LAPACK band matrix structure)

On entry, the matrix A in band storage, in rows KL+1 to 2*KL+KU+1; rows 1 to KL of the array need not be set. The j-th column of A is stored in the j-th column of the array AB as follows:

AB(KL + KU + 1 + i - j,j) = A(i,j) for max(1, j - KU) <= i <= min(m, j + KL)

This routine returns the position of AB(1,j) (fortran-1 indexing) in the above format

So to address the (i,j) position, you use the following indexing:

double *colP_j = matrix.ptrColumn(j);
double a_i_j = colP_j[kl + ku + i - j];
Parameters
jValue of the column
Returns
a pointer to the top of the column

Implements GeneralMatrix.

Definition at line 422 of file BandMatrix.cpp.

References BandMatrix::m_colPtrs.

◆ colPts()

doublereal *const * colPts ( )
virtual

Return a vector of const pointers to the columns.

Note the value of the pointers are protected by their being const. However, the value of the matrix is open to being changed.

Returns
a vector of pointers to the top of the columns of the matrices.

Implements GeneralMatrix.

Definition at line 427 of file BandMatrix.cpp.

References BandMatrix::m_colPtrs.

◆ checkRows()

size_t checkRows ( doublereal &  valueSmall) const
virtual

Check to see if we have any zero rows in the Jacobian.

This utility routine checks to see if any rows are zero. The smallest row is returned along with the largest coefficient in that row

Parameters
valueSmallOUTPUT value of the largest coefficient in the smallest row
Returns
index of the row that is most nearly zero

Implements GeneralMatrix.

Definition at line 378 of file BandMatrix.cpp.

References BandMatrix::_value(), BandMatrix::m_kl, BandMatrix::m_ku, BandMatrix::m_n, and Cantera::npos.

◆ checkColumns()

size_t checkColumns ( doublereal &  valueSmall) const
virtual

Check to see if we have any zero columns in the Jacobian.

This utility routine checks to see if any columns are zero. The smallest column is returned along with the largest coefficient in that column

Parameters
valueSmallOUTPUT value of the largest coefficient in the smallest column
Returns
index of the column that is most nearly zero

Implements GeneralMatrix.

Definition at line 400 of file BandMatrix.cpp.

References BandMatrix::_value(), BandMatrix::m_kl, BandMatrix::m_ku, BandMatrix::m_n, and Cantera::npos.

◆ info()

int info ( ) const
inline

LAPACK "info" flag after last factor/solve operation.

Definition at line 261 of file BandMatrix.h.

Referenced by MultiNewton::step().

Member Data Documentation

◆ data

vector_fp data
protected

◆ ludata

vector_fp ludata
protected

Factorized data.

Definition at line 268 of file BandMatrix.h.

Referenced by BandMatrix::BandMatrix(), BandMatrix::factor(), BandMatrix::rcond(), and BandMatrix::resize().

◆ m_n

size_t m_n
protected

◆ m_kl

size_t m_kl
protected

◆ m_ku

size_t m_ku
protected

◆ m_zero

doublereal m_zero
protected

value of zero

Definition at line 280 of file BandMatrix.h.

Referenced by BandMatrix::value().

◆ m_ipiv

std::unique_ptr<PivData> m_ipiv
protected

Pivot vector.

Definition at line 285 of file BandMatrix.h.

Referenced by BandMatrix::BandMatrix(), BandMatrix::factor(), BandMatrix::rcond(), and BandMatrix::resize().

◆ m_colPtrs

std::vector<doublereal*> m_colPtrs
protected

Vector of column pointers.

Definition at line 288 of file BandMatrix.h.

Referenced by BandMatrix::BandMatrix(), BandMatrix::colPts(), BandMatrix::ptrColumn(), and BandMatrix::resize().

◆ m_lu_col_ptrs

std::vector<double*> m_lu_col_ptrs
protected

Definition at line 289 of file BandMatrix.h.

◆ iwork_

vector_int iwork_
protected

Extra work array needed - size = n.

Definition at line 292 of file BandMatrix.h.

Referenced by BandMatrix::rcond().

◆ work_

vector_fp work_
protected

Extra dp work array needed - size = 3n.

Definition at line 295 of file BandMatrix.h.

Referenced by BandMatrix::rcond().

◆ m_info

int m_info
protected

Definition at line 297 of file BandMatrix.h.


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