Cantera
2.0
|
Class MultiJac evaluates the Jacobian of a system of equations defined by a residual function supplied by an instance of class 'OneDim. More...
#include <MultiJac.h>
Public Member Functions | |
MultiJac (OneDim &r) | |
Constructor. | |
virtual | ~MultiJac () |
Destructor. Does nothing. | |
void | eval (doublereal *x0, doublereal *resid0, double rdt) |
Evaluate the Jacobian at x0. | |
doublereal | elapsedTime () const |
Elapsed CPU time spent computing the Jacobian. | |
int | nEvals () const |
Number of Jacobian evaluations. | |
int | age () const |
Number of times 'incrementAge' has been called since the last evaluation. | |
void | incrementAge () |
Increment the Jacobian age. | |
void | updateTransient (doublereal rdt, integer *mask) |
void | setAge (int age) |
Set the age. | |
vector_int & | transientMask () |
void | incrementDiagonal (int j, doublereal d) |
void | resize (size_t n, size_t kl, size_t ku, doublereal v=0.0) |
Resize the matrix problem. | |
void | bfill (doublereal v=0.0) |
Fill or zero the matrix. | |
doublereal & | operator() (size_t i, size_t j) |
Index into the (i,j) element. | |
doublereal | operator() (size_t i, size_t j) const |
Constant index into the (i,j) element. | |
doublereal & | value (size_t i, size_t j) |
Return a changeable reference to element (i,j). | |
doublereal | value (size_t i, size_t j) const |
Return the value of element (i,j). | |
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. | |
doublereal | _value (size_t i, size_t j) const |
Return the value of the (i,j) element for (i,j) within the bandwidth. | |
virtual size_t | nRows () const |
Returns the number of rows. | |
virtual size_t | nRowsAndStruct (size_t *const iStruct=0) const |
Return the size and structure of the matrix. | |
size_t | nColumns () const |
Number of columns. | |
size_t | nSubDiagonals () const |
Number of subdiagonals. | |
size_t | nSuperDiagonals () const |
Number of superdiagonals. | |
size_t | ldim () const |
Return the number of rows of storage needed for the band storage. | |
vector_int & | ipiv () |
Return a reference to the pivot vector. | |
virtual void | mult (const doublereal *b, doublereal *prod) const |
Multiply A*b and write result to prod. | |
virtual void | leftMult (const doublereal *const b, doublereal *const prod) const |
Multiply b*A and write result to prod. | |
int | factor () |
Perform an LU decomposition, the LAPACK routine DGBTRF is used. | |
int | solve (const doublereal *const b, doublereal *const x) |
Solve the matrix problem Ax = b. | |
int | solve (doublereal *b) |
Solve the matrix problem Ax = b. | |
virtual vector_fp::iterator | begin () |
Returns an iterator for the start of the band storage data. | |
vector_fp::const_iterator | begin () const |
Returns a const iterator for the start of the band storage data. | |
vector_fp::iterator | end () |
Returns an iterator for the end of the band storage data. | |
vector_fp::const_iterator | end () const |
Returns a const iterator for the end of the band storage data. | |
virtual void | zero () |
Zero the matrix. | |
virtual int | factorQR () |
Factors the A matrix using the QR algorithm, overwriting A. | |
virtual doublereal | rcondQR () |
Returns an estimate of the inverse of the condition number for the matrix. | |
virtual doublereal | rcond (doublereal a1norm) |
Returns an estimate of the inverse of the condition number for the matrix. | |
virtual void | useFactorAlgorithm (int fAlgorithm) |
Change the way the matrix is factored. | |
virtual int | factorAlgorithm () const |
Returns the factor algorithm used. | |
virtual doublereal | oneNorm () const |
Returns the one norm of the matrix. | |
virtual GeneralMatrix * | duplMyselfAsGeneralMatrix () const |
Duplicate this object as a GeneralMatrix pointer. | |
virtual bool | factored () const |
Report whether the current matrix has been factored. | |
virtual doublereal * | ptrColumn (size_t j) |
Return a pointer to the top of column j, column values are assumed to be contiguous in memory. | |
virtual doublereal *const * | colPts () |
Return a vector of const pointers to the columns. | |
virtual void | copyData (const GeneralMatrix &y) |
Copy the data from one array into another without doing any checking. | |
virtual void | clearFactorFlag () |
Clear the factored flag. | |
virtual size_t | checkRows (doublereal &valueSmall) const |
Check to see if we have any zero rows in the jacobian. | |
virtual size_t | checkColumns (doublereal &valueSmall) const |
Check to see if we have any zero columns in the jacobian. | |
Public Attributes | |
int | matrixType_ |
Matrix type. | |
Protected Attributes | |
OneDim * | m_resid |
Residual evaluator for this jacobian. | |
vector_fp | m_r1 |
doublereal | m_rtol |
doublereal | m_atol |
doublereal | m_elapsed |
vector_fp | m_ssdiag |
vector_int | m_mask |
int | m_nevals |
int | m_age |
size_t | m_size |
size_t | m_points |
vector_fp | data |
Matrix data. | |
vector_fp | ludata |
Factorized data. | |
bool | m_factored |
Boolean indicating whether the matrix is factored. | |
size_t | m_n |
Number of rows and columns of the matrix. | |
size_t | m_kl |
Number of subdiagonals of the matrix. | |
size_t | m_ku |
Number of super diagonals of the matrix. | |
doublereal | m_zero |
value of zero | |
vector_int | m_ipiv |
Pivot vector. | |
std::vector< doublereal * > | m_colPtrs |
Vector of column pointers. | |
vector_int | iwork_ |
Extra work array needed - size = n. | |
vector_fp | work_ |
Extra dp work array needed - size = 3n. | |
Class MultiJac evaluates the Jacobian of a system of equations defined by a residual function supplied by an instance of class 'OneDim.
' The residual function may consist of several linked 1D domains, with different variables in each domain.
Definition at line 25 of file MultiJac.h.
Constructor.
Definition at line 17 of file MultiJac.cpp.
References MultiJac::m_resid, OneDim::points(), and OneDim::size().
|
inlinevirtual |
Destructor. Does nothing.
Definition at line 36 of file MultiJac.h.
void eval | ( | doublereal * | x0, |
doublereal * | resid0, | ||
double | rdt | ||
) |
Evaluate the Jacobian at x0.
The unperturbed residual function is resid0, which must be supplied on input. The third parameter 'rdt' is the reciprocal of the time step. If zero, the steady-state Jacobian is evaluated.
The array of residual values at x0 is supplied as an input.
Definition at line 54 of file MultiJac.cpp.
References BandMatrix::bfill(), DATA_PTR, OneDim::eval(), OneDim::loc(), MultiJac::m_resid, Cantera::npos, OneDim::nVars(), and BandMatrix::value().
Referenced by MultiNewton::solve(), and OneDim::solve().
|
inline |
Elapsed CPU time spent computing the Jacobian.
Definition at line 49 of file MultiJac.h.
Referenced by OneDim::saveStats().
|
inline |
Number of Jacobian evaluations.
Definition at line 54 of file MultiJac.h.
Referenced by MultiNewton::dampStep(), OneDim::saveStats(), and MultiNewton::solve().
|
inline |
Number of times 'incrementAge' has been called since the last evaluation.
Definition at line 62 of file MultiJac.h.
Referenced by MultiNewton::dampStep(), MultiJac::setAge(), and MultiNewton::solve().
|
inline |
Increment the Jacobian age.
Definition at line 69 of file MultiJac.h.
Referenced by MultiNewton::solve().
|
inline |
Set the age.
Definition at line 76 of file MultiJac.h.
References MultiJac::age().
Referenced by Domain1D::needJacUpdate().
|
inherited |
Resize the matrix problem.
All data is lost
n | size of the square matrix |
kl | band size on the lower portion of the matrix |
ku | band size on the upper portion of the matrix |
v | initial value of all matrix components. |
Definition at line 105 of file BandMatrix.cpp.
References BandMatrix::data, BandMatrix::ludata, BandMatrix::m_colPtrs, BandMatrix::m_factored, BandMatrix::m_ipiv, BandMatrix::m_kl, BandMatrix::m_ku, and BandMatrix::m_n.
|
inherited |
Fill or zero the matrix.
v | Fill value, defaults to zero. |
Definition at line 122 of file BandMatrix.cpp.
References BandMatrix::data, and BandMatrix::m_factored.
Referenced by MultiJac::eval().
|
virtualinherited |
Index into the (i,j) element.
i | row |
j | column |
Returns a changeable reference to the matrix entry
Implements GeneralMatrix.
Definition at line 134 of file BandMatrix.cpp.
References BandMatrix::value().
|
virtualinherited |
Constant index into the (i,j) element.
i | row |
j | column |
Returns an unchangeable reference to the matrix entry
Implements GeneralMatrix.
Definition at line 139 of file BandMatrix.cpp.
References BandMatrix::value().
|
inherited |
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.
i | row |
j | column |
Definition at line 144 of file BandMatrix.cpp.
References BandMatrix::data, BandMatrix::index(), BandMatrix::m_factored, BandMatrix::m_kl, and BandMatrix::m_zero.
Referenced by BandMatrix::checkColumns(), BandMatrix::checkRows(), MultiJac::eval(), BandMatrix::oneNorm(), and BandMatrix::operator()().
|
inherited |
Return the value of element (i,j).
This method does not alter the array.
i | row |
j | column |
Definition at line 153 of file BandMatrix.cpp.
References BandMatrix::data, BandMatrix::index(), and BandMatrix::m_kl.
|
inherited |
Returns the location in the internal 1D array corresponding to the (i,j) element in the banded array.
i | row |
j | column |
Definition at line 161 of file BandMatrix.cpp.
References BandMatrix::m_kl, and BandMatrix::m_ku.
Referenced by BandMatrix::_value(), and BandMatrix::value().
|
inherited |
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.
i | row |
j | column |
Definition at line 169 of file BandMatrix.cpp.
References BandMatrix::data, and BandMatrix::index().
Referenced by BandMatrix::leftMult(), and BandMatrix::mult().
|
virtualinherited |
Returns the number of rows.
Implements GeneralMatrix.
Definition at line 175 of file BandMatrix.cpp.
References BandMatrix::m_n.
Referenced by BandMatrix::factor(), BandMatrix::mult(), and Cantera::operator<<().
|
virtualinherited |
Return the size and structure of the matrix.
This is inherited from GeneralMatrix
iStruct | OUTPUT Pointer to a vector of ints that describe the structure of the matrix. istruct[0] = kl istruct[1] = ku |
Implements GeneralMatrix.
Definition at line 181 of file BandMatrix.cpp.
References BandMatrix::m_kl, BandMatrix::m_ku, and BandMatrix::m_n.
|
inherited |
Number of columns.
Definition at line 191 of file BandMatrix.cpp.
References BandMatrix::m_n.
Referenced by BandMatrix::factor(), BandMatrix::leftMult(), Cantera::operator<<(), and BandMatrix::solve().
|
inherited |
Number of subdiagonals.
Definition at line 197 of file BandMatrix.cpp.
References BandMatrix::m_kl.
Referenced by BandMatrix::factor(), and BandMatrix::solve().
|
inherited |
Number of superdiagonals.
Definition at line 203 of file BandMatrix.cpp.
References BandMatrix::m_ku.
Referenced by BandMatrix::factor(), and BandMatrix::solve().
|
inherited |
Return the number of rows of storage needed for the band storage.
Definition at line 208 of file BandMatrix.cpp.
References BandMatrix::m_kl, and BandMatrix::m_ku.
Referenced by BandMatrix::factor(), and BandMatrix::solve().
|
inherited |
Return a reference to the pivot vector.
Definition at line 213 of file BandMatrix.cpp.
References BandMatrix::m_ipiv.
Referenced by BandMatrix::factor(), and BandMatrix::solve().
|
virtualinherited |
Multiply A*b and write result to prod.
b | Vector to do the rh multiplication |
prod | OUTPUT vector to receive the result |
Implements GeneralMatrix.
Definition at line 221 of file BandMatrix.cpp.
References BandMatrix::_value(), BandMatrix::m_kl, BandMatrix::m_ku, BandMatrix::m_n, and BandMatrix::nRows().
|
virtualinherited |
Multiply b*A and write result to prod.
b | Vector to do the lh multiplication |
prod | OUTPUT vector to receive the result |
Implements GeneralMatrix.
Definition at line 241 of file BandMatrix.cpp.
References BandMatrix::_value(), BandMatrix::m_kl, BandMatrix::m_ku, BandMatrix::m_n, and BandMatrix::nColumns().
|
virtualinherited |
Perform an LU decomposition, the LAPACK routine DGBTRF is used.
The factorization is saved in ludata.
Implements GeneralMatrix.
Definition at line 263 of file BandMatrix.cpp.
References BandMatrix::data, DATA_PTR, BandMatrix::ipiv(), BandMatrix::ldim(), BandMatrix::ludata, BandMatrix::m_factored, BandMatrix::nColumns(), BandMatrix::nRows(), BandMatrix::nSubDiagonals(), and BandMatrix::nSuperDiagonals().
Referenced by BandMatrix::factorQR(), and BandMatrix::solve().
|
inherited |
Solve the matrix problem Ax = b.
b | INPUT rhs of the problem |
x | OUTPUT solution to the problem |
Definition at line 282 of file BandMatrix.cpp.
References BandMatrix::m_n.
Referenced by MultiNewton::step().
|
virtualinherited |
Solve the matrix problem Ax = b.
b | INPUT rhs of the problem OUTPUT solution to the problem |
Implements GeneralMatrix.
Definition at line 288 of file BandMatrix.cpp.
References DATA_PTR, BandMatrix::factor(), BandMatrix::ipiv(), BandMatrix::ldim(), BandMatrix::ludata, BandMatrix::m_factored, BandMatrix::nColumns(), BandMatrix::nSubDiagonals(), and BandMatrix::nSuperDiagonals().
|
virtualinherited |
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 308 of file BandMatrix.cpp.
References BandMatrix::data, and BandMatrix::m_factored.
|
virtualinherited |
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 320 of file BandMatrix.cpp.
References BandMatrix::data.
|
inherited |
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 314 of file BandMatrix.cpp.
References BandMatrix::data, and BandMatrix::m_factored.
|
inherited |
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 325 of file BandMatrix.cpp.
References BandMatrix::data.
|
virtualinherited |
Zero the matrix.
Implements GeneralMatrix.
Definition at line 128 of file BandMatrix.cpp.
References BandMatrix::data, and BandMatrix::m_factored.
|
virtualinherited |
Factors the A matrix using the QR algorithm, overwriting A.
we set m_factored to 2 to indicate the matrix is now QR factored
Implements GeneralMatrix.
Definition at line 354 of file BandMatrix.cpp.
References BandMatrix::factor().
|
virtualinherited |
Returns an estimate of the inverse of the condition number for the matrix.
The matrix must have been previously factored using the QR algorithm
Implements GeneralMatrix.
Definition at line 367 of file BandMatrix.cpp.
References BandMatrix::oneNorm(), and BandMatrix::rcond().
|
virtualinherited |
Returns an estimate of the inverse of the condition number for the matrix.
The matrix must have been previously factored using the LU algorithm
a1norm | Norm of the matrix |
Implements GeneralMatrix.
Definition at line 381 of file BandMatrix.cpp.
References DATA_PTR, Cantera::int2str(), BandMatrix::iwork_, BandMatrix::ludata, BandMatrix::m_factored, BandMatrix::m_ipiv, BandMatrix::m_kl, BandMatrix::m_ku, BandMatrix::m_n, BandMatrix::work_, and Cantera::writelogf().
Referenced by BandMatrix::rcondQR().
|
virtualinherited |
Change the way the matrix is factored.
fAlgorithm | integer 0 LU factorization 1 QR factorization |
Implements GeneralMatrix.
Definition at line 418 of file BandMatrix.cpp.
|
virtualinherited |
Returns the factor algorithm used.
0 LU decomposition 1 QR decomposition
This routine will always return 0
Implements GeneralMatrix.
Definition at line 423 of file BandMatrix.cpp.
|
virtualinherited |
Returns the one norm of the matrix.
Implements GeneralMatrix.
Definition at line 429 of file BandMatrix.cpp.
References BandMatrix::m_colPtrs, BandMatrix::m_kl, BandMatrix::m_ku, BandMatrix::m_n, and BandMatrix::value().
Referenced by BandMatrix::rcondQR().
|
virtualinherited |
Duplicate this object as a GeneralMatrix pointer.
Implements GeneralMatrix.
Definition at line 499 of file BandMatrix.cpp.
References BandMatrix::BandMatrix().
|
virtualinherited |
Report whether the current matrix has been factored.
Implements GeneralMatrix.
Definition at line 505 of file BandMatrix.cpp.
References BandMatrix::m_factored.
|
virtualinherited |
Return a pointer to the top of column j, column values are assumed to be contiguous in memory.
The LAPACK bandstructure has column values which are contiguous in memory:
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];
j | Value of the column |
Implements GeneralMatrix.
Definition at line 516 of file BandMatrix.cpp.
References BandMatrix::m_colPtrs.
|
virtualinherited |
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.
Implements GeneralMatrix.
Definition at line 529 of file BandMatrix.cpp.
References BandMatrix::m_colPtrs.
|
virtualinherited |
Copy the data from one array into another without doing any checking.
This differs from the assignment operator as no resizing is done and memcpy() is used.
y | Array to be copied |
Implements GeneralMatrix.
Definition at line 539 of file BandMatrix.cpp.
References BandMatrix::data, DATA_PTR, BandMatrix::m_factored, BandMatrix::m_kl, BandMatrix::m_ku, BandMatrix::m_n, and GeneralMatrix::ptrColumn().
|
virtualinherited |
Clear the factored flag.
Implements GeneralMatrix.
Definition at line 550 of file BandMatrix.cpp.
References BandMatrix::m_factored.
|
virtualinherited |
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
valueSmall | OUTPUT value of the largest coefficient in the smallest row |
Implements GeneralMatrix.
Definition at line 447 of file BandMatrix.cpp.
References BandMatrix::m_kl, BandMatrix::m_ku, BandMatrix::m_n, Cantera::npos, and BandMatrix::value().
|
virtualinherited |
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
valueSmall | OUTPUT value of the largest coefficient in the smallest column |
Implements GeneralMatrix.
Definition at line 473 of file BandMatrix.cpp.
References BandMatrix::m_kl, BandMatrix::m_ku, BandMatrix::m_n, Cantera::npos, and BandMatrix::value().
|
protected |
Residual evaluator for this jacobian.
This is a pointer to the residual evaluator. This object isn't owned by this jacobian object.
Definition at line 93 of file MultiJac.h.
Referenced by MultiJac::eval(), and MultiJac::MultiJac().
|
protectedinherited |
Matrix data.
Definition at line 388 of file BandMatrix.h.
Referenced by BandMatrix::_value(), BandMatrix::BandMatrix(), BandMatrix::begin(), BandMatrix::bfill(), BandMatrix::copyData(), BandMatrix::end(), BandMatrix::factor(), BandMatrix::operator=(), BandMatrix::resize(), BandMatrix::value(), and BandMatrix::zero().
|
protectedinherited |
Factorized data.
Definition at line 391 of file BandMatrix.h.
Referenced by BandMatrix::BandMatrix(), BandMatrix::factor(), BandMatrix::operator=(), BandMatrix::rcond(), BandMatrix::resize(), and BandMatrix::solve().
|
protectedinherited |
Boolean indicating whether the matrix is factored.
Definition at line 394 of file BandMatrix.h.
Referenced by BandMatrix::BandMatrix(), BandMatrix::begin(), BandMatrix::bfill(), BandMatrix::clearFactorFlag(), BandMatrix::copyData(), BandMatrix::end(), BandMatrix::factor(), BandMatrix::factored(), BandMatrix::operator=(), BandMatrix::rcond(), BandMatrix::resize(), BandMatrix::solve(), BandMatrix::value(), and BandMatrix::zero().
|
protectedinherited |
Number of rows and columns of the matrix.
Definition at line 397 of file BandMatrix.h.
Referenced by BandMatrix::BandMatrix(), BandMatrix::checkColumns(), BandMatrix::checkRows(), BandMatrix::copyData(), BandMatrix::leftMult(), BandMatrix::mult(), BandMatrix::nColumns(), BandMatrix::nRows(), BandMatrix::nRowsAndStruct(), BandMatrix::oneNorm(), BandMatrix::operator=(), BandMatrix::rcond(), BandMatrix::resize(), and BandMatrix::solve().
|
protectedinherited |
Number of subdiagonals of the matrix.
Definition at line 400 of file BandMatrix.h.
Referenced by BandMatrix::BandMatrix(), BandMatrix::checkColumns(), BandMatrix::checkRows(), BandMatrix::copyData(), BandMatrix::index(), BandMatrix::ldim(), BandMatrix::leftMult(), BandMatrix::mult(), BandMatrix::nRowsAndStruct(), BandMatrix::nSubDiagonals(), BandMatrix::oneNorm(), BandMatrix::operator=(), BandMatrix::rcond(), BandMatrix::resize(), and BandMatrix::value().
|
protectedinherited |
Number of super diagonals of the matrix.
Definition at line 403 of file BandMatrix.h.
Referenced by BandMatrix::BandMatrix(), BandMatrix::checkColumns(), BandMatrix::checkRows(), BandMatrix::copyData(), BandMatrix::index(), BandMatrix::ldim(), BandMatrix::leftMult(), BandMatrix::mult(), BandMatrix::nRowsAndStruct(), BandMatrix::nSuperDiagonals(), BandMatrix::oneNorm(), BandMatrix::operator=(), BandMatrix::rcond(), and BandMatrix::resize().
|
protectedinherited |
|
protectedinherited |
Pivot vector.
Definition at line 409 of file BandMatrix.h.
Referenced by BandMatrix::BandMatrix(), BandMatrix::ipiv(), BandMatrix::operator=(), BandMatrix::rcond(), and BandMatrix::resize().
|
protectedinherited |
Vector of column pointers.
Definition at line 412 of file BandMatrix.h.
Referenced by BandMatrix::BandMatrix(), BandMatrix::colPts(), BandMatrix::oneNorm(), BandMatrix::operator=(), BandMatrix::ptrColumn(), and BandMatrix::resize().
|
protectedinherited |
Extra work array needed - size = n.
Definition at line 415 of file BandMatrix.h.
Referenced by BandMatrix::rcond().
|
protectedinherited |
Extra dp work array needed - size = 3n.
Definition at line 418 of file BandMatrix.h.
Referenced by BandMatrix::rcond().
|
inherited |
Matrix type.
0 Square 1 Banded
Definition at line 236 of file GeneralMatrix.h.
Referenced by NonlinearSolver::beuler_jac(), NonlinearSolver::doAffineNewtonSolve(), GeneralMatrix::operator=(), and NonlinearSolver::scaleMatrix().