H2Lib
3.0
|
Triangular and orthogonal factorizations. More...
Functions | |
void | diagsolve_amatrix_avector (bool atrans, pcamatrix a, pavector x) |
Solve and variants for a diagonal matrix . More... | |
void | diagsolve_amatrix (bool atrans, pcamatrix a, bool xtrans, pamatrix x) |
Solve and variants for a diagonal matrix . More... | |
void | diageval_amatrix_avector (bool atrans, pcamatrix a, pavector x) |
Evaluate and variants for a diagonal matrix . More... | |
void | diageval_amatrix (bool atrans, pcamatrix a, bool xtrans, pamatrix x) |
Evaluate and variants for a diagonal matrix . More... | |
void | diageval_realavector_amatrix (field alpha, bool atrans, pcrealavector a, bool xtrans, pamatrix x) |
Evaluate and variants for a diagonal matrix represented by a vector of reals. More... | |
void | triangularsolve_amatrix_avector (bool alower, bool aunit, bool atrans, pcamatrix a, pavector x) |
Solve for a triangular matrix . More... | |
void | triangularsolve_amatrix (bool alower, bool aunit, bool atrans, pcamatrix a, bool xtrans, pamatrix x) |
Solve and variants for a triangular matrix . More... | |
void | triangulareval_amatrix_avector (bool alower, bool aunit, bool atrans, pcamatrix a, pavector x) |
Evaluate for a triangular matrix . More... | |
void | triangulareval_amatrix (bool alower, bool aunit, bool atrans, pcamatrix a, bool xtrans, pamatrix x) |
Evaluate and variants for a triangular matrix . More... | |
void | triangularaddmul_amatrix (field alpha, bool alower, bool atrans, pcamatrix a, bool blower, bool btrans, pcamatrix b, pamatrix c) |
Add the product of two triangular matrices and to a matrix , . More... | |
void | copy_lower_amatrix (pcamatrix a, bool aunit, pamatrix b) |
Copy the lower triangular part of a matrix. More... | |
void | copy_upper_amatrix (pcamatrix a, bool aunit, pamatrix b) |
Copy the upper triangular part of a matrix. More... | |
uint | lrdecomp_amatrix (pamatrix a) |
Compute the LR decomposition of a matrix. More... | |
uint | lrdecomp_blocks_amatrix (pamatrix a, uint blocksize) |
Compute the LR decomposition of a matrix using a block-based algorithm. More... | |
void | lrdecomp_tasks_amatrix (pamatrix a, uint blocksize) |
Compute the LR decomposition of a matrix using a block-based algorithm. More... | |
void | lrsolve_n_amatrix_avector (pcamatrix a, pavector x) |
Solve the linear system using a LR factorization. More... | |
void | lrsolve_t_amatrix_avector (pcamatrix a, pavector x) |
Solve the linear system using a LR factorization. More... | |
void | lrsolve_amatrix_avector (bool atrans, pcamatrix a, pavector x) |
Solve the linear system or using a LR factorization. More... | |
void | lrsolve_amatrix (pcamatrix a, pamatrix x) |
Solve the linear system using a LR factorization. More... | |
void | lreval_n_amatrix_avector (pcamatrix a, pavector x) |
Evaluate using a LR factorization. More... | |
void | lreval_t_amatrix_avector (pcamatrix a, pavector x) |
Evaluate using a LR factorization. More... | |
void | lreval_amatrix_avector (bool atrans, pcamatrix a, pavector x) |
Evaluate or using a LR factorization. More... | |
uint | choldecomp_amatrix (pamatrix a) |
Compute the Cholesky decomposition of a self-adjoint positive definite matrix. More... | |
uint | choldecomp_blocks_amatrix (pamatrix a, uint blocksize) |
Compute the Cholesky decomposition of a self-adjoint positive definite matrix using a block-based algorithm. More... | |
void | choldecomp_tasks_amatrix (pamatrix a, uint blocksize) |
Compute the Cholesky decomposition of a self-adjoint positive definite matrix using a block-based algorithm. More... | |
void | cholsolve_amatrix_avector (pcamatrix a, pavector x) |
Solve the linear system using a Cholesky factorization. More... | |
void | cholsolve_amatrix (pcamatrix a, pamatrix x) |
Solve the linear system using a Cholesky factorization. More... | |
void | choleval_amatrix_avector (pcamatrix a, pavector x) |
Evaluate using a Cholesky factorization. More... | |
uint | ldltdecomp_amatrix (pamatrix a) |
Compute the LDLT decomposition of a self-adjoint matrix. More... | |
void | ldltsolve_amatrix_avector (pcamatrix a, pavector x) |
Solve the linear system using a LDLT factorization. More... | |
void | ldltsolve_amatrix (pcamatrix a, pamatrix x) |
Solve the linear system using a LDLT factorization. More... | |
void | qrdecomp_amatrix (pamatrix a, pavector tau) |
Compute the QR decomposition of a matrix. More... | |
uint | qrdecomp_pivot_amatrix (pamatrix a, pavector tau, uint *colpiv) |
Compute the QR decomposition of a matrix with column pivoting. More... | |
uint | qrdecomp_rank_amatrix (pamatrix a, pavector tau, pctruncmode tm, real eps, uint *colpiv) |
Compute the QR decomposition of a matrix with column pivoting, exit early if the remainder becomes small enough. More... | |
void | qreval_amatrix_avector (bool qtrans, pcamatrix a, pcavector tau, pavector x) |
Evaluate and variants for the matrix of a QR decomposition. More... | |
void | qreval_amatrix (bool qtrans, pcamatrix a, pcavector tau, pamatrix x) |
Evaluate and variants for the matrix of a QR decomposition. More... | |
void | qrsolve_amatrix_avector (pcamatrix a, pcavector tau, pavector x) |
Solve the linear system using a QR factorization. More... | |
void | qrinvert_amatrix (pamatrix a) |
Compute the inverse of a matrix by QR decomposition. More... | |
void | qrexpand_amatrix (pcamatrix a, pcavector tau, pamatrix q) |
Compute the factor of a QR factorization. More... | |
Triangular and orthogonal factorizations.
The factorizations module covers triangular factorizations like the LR factorization (without pivoting), the Cholesky factorization and the LDL^* factorization, and the orthogonal QR factorization.
Compute the Cholesky decomposition of a self-adjoint positive definite matrix.
a | Original matrix . Lower part gets overwritten by . |
Compute the Cholesky decomposition of a self-adjoint positive definite matrix using a block-based algorithm.
a | Original matrix . Lower part gets overwritten by . |
blocksize | Minimal block size. |
Compute the Cholesky decomposition of a self-adjoint positive definite matrix using a block-based algorithm.
If OpenMP is enabled and supports tasks, the function will try to utilize multiple threads.
a | Original matrix . Lower part gets overwritten by . |
blocksize | Minimal block size. |
Evaluate using a Cholesky factorization.
a | Cholesky factorization of as provided by choldecomp_amatrix. |
x | Right-hand side vector , that gets overwritten by the result. |
Solve the linear system using a Cholesky factorization.
a | Cholesky factorization of as provided by choldecomp_amatrix. |
x | Right-hand side , gets overwritten by . |
Solve the linear system using a Cholesky factorization.
a | Cholesky factorization of as provided by choldecomp_amatrix. |
x | Right-hand side , gets overwritten by . |
Copy the lower triangular part of a matrix.
a | Source matrix. |
aunit | Set if has unit diagonal. |
b | Target matrix. |
Copy the upper triangular part of a matrix.
a | Source matrix. |
aunit | Set if has unit diagonal. |
b | Target matrix. |
Evaluate and variants for a diagonal matrix .
atrans | Set if has to be transposed. |
a | Diagonal matrix. |
xtrans | Set if and have to be transposed. |
x | Right-hand side , gets overwritten by solution . |
Evaluate and variants for a diagonal matrix .
atrans | Set if has to be transposed. |
a | Diagonal matrix. |
x | Right-hand side , gets overwritten by solution . |
void diageval_realavector_amatrix | ( | field | alpha, |
bool | atrans, | ||
pcrealavector | a, | ||
bool | xtrans, | ||
pamatrix | x | ||
) |
Evaluate and variants for a diagonal matrix represented by a vector of reals.
alpha | Scaling factor for . |
atrans | Set if has to be transposed. |
a | Vector containing the diagonal of . |
xtrans | Set if and have to be transposed. |
x | Right-hand side , gets overwritten by solution . |
Solve and variants for a diagonal matrix .
atrans | Set if has to be transposed. |
a | Diagonal matrix. |
xtrans | Set if and have to be transposed. |
x | Right-hand side , gets overwritten by solution . |
Solve and variants for a diagonal matrix .
atrans | Set if adjoint matrix should be used. |
a | Diagonal matrix. |
x | Right-hand side, gets overwritten by solution. |
Compute the LDLT decomposition of a self-adjoint matrix.
a | Original matrix . Diagonal part gets overwritten by , strictly lower part by . has unit diagonal. |
Solve the linear system using a LDLT factorization.
a | LDLT factorization of as provided by ldltdecomp_amatrix. |
x | Right-hand side , gets overwritten by . |
Solve the linear system using a LDLT factorization.
a | LDLT factorization of as provided by ldltdecomp_amatrix. |
x | Right-hand side , gets overwritten by . |
Compute the LR decomposition of a matrix.
a | Original matrix . Upper part gets overwritten by , strictly lower part by . has unit diagonal. |
Compute the LR decomposition of a matrix using a block-based algorithm.
a | Original matrix . Lower part gets overwritten by . |
blocksize | Minimal block size. |
Compute the LR decomposition of a matrix using a block-based algorithm.
If OpenMP is enabled and supports tasks, the function will try to utilize multiple threads.
a | Original matrix . Lower part gets overwritten by . |
blocksize | Minimal block size. |
Evaluate or using a LR factorization.
atrans | Set if instead of is to be used. |
a | LR factorization of , as provided by lrdecomp_amatrix. |
x | Right-hand side vector , that gets overwritten by the result. |
Evaluate using a LR factorization.
a | LR factorization of , as provided by lrdecomp_amatrix. |
x | Right-hand side vector , that gets overwritten by the result. |
Evaluate using a LR factorization.
a | LR factorization of , as provided by lrdecomp_amatrix. |
x | Right-hand side vector , that gets overwritten by the result. |
Solve the linear system using a LR factorization.
a | LR factorization of , as provided by lrdecomp_amatrix. |
x | Right-hand side , gets overwritten by . |
Solve the linear system or using a LR factorization.
atrans | Set if instead of is to be used. |
a | LR factorization of , as provided by lrdecomp_amatrix. |
x | Right-hand side , gets overwritten by . |
Solve the linear system using a LR factorization.
a | LR factorization of , as provided by lrdecomp_amatrix. |
x | Right-hand side , gets overwritten by . |
Solve the linear system using a LR factorization.
a | LR factorization of , as provided by lrdecomp_amatrix. |
x | Right-hand side , gets overwritten by . |
Compute the QR decomposition of a matrix.
a | Original matrix . Upper triangular part gets overwritten by , strictly lower triangular part by Householder vectors. |
tau | Scaling factors of Householder reflections . |
Compute the QR decomposition of a matrix with column pivoting.
a | Original matrix . Upper triangular part gets overwritten by , strictly lower triangular part by Householder vectors. |
tau | Scaling factors of Householder reflections. |
colpiv | If not null, will be filled with column pivots. |
Compute the QR decomposition of a matrix with column pivoting, exit early if the remainder becomes small enough.
a | Original matrix . Upper triangular part gets overwritten by , strictly lower triangular part by Householder vectors. |
tau | Scaling factors of Householder reflections. |
tm | Truncation mode. |
eps | Truncation accuracy. |
colpiv | If not null, will be filled with column pivots. |
Evaluate and variants for the matrix of a QR decomposition.
qtrans | Set if is to be used instead of . |
a | Householder vectors as provided by qrdecomp_amatrix. |
tau | Scaling factors as provided by qrdecomp_amatrix. |
x | Right-hand side matrix , gets overwritten by . |
Evaluate and variants for the matrix of a QR decomposition.
qtrans | Set if is to be used instead of . |
a | Householder vectors as provided by qrdecomp_amatrix. |
tau | Scaling factors as provided by qrdecomp_amatrix. |
x | Right-hand side vector , gets overwritten by . |
Compute the factor of a QR factorization.
a | Householder vectors as provided by qrdecomp_amatrix. |
tau | Scaling factors as provided by qrdecomp_amatrix. |
q | Matrix, gets overwritten by . |
void qrinvert_amatrix | ( | pamatrix | a | ) |
Compute the inverse of a matrix by QR decomposition.
a | Matrix , gets overwritten by . |
Solve the linear system using a QR factorization.
a | QR factorization of as provided by qrdecomp_amatrix. |
tau | Scaling factors as provided by qrdecomp_amatrix. |
x | Right-hand side vector , gets overwritten by . |
void triangularaddmul_amatrix | ( | field | alpha, |
bool | alower, | ||
bool | atrans, | ||
pcamatrix | a, | ||
bool | blower, | ||
bool | btrans, | ||
pcamatrix | b, | ||
pamatrix | c | ||
) |
Add the product of two triangular matrices and to a matrix , .
alpha | Scaling factor . |
alower | Set if is lower triangular, otherwise it is assumed to be upper triangular. |
atrans | Set if is to be used instead of . |
a | Left triangular matrix . |
blower | Set if is lower triangular, otherwise if is assumed to be upper triangular. |
btrans | Set if is to be used instead of . |
b | Right triangular matrix . |
c | Target matrix . |
void triangulareval_amatrix | ( | bool | alower, |
bool | aunit, | ||
bool | atrans, | ||
pcamatrix | a, | ||
bool | xtrans, | ||
pamatrix | x | ||
) |
Evaluate and variants for a triangular matrix .
alower | Set if is lower triangular, otherwise it is assumed to be upper triangular. |
atrans | Set if is to be used instead of . |
aunit | Set if has unit diagonal. |
a | Triangular matrix . |
xtrans | Set if and are to be used instead of and . |
x | Right-hand side matrix , gets overwritten by . |
void triangulareval_amatrix_avector | ( | bool | alower, |
bool | aunit, | ||
bool | atrans, | ||
pcamatrix | a, | ||
pavector | x | ||
) |
Evaluate for a triangular matrix .
alower | Set if is lower triangular, otherwise it is assumed to be upper triangular. |
aunit | Set if has unit diagonal. |
atrans | |
a | Triangular matrix . |
x | Right-hand side vector , gets overwritten by . |
void triangularsolve_amatrix | ( | bool | alower, |
bool | aunit, | ||
bool | atrans, | ||
pcamatrix | a, | ||
bool | xtrans, | ||
pamatrix | x | ||
) |
Solve and variants for a triangular matrix .
alower | Set if is lower triangular, otherwise it is assumed to be upper triangular. |
atrans | Set if is to be used instead of . |
aunit | Set if has unit diagonal. |
a | Triangular matrix . |
x | Right-hand side , gets overwritten by solution . |
xtrans | Set if has to be used instead of . |
void triangularsolve_amatrix_avector | ( | bool | alower, |
bool | aunit, | ||
bool | atrans, | ||
pcamatrix | a, | ||
pavector | x | ||
) |
Solve for a triangular matrix .
alower | Set if is lower triangular, otherwise it is assumed to be upper triangular. |
atrans | Set if is to be used instead of . |
aunit | Set of has unit diagonal. |
a | Triangular matrix . |
x | Right-hand side , gets overwritten by solution |