H2Lib  3.0
Functions
factorizations

Triangular and orthogonal factorizations. More...

Functions

void diagsolve_amatrix_avector (bool atrans, pcamatrix a, pavector x)
 Solve $A x = b$ and variants for a diagonal matrix $A$. More...
 
void diagsolve_amatrix (bool atrans, pcamatrix a, bool xtrans, pamatrix x)
 Solve $A X = B$ and variants for a diagonal matrix $A$. More...
 
void diageval_amatrix_avector (bool atrans, pcamatrix a, pavector x)
 Evaluate $b = A x$ and variants for a diagonal matrix $A$. More...
 
void diageval_amatrix (bool atrans, pcamatrix a, bool xtrans, pamatrix x)
 Evaluate $B = A X$ and variants for a diagonal matrix $A$. More...
 
void diageval_realavector_amatrix (field alpha, bool atrans, pcrealavector a, bool xtrans, pamatrix x)
 Evaluate $B = A X$ and variants for a diagonal matrix $A$ represented by a vector of reals. More...
 
void triangularsolve_amatrix_avector (bool alower, bool aunit, bool atrans, pcamatrix a, pavector x)
 Solve $A x = b$ for a triangular matrix $A$. More...
 
void triangularsolve_amatrix (bool alower, bool aunit, bool atrans, pcamatrix a, bool xtrans, pamatrix x)
 Solve $A X = B$ and variants for a triangular matrix $A$. More...
 
void triangulareval_amatrix_avector (bool alower, bool aunit, bool atrans, pcamatrix a, pavector x)
 Evaluate $b = A x$ for a triangular matrix $A$. More...
 
void triangulareval_amatrix (bool alower, bool aunit, bool atrans, pcamatrix a, bool xtrans, pamatrix x)
 Evaluate $B = A X$ and variants for a triangular matrix $A$. 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 $A$ and $B$ to a matrix $C$, $C \gets C + \alpha A B$. 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 $A=LR$ of a matrix. More...
 
uint lrdecomp_blocks_amatrix (pamatrix a, uint blocksize)
 Compute the LR decomposition $A=LR$ of a matrix using a block-based algorithm. More...
 
void lrdecomp_tasks_amatrix (pamatrix a, uint blocksize)
 Compute the LR decomposition $A=LR$ of a matrix using a block-based algorithm. More...
 
void lrsolve_n_amatrix_avector (pcamatrix a, pavector x)
 Solve the linear system $A x = b$ using a LR factorization. More...
 
void lrsolve_t_amatrix_avector (pcamatrix a, pavector x)
 Solve the linear system $A^* x = b$ using a LR factorization. More...
 
void lrsolve_amatrix_avector (bool atrans, pcamatrix a, pavector x)
 Solve the linear system $A x = b$ or $A^* x = b$ using a LR factorization. More...
 
void lrsolve_amatrix (pcamatrix a, pamatrix x)
 Solve the linear system $A X = B$ using a LR factorization. More...
 
void lreval_n_amatrix_avector (pcamatrix a, pavector x)
 Evaluate $x \gets A x$ using a LR factorization. More...
 
void lreval_t_amatrix_avector (pcamatrix a, pavector x)
 Evaluate $x \gets A^* x$ using a LR factorization. More...
 
void lreval_amatrix_avector (bool atrans, pcamatrix a, pavector x)
 Evaluate $x \gets A x$ or $x \gets A^* x$ using a LR factorization. More...
 
uint choldecomp_amatrix (pamatrix a)
 Compute the Cholesky decomposition $A=LL^*$ of a self-adjoint positive definite matrix. More...
 
uint choldecomp_blocks_amatrix (pamatrix a, uint blocksize)
 Compute the Cholesky decomposition $A=LL^*$ 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 $A=LL^*$ 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 $A x = b$ using a Cholesky factorization. More...
 
void cholsolve_amatrix (pcamatrix a, pamatrix x)
 Solve the linear system $A X = B$ using a Cholesky factorization. More...
 
