H2Lib  3.0
Functions
harith.h File Reference
#include "hmatrix.h"
#include "truncation.h"

Go to the source code of this file.

Functions

void trunc_rkmatrix (pctruncmode tm, real eps, prkmatrix r)
 Truncate an rkmatrix, $A \gets \operatorname{trunc}(A,\epsilon)$. More...
 
void add_amatrix_rkmatrix (field alpha, bool atrans, pcamatrix a, pctruncmode tm, real eps, prkmatrix b)
 Truncated addition of a matrix to a low-rank matrices in rkmatrix representation, $B \gets \operatorname{trunc}(B+\alpha A,\epsilon)$. More...
 
void add_rkmatrix_amatrix (field alpha, bool atrans, pcrkmatrix a, pamatrix b)
 Truncated addition of a full matrix to a low-rank matrix. More...
 
void add_hmatrix_amatrix (field alpha, bool atrans, pchmatrix a, pamatrix b)
 Addition of a hierarchical matrix to a matrix, $B \gets B+\alpha A$. More...
 
void add_rkmatrix (field alpha, pcrkmatrix src, pctruncmode tm, real eps, prkmatrix trg)
 Truncated addition of two low-rank matrices in rkmatrix representation, $A \gets \operatorname{trunc}(A+\alpha B,\epsilon)$. More...
 
void add_amatrix_hmatrix (field alpha, bool atrans, pcamatrix a, pctruncmode tm, real eps, phmatrix b)
 Truncated addition of two low-rank matrices in rkmatrix. More...
 
void add_lower_amatrix_hmatrix (field alpha, bool atrans, pcamatrix a, pctruncmode tm, real eps, phmatrix b)
 Locally truncated addition of a matrix in amatrix representation to the lower triangular part of a hierarchical matrix in hmatrix representation, $\operatorname{lower}(B) \gets \operatorname{blocktrunc}(\operatorname{lower}(B + \alpha A),\epsilon)$. More...
 
void add_rkmatrix_hmatrix (field alpha, pcrkmatrix a, pctruncmode tm, real eps, phmatrix b)
 Locally truncated addition of a low-rank matrix in rkmatrix representation to a hierarchical matrix in hmatrix representation, $B \gets \operatorname{blocktrunc}(B + \alpha A,\epsilon)$. More...
 
void add_lower_rkmatrix_hmatrix (field alpha, pcrkmatrix a, pctruncmode tm, real eps, phmatrix b)
 Locally truncated addition of a low-rank matrix in rkmatrix representation to the lower triangular part of a hierarchical matrix in hmatrix representation, $\operatorname{lower}(B) \gets \operatorname{blocktrunc}(\operatorname{lower}(B + \alpha A),\epsilon)$. More...
 
void add_hmatrix (field alpha, pchmatrix a, pctruncmode tm, real eps, phmatrix b)
 Locally truncated addition of a hierarchical matrix to a hierarchical matrix representation, $B \gets \operatorname{blocktrunc}(B + \alpha A,\epsilon)$. More...
 
phmatrix split_sub_amatrix (pcamatrix f, pccluster rc, pccluster cc, bool rsplit, bool csplit)
 Split an amatrix into submatrices. More...
 
phmatrix split_rkmatrix (pcrkmatrix r, pccluster rc, pccluster cc, bool rsplit, bool csplit, bool copy)
 Split an rkmatrix into submatrices. More...
 
phmatrix split_sub_rkmatrix (pcrkmatrix r, pccluster rc, pccluster cc, bool rsplit, bool csplit)
 Split an rkmatrix into submatrices referencing the original matrix. More...
 
void merge_rkmatrix (bool colmerge, pcrkmatrix src, pctruncmode tm, real eps, prkmatrix trg)
 Split a leaf hmatrix into submatrices. More...
 
prkmatrix merge_hmatrix_rkmatrix (pchmatrix s, pctruncmode tm, real eps)
 Merge a depth-one hmatrix containing only rkmatrix submatrices into a new rkmatrix, $A \gets \operatorname{succtrunc}\left(\begin{pmatrix} B_{11} & \ldots & B_{1m}\\ \vdots & \ddots & \vdots\\ B_{n1} & \ldots & B_{nm} \end{pmatrix},\epsilon\right)$. More...
 
void addmul_rkmatrix_amatrix_amatrix (field alpha, bool atrans, pcrkmatrix a, bool btrans, pcamatrix b, bool ctrans, pamatrix c)
 Multiply a matrix by an rkmatrix, $C \gets C + \alpha A B$. More...
 
void addmul_hmatrix_amatrix_amatrix (field alpha, bool atrans, pchmatrix a, bool btrans, pcamatrix bp, bool ctrans, pamatrix cp)
 Multiply a matrix by a hmatrix, $C \gets C + \alpha A B$. More...
 
void addmul_hmatrix (field alpha, bool xtrans, pchmatrix x, bool ytrans, pchmatrix y, pctruncmode tm, real eps, phmatrix z)
 Multiply two H-matrices, $Z \gets \operatorname{succtrunc}(Z + \alpha X Y,\epsilon)$. More...
 
