Cantera  2.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 /*
8  * Copyright 2004 Sandia Corporation. Under the terms of Contract
9  * DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government
10  * retains certain rights in this software.
11  * See file License.txt for licensing information.
12  */
13 
14 #ifndef CT_GENERALMATRIX_H
15 #define CT_GENERALMATRIX_H
16 
17 #include "cantera/base/ct_defs.h"
18 
19 namespace Cantera
20 {
21 
22 //! Generic matrix
24 {
25 public:
26  //! Base Constructor
27  /*!
28  * @param matType Matrix type
29  * 0 full
30  * 1 banded
31  */
32  GeneralMatrix(int matType);
33 
34  //! Copy Constructor
35  /*!
36  * @param right Object to be copied
37  */
38  GeneralMatrix(const GeneralMatrix& right);
39 
40  //! Assignment operator
41  /*!
42  * @param right Object to be copied
43  */
44  GeneralMatrix& operator=(const GeneralMatrix& right);
45 
46  //! Destructor. Does nothing.
47  virtual ~GeneralMatrix();
48 
49  //! Duplicator member function
50  /*!
51  * This function will duplicate the matrix given a generic GeneralMatrix pointer
52  *
53  * @return Returns a pointer to the malloced object
54  */
55  virtual GeneralMatrix* duplMyselfAsGeneralMatrix() const = 0;
56 
57  //! Zero the matrix elements
58  virtual void zero() = 0;
59 
60  //! Multiply A*b and write result to prod.
61  /*!
62  * @param b Vector to do the rh multiplication
63  * @param prod OUTPUT vector to receive the result
64  */
65  virtual void mult(const doublereal* b, doublereal* prod) const = 0;
66 
67  //! Multiply b*A and write result to prod.
68  /*!
69  * @param b Vector to do the lh multiplication
70  * @param prod OUTPUT vector to receive the result
71  */
72  virtual void leftMult(const doublereal* const b, doublereal* const prod) const = 0;
73 
74  //! Factors the A matrix, overwriting A.
75  /*
76  * We flip m_factored boolean to indicate that the matrix is now A-1.
77  */
78  virtual int factor() = 0;
79 
80  //! Factors the A matrix using the QR algorithm, overwriting A
81  /*!
82  * we set m_factored to 2 to indicate the matrix is now QR factored
83  *
84  * @return Returns the info variable from lapack
85  */
86  virtual int factorQR() = 0;
87 
88  //! Returns an estimate of the inverse of the condition number for the matrix
89  /*!
90  * The matrix must have been previously factored using the QR algorithm
91  *
92  * @return returns the inverse of the condition number
93  */
94  virtual doublereal rcondQR() = 0;
95 
96  //! Returns an estimate of the inverse of the condition number for the matrix
97  /*!
98  * The matrix must have been previously factored using the LU algorithm
99  *
100  * @param a1norm Norm of the matrix
101  *
102  * @return returns the inverse of the condition number
103  */
104  virtual doublereal rcond(doublereal a1norm) = 0;
105 
106  //! Change the way the matrix is factored
107  /*!
108  * @param fAlgorithm integer
109  * 0 LU factorization
110  * 1 QR factorization
111  */
112  virtual void useFactorAlgorithm(int fAlgorithm) = 0;
113 
114  //! Return the factor algorithm used
115  /*!
116  *
117  */
118  virtual int factorAlgorithm() const = 0;
119 
120  //! Calculate the one norm of the matrix
121  /*!
122  * Returns the one norm of the matrix
123  */
124  virtual doublereal oneNorm() const = 0;
125 
126 
127  //! Return the number of rows in the matrix
128  virtual size_t nRows() const = 0;
129 
130 
131  //! Return the size and structure of the matrix
132  /*!
133  * This is inherited from GeneralMatrix
134  *
135  * @param iStruct OUTPUT Pointer to a vector of ints that describe the structure of the matrix.
136  *
137  * @return returns the number of rows and columns in the matrix.
138  */
139  virtual size_t nRowsAndStruct(size_t* const iStruct = 0) const = 0;
140 
141  //! clear the factored flag
142  virtual void clearFactorFlag() = 0;
143 
144  //! Solves the Ax = b system returning x in the b spot.
145  /*!
146  * @param b Vector for the rhs of the equation system
147  */
148  virtual int solve(doublereal* b) = 0;
149 
150  //! true if the current factorization is up to date with the matrix
151  virtual bool factored() const = 0;
152 
153  //! Return a pointer to the top of column j, columns are assumed to be contiguous in memory
154  /*!
155  * @param j Value of the column
156  *
157  * @return Returns a pointer to the top of the column
158  */
159  virtual doublereal* ptrColumn(size_t j) = 0;
160 
161  //! Index into the (i,j) element
162  /*!
163  * @param i row
164  * @param j column
165  *
166  * Returns a changeable reference to the matrix entry
167  */
168  virtual doublereal& operator()(size_t i, size_t j) = 0;
169 
170 
171  //! Constant Index into the (i,j) element
172  /*!
173  * @param i row
174  * @param j column
175  *
176  * Returns an unchangeable reference to the matrix entry
177  */
178  virtual doublereal operator()(size_t i, size_t j) const = 0;
179 
180  //! Copy the data from one array into another without doing any checking
181  /*!
182  * This differs from the assignment operator as no resizing is done and memcpy() is used.
183  * @param y Array to be copied
184  */
185  virtual void copyData(const GeneralMatrix& y) = 0;
186 
187  //! Return an iterator pointing to the first element
188  /*!
189  * We might drop this later
190  */
191  virtual vector_fp::iterator begin() = 0;
192 
193  //! Return a const iterator pointing to the first element
194  /*!
195  * We might drop this later
196  */
197  virtual vector_fp::const_iterator begin() const = 0;
198 
199  //! Return a vector of const pointers to the columns
200  /*!
201  * Note the value of the pointers are protected by their being const.
202  * However, the value of the matrix is open to being changed.
203  *
204  * @return returns a vector of pointers to the top of the columns
205  * of the matrices.
206  */
207  virtual doublereal* const* colPts() = 0;
208 
209  //! Check to see if we have any zero rows in the jacobian
210  /*!
211  * This utility routine checks to see if any rows are zero.
212  * The smallest row is returned along with the largest coefficient in that row
213  *
214  * @param valueSmall OUTPUT value of the largest coefficient in the smallest row
215  *
216  * @return index of the row that is most nearly zero
217  */
218  virtual size_t checkRows(doublereal& valueSmall) const = 0;
219 
220  //! Check to see if we have any zero columns in the jacobian
221  /*!
222  * This utility routine checks to see if any columns are zero.
223  * The smallest column is returned along with the largest coefficient in that column
224  *
225  * @param valueSmall OUTPUT value of the largest coefficient in the smallest column
226  *
227  * @return index of the column that is most nearly zero
228  */
229  virtual size_t checkColumns(doublereal& valueSmall) const = 0;
230 
231  //! Matrix type
232  /*!
233  * 0 Square
234  * 1 Banded
235  */
237 
238 };
239 }
240 #endif