void choleval_amatrix_avector (pcamatrix a, pavector x)
 Evaluate $x \gets A x$ using a Cholesky factorization. More...
 
uint ldltdecomp_amatrix (pamatrix a)
 Compute the LDLT decomposition $A=LDL^*$ of a self-adjoint matrix. More...
 
void ldltsolve_amatrix_avector (pcamatrix a, pavector x)
 Solve the linear system $A x = b$ using a LDLT factorization. More...
 
void ldltsolve_amatrix (pcamatrix a, pamatrix x)
 Solve the linear system $A X = B$ using a LDLT factorization. More...
 
void qrdecomp_amatrix (pamatrix a, pavector tau)
 Compute the QR decomposition $A=QR$ of a matrix. More...
 
uint qrdecomp_pivot_amatrix (pamatrix a, pavector tau, uint *colpiv)
 Compute the QR decomposition $A=QR$ 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 $A=QR$ 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 $b=Qx$ and variants for the matrix $Q$ of a QR decomposition. More...
 
void qreval_amatrix (bool qtrans, pcamatrix a, pcavector tau, pamatrix x)
 Evaluate $B=QX$ and variants for the matrix $Q$ of a QR decomposition. More...
 
void qrsolve_amatrix_avector (pcamatrix a, pcavector tau, pavector x)
 Solve the linear system $A x = b$ using a QR factorization. More...
 
void qrinvert_amatrix (pamatrix a)
 Compute the inverse $A^{-1}$ of a matrix by QR decomposition. More...
 
void qrexpand_amatrix (pcamatrix a, pcavector tau, pamatrix q)
 Compute the factor $Q$ of a QR factorization. More...
 

Detailed Description

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.

Function Documentation

uint choldecomp_amatrix ( pamatrix  a)

Compute the Cholesky decomposition $A=LL^*$ of a self-adjoint positive definite matrix.

Parameters
aOriginal matrix $A$. Lower part gets overwritten by $L$.
Returns
Zero if successful, number of first failed step otherwise.
uint choldecomp_blocks_amatrix ( pamatrix  a,
uint  blocksize 
)

Compute the Cholesky decomposition $A=LL^*$ of a self-adjoint positive definite matrix using a block-based algorithm.

Parameters
aOriginal matrix $A$. Lower part gets overwritten by $L$.
blocksizeMinimal block size.
Returns
Zero if successful, number of first failed step otherwise.
void choldecomp_tasks_amatrix ( pamatrix  a,
uint  blocksize 
)

Compute the Cholesky decomposition $A=LL^*$ 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.

Parameters
aOriginal matrix $A$. Lower part gets overwritten by $L$.
blocksizeMinimal block size.
void choleval_amatrix_avector ( pcamatrix  a,
pavector  x 
)

Evaluate $x \gets A x$ using a Cholesky factorization.

Parameters
aCholesky factorization of $A$ as provided by choldecomp_amatrix.
xRight-hand side vector $x$, that gets overwritten by the result.
void cholsolve_amatrix ( pcamatrix  a,
pamatrix  x 
)

Solve the linear system $A X = B$ using a Cholesky factorization.

Parameters
aCholesky factorization of $A$ as provided by choldecomp_amatrix.
xRight-hand side $B$, gets overwritten by $X$.
void cholsolve_amatrix_avector ( pcamatrix  a,
pavector  x 
)

Solve the linear system $A x = b$ using a Cholesky factorization.

Parameters
aCholesky factorization of $A$ as provided by choldecomp_amatrix.
xRight-hand side $b$, gets overwritten by $x$.
void copy_lower_amatrix ( pcamatrix  a,
bool  aunit,
pamatrix  b 
)

Copy the lower triangular part of a matrix.

Parameters
aSource matrix.
aunitSet if $A$ has unit diagonal.
bTarget matrix.
void copy_upper_amatrix ( pcamatrix  a,
bool  aunit,
pamatrix  b 
)

Copy the upper triangular part of a matrix.

Parameters
aSource matrix.
aunitSet if $A$ has unit diagonal.
bTarget matrix.
void diageval_amatrix ( bool  atrans,
pcamatrix  a,
bool  xtrans,
pamatrix  x 
)

