25 #ifndef UQ_GSL_MATRIX_H
26 #define UQ_GSL_MATRIX_H
32 #include <queso/Matrix.h>
33 #include <gsl/gsl_matrix.h>
34 #include <gsl/gsl_permutation.h>
35 #include <queso/GslVector.h>
105 double&
operator()(
unsigned int i,
unsigned int j);
109 const double&
operator()(
unsigned int i,
unsigned int j)
const;
130 unsigned int rank (
double absoluteZeroThreshold,
double relativeZeroThreshold)
const;
250 void cwSet (
double value);
253 void cwSet (
unsigned int rowId,
unsigned int colId,
const GslMatrix& mat);
260 void zeroLower (
bool includeDiagonal =
false);
265 void zeroUpper (
bool includeDiagonal =
false);
279 bool checkForExactNumRowsMatching,
280 bool checkForExactNumColsMatching);
285 const std::vector<const GslMatrix* >& matrices,
286 bool checkForExactNumRowsMatching,
287 bool checkForExactNumColsMatching);
292 const std::vector< GslMatrix* >& matrices,
293 bool checkForExactNumRowsMatching,
294 bool checkForExactNumColsMatching);
299 const std::vector<const GslMatrix* >& matrices,
300 bool checkForExactNumRowsMatching,
301 bool checkForExactNumColsMatching);
306 const std::vector< GslMatrix* >& matrices,
307 bool checkForExactNumRowsMatching,
308 bool checkForExactNumColsMatching);
313 const std::vector<const GslMatrix* >& matrices,
314 bool checkForExactNumRowsMatching,
315 bool checkForExactNumColsMatching);
320 const std::vector< GslMatrix* >& matrices,
321 bool checkForExactNumRowsMatching,
322 bool checkForExactNumColsMatching);
329 bool checkForExactNumRowsMatching,
330 bool checkForExactNumColsMatching);
337 bool checkForExactNumRowsMatching,
338 bool checkForExactNumColsMatching);
358 void print (std::ostream& os)
const;
362 const std::string& fileName,
363 const std::string& fileType,
364 const std::set<unsigned int>& allowedSubEnvIds)
const;
368 const std::string& fileType,
369 const std::set<unsigned int>& allowedSubEnvIds);
460 #endif // UQ_GSL_MATRIX_H
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.
GslMatrix & operator/=(double a)
Stores in this the coordinate-wise division of this by a.
void setColumn(const unsigned int column_num, const GslVector &column)
This function copies vector column into the column_num-th column of this matrix.
GslMatrix & operator+=(const GslMatrix &rhs)
Stores in this the coordinate-wise addition of this and rhs.
double max() const
Returns the maximum element value of the matrix.
Class for matrix operations using GSL library.
GslMatrix operator*(double a, const GslMatrix &mat)
gsl_matrix * m_LU
GSL matrix for the LU decomposition of m_mat.
GslMatrix operator-(const GslMatrix &m1, const GslMatrix &m2)
GslMatrix operator+(const GslMatrix &m1, const GslMatrix &m2)
GslMatrix rightDiagScaling(const GslMatrix &mat, const GslVector &vec)
void zeroUpper(bool includeDiagonal=false)
This function sets all the entries above the main diagonal of this matrix to zero.
GslMatrix matrixProduct(const GslVector &v1, const GslVector &v2)
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...
unsigned int numRowsLocal() const
Returns the local row dimension of this matrix.
std::ostream & operator<<(std::ostream &os, const BaseEnvironment &obj)
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.
void eigen(GslVector &eigenValues, GslMatrix *eigenVectors) const
This function computes the eigenvalues of a real symmetric matrix.
GslMatrix leftDiagScaling(const GslVector &vec, const GslMatrix &mat)
GslVector invertMultiply(const GslVector &b) const
This function calculates the inverse of this matrix and multiplies it with vector b...
double & operator()(unsigned int i, unsigned int j)
Element access method (non-const).
double normFrob() const
Returns the Frobenius norm of this matrix.
GslVector * m_svdSvec
m_svdSvec stores the diagonal of the N-by-N diagonal matrix of singular values S after the singular v...
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...
GslMatrix * m_svdUmat
m_svdUmat stores the M-by-N orthogonal matrix U after the singular value decomposition of a matrix...
void filterLargeValues(double thresholdValue)
This function sets to zero (filters) all entries of this matrix which are greater than thresholdValue...
This (virtual) class sets up the environment underlying the use of the QUESO library by an executable...
void print(std::ostream &os) const
Print method. Defines the behavior of the ostream << operator inherited from the Object class...
void mpiSum(const MpiComm &comm, GslMatrix &M_global) const
A class for partitioning vectors and matrices.
GslVector invertMultiplyForceLU(const GslVector &b) const
This function calculates the inverse of this matrix and multiplies it with vector b...
double m_lnDeterminant
The natural logarithm of the determinant of this matrix.
GslMatrix & operator-=(const GslMatrix &rhs)
Stores in this the coordinate-wise subtraction of this by rhs.
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...
void zeroLower(bool includeDiagonal=false)
This function sets all the entries bellow the main diagonal of this matrix to zero.
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.
GslMatrix transpose() const
This function calculated the transpose of this matrix (square).
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).
Class for vector operations using GSL library.
GslMatrix * m_svdVTmat
m_svdVmatT stores the transpose of N-by-N orthogonal square matrix V, namely V^T, after the singular ...
void cwSet(double value)
Component-wise set all values to this with value.
GslMatrix & operator=(const GslMatrix &rhs)
Copies values from matrix rhs to this.
void largestEigen(double &eigenValue, GslVector &eigenVector) const
This function finds largest eigenvalue, namely eigenValue, of this matrix and its corresponding eigen...
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.
GslMatrix()
Default Constructor.
Map * m_svdColMap
Mapping for matrices involved in the singular value decomposition (svd) routine.
const GslMatrix & svdMatV() const
This function calls private member GslMatrix::internalSvd() to set a N-by-N orthogonal square matrix ...
void setRow(const unsigned int row_num, const GslVector &row)
This function copies vector row into the row_num-th column of this matrix.
The QUESO MPI Communicator Class.
gsl_matrix * data()
Returns this matrix.
GslMatrix inverse() const
This function calculated the inverse of this matrix (square).
unsigned int rank(double absoluteZeroThreshold, double relativeZeroThreshold) const
This function returns the number of singular values of this matrix (rank).
Class for matrix operations (virtual).
void copy(const GslMatrix &src)
In this function this matrix receives a copy of matrix src.
void smallestEigen(double &eigenValue, GslVector &eigenVector) const
This function finds smallest eigenvalue, namely eigenValue, of this matrix and its corresponding eige...
GslVector multiply(const GslVector &x) const
This function multiplies this matrix by vector x and returns the resulting vector.
void subReadContents(const std::string &fileName, const std::string &fileType, const std::set< unsigned int > &allowedSubEnvIds)
Read contents of subenvironment from file fileName.
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.
GslMatrix & operator*=(double a)
Stores in this the coordinate-wise multiplication of this and a.
void filterSmallValues(double thresholdValue)
This function sets to zero (filters) all entries of this matrix which are smaller than thresholdValue...
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...
gsl_permutation * m_permutation
The permutation matrix of a LU decomposition.
double normMax() const
Returns the Frobenius norm of this matrix.
void matlabLinearInterpExtrap(const GslVector &x1Vec, const GslMatrix &y1Mat, const GslVector &x2Vec)
gsl_matrix * m_mat
GSL matrix, also referred to as this matrix.
void resetLU()
In this function resets the LU decomposition of this matrix, as well as deletes the private member po...
int m_signum
m_signum stores the sign of the permutation of the LU decomposition PA = LU.
unsigned int numRowsGlobal() const
Returns the global row dimension of this matrix.
int internalSvd() const
This function factorizes the M-by-N matrix A into the singular value decomposition A = U S V^T for M ...
double lnDeterminant() const
Calculates the ln(determinant) of this matrix.
double determinant() const
Calculates the determinant of this matrix.
bool m_isSingular
Indicates whether or not this matrix is singular.
double m_determinant
The determinant of this matrix.
GslMatrix * m_svdVmat
m_svdVmat stores the N-by-N orthogonal square matrix V after the singular value decomposition of a ma...
unsigned int numCols() const
Returns the column dimension of this matrix.
GslMatrix * m_inverse
Inverse matrix of this.
void cwExtract(unsigned int rowId, unsigned int colId, GslMatrix &mat) const
const GslMatrix & svdMatU() const
This function calls private member GslMatrix::internalSvd() to set a M-by-N orthogonal matrix U of th...
const BaseEnvironment & env() const
int chol()
Computes Cholesky factorization of a real symmetric positive definite matrix this.