25 #ifndef UQ_GSL_MATRIX_H
26 #define UQ_GSL_MATRIX_H
32 #include <queso/Defines.h>
33 #include <queso/GslVector.h>
34 #include <queso/Matrix.h>
36 #include <gsl/gsl_matrix.h>
37 #include <gsl/gsl_permutation.h>
103 double&
operator()(
unsigned int i,
unsigned int j)
111 return *gsl_matrix_ptr(
m_mat,i,j);
115 const double&
operator()(
unsigned int i,
unsigned int j)
const
119 return *gsl_matrix_ptr(
m_mat,i,j);
141 unsigned int rank (
double absoluteZeroThreshold,
double relativeZeroThreshold)
const;
261 void cwSet (
double value);
264 void cwSet (
unsigned int rowId,
unsigned int colId,
const GslMatrix& mat);
271 void zeroLower (
bool includeDiagonal =
false);
276 void zeroUpper (
bool includeDiagonal =
false);
290 bool checkForExactNumRowsMatching,
291 bool checkForExactNumColsMatching);
296 const std::vector<const GslMatrix* >& matrices,
297 bool checkForExactNumRowsMatching,
298 bool checkForExactNumColsMatching);
303 const std::vector< GslMatrix* >& matrices,
304 bool checkForExactNumRowsMatching,
305 bool checkForExactNumColsMatching);
310 const std::vector<const GslMatrix* >& matrices,
311 bool checkForExactNumRowsMatching,
312 bool checkForExactNumColsMatching);
317 const std::vector< GslMatrix* >& matrices,
318 bool checkForExactNumRowsMatching,
319 bool checkForExactNumColsMatching);
324 const std::vector<const GslMatrix* >& matrices,
325 bool checkForExactNumRowsMatching,
326 bool checkForExactNumColsMatching);
331 const std::vector< GslMatrix* >& matrices,
332 bool checkForExactNumRowsMatching,
333 bool checkForExactNumColsMatching);
340 bool checkForExactNumRowsMatching,
341 bool checkForExactNumColsMatching);
348 bool checkForExactNumRowsMatching,
349 bool checkForExactNumColsMatching);
369 void print (std::ostream& os)
const;
373 const std::string& fileName,
374 const std::string& fileType,
375 const std::set<unsigned int>& allowedSubEnvIds)
const;
379 const std::string& fileType,
380 const std::set<unsigned int>& allowedSubEnvIds);
475 #endif // UQ_GSL_MATRIX_H
GslMatrix transpose() const
This function calculated the transpose of this matrix (square).
void cwSet(double value)
Component-wise set all values to this with value.
int chol()
Computes Cholesky factorization of a real symmetric positive definite matrix this.
const GslMatrix & svdMatU() const
This function calls private member GslMatrix::internalSvd() to set a M-by-N orthogonal matrix U of th...
double normFrob() const
Returns the Frobenius norm of this matrix.
void mpiSum(const MpiComm &comm, GslMatrix &M_global) const
void setRow(const unsigned int row_num, const GslVector &row)
This function copies vector row into the row_num-th column of this matrix.
GslMatrix & operator-=(const GslMatrix &rhs)
Stores in this the coordinate-wise subtraction of this by rhs.
GslMatrix & operator*=(double a)
Stores in this the coordinate-wise multiplication of this and a.
GslMatrix()
Default Constructor.
void setColumn(const unsigned int column_num, const GslVector &column)
This function copies vector column into the column_num-th column of this matrix.
Map * m_svdColMap
Mapping for matrices involved in the singular value decomposition (svd) routine.
void zeroUpper(bool includeDiagonal=false)
This function sets all the entries above the main diagonal of this matrix to zero.
double m_determinant
The determinant of this matrix.
Class for matrix operations (virtual).
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.
bool m_isSingular
Indicates whether or not this matrix is singular.
A class for partitioning vectors and matrices.
void cwExtract(unsigned int rowId, unsigned int colId, GslMatrix &mat) const
#define queso_require_less_msg(expr1, expr2, msg)
unsigned int numCols() const
Returns the column dimension of this matrix.
int m_signum
m_signum stores the sign of the permutation of the LU decomposition PA = LU.
void zeroLower(bool includeDiagonal=false)
This function sets all the entries bellow the main diagonal of this matrix to zero.
GslMatrix matrixProduct(const GslVector &v1, const GslVector &v2)
unsigned int numRowsLocal() const
Returns the local row dimension of this matrix.
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.
GslMatrix operator*(double a, const GslMatrix &mat)
void print(std::ostream &os) const
Print method. Defines the behavior of the ostream << operator inherited from the Object class...
double m_lnDeterminant
The natural logarithm of the determinant of this matrix.
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...
double lnDeterminant() const
Calculates the ln(determinant) 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 ...
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...
GslMatrix * m_svdVmat
m_svdVmat stores the N-by-N orthogonal square matrix V after the singular value decomposition of a ma...
void matlabLinearInterpExtrap(const GslVector &x1Vec, const GslMatrix &y1Mat, const GslVector &x2Vec)
GslVector * m_svdSvec
m_svdSvec stores the diagonal of the N-by-N diagonal matrix of singular values S after the singular v...
GslMatrix leftDiagScaling(const GslVector &vec, const GslMatrix &mat)
std::ostream & operator<<(std::ostream &os, const BaseEnvironment &obj)
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 rightDiagScaling(const GslMatrix &mat, const GslVector &vec)
gsl_matrix * m_mat
GSL matrix, also referred to as this matrix.
Class for matrix operations using GSL library.
gsl_matrix * m_LU
GSL matrix for the LU decomposition of m_mat.
double & operator()(unsigned int i, unsigned int j)
Element access method (non-const).
const BaseEnvironment & env() const
double max() const
Returns the maximum element value of the matrix.
This (virtual) class sets up the environment underlying the use of the QUESO library by an executable...
void filterLargeValues(double thresholdValue)
This function sets to zero (filters) all entries of this matrix which are greater than thresholdValue...
gsl_permutation * m_permutation
The permutation matrix of a LU decomposition.
GslMatrix * m_inverse
Inverse matrix of this.
void smallestEigen(double &eigenValue, GslVector &eigenVector) const
This function finds smallest eigenvalue, namely eigenValue, of this matrix and its corresponding eige...
GslVector invertMultiplyForceLU(const GslVector &b) const
This function calculates the inverse of this matrix and multiplies it with vector b...
unsigned int numRowsGlobal() const
Returns the global row dimension of this matrix.
double normMax() const
Returns the Frobenius norm of this matrix.
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.
gsl_matrix * data()
Returns this matrix.
GslMatrix operator-(const GslMatrix &m1, const GslMatrix &m2)
GslMatrix inverse() const
This function calculated the inverse of this matrix (square).
void subReadContents(const std::string &fileName, const std::string &fileType, const std::set< unsigned int > &allowedSubEnvIds)
Read contents of subenvironment from file fileName.
void filterSmallValues(double thresholdValue)
This function sets to zero (filters) all entries of this matrix which are smaller than thresholdValue...
GslMatrix & operator+=(const GslMatrix &rhs)
Stores in this the coordinate-wise addition of this and rhs.
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...
GslVector multiply(const GslVector &x) const
This function multiplies this matrix by vector x and returns the resulting vector.
GslMatrix operator+(const GslMatrix &m1, const GslMatrix &m2)
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.
const double & operator()(unsigned int i, unsigned int j) const
Element access method (const).
GslMatrix * m_svdUmat
m_svdUmat stores the M-by-N orthogonal matrix U after the singular value decomposition of a matrix...
Class for vector operations using GSL library.
const GslMatrix & svdMatV() const
This function calls private member GslMatrix::internalSvd() to set a N-by-N orthogonal square matrix ...
The QUESO MPI Communicator Class.
void copy(const GslMatrix &src)
In this function this matrix receives a copy of matrix src.
GslMatrix & operator=(const GslMatrix &rhs)
Copies values from matrix rhs to this.
void resetLU()
In this function resets the LU decomposition of this matrix, as well as deletes the private member po...
unsigned int rank(double absoluteZeroThreshold, double relativeZeroThreshold) const
This function returns the number of singular values of this matrix (rank).
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...
GslVector invertMultiply(const GslVector &b) const
This function calculates the inverse of this matrix and multiplies it with vector b...
GslMatrix & operator/=(double a)
Stores in this the coordinate-wise division of this by a.
GslMatrix * m_svdVTmat
m_svdVmatT stores the transpose of N-by-N orthogonal square matrix V, namely V^T, after the singular ...
void largestEigen(double &eigenValue, GslVector &eigenVector) const
This function finds largest eigenvalue, namely eigenValue, of this matrix and its corresponding eigen...
void eigen(GslVector &eigenValues, GslMatrix *eigenVectors) const
This function computes the eigenvalues of a real symmetric matrix.
double determinant() const
Calculates the determinant of this matrix.
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).