Evaluate $B = A X$ and variants for a diagonal matrix $A$.

Parameters
atransSet if $A$ has to be transposed.
aDiagonal matrix.
xtransSet if $X$ and $B$ have to be transposed.
xRight-hand side $B$, gets overwritten by solution $X$.
void diageval_amatrix_avector ( bool  atrans,
pcamatrix  a,
pavector  x 
)

Evaluate $b = A x$ and variants for a diagonal matrix $A$.

Parameters
atransSet if $A$ has to be transposed.
aDiagonal matrix.
xRight-hand side $b$, gets overwritten by solution $x$.
void diageval_realavector_amatrix ( field  alpha,
bool  atrans,
pcrealavector  a,
bool  xtrans,
pamatrix  x 
)

Evaluate $B = A X$ and variants for a diagonal matrix $A$ represented by a vector of reals.

Parameters
alphaScaling factor for $A$.
atransSet if $A$ has to be transposed.
aVector containing the diagonal of $A$.
xtransSet if $X$ and $B$ have to be transposed.
xRight-hand side $B$, gets overwritten by solution $X$.
void diagsolve_amatrix ( bool  atrans,
pcamatrix  a,
bool  xtrans,
pamatrix  x 
)

Solve $A X = B$ and variants for a diagonal matrix $A$.

Parameters
atransSet if $A$ has to be transposed.
aDiagonal matrix.
xtransSet if $X$ and $B$ have to be transposed.
xRight-hand side $B$, gets overwritten by solution $X$.
void diagsolve_amatrix_avector ( bool  atrans,
pcamatrix  a,
pavector  x 
)

Solve $A x = b$ and variants for a diagonal matrix $A$.

Parameters
atransSet if adjoint matrix should be used.
aDiagonal matrix.
xRight-hand side, gets overwritten by solution.
uint ldltdecomp_amatrix ( pamatrix  a)

Compute the LDLT decomposition $A=LDL^*$ of a self-adjoint matrix.

Parameters
aOriginal matrix $A$. Diagonal part gets overwritten by $D$, strictly lower part by $L$. $L$ has unit diagonal.
Returns
Zero if successful, number of first failed step otherwise.
void ldltsolve_amatrix ( pcamatrix  a,
pamatrix  x 
)

Solve the linear system $A X = B$ using a LDLT factorization.

Parameters
aLDLT factorization of $A$ as provided by ldltdecomp_amatrix.
xRight-hand side $B$, gets overwritten by $X$.
void ldltsolve_amatrix_avector ( pcamatrix  a,
pavector  x 
)

Solve the linear system $A x = b$ using a LDLT factorization.

Parameters
aLDLT factorization of $A$ as provided by ldltdecomp_amatrix.
xRight-hand side $b$, gets overwritten by $x$.
uint lrdecomp_amatrix ( pamatrix  a)

Compute the LR decomposition $A=LR$ of a matrix.

Parameters
aOriginal matrix $A$. Upper part gets overwritten by $R$, strictly lower part by $L$. $L$ has unit diagonal.
Returns
Zero if successful, number of first failed step otherwise.
uint lrdecomp_blocks_amatrix ( pamatrix  a,
uint  blocksize 
)

Compute the LR decomposition $A=LR$ of a matrix using a block-based algorithm.

Parameters
aOriginal matrix $A$. Lower part gets overwritten by $L$.
blocksizeMinimal block size.
Returns
Zero if successful, number of first failed step otherwise.
void lrdecomp_tasks_amatrix ( pamatrix  a,
uint  blocksize 
)

Compute the LR decomposition $A=LR$ of a matrix using a block-based algorithm.

If OpenMP is enabled and supports tasks, the function will try to utilize multiple threads.

Parameters
aOriginal matrix $A$. Lower part gets overwritten by $L$.
blocksizeMinimal block size.
void lreval_amatrix_avector ( bool  atrans,
pcamatrix  a,
pavector  x 
)

Evaluate $x \gets A x$ or $x \gets A^* x$ using a LR factorization.

