Cantera  2.4.0
GeneralMatrix.h
Go to the documentation of this file.
1 /**
2  * @file GeneralMatrix.h
3  * Declarations for the class GeneralMatrix which is a virtual base class for matrices handled by solvers
4  * (see class \ref numerics and \link Cantera::GeneralMatrix GeneralMatrix\endlink).
5  */
6 
7 // This file is part of Cantera. See License.txt in the top-level directory or
8 // at http://www.cantera.org/license.txt for license and copyright information.
9 
10 #ifndef CT_GENERALMATRIX_H
11 #define CT_GENERALMATRIX_H
12 
13 #include "cantera/base/ct_defs.h"
15 #include "cantera/base/global.h"
16 
17 namespace Cantera
18 {
19 
20 //! Generic matrix
22 {
23 public:
24  //! Base Constructor
25  GeneralMatrix() : m_factored(false) {}
26 
27  virtual ~GeneralMatrix() {}
28 
29  //! Zero the matrix elements
30  virtual void zero() = 0;
31 
32  //! Multiply A*b and write result to prod.
33  /*!
34  * @param b Vector to do the rh multiplication
35  * @param prod OUTPUT vector to receive the result
36  */
37  virtual void mult(const doublereal* b, doublereal* prod) const = 0;
38 
39  //! Multiply b*A and write result to prod.
40  /*!
41  * @param b Vector to do the lh multiplication
42  * @param prod OUTPUT vector to receive the result
43  */
44  virtual void leftMult(const doublereal* const b, doublereal* const prod) const = 0;
45 
46  //! Factors the A matrix, overwriting A.
47  /*!
48  * We flip m_factored boolean to indicate that the matrix is now A-1.
49  */
50  virtual int factor() = 0;
51 
52  //! Factors the A matrix using the QR algorithm, overwriting A
53  /*!
54  * we set m_factored to 2 to indicate the matrix is now QR factored
55  *
56  * @returns the info variable from LAPACK
57  */
58  virtual int factorQR() {
59  throw NotImplementedError("GeneralMatrix::factorQR");
60  }
61 
62  //! Returns an estimate of the inverse of the condition number for the matrix
63  /*!
64  * The matrix must have been previously factored using the QR algorithm
65  *
66  * @returns the inverse of the condition number
67  */
68  virtual doublereal rcondQR() {
69  throw NotImplementedError("GeneralMatrix::rcondQR");
70  }
71 
72  //! Returns an estimate of the inverse of the condition number for the matrix
73  /*!
74  * The matrix must have been previously factored using the LU algorithm
75  *
76  * @param a1norm Norm of the matrix
77  * @returns the inverse of the condition number
78  */
79  virtual doublereal rcond(doublereal a1norm) = 0;
80 
81  //! Change the way the matrix is factored
82  /*!
83  * @param fAlgorithm integer
84  * 0 LU factorization
85  * 1 QR factorization
86  */
87  virtual void useFactorAlgorithm(int fAlgorithm) {
88  throw NotImplementedError("GeneralMatrix::useFactorAlgorithm");
89  };
90 
91  //! Return the factor algorithm used
92  virtual int factorAlgorithm() const = 0;
93 
94  //! Calculate the one norm of the matrix
95  virtual doublereal oneNorm() const = 0;
96 
97  //! Return the number of rows in the matrix
98  virtual size_t nRows() const = 0;
99 
100  //! clear the factored flag
101  virtual void clearFactorFlag() {
102  m_factored = 0;
103  };
104 
105  //! Solves the Ax = b system returning x in the b spot.
106  /*!
107  * @param b Vector for the RHS of the equation system
108  * @param nrhs Number of right-hand sides to solve, default 1
109  * @param ldb Leading dimension of the right-hand side array. Defaults to
110  * nRows()
111  */
112  virtual int solve(doublereal* b, size_t nrhs=1, size_t ldb=0) = 0;
113 
114  //! true if the current factorization is up to date with the matrix
115  virtual bool factored() const {
116  return (m_factored != 0);
117  }
118 
119  //! Return a pointer to the top of column j, columns are assumed to be
120  //! contiguous in memory
121  /*!
122  * @param j Value of the column
123  * @returns a pointer to the top of the column
124  */
125  virtual doublereal* ptrColumn(size_t j) = 0;
126 
127  //! Index into the (i,j) element
128  /*!
129  * @param i row
130  * @param j column
131  * @returns a changeable reference to the matrix entry
132  */
133  virtual doublereal& operator()(size_t i, size_t j) = 0;
134 
135  //! Constant Index into the (i,j) element
136  /*!
137  * @param i row
138  * @param j column
139  * @returns an unchangeable reference to the matrix entry
140  */
141  virtual doublereal operator()(size_t i, size_t j) const = 0;
142 
143  //! Return an iterator pointing to the first element
144  /*!
145  * We might drop this later
146  */
147  virtual vector_fp::iterator begin() = 0;
148 
149  //! Return a const iterator pointing to the first element
150  /*!
151  * We might drop this later
152  */
153  virtual vector_fp::const_iterator begin() const = 0;
154 
155  //! Return a vector of const pointers to the columns
156  /*!
157  * Note the value of the pointers are protected by their being const.
158  * However, the value of the matrix is open to being changed.
159  *
160  * @returns a vector of pointers to the top of the columns of the matrices.
161  */
162  virtual doublereal* const* colPts() = 0;
163 
164  //! Check to see if we have any zero rows in the Jacobian
165  /*!
166  * This utility routine checks to see if any rows are zero. The smallest row
167  * is returned along with the largest coefficient in that row
168  *
169  * @param valueSmall OUTPUT value of the largest coefficient in the smallest row
170  * @return index of the row that is most nearly zero
171  */
172  virtual size_t checkRows(doublereal& valueSmall) const = 0;
173 
174  //! Check to see if we have any zero columns in the Jacobian
175  /*!
176  * This utility routine checks to see if any columns are zero. The smallest
177  * column is returned along with the largest coefficient in that column
178  *
179  * @param valueSmall OUTPUT value of the largest coefficient in the smallest column
180  * @return index of the column that is most nearly zero
181  */
182  virtual size_t checkColumns(doublereal& valueSmall) const = 0;
183 
184 protected:
185  //! Indicates whether the matrix is factored. 0 for unfactored; Non-zero
186  //! values indicate a particular factorization (LU=1, QR=2).
188 };
189 
190 }
191 #endif
virtual doublereal * ptrColumn(size_t j)=0
Return a pointer to the top of column j, columns are assumed to be contiguous in memory.
An error indicating that an unimplemented function has been called.
Definition: ctexceptions.h:187
int m_factored
Indicates whether the matrix is factored.
virtual doublereal oneNorm() const =0
Calculate the one norm of the matrix.
virtual void clearFactorFlag()
clear the factored flag
virtual size_t nRows() const =0
Return the number of rows in the matrix.
virtual doublereal rcond(doublereal a1norm)=0
Returns an estimate of the inverse of the condition number for the matrix.
GeneralMatrix()
Base Constructor.
Definition: GeneralMatrix.h:25
virtual vector_fp::iterator begin()=0
Return an iterator pointing to the first element.
This file contains definitions of terms that are used in internal routines and are unlikely to need m...
Generic matrix.
Definition: GeneralMatrix.h:21
virtual int factorQR()
Factors the A matrix using the QR algorithm, overwriting A.
Definition: GeneralMatrix.h:58
virtual int factorAlgorithm() const =0
Return the factor algorithm used.
This file contains definitions for utility functions and text for modules, inputfiles, logs, textlogs, (see Input File Handling, Diagnostic Output, and Writing messages to the screen).
virtual int factor()=0
Factors the A matrix, overwriting A.
virtual void zero()=0
Zero the matrix elements.
virtual doublereal & operator()(size_t i, size_t j)=0
Index into the (i,j) element.
virtual bool factored() const
true if the current factorization is up to date with the matrix
virtual size_t checkRows(doublereal &valueSmall) const =0
Check to see if we have any zero rows in the Jacobian.
virtual void leftMult(const doublereal *const b, doublereal *const prod) const =0
Multiply b*A and write result to prod.
virtual void mult(const doublereal *b, doublereal *prod) const =0
Multiply A*b and write result to prod.
virtual int solve(doublereal *b, size_t nrhs=1, size_t ldb=0)=0
Solves the Ax = b system returning x in the b spot.
virtual size_t checkColumns(doublereal &valueSmall) const =0
Check to see if we have any zero columns in the Jacobian.
virtual doublereal *const * colPts()=0
Return a vector of const pointers to the columns.
Namespace for the Cantera kernel.
Definition: AnyMap.cpp:8
virtual doublereal rcondQR()
Returns an estimate of the inverse of the condition number for the matrix.
Definition: GeneralMatrix.h:68
Definitions for the classes that are thrown when Cantera experiences an error condition (also contain...
virtual void useFactorAlgorithm(int fAlgorithm)
Change the way the matrix is factored.
Definition: GeneralMatrix.h:87