To handle dense matrices. More...
#include <DMatrix.h>

Public Member Functions | |
DMatrix () | |
Default constructor. More... | |
DMatrix (size_t nr) | |
Constructor for a matrix with nr rows and nr columns. More... | |
DMatrix (size_t nr, size_t nc) | |
Constructor for a matrix with nr rows and nc columns. More... | |
DMatrix (Vect< T_ > &v) | |
Constructor that uses a Vect instance. The class uses the memory space occupied by this vector. More... | |
DMatrix (const DMatrix< T_ > &m) | |
Copy Constructor. More... | |
DMatrix (Mesh &mesh, size_t dof=0, int is_diagonal=false) | |
Constructor using mesh to initialize structure of matrix. More... | |
~DMatrix () | |
Destructor. | |
void | setDiag () |
Store diagonal entries in a separate internal vector. | |
void | setDiag (const T_ &a) |
Set matrix as diagonal and assign its diagonal entries as a constant. More... | |
void | setDiag (const vector< T_ > &d) |
Set matrix as diagonal and assign its diagonal entries. More... | |
void | setSize (size_t size) |
Set size (number of rows) of matrix. More... | |
void | setSize (size_t nr, size_t nc) |
Set size (number of rows and columns) of matrix. More... | |
void | getColumn (size_t j, Vect< T_ > &v) const |
Get j -th column vector. More... | |
Vect< T_ > | getColumn (size_t j) const |
Get j -th column vector. More... | |
void | getRow (size_t i, Vect< T_ > &v) const |
Get i -th row vector. More... | |
Vect< T_ > | getRow (size_t i) const |
Get i -th row vector. More... | |
void | set (size_t i, size_t j, const T_ &val) |
Assign a constant value to an entry of the matrix. More... | |
void | reset () |
Set matrix to 0 and reset factorization parameter. More... | |
void | setRow (size_t i, const Vect< T_ > &v) |
Copy a given vector to a prescribed row in the matrix. More... | |
void | setColumn (size_t j, const Vect< T_ > &v) |
Copy a given vector to a prescribed column in the matrix. More... | |
void | MultAdd (T_ a, const Vect< T_ > &x, Vect< T_ > &y) const |
Multiply matrix by vector a*x and add result to y . More... | |
void | MultAdd (const Vect< T_ > &x, Vect< T_ > &y) const |
Multiply matrix by vector x and add result to y . More... | |
void | Mult (const Vect< T_ > &x, Vect< T_ > &y) const |
Multiply matrix by vector x and save result in y . More... | |
void | TMult (const Vect< T_ > &x, Vect< T_ > &y) const |
Multiply transpose of matrix by vector x and add result in y . More... | |
void | add (size_t i, size_t j, const T_ &val) |
Add constant val to entry (i,j) of the matrix. More... | |
void | Axpy (T_ a, const DMatrix< T_ > &m) |
Add to matrix the product of a matrix by a scalar. More... | |
void | Axpy (T_ a, const Matrix< T_ > *m) |
Add to matrix the product of a matrix by a scalar. More... | |
int | setQR () |
Construct a QR factorization of the matrix. More... | |
int | setTransQR () |
Construct a QR factorization of the transpose of the matrix. More... | |
int | solveQR (const Vect< T_ > &b, Vect< T_ > &x) |
Solve a linear system by QR decomposition. More... | |
int | solveTransQR (const Vect< T_ > &b, Vect< T_ > &x) |
Solve a transpose linear system by QR decomposition. More... | |
T_ | operator() (size_t i, size_t j) const |
Operator () (Constant version). Return a(i,j) More... | |
T_ & | operator() (size_t i, size_t j) |
Operator () (Non constant version). Return a(i,j) More... | |
int | setLU () |
Factorize the matrix (LU factorization) More... | |
int | setTransLU () |
Factorize the transpose of the matrix (LU factorization) More... | |
int | solve (Vect< T_ > &b, bool fact=true) |
Solve linear system. More... | |
int | solveTrans (Vect< T_ > &b, bool fact=true) |
Solve the transpose linear system. More... | |
int | solve (const Vect< T_ > &b, Vect< T_ > &x, bool fact=true) |
Solve linear system. More... | |
int | solveTrans (const Vect< T_ > &b, Vect< T_ > &x, bool fact=true) |
Solve the transpose linear system. More... | |
DMatrix & | operator= (DMatrix< T_ > &m) |
Operator = More... | |
DMatrix & | operator+= (const DMatrix< T_ > &m) |
Operator +=. More... | |
DMatrix & | operator-= (const DMatrix< T_ > &m) |
Operator -=. More... | |
DMatrix & | operator= (const T_ &x) |
Operator = More... | |
DMatrix & | operator*= (const T_ &x) |
Operator *= More... | |
DMatrix & | operator+= (const T_ &x) |
Operator += More... | |
DMatrix & | operator-= (const T_ &x) |
Operator -= More... | |
const T_ * | getArray () |
Return matrix as C-Array. More... | |
T_ | get (size_t i, size_t j) const |
Return entry (i,j) of matrix. | |
void | add (size_t i, const T_ &val) |
Add val to entry i . | |
![]() | |
Matrix () | |
Default constructor. More... | |
Matrix (const Matrix< T_ > &m) | |
Copy Constructor. | |
virtual | ~Matrix () |
Destructor. | |
virtual void | reset () |
Set matrix to 0 and reset factorization parameter. More... | |
size_t | getNbRows () const |
Return number of rows. | |
size_t | getNbColumns () const |
Return number of columns. | |
void | setPenal (real_t p) |
Set Penalty Parameter (For boundary condition prescription). | |
void | setDiagonal () |
Set the matrix as diagonal. | |
T_ | getDiag (size_t k) const |
Return k -th diagonal entry of matrix. More... | |
size_t | size () const |
Return matrix dimension (Number of rows and columns). | |
virtual void | MultAdd (const Vect< T_ > &x, Vect< T_ > &y) const =0 |
Multiply matrix by vector x and add to y | |
virtual void | MultAdd (T_ a, const Vect< T_ > &x, Vect< T_ > &y) const =0 |
Multiply matrix by vector a*x and add to y | |
virtual void | Mult (const Vect< T_ > &x, Vect< T_ > &y) const =0 |
Multiply matrix by vector x and save in y | |
virtual void | TMult (const Vect< T_ > &v, Vect< T_ > &w) const =0 |
Multiply transpose of matrix by vector x and save in y | |
virtual void | Axpy (T_ a, const Matrix< T_ > *x)=0 |
Add to matrix the product of a matrix by a scalar. More... | |
void | setDiagonal (Mesh &mesh) |
Initialize matrix storage in the case where only diagonal terms are stored. More... | |
virtual void | clear () |
brief Set all matrix entries to zero | |
void | Assembly (const Element &el, T_ *a) |
Assembly of element matrix into global matrix. More... | |
void | Assembly (const Side &sd, T_ *a) |
Assembly of side matrix into global matrix. More... | |
void | Prescribe (Vect< T_ > &b, const Vect< T_ > &u, int flag=0) |
Impose by a penalty method an essential boundary condition, using the Mesh instance provided by the constructor. More... | |
void | Prescribe (int dof, int code, Vect< T_ > &b, const Vect< T_ > &u, int flag=0) |
Impose by a penalty method an essential boundary condition to a given degree of freedom for a given code. More... | |
void | Prescribe (Vect< T_ > &b, int flag=0) |
Impose by a penalty method a homegeneous (=0) essential boundary condition. More... | |
void | Prescribe (size_t dof, Vect< T_ > &b, const Vect< T_ > &u, int flag=0) |
Impose by a penalty method an essential boundary condition when only one DOF is treated. More... | |
void | PrescribeSide () |
Impose by a penalty method an essential boundary condition when DOFs are supported by sides. More... | |
virtual void | add (size_t i, size_t j, const T_ &val)=0 |
Add val to entry (i,j) . | |
virtual void | add (size_t i, const T_ &val)=0 |
Add val to entry i . | |
virtual int | Factor ()=0 |
Factorize matrix. Available only if the storage class enables it. | |
virtual int | solve (Vect< T_ > &b, bool fact=true)=0 |
Solve the linear system. More... | |
virtual int | solve (const Vect< T_ > &b, Vect< T_ > &x, bool fact=true)=0 |
Solve the linear system. More... | |
int | FactorAndSolve (Vect< T_ > &b) |
Factorize matrix and solve the linear system. More... | |
int | FactorAndSolve (const Vect< T_ > &b, Vect< T_ > &x) |
Factorize matrix and solve the linear system. More... | |
size_t | getLength () const |
Return number of stored terms in matrix. | |
int | isDiagonal () const |
Say if matrix is diagonal or not. | |
int | isFactorized () const |
Say if matrix is factorized or not. More... | |
virtual size_t | getColInd (size_t i) const |
Return Column index for column i (See the description for class SpMatrix). | |
virtual size_t | getRowPtr (size_t i) const |
Return Row pointer for row i (See the description for class SpMatrix). | |
virtual void | set (size_t i, size_t j, const T_ &val)=0 |
Assign a value to an entry of the matrix. More... | |
virtual T_ & | operator() (size_t i, size_t j)=0 |
Operator () (Non constant version). More... | |
virtual T_ | operator() (size_t i, size_t j) const =0 |
Operator () (Non constant version). More... | |
T_ | operator() (size_t i) const |
Operator () with one argument (Constant version). More... | |
T_ & | operator() (size_t i) |
Operator () with one argument (Non Constant version). More... | |
T_ & | operator[] (size_t k) |
Operator [] (Non constant version). More... | |
T_ | operator[] (size_t k) const |
Operator [] (Constant version). More... | |
Matrix & | operator= (Matrix< T_ > &m) |
Operator =. More... | |
Matrix & | operator+= (const Matrix< T_ > &m) |
Operator +=. More... | |
Matrix & | operator-= (const Matrix< T_ > &m) |
Operator -=. More... | |
Matrix & | operator= (const T_ &x) |
Operator =. More... | |
Matrix & | operator*= (const T_ &x) |
Operator *=. More... | |
Matrix & | operator+= (const T_ &x) |
Operator +=. More... | |
Matrix & | operator-= (const T_ &x) |
Operator -=. More... | |
virtual T_ | get (size_t i, size_t j) const =0 |
Return entry (i,j) of matrix if this one is stored, 0 else. | |
Detailed Description
class OFELI::DMatrix< T_ >
To handle dense matrices.
This class enables storing and manipulating general dense matrices. Matrices can be square or rectangle ones.
- Template Parameters
-
T_ Data type (double, float, complex<double>, ...)
- Copyright
- GNU Lesser Public License
Constructor & Destructor Documentation
◆ DMatrix() [1/6]
DMatrix | ( | ) |
Default constructor.
Initializes a zero-dimension matrix.
◆ DMatrix() [2/6]
DMatrix | ( | size_t | nr | ) |
Constructor for a matrix with nr
rows and nr
columns.
Matrix entries are set to 0
.
◆ DMatrix() [3/6]
DMatrix | ( | size_t | nr, |
size_t | nc | ||
) |
Constructor for a matrix with nr
rows and nc
columns.
Matrix entries are set to 0.
◆ DMatrix() [4/6]
Constructor that uses a Vect instance. The class uses the memory space occupied by this vector.
- Parameters
-
[in] v Vector to copy
◆ DMatrix() [5/6]
◆ DMatrix() [6/6]
Constructor using mesh to initialize structure of matrix.
- Parameters
-
[in] mesh Mesh instance for which matrix graph is determined. [in] dof Option parameter, with default value 0
.
dof=1
means that only one degree of freedom for each node (or element or side) is taken to determine matrix structure. The valuedof=0
means that matrix structure is determined using all DOFs.[in] is_diagonal Boolean argument to say is the matrix is actually a diagonal matrix or not.
Member Function Documentation
◆ add()
|
virtual |
Add constant val
to entry (i,j)
of the matrix.
- Parameters
-
[in] i row index [in] j column index [in] val Constant to add
Implements Matrix< T_ >.
◆ Axpy() [1/2]
void Axpy | ( | T_ | a, |
const DMatrix< T_ > & | m | ||
) |
Add to matrix the product of a matrix by a scalar.
- Parameters
-
[in] a Scalar to premultiply [in] m Matrix by which a
is multiplied. The result is added to current instance
◆ Axpy() [2/2]
|
virtual |
Add to matrix the product of a matrix by a scalar.
- Parameters
-
[in] a Scalar to premultiply [in] m Matrix by which a
is multiplied. The result is added to current instance
Implements Matrix< T_ >.
◆ getArray()
const T_ * getArray | ( | ) |
Return matrix as C-Array.
Matrix is stored row by row.
◆ getColumn() [1/2]
Vect< T_ > getColumn | ( | size_t | j | ) | const |
Get j
-th column vector.
- Parameters
-
[in] j Index of column to extract
- Returns
- Vect instance where the column is stored
- Remarks
- Vector v does not need to be sized before. It is resized in the function
◆ getColumn() [2/2]
void getColumn | ( | size_t | j, |
Vect< T_ > & | v | ||
) | const |
Get j
-th column vector.
- Parameters
-
[in] j Index of column to extract [out] v Reference to Vect instance where the column is stored
- Remarks
- Vector v does not need to be sized before. It is resized in the function
◆ getRow() [1/2]
Vect< T_ > getRow | ( | size_t | i | ) | const |
Get i
-th row vector.
- Parameters
-
[in] i Index of row to extract
- Returns
- Vect instance where the row is stored
- Remarks
- Vector v does not need to be sized before. It is resized in the function
◆ getRow() [2/2]
void getRow | ( | size_t | i, |
Vect< T_ > & | v | ||
) | const |
Get i
-th row vector.
- Parameters
-
[in] i Index of row to extract [out] v Reference to Vect instance where the row is stored
- Remarks
- Vector v does not need to be sized before. It is resized in the function
◆ Mult()
Multiply matrix by vector x
and save result in y
.
- Parameters
-
[in] x Vector to add to y
[out] y Result.
Implements Matrix< T_ >.
◆ MultAdd() [1/2]
Multiply matrix by vector x
and add result to y
.
- Parameters
-
[in] x Vector to add to y
[in,out] y on input, vector to add to. On output, result.
Implements Matrix< T_ >.
◆ MultAdd() [2/2]
Multiply matrix by vector a*x
and add result to y
.
- Parameters
-
[in] a constant to multiply by [in] x Vector to multiply by a
[in,out] y on input, vector to add to. On output, result.
Implements Matrix< T_ >.
◆ operator()() [1/2]
|
virtual |
Operator () (Non constant version). Return a(i,j)
- Parameters
-
[in] i row index [in] j column index
Implements Matrix< T_ >.
◆ operator()() [2/2]
|
virtual |
Operator () (Constant version). Return a(i,j)
- Parameters
-
[in] i row index [in] j column index
Implements Matrix< T_ >.
◆ operator*=()
DMatrix & operator*= | ( | const T_ & | x | ) |
Operator *=
Premultiply matrix entries by constant value x
.
◆ operator+=() [1/2]
Operator +=.
Add matrix m
to current matrix instance.
◆ operator+=() [2/2]
DMatrix & operator+= | ( | const T_ & | x | ) |
Operator +=
Add constant value x
to matrix entries
◆ operator-=() [1/2]
Operator -=.
Subtract matrix m
from current matrix instance.
◆ operator-=() [2/2]
DMatrix & operator-= | ( | const T_ & | x | ) |
Operator -=
Subtract constant value x
from matrix entries.
◆ operator=() [1/2]
DMatrix & operator= | ( | const T_ & | x | ) |
Operator =
Assign matrix to identity times x
◆ operator=() [2/2]
◆ reset()
|
virtual |
Set matrix to 0 and reset factorization parameter.
- Warning
- This function must be used if after a factorization, the matrix has modified
Reimplemented from Matrix< T_ >.
◆ set()
|
virtual |
Assign a constant value to an entry of the matrix.
- Parameters
-
[in] i row index of matrix [in] j column index of matrix [in] val Value to assign to a(i,j)
.
Implements Matrix< T_ >.
◆ setColumn()
void setColumn | ( | size_t | j, |
const Vect< T_ > & | v | ||
) |
Copy a given vector to a prescribed column in the matrix.
- Parameters
-
[in] j column index to be assigned [in] v Vect instance to copy
◆ setDiag() [1/2]
void setDiag | ( | const T_ & | a | ) |
Set matrix as diagonal and assign its diagonal entries as a constant.
- Parameters
-
[in] a Value to assign to all diagonal entries
◆ setDiag() [2/2]
void setDiag | ( | const vector< T_ > & | d | ) |
Set matrix as diagonal and assign its diagonal entries.
- Parameters
-
[in] d Vector entries to assign to matrix diagonal entries
◆ setLU()
int setLU | ( | ) |
Factorize the matrix (LU factorization)
LU factorization of the matrix is realized. Note that since this is an in place factorization, the contents of the matrix are modified.
- Returns
-
0
if factorization was normally performed, -
n
if then
-th pivot is null.
-
- Remarks
- A flag in this class indicates after factorization that this one has been realized, so that, if the member function solve is called after this no further factorization is done.
◆ setQR()
int setQR | ( | ) |
Construct a QR factorization of the matrix.
This function constructs the QR decomposition using the Householder method. The upper triangular matrix R is returned in the upper triangle of the current matrix, except for the diagonal elements of R which are stored in an internal vector. The orthogonal matrix Q is represented as a product of n-1 Householder matrices Q1 . . . Qn-1, where Qj = 1 - uj.uj /cj . The i-th component of uj is zero for i = 1, ..., j-1 while the nonzero components are returned in a[i][j] for i = j, ..., n.
- Returns
0
if the decomposition was successful,k
is thek-th
row is singular
- Remarks
- The matrix can be square or rectangle
◆ setRow()
void setRow | ( | size_t | i, |
const Vect< T_ > & | v | ||
) |
Copy a given vector to a prescribed row in the matrix.
- Parameters
-
[in] i row index to be assigned [in] v Vect instance to copy
◆ setSize() [1/2]
void setSize | ( | size_t | nr, |
size_t | nc | ||
) |
Set size (number of rows and columns) of matrix.
- Parameters
-
[in] nr Number of rows. [in] nc Number of columns.
◆ setSize() [2/2]
void setSize | ( | size_t | size | ) |
Set size (number of rows) of matrix.
- Parameters
-
[in] size Number of rows and columns.
◆ setTransLU()
int setTransLU | ( | ) |
Factorize the transpose of the matrix (LU factorization)
LU factorization of the transpose of the matrix is realized. Note that since this is an in place factorization, the contents of the matrix are modified.
- Returns
-
0
if factorization was normally performed, -
n
if then
-th pivot is null.
-
- Remarks
- A flag in this class indicates after factorization that this one has been realized, so that, if the member function solve is called after this no further factorization is done.
◆ setTransQR()
int setTransQR | ( | ) |
Construct a QR factorization of the transpose of the matrix.
This function constructs the QR decomposition using the Householder method. The upper triangular matrix R is returned in the upper triangle of the current matrix, except for the diagonal elements of R which are stored in an internal vector. The orthogonal matrix Q is represented as a product of n-1 Householder matrices Q1 . . . Qn-1, where Qj = 1 - uj.uj /cj . The i-th component of uj is zero for i = 1, ..., j-1 while the nonzero components are returned in a[i][j] for i = j, ..., n.
- Returns
0
if the decomposition was successful,k
is thek-th
row is singular
- Remarks
- The matrix can be square or rectangle
◆ solve() [1/2]
Solve linear system.
The linear system having the current instance as a matrix is solved by using the LU decomposition. Solution is thus realized after a factorization step and a forward/backward substitution step. The factorization step is realized only if this was not already done.
Note that this function modifies the matrix contents is a factorization is performed. Naturally, if the the matrix has been modified after using this function, the user has to refactorize it using the function setLU. This is because the class has no non-expensive way to detect if the matrix has been modified. The function setLU realizes the factorization step only.
- Parameters
-
[in] b Vect instance that contains right-hand side. [out] x Vect instance that contains solution [in] fact Set true if matrix is to be factorized (Default value), false if not
- Returns
-
0
if solution was normally performed, -
n
if then
-th pivot is null.
-
Implements Matrix< T_ >.
◆ solve() [2/2]
|
virtual |
Solve linear system.
The linear system having the current instance as a matrix is solved by using the LU decomposition. Solution is thus realized after a factorization step and a forward/backward substitution step. The factorization step is realized only if this was not already done.
Note that this function modifies the matrix contents is a factorization is performed. Naturally, if the the matrix has been modified after using this function, the user has to refactorize it using the function setLU. This is because the class has no non-expensive way to detect if the matrix has been modified. The function setLU realizes the factorization step only.
- Parameters
-
[in,out] b Vect instance that contains right-hand side on input and solution on output. [in] fact Set true if matrix is to be factorized (Default value), false if not
- Returns
-
0
if solution was normally performed, -
n
if then
-th pivot is null.
-
Implements Matrix< T_ >.
◆ solveQR()
Solve a linear system by QR decomposition.
This function constructs the QR decomposition, if this was not already done by using the member function QR and solves the linear system
- Parameters
-
[in] b Right-hand side vector [out] x Solution vector. Must have been sized before using this function.
- Returns
- The same value as returned by the function QR
◆ solveTrans() [1/2]
Solve the transpose linear system.
The linear system having the current instance as a transpose matrix is solved by using the LU decomposition. Solution is thus realized after a factorization step and a forward/backward substitution step. The factorization step is realized only if this was not already done.
Note that this function modifies the matrix contents is a factorization is performed. Naturally, if the the matrix has been modified after using this function, the user has to refactorize it using the function setLU. This is because the class has no non-expensive way to detect if the matrix has been modified. The function setLU realizes the factorization step only.
- Parameters
-
[in] b Vect instance that contains right-hand side. [out] x Vect instance that contains solution [in] fact Set true if matrix is to be factorized (Default value), false if not
- Returns
-
0
if solution was normally performed, -
n
if then
-th pivot is null.
-
◆ solveTrans() [2/2]
int solveTrans | ( | Vect< T_ > & | b, |
bool | fact = true |
||
) |
Solve the transpose linear system.
The linear system having the current instance as a transpose matrix is solved by using the LU decomposition. Solution is thus realized after a factorization step and a forward/backward substitution step. The factorization step is realized only if this was not already done.
Note that this function modifies the matrix contents is a factorization is performed. Naturally, if the the matrix has been modified after using this function, the user has to refactorize it using the function setLU. This is because the class has no non-expensive way to detect if the matrix has been modified. The function setLU realizes the factorization step only.
- Parameters
-
[in,out] b Vect instance that contains right-hand side on input and solution on output. [in] fact Set true if matrix is to be factorized (Default value), false if not
- Returns
-
0
if solution was normally performed, -
n
if then
-th pivot is null.
-
◆ solveTransQR()
Solve a transpose linear system by QR decomposition.
This function constructs the QR decomposition, if this was not already done by using the member function QR and solves the linear system
- Parameters
-
[in] b Right-hand side vector [out] x Solution vector. Must have been sized before using this function.
- Returns
- The same value as returned by the function QR
◆ TMult()
Multiply transpose of matrix by vector x
and add result in y
.
- Parameters
-
[in] x Vector to add to y
[in,out] y on input, vector to add to. On output, result.
Implements Matrix< T_ >.