Parameters
atransSet if $A^*$ instead of $A$ is to be used.
aLR factorization of $A$, as provided by lrdecomp_amatrix.
xRight-hand side vector $x$, that gets overwritten by the result.
void lreval_n_amatrix_avector ( pcamatrix  a,
pavector  x 
)

Evaluate $x \gets A x$ using a LR factorization.

Parameters
aLR factorization of $A$, as provided by lrdecomp_amatrix.
xRight-hand side vector $x$, that gets overwritten by the result.
void lreval_t_amatrix_avector ( pcamatrix  a,
pavector  x 
)

Evaluate $x \gets A^* x$ using a LR factorization.

Parameters
aLR factorization of $A$, as provided by lrdecomp_amatrix.
xRight-hand side vector $x$, that gets overwritten by the result.
void lrsolve_amatrix ( pcamatrix  a,
pamatrix  x 
)

Solve the linear system $A X = B$ using a LR factorization.

Parameters
aLR factorization of $A$, as provided by lrdecomp_amatrix.
xRight-hand side $B$, gets overwritten by $X$.
void lrsolve_amatrix_avector ( bool  atrans,
pcamatrix  a,
pavector  x 
)

Solve the linear system $A x = b$ or $A^* x = b$ using a LR factorization.

Parameters
atransSet if $A^*$ instead of $A$ is to be used.
aLR factorization of $A$, as provided by lrdecomp_amatrix.
xRight-hand side $b$, gets overwritten by $x$.
void lrsolve_n_amatrix_avector ( pcamatrix  a,
pavector  x 
)

Solve the linear system $A x = b$ using a LR factorization.

Parameters
aLR factorization of $A$, as provided by lrdecomp_amatrix.
xRight-hand side $b$, gets overwritten by $x$.
void lrsolve_t_amatrix_avector ( pcamatrix  a,
pavector  x 
)

Solve the linear system $A^* x = b$ using a LR factorization.

Parameters
aLR factorization of $A$, as provided by lrdecomp_amatrix.
xRight-hand side $b$, gets overwritten by $x$.
void qrdecomp_amatrix ( pamatrix  a,
pavector  tau 
)

Compute the QR decomposition $A=QR$ of a matrix.

Parameters
aOriginal matrix $A$. Upper triangular part gets overwritten by $R$, strictly lower triangular part by Householder vectors.
tauScaling factors of Householder reflections .
uint qrdecomp_pivot_amatrix ( pamatrix  a,
pavector  tau,
uint colpiv 
)

Compute the QR decomposition $A=QR$ of a matrix with column pivoting.

Parameters
aOriginal matrix $A$. Upper triangular part gets overwritten by $R$, strictly lower triangular part by Householder vectors.
tauScaling factors of Householder reflections.
colpivIf not null, will be filled with column pivots.
Returns
Number of elementary reflections.
uint qrdecomp_rank_amatrix ( pamatrix  a,
pavector  tau,
pctruncmode  tm,
real  eps,
uint colpiv 
)

Compute the QR decomposition $A=QR$ of a matrix with column pivoting, exit early if the remainder becomes small enough.

Parameters
aOriginal matrix $A$. Upper triangular part gets overwritten by $R$, strictly lower triangular part by Householder vectors.
tauScaling factors of Householder reflections.
tmTruncation mode.
epsTruncation accuracy.
colpivIf not null, will be filled with column pivots.
Returns
Number of elementary reflections.
void qreval_amatrix ( bool  qtrans,
pcamatrix  a,
pcavector  tau,
pamatrix  x 
)

Evaluate $B=QX$ and variants for the matrix $Q$ of a QR decomposition.

Parameters
qtransSet if $Q^*$ is to be used instead of $Q$.
aHouseholder vectors as provided by qrdecomp_amatrix.
tauScaling factors as provided by qrdecomp_amatrix.
xRight-hand side matrix $X$, gets overwritten by $B$.
void qreval_amatrix_avector ( bool  qtrans,
pcamatrix  a,
pcavector  tau,
pavector  x 
)

