queso-0.57.1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
GslMatrix.h
Go to the documentation of this file.
1 //-----------------------------------------------------------------------bl-
2 //--------------------------------------------------------------------------
3 //
4 // QUESO - a library to support the Quantification of Uncertainty
5 // for Estimation, Simulation and Optimization
6 //
7 // Copyright (C) 2008-2017 The PECOS Development Team
8 //
9 // This library is free software; you can redistribute it and/or
10 // modify it under the terms of the Version 2.1 GNU Lesser General
11 // Public License as published by the Free Software Foundation.
12 //
13 // This library is distributed in the hope that it will be useful,
14 // but WITHOUT ANY WARRANTY; without even the implied warranty of
15 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
16 // Lesser General Public License for more details.
17 //
18 // You should have received a copy of the GNU Lesser General Public
19 // License along with this library; if not, write to the Free Software
20 // Foundation, Inc. 51 Franklin Street, Fifth Floor,
21 // Boston, MA 02110-1301 USA
22 //
23 //-----------------------------------------------------------------------el-
24 
25 #ifndef UQ_GSL_MATRIX_H
26 #define UQ_GSL_MATRIX_H
27 
32 #include <queso/Defines.h>
33 #include <queso/GslVector.h>
34 #include <queso/Matrix.h>
35 
36 #include <gsl/gsl_matrix.h>
37 #include <gsl/gsl_permutation.h>
38 
39 namespace QUESO {
40 
49 class GslMatrix : public Matrix
50 {
51 public:
53 
54 
57  const Map& map,
58  unsigned int numCols);
59 
61  GslMatrix(const BaseEnvironment& env,
62  const Map& map,
63  double diagValue); // MATLAB eye
64 
66  GslMatrix(const GslVector& v,
67  double diagValue); // MATLAB eye
68 
70 
71  GslMatrix(const GslVector& v); // MATLAB diag
72 
74 
75  GslMatrix(const GslMatrix& B);
76 
78  ~GslMatrix();
80 
81 
83 
84  GslMatrix& operator= (const GslMatrix& rhs);
86 
88  GslMatrix& operator*=(double a);
89 
91  GslMatrix& operator/=(double a);
92 
94  GslMatrix& operator+=(const GslMatrix& rhs);
95 
97  GslMatrix& operator-=(const GslMatrix& rhs);
99 
100 
102 
103  double& operator()(unsigned int i, unsigned int j)
105  {
106  // Changing the matrix changes its decomposition(s). We'll
107  // invalidate for now; recalculate later if needed.
108  this->resetLU();
109  queso_require_less_msg(i, m_mat->size1, "i is too large");
110  queso_require_less_msg(j, m_mat->size2, "j is too large");
111  return *gsl_matrix_ptr(m_mat,i,j);
112  }
113 
115  const double& operator()(unsigned int i, unsigned int j) const
116  {
117  queso_require_less_msg(i, m_mat->size1, "i is too large");
118  queso_require_less_msg(j, m_mat->size2, "j is too large");
119  return *gsl_matrix_ptr(m_mat,i,j);
120  }
121 
122 
124 
126 
127  unsigned int numRowsLocal () const;
129 
131  unsigned int numRowsGlobal () const;
132 
134  unsigned int numCols () const;
135 
137  double max () const;
138 
140 
141  unsigned int rank (double absoluteZeroThreshold, double relativeZeroThreshold) const;
142 
144  GslMatrix transpose () const;
145 
147  GslMatrix inverse () const;
148 
149 
151  double determinant () const;
152 
154  double lnDeterminant () const;
155 
157 
159 
160  double normFrob () const;
162 
164  double normMax () const;
166 
167 
168 
169 
170 
172 
173 
175 
176  int chol ();
177 
179  int svd (GslMatrix& matU, GslVector& vecS, GslMatrix& matVt) const;
180 
182 
183  const GslMatrix& svdMatU () const;
184 
186 
187  const GslMatrix& svdMatV () const;
188 
190  int svdSolve (const GslVector& rhsVec, GslVector& solVec) const;
191 
193  int svdSolve (const GslMatrix& rhsMat, GslMatrix& solMat) const;
194 
195 
196 
198  GslVector multiply (const GslVector& x) const;
199 
201  //the preallocated vector \c y.
202  void multiply (const GslVector& x, GslVector& y) const;
203 
205  GslMatrix multiply (const GslMatrix& X) const;
206 
208  //the preallocated zeroed matrix \c Y.
209  void multiply (const GslMatrix& X, GslMatrix& Y) const;
210 
212 
213  GslVector invertMultiply (const GslVector& b) const;
214 
216 
218  void invertMultiply (const GslVector& b, GslVector& x) const;
219 
221 
222  GslMatrix invertMultiply (const GslMatrix& B) const;
223 
225 
227  void invertMultiply (const GslMatrix& B, GslMatrix& X) const;
228 
230 
232  GslVector invertMultiplyForceLU (const GslVector& b) const;
233 
235 
236  void invertMultiplyForceLU (const GslVector& b, GslVector& x) const;
237 
239  void getColumn (const unsigned int column_num, GslVector& column) const;
240 
242  GslVector getColumn (const unsigned int column_num) const;
243 
245  void setColumn (const unsigned int column_num, const GslVector& column);
246 
248  void getRow (const unsigned int row_num, GslVector& row) const;
249 
251  GslVector getRow (const unsigned int row_num) const;
252 
254  void setRow (const unsigned int row_num, const GslVector& row);
255 
257  void eigen (GslVector& eigenValues, GslMatrix* eigenVectors) const;
258 
260  void largestEigen (double& eigenValue, GslVector& eigenVector) const;
261 
263  void smallestEigen (double& eigenValue, GslVector& eigenVector) const;
264 
265 
267 
269 
270 
272  void cwSet (double value);
273 
275  void cwSet (unsigned int rowId, unsigned int colId, const GslMatrix& mat);
276 
277  void cwExtract (unsigned int rowId, unsigned int colId, GslMatrix& mat) const;
278 
280 
282  void zeroLower (bool includeDiagonal = false);
283 
285 
287  void zeroUpper (bool includeDiagonal = false);
288 
290 
291  void filterSmallValues (double thresholdValue);
292 
294 
295  void filterLargeValues (double thresholdValue);
296 
298  void fillWithTranspose (unsigned int rowId,
299  unsigned int colId,
300  const GslMatrix& mat,
301  bool checkForExactNumRowsMatching,
302  bool checkForExactNumColsMatching);
303 
305  void fillWithBlocksDiagonally (unsigned int rowId,
306  unsigned int colId,
307  const std::vector<const GslMatrix* >& matrices,
308  bool checkForExactNumRowsMatching,
309  bool checkForExactNumColsMatching);
310 
312  void fillWithBlocksDiagonally (unsigned int rowId,
313  unsigned int colId,
314  const std::vector< GslMatrix* >& matrices,
315  bool checkForExactNumRowsMatching,
316  bool checkForExactNumColsMatching);
317 
319  void fillWithBlocksHorizontally(unsigned int rowId,
320  unsigned int colId,
321  const std::vector<const GslMatrix* >& matrices,
322  bool checkForExactNumRowsMatching,
323  bool checkForExactNumColsMatching);
324 
326  void fillWithBlocksHorizontally(unsigned int rowId,
327  unsigned int colId,
328  const std::vector< GslMatrix* >& matrices,
329  bool checkForExactNumRowsMatching,
330  bool checkForExactNumColsMatching);
331 
333  void fillWithBlocksVertically (unsigned int rowId,
334  unsigned int colId,
335  const std::vector<const GslMatrix* >& matrices,
336  bool checkForExactNumRowsMatching,
337  bool checkForExactNumColsMatching);
338 
340  void fillWithBlocksVertically (unsigned int rowId,
341  unsigned int colId,
342  const std::vector< GslMatrix* >& matrices,
343  bool checkForExactNumRowsMatching,
344  bool checkForExactNumColsMatching);
345 
347  void fillWithTensorProduct (unsigned int rowId,
348  unsigned int colId,
349  const GslMatrix& mat1,
350  const GslMatrix& mat2,
351  bool checkForExactNumRowsMatching,
352  bool checkForExactNumColsMatching);
353 
355  void fillWithTensorProduct (unsigned int rowId,
356  unsigned int colId,
357  const GslMatrix& mat1,
358  const GslVector& vec2,
359  bool checkForExactNumRowsMatching,
360  bool checkForExactNumColsMatching);
362 
363 
365 
366  gsl_matrix* data ();
368 
369  void mpiSum (const MpiComm& comm, GslMatrix& M_global) const;
370 
371  void matlabLinearInterpExtrap (const GslVector& x1Vec, const GslMatrix& y1Mat, const GslVector& x2Vec);
373 
374 
375 
377 
378 
380  void print (std::ostream& os) const;
381 
383  void subWriteContents (const std::string& varNamePrefix,
384  const std::string& fileName,
385  const std::string& fileType,
386  const std::set<unsigned int>& allowedSubEnvIds) const;
387 
389  void subReadContents (const std::string& fileName,
390  const std::string& fileType,
391  const std::set<unsigned int>& allowedSubEnvIds);
393 
394 private:
396 
397  GslMatrix();
398 
400  void copy (const GslMatrix& src);
401 
403  void resetLU ();
404 
406  int internalSvd () const;
407 
409  gsl_matrix* m_mat;
410 
412  mutable gsl_matrix* m_LU;
413 
416 
418  mutable Map* m_svdColMap;
419 
421 
425 
427 
432 
434 
438 
440 
444 
446  mutable double m_determinant;
447 
449  mutable double m_lnDeterminant;
450 
452 
455  mutable gsl_permutation* m_permutation;
456 
458 
461  mutable int m_signum;
462 
464  mutable bool m_isSingular;
465 };
466 
467 GslMatrix operator* (double a, const GslMatrix& mat);
468 GslVector operator* (const GslMatrix& mat, const GslVector& vec);
469 GslMatrix operator* (const GslMatrix& m1, const GslMatrix& m2 );
470 GslMatrix operator+ (const GslMatrix& m1, const GslMatrix& m2 );
471 GslMatrix operator- (const GslMatrix& m1, const GslMatrix& m2 );
472 GslMatrix matrixProduct (const GslVector& v1, const GslVector& v2 );
473 
474 // Row \c i of the returned matrix is equal to row \c i of \c mat multiplied by element \c i of \c v.
475 GslMatrix leftDiagScaling (const GslVector& vec, const GslMatrix& mat);
476 
477 // Column \c j of the returned matrix is equal to column \c j of \c mat multiplied by element \c j of \c v.
478 GslMatrix rightDiagScaling(const GslMatrix& mat, const GslVector& vec);
479 std::ostream& operator<< (std::ostream& os, const GslMatrix& obj);
480 
481 } // End namespace QUESO
482 
483 #endif // UQ_GSL_MATRIX_H
gsl_permutation * m_permutation
The permutation matrix of a LU decomposition.
Definition: GslMatrix.h:455
GslMatrix * m_svdVmat
m_svdVmat stores the N-by-N orthogonal square matrix V after the singular value decomposition of a ma...
Definition: GslMatrix.h:437
GslVector invertMultiply(const GslVector &b) const
This function calculates the inverse of this matrix and multiplies it with vector b...
Definition: GslMatrix.C:1222
unsigned int numRowsGlobal() const
Returns the global row dimension of this matrix.
Definition: GslMatrix.C:281
GslVector invertMultiplyForceLU(const GslVector &b) const
This function calculates the inverse of this matrix and multiplies it with vector b...
Definition: GslMatrix.C:1373
void fillWithBlocksDiagonally(unsigned int rowId, unsigned int colId, const std::vector< const GslMatrix * > &matrices, bool checkForExactNumRowsMatching, bool checkForExactNumColsMatching)
This function fills this matrix diagonally with const block matrices.
Definition: GslMatrix.C:753
A class for partitioning vectors and matrices.
Definition: Map.h:49
const BaseEnvironment & env() const
Definition: Matrix.C:47
void matlabLinearInterpExtrap(const GslVector &x1Vec, const GslMatrix &y1Mat, const GslVector &x2Vec)
Definition: GslMatrix.C:1741
This (virtual) class sets up the environment underlying the use of the QUESO library by an executable...
Definition: Environment.h:198
GslMatrix & operator/=(double a)
Stores in this the coordinate-wise division of this by a.
Definition: GslMatrix.C:187
int m_signum
m_signum stores the sign of the permutation of the LU decomposition PA = LU.
Definition: GslMatrix.h:461
void zeroLower(bool includeDiagonal=false)
This function sets all the entries bellow the main diagonal of this matrix to zero.
Definition: GslMatrix.C:591
GslVector multiply(const GslVector &x) const
This function multiplies this matrix by vector x and returns the resulting vector.
Definition: GslMatrix.C:1155
unsigned int numRowsLocal() const
Returns the local row dimension of this matrix.
Definition: GslMatrix.C:275
void subReadContents(const std::string &fileName, const std::string &fileType, const std::set< unsigned int > &allowedSubEnvIds)
Read contents of subenvironment from file fileName.
Definition: GslMatrix.C:1841
GslMatrix * m_svdUmat
m_svdUmat stores the M-by-N orthogonal matrix U after the singular value decomposition of a matrix...
Definition: GslMatrix.h:424
GslMatrix()
Default Constructor.
std::ostream & operator<<(std::ostream &os, const SequenceStatisticalOptions &obj)
void filterSmallValues(double thresholdValue)
This function sets to zero (filters) all entries of this matrix which are smaller than thresholdValue...
Definition: GslMatrix.C:647
const GslMatrix & svdMatV() const
This function calls private member GslMatrix::internalSvd() to set a N-by-N orthogonal square matrix ...
Definition: GslMatrix.C:522
GslMatrix & operator=(const GslMatrix &rhs)
Copies values from matrix rhs to this.
Definition: GslMatrix.C:169
void filterLargeValues(double thresholdValue)
This function sets to zero (filters) all entries of this matrix which are greater than thresholdValue...
Definition: GslMatrix.C:670
void setRow(const unsigned int row_num, const GslVector &row)
This function copies vector row into the row_num-th column of this matrix.
Definition: GslMatrix.C:1664
gsl_matrix * data()
Returns this matrix.
Definition: GslMatrix.C:1766
GslMatrix operator-(const GslMatrix &m1, const GslMatrix &m2)
Definition: GslMatrix.C:2029
void resetLU()
In this function resets the LU decomposition of this matrix, as well as deletes the private member po...
Definition: GslMatrix.C:232
void smallestEigen(double &eigenValue, GslVector &eigenVector) const
This function finds smallest eigenvalue, namely eigenValue, of this matrix and its corresponding eige...
Definition: GslMatrix.C:1515
GslMatrix inverse() const
This function calculated the inverse of this matrix (square).
Definition: GslMatrix.C:713
Class for matrix operations using GSL library.
Definition: GslMatrix.h:49
bool m_isSingular
Indicates whether or not this matrix is singular.
Definition: GslMatrix.h:464
double determinant() const
Calculates the determinant of this matrix.
Definition: GslMatrix.C:1038
int svdSolve(const GslVector &rhsVec, GslVector &solVec) const
This function solves the system A x = b using the singular value decomposition (U, S, V) of A which must have been computed previously with GslMatrix::svd (x=solVec, b=rhsVec).
Definition: GslMatrix.C:456
GslMatrix * m_svdVTmat
m_svdVmatT stores the transpose of N-by-N orthogonal square matrix V, namely V^T, after the singular ...
Definition: GslMatrix.h:443
int chol()
Computes Cholesky factorization of a real symmetric positive definite matrix this.
Definition: GslMatrix.C:408
void getRow(const unsigned int row_num, GslVector &row) const
This function gets the row_num-th column of this matrix and stores it into vector row...
Definition: GslMatrix.C:1612
double & operator()(unsigned int i, unsigned int j)
Element access method (non-const).
Definition: GslMatrix.h:104
void largestEigen(double &eigenValue, GslVector &eigenVector) const
This function finds largest eigenvalue, namely eigenValue, of this matrix and its corresponding eigen...
Definition: GslMatrix.C:1446
Map * m_svdColMap
Mapping for matrices involved in the singular value decomposition (svd) routine.
Definition: GslMatrix.h:418
double m_lnDeterminant
The natural logarithm of the determinant of this matrix.
Definition: GslMatrix.h:449
unsigned int rank(double absoluteZeroThreshold, double relativeZeroThreshold) const
This function returns the number of singular values of this matrix (rank).
Definition: GslMatrix.C:1120
void getColumn(const unsigned int column_num, GslVector &column) const
This function gets the column_num-th column of this matrix and stores it into vector column...
Definition: GslMatrix.C:1586
Class for matrix operations (virtual).
Definition: Matrix.h:46
void fillWithBlocksHorizontally(unsigned int rowId, unsigned int colId, const std::vector< const GslMatrix * > &matrices, bool checkForExactNumRowsMatching, bool checkForExactNumColsMatching)
This function fills this matrix horizontally with const block matrices.
Definition: GslMatrix.C:825
void copy(const GslMatrix &src)
In this function this matrix receives a copy of matrix src.
Definition: GslMatrix.C:220
double normMax() const
Returns the Frobenius norm of this matrix.
Definition: GslMatrix.C:311
GslMatrix operator*(double a, const GslMatrix &mat)
Definition: GslMatrix.C:1982
void eigen(GslVector &eigenValues, GslMatrix *eigenVectors) const
This function computes the eigenvalues of a real symmetric matrix.
Definition: GslMatrix.C:1420
GslMatrix transpose() const
This function calculated the transpose of this matrix (square).
Definition: GslMatrix.C:693
Class for vector operations using GSL library.
Definition: GslVector.h:49
void cwSet(double value)
Component-wise set all values to this with value.
Definition: GslMatrix.C:347
const GslMatrix & svdMatU() const
This function calls private member GslMatrix::internalSvd() to set a M-by-N orthogonal matrix U of th...
Definition: GslMatrix.C:512
void mpiSum(const MpiComm &comm, GslMatrix &M_global) const
Definition: GslMatrix.C:1698
void print(std::ostream &os) const
Print method. Defines the behavior of the ostream &lt;&lt; operator inherited from the Object class...
Definition: GslMatrix.C:1772
GslMatrix & operator*=(double a)
Stores in this the coordinate-wise multiplication of this and a.
Definition: GslMatrix.C:177
~GslMatrix()
Destructor.
Definition: GslMatrix.C:161
void zeroUpper(bool includeDiagonal=false)
This function sets all the entries above the main diagonal of this matrix to zero.
Definition: GslMatrix.C:619
GslMatrix & operator-=(const GslMatrix &rhs)
Stores in this the coordinate-wise subtraction of this by rhs.
Definition: GslMatrix.C:207
GslMatrix matrixProduct(const GslVector &v1, const GslVector &v2)
Definition: GslMatrix.C:2036
void subWriteContents(const std::string &varNamePrefix, const std::string &fileName, const std::string &fileType, const std::set< unsigned int > &allowedSubEnvIds) const
Write contents of subenvironment in file fileName.
Definition: GslMatrix.C:1801
gsl_matrix * m_mat
GSL matrix, also referred to as this matrix.
Definition: GslMatrix.h:409
void fillWithTensorProduct(unsigned int rowId, unsigned int colId, const GslMatrix &mat1, const GslMatrix &mat2, bool checkForExactNumRowsMatching, bool checkForExactNumColsMatching)
This function calculates the tensor product of matrices mat1 and mat2 and stores it in this matrix...
Definition: GslMatrix.C:953
double normFrob() const
Returns the Frobenius norm of this matrix.
Definition: GslMatrix.C:293
GslVector * m_svdSvec
m_svdSvec stores the diagonal of the N-by-N diagonal matrix of singular values S after the singular v...
Definition: GslMatrix.h:431
The QUESO MPI Communicator Class.
Definition: MpiComm.h:203
double m_determinant
The determinant of this matrix.
Definition: GslMatrix.h:446
double lnDeterminant() const
Calculates the ln(determinant) of this matrix.
Definition: GslMatrix.C:1079
void fillWithTranspose(unsigned int rowId, unsigned int colId, const GslMatrix &mat, bool checkForExactNumRowsMatching, bool checkForExactNumColsMatching)
This function stores the transpose of this matrix into this matrix.
Definition: GslMatrix.C:1014
void cwExtract(unsigned int rowId, unsigned int colId, GslMatrix &mat) const
Definition: GslMatrix.C:385
unsigned int numCols() const
Returns the column dimension of this matrix.
Definition: GslMatrix.C:287
int svd(GslMatrix &matU, GslVector &vecS, GslMatrix &matVt) const
Checks for the dimension of this matrix, matU, VecS and matVt, and calls the protected routine intern...
Definition: GslMatrix.C:435
GslMatrix leftDiagScaling(const GslVector &vec, const GslMatrix &mat)
Definition: GslMatrix.C:2052
gsl_matrix * m_LU
GSL matrix for the LU decomposition of m_mat.
Definition: GslMatrix.h:412
void setColumn(const unsigned int column_num, const GslVector &column)
This function copies vector column into the column_num-th column of this matrix.
Definition: GslMatrix.C:1681
const Map & map() const
Definition: Matrix.C:54
int internalSvd() const
This function factorizes the M-by-N matrix A into the singular value decomposition A = U S V^T for M ...
Definition: GslMatrix.C:532
GslMatrix * m_inverse
Inverse matrix of this.
Definition: GslMatrix.h:415
GslMatrix rightDiagScaling(const GslMatrix &mat, const GslVector &vec)
Definition: GslMatrix.C:2073
GslMatrix & operator+=(const GslMatrix &rhs)
Stores in this the coordinate-wise addition of this and rhs.
Definition: GslMatrix.C:196
double max() const
Returns the maximum element value of the matrix.
Definition: GslMatrix.C:329
void fillWithBlocksVertically(unsigned int rowId, unsigned int colId, const std::vector< const GslMatrix * > &matrices, bool checkForExactNumRowsMatching, bool checkForExactNumColsMatching)
This function fills this matrix vertically with const block matrices.
Definition: GslMatrix.C:889
GslMatrix operator+(const GslMatrix &m1, const GslMatrix &m2)
Definition: GslMatrix.C:2022
const double & operator()(unsigned int i, unsigned int j) const
Element access method (const).
Definition: GslMatrix.h:115

Generated on Tue Jun 5 2018 19:48:54 for queso-0.57.1 by  doxygen 1.8.5