void addmul_lower_hmatrix (field alpha, bool xtrans, pchmatrix x, bool ytrans, pchmatrix y, pctruncmode tm, real eps, phmatrix z)
 Multiply two H-matrices, computing only the lower triangular part of the result, $\operatorname{lower}(Z) \gets \operatorname{succtrunc}(\operatorname{lower}(Z + \alpha X Y),\epsilon)$. More...
 
void invert_hmatrix (phmatrix a, phmatrix work, pctruncmode tm, real eps)
 Invert an H-matrix, $A \gets \operatorname{succtrunc}(A^{-1},\epsilon)$. More...
 
phmatrix clone_lower_hmatrix (bool aunit, pchmatrix a)
 Create a duplicate of the lower triangular part of an H-matrix. More...
 
phmatrix clone_upper_hmatrix (bool aunit, pchmatrix a)
 Create a duplicate of the upper triangular part of an H-matrix. More...
 
void triangularinvmul_hmatrix_avector (bool alower, bool aunit, bool atrans, pchmatrix a, pavector xp)
 Solve a triangular system, $x \gets R^{-1} x$ or $x \gets R^{-*} x$. More...
 
void triangularsolve_hmatrix_avector (bool alower, bool aunit, bool atrans, pchmatrix a, pavector x)
 Solve a triangular system $x \gets A^{-1} x$ or $x \gets (A^*)^{-1} x$. More...
 
void triangularinvmul_hmatrix_amatrix (bool alower, bool aunit, bool atrans, pchmatrix a, bool xtrans, pamatrix xp)
 Solve a triangular system $X \gets A^{-1} X$ or $X \gets (A^*)^{-1} X$. More...
 
void triangularinvmul_hmatrix (bool alower, bool aunit, bool atrans, pchmatrix a, pctruncmode tm, real eps, bool xtrans, phmatrix xp)
 Solve a triangular system, $X \gets \operatorname{succtrunc}(A^{-1} X, \epsilon)$ or $X \gets \operatorname{succtrunc}(A^{-*} X, \epsilon)$. More...
 
void triangularinvmul_amatrix_hmatrix (bool alower, bool aunit, bool atrans, pcamatrix a, pctruncmode tm, real eps, bool xtrans, phmatrix xp)
 Solve a triangular system, $X \gets \operatorname{succtrunc}(A^{-1} X, \epsilon)$ or $X \gets \operatorname{succtrunc}(A^{-*} X, \epsilon)$. More...
 
void triangularmul_hmatrix_avector (bool alower, bool aunit, bool atrans, pchmatrix a, pavector xp)
 Evaluate a triangular system $x \gets A x$ or $x \gets A^* x$. More...
 
void triangulareval_hmatrix_avector (bool alower, bool aunit, bool atrans, pchmatrix a, pavector x)
 Evaluate a triangular system $x \gets A x$ or $x \gets A^* x$. More...
 
void triangularmul_hmatrix_amatrix (bool alower, bool aunit, bool atrans, pchmatrix a, bool xtrans, pamatrix xp)
 Evaluate a triangular system $X \gets A X$ or $X \gets A^* X$. More...
 
void triangularmul_hmatrix (bool alower, bool aunit, bool atrans, pchmatrix a, pctruncmode tm, real eps, bool xtrans, phmatrix xp)
 Evaluate a triangular system, $X \gets \operatorname{succtrunc}(A X, \epsilon)$ or $X \gets \operatorname{succtrunc}(A^{*} X, \epsilon)$. More...
 
void lrdecomp_hmatrix (phmatrix a, pctruncmode tm, real eps)
 Compute the LR factorization, $A \approx L R$. More...
 
void lrsolve_n_hmatrix_avector (pchmatrix a, pavector x)
 Solve the linear systems $A x = b$ using the LR factorization provided by lrdecomp_hmatrix. More...
 
void lrsolve_t_hmatrix_avector (pchmatrix a, pavector x)
 Solve the linear systems $A^* x = b$ using the LR factorization provided by lrdecomp_hmatrix. More...
 
void lrsolve_hmatrix_avector (bool atrans, pchmatrix a, pavector x)
 Solve the linear systems $A x = b$ or $A^* x = b$ using the LR factorization provided by lrdecomp_hmatrix. More...
 
void lreval_n_hmatrix_avector (pchmatrix a, pavector x)
 Evaluate the linear system $x \gets A x$ using the LR factorization provided by lrdecomp_hmatrix. More...
 
void lreval_t_hmatrix_avector (pchmatrix a, pavector x)
 Evaluate the linear system $x \gets A^* x$ using the LR factorization provided by lrdecomp_hmatrix. More...
 
void lreval_hmatrix_avector (bool atrans, pchmatrix a, pavector x)
 Evaluate the linear system $x \gets A x$ or $x \gets A^* x$ using the LR factorization provided by lrdecomp_hmatrix. More...
 
void choldecomp_hmatrix (phmatrix a, pctruncmode tm, real eps)
 Compute the Cholesky factorization, $A \approx L L^*$. More...
 
void cholsolve_hmatrix_avector (pchmatrix a, pavector x)
 Solve the linear system $A x = b$ using the Cholesky factorization provided by choldecomp_hmatrix. More...
 
void choleval_hmatrix_avector (pchmatrix a, pavector x)
 Evaluate the linear system $x \gets A x$ using the Cholesky factorization provided by choldecomp_hmatrix. More...
 

Detailed Description

Author
Steffen Börm