Evaluate $b=Qx$ and variants for the matrix $Q$ of a QR decomposition.

Parameters
qtransSet if $Q^*$ is to be used instead of $Q$.
aHouseholder vectors as provided by qrdecomp_amatrix.
tauScaling factors as provided by qrdecomp_amatrix.
xRight-hand side vector $x$, gets overwritten by $b$.
void qrexpand_amatrix ( pcamatrix  a,
pcavector  tau,
pamatrix  q 
)

Compute the factor $Q$ of a QR factorization.

Parameters
aHouseholder vectors as provided by qrdecomp_amatrix.
tauScaling factors as provided by qrdecomp_amatrix.
qMatrix, gets overwritten by $Q$.
void qrinvert_amatrix ( pamatrix  a)

Compute the inverse $A^{-1}$ of a matrix by QR decomposition.

Parameters
aMatrix $A$, gets overwritten by $A^{-1}$.
void qrsolve_amatrix_avector ( pcamatrix  a,
pcavector  tau,
pavector  x 
)

Solve the linear system $A x = b$ using a QR factorization.

Parameters
aQR factorization of $A$ as provided by qrdecomp_amatrix.
tauScaling factors as provided by qrdecomp_amatrix.
xRight-hand side vector $b$, gets overwritten by $x$.
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 $A$ and $B$ to a matrix $C$, $C \gets C + \alpha A B$.

Parameters
alphaScaling factor $\alpha$.
alowerSet if $A$ is lower triangular, otherwise it is assumed to be upper triangular.
atransSet if $A^*$ is to be used instead of $A$.
aLeft triangular matrix $A$.
blowerSet if $B$ is lower triangular, otherwise if is assumed to be upper triangular.
btransSet if $B^*$ is to be used instead of $B$.
bRight triangular matrix $B$.
cTarget matrix $C$.
void triangulareval_amatrix ( bool  alower,
bool  aunit,
bool  atrans,
pcamatrix  a,
bool  xtrans,
pamatrix  x 
)

Evaluate $B = A X$ and variants for a triangular matrix $A$.

Parameters
alowerSet if $A$ is lower triangular, otherwise it is assumed to be upper triangular.
atransSet if $A^*$ is to be used instead of $A$.
aunitSet if $A$ has unit diagonal.
aTriangular matrix $A$.
xtransSet if $X^*$ and $B^*$ are to be used instead of $X$ and $B$.
xRight-hand side matrix $X$, gets overwritten by $B$.
void triangulareval_amatrix_avector ( bool  alower,
bool  aunit,
bool  atrans,
pcamatrix  a,
pavector  x 
)

Evaluate $b = A x$ for a triangular matrix $A$.

Parameters
alowerSet if $A$ is lower triangular, otherwise it is assumed to be upper triangular.
aunitSet if $A$ has unit diagonal.
atrans
aTriangular matrix $A$.
xRight-hand side vector $x$, gets overwritten by $b$.
void triangularsolve_amatrix ( bool  alower,
bool  aunit,
bool  atrans,
pcamatrix  a,
bool  xtrans,
pamatrix  x 
)

Solve $A X = B$ and variants for a triangular matrix $A$.

Parameters
alowerSet if $A$ is lower triangular, otherwise it is assumed to be upper triangular.
atransSet if $A^*$ is to be used instead of $A$.
aunitSet if $A$ has unit diagonal.
aTriangular matrix $A$.
xRight-hand side $B$, gets overwritten by solution $X$.
xtransSet if $X^*$ has to be used instead of $X$.
void triangularsolve_amatrix_avector ( bool  alower,
bool  aunit,
bool  atrans,
pcamatrix  a,
pavector  x 
)

Solve $A x = b$ for a triangular matrix $A$.

Parameters
alowerSet if $A$ is lower triangular, otherwise it is assumed to be upper triangular.
atransSet if $A^*$ is to be used instead of $A$.
aunitSet of $A$ has unit diagonal.
aTriangular matrix $A$.
xRight-hand side $b$, gets overwritten by solution $x$