![]() |
H2Lib
3.0
|
Algebraic operations with hierarchical matrices. More...
Data Structures | |
struct | _haccum |
Accumulator for H-matrix products. More... | |
Typedefs | |
typedef struct _haccum | haccum |
Accumulator for H-matrix products. | |
typedef haccum * | phaccum |
Pointer to haccum object. | |
typedef const haccum * | pchaccum |
Pointer to constant haccum object. | |
typedef struct _hprodentry | hprodentry |
Unhandled H-matrix products. | |
typedef hprodentry * | phprodentry |
Points to hprodentry object. | |
Functions | |
void | trunc_rkmatrix (pctruncmode tm, real eps, prkmatrix r) |
Truncate an rkmatrix, ![]() | |
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, ![]() | |
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, ![]() | |
void | add_rkmatrix (field alpha, pcrkmatrix src, pctruncmode tm, real eps, prkmatrix trg) |
Truncated addition of two low-rank matrices in rkmatrix representation, ![]() | |
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, ![]() | |
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, ![]() | |
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, ![]() | |
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, ![]() | |
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, ![]() | |
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, ![]() | |
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, ![]() | |
void | addmul_hmatrix (field alpha, bool xtrans, pchmatrix x, bool ytrans, pchmatrix y, pctruncmode tm, real eps, phmatrix z) |
Multiply two H-matrices, ![]() | |
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, ![]() | |
void | invert_hmatrix (phmatrix a, phmatrix work, pctruncmode tm, real eps) |
Invert an H-matrix, ![]() | |
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, ![]() ![]() | |
void | triangularsolve_hmatrix_avector (bool alower, bool aunit, bool atrans, pchmatrix a, pavector x) |
Solve a triangular system ![]() ![]() | |
void | triangularinvmul_hmatrix_amatrix (bool alower, bool aunit, bool atrans, pchmatrix a, bool xtrans, pamatrix xp) |
Solve a triangular system ![]() ![]() | |
void | triangularinvmul_hmatrix (bool alower, bool aunit, bool atrans, pchmatrix a, pctruncmode tm, real eps, bool xtrans, phmatrix xp) |
Solve a triangular system, ![]() ![]() | |
void | triangularinvmul_amatrix_hmatrix (bool alower, bool aunit, bool atrans, pcamatrix a, pctruncmode tm, real eps, bool xtrans, phmatrix xp) |
Solve a triangular system, ![]() ![]() | |
void | triangularmul_hmatrix_avector (bool alower, bool aunit, bool atrans, pchmatrix a, pavector xp) |
Evaluate a triangular system ![]() ![]() | |
void | triangulareval_hmatrix_avector (bool alower, bool aunit, bool atrans, pchmatrix a, pavector x) |
Evaluate a triangular system ![]() ![]() | |
void | triangularmul_hmatrix_amatrix (bool alower, bool aunit, bool atrans, pchmatrix a, bool xtrans, pamatrix xp) |
Evaluate a triangular system ![]() ![]() | |
void | triangularmul_hmatrix (bool alower, bool aunit, bool atrans, pchmatrix a, pctruncmode tm, real eps, bool xtrans, phmatrix xp) |
Evaluate a triangular system, ![]() ![]() | |
void | lrdecomp_hmatrix (phmatrix a, pctruncmode tm, real eps) |
Compute the LR factorization, ![]() | |
void | lrsolve_n_hmatrix_avector (pchmatrix a, pavector x) |
Solve the linear systems ![]() | |
void | lrsolve_t_hmatrix_avector (pchmatrix a, pavector x) |
Solve the linear systems ![]() | |
void | lrsolve_hmatrix_avector (bool atrans, pchmatrix a, pavector x) |
Solve the linear systems ![]() ![]() | |
void | lreval_n_hmatrix_avector (pchmatrix a, pavector x) |
Evaluate the linear system ![]() | |
void | lreval_t_hmatrix_avector (pchmatrix a, pavector x) |
Evaluate the linear system ![]() | |
void | lreval_hmatrix_avector (bool atrans, pchmatrix a, pavector x) |
Evaluate the linear system ![]() ![]() | |
void | choldecomp_hmatrix (phmatrix a, pctruncmode tm, real eps) |
Compute the Cholesky factorization, ![]() | |
void | cholsolve_hmatrix_avector (pchmatrix a, pavector x) |
Solve the linear system ![]() | |
void | choleval_hmatrix_avector (pchmatrix a, pavector x) |
Evaluate the linear system ![]() | |
phaccum | new_haccum (phmatrix z, pctruncmode tm, real eps) |
Create a new H-matrix accumulator. More... | |
void | del_haccum (phaccum ha) |
Delete an accumulator. More... | |
void | addproduct_haccum (field alpha, bool xtrans, pchmatrix x, bool ytrans, pchmatrix y, phaccum za) |
Add a product ![]() | |
phaccum * | split_haccum (phmatrix z, phaccum za) |
Create H-matrix accumulators for submatrices. More... | |
void | flush_haccum (phaccum ha) |
Flush an H-matrix accumulator, i.e., add all accumulated products to the target matrix and clear the accumulator. More... | |
void | add_amatrix_destructive_rkmatrix (field alpha, bool atrans, pamatrix a, pctruncmode tm, real eps, prkmatrix b) |
Truncated addition of a matrix to a low-rank matrices in rkmatrix representation, ![]() | |
void | add_amatrix_destructive_hmatrix (field alpha, bool atrans, pamatrix a, pctruncmode tm, real eps, phmatrix b) |
Locally truncated addition of a matrix in amatrix representation to a hierarchical matrix in hmatrix representation, ![]() | |
void | addmul2_hmatrix (field alpha, bool xtrans, pchmatrix x, bool ytrans, pchmatrix y, pctruncmode tm, real eps, phmatrix z) |
Multiply two H-matrices using accumulators, ![]() | |
void | lowersolve_haccum (bool aunit, bool atrans, pchmatrix a, bool xtrans, phaccum xa) |
Solve a lower triangular system using accumulators, ![]() ![]() | |
void | lowersolve2_hmatrix_hmatrix (bool aunit, bool atrans, pchmatrix a, pctruncmode tm, real eps, bool xtrans, phmatrix x) |
Solve a lower triangular system using accumulators, ![]() ![]() | |
void | uppersolve_haccum (bool aunit, bool atrans, pchmatrix a, bool xtrans, phaccum xa) |
Solve an upper triangular system using accumulators, ![]() ![]() | |
void | uppersolve2_hmatrix_hmatrix (bool aunit, bool atrans, pchmatrix a, pctruncmode tm, real eps, bool xtrans, phmatrix x) |
Solve an upper triangular system using accumulators, ![]() ![]() | |
void | lrdecomp_haccum (phaccum aa) |
Compute the LR factorization using accumulators, ![]() | |
void | lrdecomp2_hmatrix (phmatrix a, pctruncmode tm, real eps) |
Compute the LR factorization using accumulators, ![]() | |
void | choldecomp_haccum (phaccum aa) |
Compute the Cholesky factorization using accumulators, ![]() | |
void | choldecomp2_hmatrix (phmatrix a, pctruncmode tm, real eps) |
Compute the Cholesky factorization using accumulators, ![]() | |
Algebraic operations with hierarchical matrices.
void add_amatrix_destructive_hmatrix | ( | field | alpha, |
bool | atrans, | ||
pamatrix | a, | ||
pctruncmode | tm, | ||
real | eps, | ||
phmatrix | b | ||
) |
Locally truncated addition of a matrix in amatrix representation to a hierarchical matrix in hmatrix representation, .
alpha | Scaling factor ![]() |
atrans | Set if ![]() ![]() |
a | Matrix ![]() |
tm | Truncation mode. |
eps | Truncation accuracy ![]() |
b | Target hmatrix ![]() |
void add_amatrix_destructive_rkmatrix | ( | field | alpha, |
bool | atrans, | ||
pamatrix | a, | ||
pctruncmode | tm, | ||
real | eps, | ||
prkmatrix | b | ||
) |
Truncated addition of a matrix to a low-rank matrices in rkmatrix representation, .
alpha | Scaling factor ![]() |
atrans | Set if ![]() ![]() |
a | Matrix ![]() |
tm | Truncation mode. |
eps | Truncation accuracy ![]() |
b | Target matrix ![]() |
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.
*Locally truncated addition of a matrix in amatrix representation to a hierarchical matrix in hmatrix representation, .
alpha | Scaling factor ![]() |
atrans | Set if ![]() ![]() |
a | Matrix ![]() |
tm | Truncation mode. |
eps | Truncation accuracy ![]() |
b | Target hmatrix ![]() |
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, .
alpha | Scaling factor ![]() |
atrans | Set if ![]() ![]() |
a | Matrix ![]() |
tm | Truncation mode. |
eps | Truncation accuracy ![]() |
b | Target matrix ![]() |
void add_hmatrix | ( | field | alpha, |
pchmatrix | a, | ||
pctruncmode | tm, | ||
real | eps, | ||
phmatrix | b | ||
) |
Addition of a hierarchical matrix to a matrix, .
alpha | Scaling factor ![]() |
atrans | Set if ![]() ![]() |
a | Matrix ![]() |
b | Target matrix ![]() |
void add_lower_rkmatrix_hmatrix | ( | field | alpha, |
pcrkmatrix | a, | ||
pctruncmode | tm, | ||
real | eps, | ||
phmatrix | b | ||
) |
void add_rkmatrix | ( | field | alpha, |
pcrkmatrix | src, | ||
pctruncmode | tm, | ||
real | eps, | ||
prkmatrix | trg | ||
) |
Truncated addition of two low-rank matrices in rkmatrix representation, .
alpha | Scaling factor ![]() |
src | Source matrix ![]() |
tm | Truncation mode. |
eps | Truncation accuracy ![]() |
trg | Target matrix ![]() |
void add_rkmatrix_amatrix | ( | field | alpha, |
bool | atrans, | ||
pcrkmatrix | a, | ||
pamatrix | b | ||
) |
Truncated addition of a full matrix to a low-rank matrix.
*Addition of a low-rank matrix to a matrix, .
alpha | Scaling factor ![]() |
atrans | Set if ![]() ![]() |
a | Matrix ![]() |
b | Target matrix ![]() |
void add_rkmatrix_hmatrix | ( | field | alpha, |
pcrkmatrix | a, | ||
pctruncmode | tm, | ||
real | eps, | ||
phmatrix | b | ||
) |
void addmul2_hmatrix | ( | field | alpha, |
bool | xtrans, | ||
pchmatrix | x, | ||
bool | ytrans, | ||
pchmatrix | y, | ||
pctruncmode | tm, | ||
real | eps, | ||
phmatrix | z | ||
) |
Multiply two H-matrices using accumulators, .
Compared to the standard implementation addmul_hmatrix, this algorithm uses H-matrix accumulators to reduce the number of singular value decompositions.
alpha | Scaling factor ![]() |
xtrans | Set if ![]() ![]() |
x | Hierarchical matrix ![]() |
ytrans | Set if ![]() ![]() |
y | Hierarchical matrix ![]() |
tm | Truncation mode. |
eps | Truncation accuracy ![]() |
z | Target matrix ![]() |
void addmul_hmatrix | ( | field | alpha, |
bool | xtrans, | ||
pchmatrix | x, | ||
bool | ytrans, | ||
pchmatrix | y, | ||
pctruncmode | tm, | ||
real | eps, | ||
phmatrix | z | ||
) |
Multiply two H-matrices, .
alpha | Scaling factor ![]() |
xtrans | Set if ![]() ![]() |
x | Hierarchical matrix ![]() |
ytrans | Set if ![]() ![]() |
y | Hierarchical matrix ![]() |
tm | Truncation mode. |
eps | Truncation accuracy ![]() |
z | Target matrix ![]() |
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, .
alpha | Scaling factor ![]() |
atrans | Set if ![]() ![]() |
a | Hierarchical matrix ![]() |
btrans | Set if ![]() ![]() |
bp | Matrix ![]() |
ctrans | Set if ![]() ![]() |
cp | Target matrix ![]() |
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, .
alpha | Scaling factor ![]() |
xtrans | Set if ![]() ![]() |
x | Hierarchical matrix ![]() |
ytrans | Set if ![]() ![]() |
y | Hierarchical matrix ![]() |
tm | Truncation mode. |
eps | Truncation accuracy ![]() |
z | Target matrix ![]() |
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, .
alpha | Scaling factor ![]() |
atrans | Set if ![]() ![]() |
a | Low-rank matrix ![]() |
btrans | Set if ![]() ![]() |
b | Matrix ![]() |
ctrans | Set if ![]() ![]() |
c | Target matrix ![]() |
void addproduct_haccum | ( | field | alpha, |
bool | xtrans, | ||
pchmatrix | x, | ||
bool | ytrans, | ||
pchmatrix | y, | ||
phaccum | za | ||
) |
Add a product to an accumulator.
alpha | Scaling factor. |
xtrans | Set if ![]() ![]() |
x | H-matrix ![]() |
ytrans | Set if ![]() ![]() |
y | H-matrix ![]() |
za | Accumulator for the result. |
void choldecomp2_hmatrix | ( | phmatrix | a, |
pctruncmode | tm, | ||
real | eps | ||
) |
Compute the Cholesky factorization using accumulators, .
The lower triangular part of is stored in the lower triangular part of the source matrix.
The strictly upper triangular part of the source matrix is not used.
a | Source matrix ![]() ![]() |
tm | Truncation mode. |
eps | Truncation accuracy. |
void choldecomp_haccum | ( | phaccum | aa | ) |
Compute the Cholesky factorization using accumulators, .
The lower triangular part of is stored in the lower triangular part of the source matrix.
The strictly upper triangular part of the source matrix is not used.
aa | Accumulator representation of ![]() |
void choldecomp_hmatrix | ( | phmatrix | a, |
pctruncmode | tm, | ||
real | eps | ||
) |
Compute the Cholesky factorization, .
The lower triangular part of is stored in the lower triangular part of the source matrix.
The strictly upper triangular part of the source matrix is not used.
a | Source matrix ![]() ![]() |
tm | Truncation mode. |
eps | Truncation accuracy. |
Evaluate the linear system using the Cholesky factorization provided by choldecomp_hmatrix.
a | Matrix containing the Cholesky factorization in the form returned by choldecomp_hmatrix. |
x | Right-hand side vector ![]() |
Solve the linear system using the Cholesky factorization provided by choldecomp_hmatrix.
a | Matrix containing the Cholesky factorization in the form returned by choldecomp_hmatrix. |
x | Right-hand side vector ![]() ![]() |
Create a duplicate of the lower triangular part of an H-matrix.
aunit | Set if the duplicate has to have unit diagonal. |
a | Source matrix ![]() |
aunit==true
. Create a duplicate of the upper triangular part of an H-matrix.
aunit | Set if the duplicate has to have unit diagonal. |
a | Source matrix ![]() |
aunit==true
. void del_haccum | ( | phaccum | ha | ) |
Delete an accumulator.
ha | Accumulator to be deleted. |
void flush_haccum | ( | phaccum | ha | ) |
Flush an H-matrix accumulator, i.e., add all accumulated products to the target matrix and clear the accumulator.
ha | H-matrix accumulator. |
void invert_hmatrix | ( | phmatrix | a, |
phmatrix | work, | ||
pctruncmode | tm, | ||
real | eps | ||
) |
Invert an H-matrix, .
a | Matrix ![]() ![]() |
work | Auxiliary storage, has to be zero and be of the same structure as a . Consider using clonestructure_hmatrix. |
tm | Truncation mode. |
eps | Truncation accuracy. |
void lowersolve2_hmatrix_hmatrix | ( | bool | aunit, |
bool | atrans, | ||
pchmatrix | a, | ||
pctruncmode | tm, | ||
real | eps, | ||
bool | xtrans, | ||
phmatrix | x | ||
) |
Solve a lower triangular system using accumulators, or
.
aunit | Set if ![]() |
atrans | Set if ![]() ![]() |
a | Matrix containing the lower triangular part of ![]() |
tm | Truncation mode. |
eps | Truncation accuracy. |
xtrans | set if ![]() ![]() |
x | Input matrix ![]() a->rc . |
Solve a lower triangular system using accumulators, or
.
aunit | Set if ![]() |
atrans | Set if ![]() ![]() |
a | Matrix containing the lower triangular part of ![]() |
xtrans | set if ![]() ![]() |
xa | Accumulator representation of ![]() |
void lrdecomp2_hmatrix | ( | phmatrix | a, |
pctruncmode | tm, | ||
real | eps | ||
) |
Compute the LR factorization using accumulators, .
The lower triangular part has unit diagonal, only its strict lower triangular part is stored in the strict lower triangular part of the source matrix.
The upper triangular part is stored in the upper triangular part of the source matrix.
a | Source matrix ![]() ![]() ![]() |
tm | Truncation mode. |
eps | Truncation accuracy. |
void lrdecomp_haccum | ( | phaccum | aa | ) |
Compute the LR factorization using accumulators, .
The lower triangular part has unit diagonal, only its strict lower triangular part is stored in the strict lower triangular part of the source matrix.
The upper triangular part is stored in the upper triangular part of the source matrix.
aa | Accumulator representation of ![]() ![]() ![]() |
void lrdecomp_hmatrix | ( | phmatrix | a, |
pctruncmode | tm, | ||
real | eps | ||
) |
Compute the LR factorization, .
The lower triangular part has unit diagonal, only its strict lower triangular part is stored in the strict lower triangular part of the source matrix.
The upper triangular part is stored in the upper triangular part of the source matrix.
a | Source matrix ![]() ![]() ![]() |
tm | Truncation mode. |
eps | Truncation accuracy. |
Evaluate the linear system or
using the LR factorization provided by lrdecomp_hmatrix.
atrans | Set if ![]() ![]() |
a | Matrix containing the LR factorization in the form returned by lrdecomp_hmatrix. |
x | Right-hand side vector ![]() |
Evaluate the linear system using the LR factorization provided by lrdecomp_hmatrix.
a | Matrix containing the LR factorization in the form returned by lrdecomp_hmatrix. |
x | Right-hand side vector ![]() |
Evaluate the linear system using the LR factorization provided by lrdecomp_hmatrix.
a | Matrix containing the LR factorization in the form returned by lrdecomp_hmatrix. |
x | Right-hand side vector ![]() |
Solve the linear systems or
using the LR factorization provided by lrdecomp_hmatrix.
atrans | Set if ![]() ![]() |
a | Matrix containing the LR factorization in the form returned by lrdecomp_hmatrix. |
x | Right-hand side vector ![]() ![]() |
Solve the linear systems using the LR factorization provided by lrdecomp_hmatrix.
a | Matrix containing the LR factorization in the form returned by lrdecomp_hmatrix. |
x | Right-hand side vector ![]() ![]() |
Solve the linear systems using the LR factorization provided by lrdecomp_hmatrix.
a | Matrix containing the LR factorization in the form returned by lrdecomp_hmatrix. |
x | Right-hand side vector ![]() ![]() |
prkmatrix merge_hmatrix_rkmatrix | ( | pchmatrix | s, |
pctruncmode | tm, | ||
real | eps | ||
) |
void merge_rkmatrix | ( | bool | colmerge, |
pcrkmatrix | src, | ||
pctruncmode | tm, | ||
real | eps, | ||
prkmatrix | trg | ||
) |
Split a leaf hmatrix into submatrices.
** Split a leaf hmatrix into submatrices referencing Block-merge two rk matrices into a new rk matrix, or
.
If colmerge
is not set, and
have to have the same number of rows, otherwise the same number of columns.
colmerge | Set to merge into a column instead of a row. |
src | Source matrix ![]() |
tm | Truncation mode. |
eps | Truncation accuracy ![]() |
trg | Target matrix ![]() |
phaccum new_haccum | ( | phmatrix | z, |
pctruncmode | tm, | ||
real | eps | ||
) |
Create a new H-matrix accumulator.
z | Target matrix. |
tm | Truncation mode. |
eps | Truncation accuracy. |
z
Create H-matrix accumulators for submatrices.
z | Target matrix, row and column cluster have to conincide with ha->z->rc and ha->z->cc . |
za | H-matrix accumulator. |
z
, numbered as in z->son
. phmatrix split_rkmatrix | ( | pcrkmatrix | r, |
pccluster | rc, | ||
pccluster | cc, | ||
bool | rsplit, | ||
bool | csplit, | ||
bool | copy | ||
) |
Split an rkmatrix into submatrices.
r | Original matrix. |
rc | Row cluster. |
cc | Column cluster. |
rsplit | Set to split in row direction. |
csplit | Set to split in column direction. |
copy | Set to copy values into new matrix, otherwise it will be zero. |
Split an rkmatrix into submatrices referencing the original matrix.
r | Original matrix. |
rc | Row cluster. |
cc | Column cluster. |
rsplit | Set to split in row direction. |
csplit | Set to split in column direction. |
void triangulareval_hmatrix_avector | ( | bool | alower, |
bool | aunit, | ||
bool | atrans, | ||
pchmatrix | a, | ||
pavector | x | ||
) |
Evaluate a triangular system or
.
alower | Set if ![]() |
aunit | Set if ![]() |
atrans | Set if ![]() ![]() |
a | Matrix containing the lower or upper triangular part of ![]() |
x | Input vector ![]() |
void triangularinvmul_amatrix_hmatrix | ( | bool | alower, |
bool | aunit, | ||
bool | atrans, | ||
pcamatrix | a, | ||
pctruncmode | tm, | ||
real | eps, | ||
bool | xtrans, | ||
phmatrix | xp | ||
) |
Solve a triangular system, or
.
alower | Set if ![]() |
aunit | Set if ![]() |
atrans | Set if ![]() ![]() |
a | Matrix containing the lower or upper triangular part of ![]() |
tm | Truncation mode. |
eps | Truncation accuracy. |
xtrans | set if ![]() ![]() |
xp | Input matrix ![]() a->rc . |
void triangularinvmul_hmatrix | ( | bool | alower, |
bool | aunit, | ||
bool | atrans, | ||
pchmatrix | a, | ||
pctruncmode | tm, | ||
real | eps, | ||
bool | xtrans, | ||
phmatrix | xp | ||
) |
Solve a triangular system, or
.
alower | Set if ![]() |
aunit | Set if ![]() |
atrans | Set if ![]() ![]() |
a | Matrix containing the triangular part of ![]() |
tm | Truncation mode. |
eps | Truncation accuracy. |
xtrans | set if ![]() ![]() |
xp | Input matrix ![]() a->rc . |
void triangularinvmul_hmatrix_amatrix | ( | bool | alower, |
bool | aunit, | ||
bool | atrans, | ||
pchmatrix | a, | ||
bool | xtrans, | ||
pamatrix | xp | ||
) |
Solve a triangular system or
.
alower | Set if ![]() |
aunit | Set if ![]() |
atrans | Set if ![]() ![]() |
a | Matrix containing the lower or upper triangular part of ![]() |
xtrans | set if ![]() ![]() |
xp | Input matrix ![]() a->rc . |
void triangularinvmul_hmatrix_avector | ( | bool | alower, |
bool | aunit, | ||
bool | atrans, | ||
pchmatrix | a, | ||
pavector | xp | ||
) |
Solve a triangular system, or
.
alower | Set if ![]() |
aunit | Set if ![]() |
atrans | Set if ![]() ![]() |
a | Matrix containing the upper triangular part of ![]() |
xp | Input vector ![]() a->cc respectively a->rc . |
void triangularmul_hmatrix | ( | bool | alower, |
bool | aunit, | ||
bool | atrans, | ||
pchmatrix | a, | ||
pctruncmode | tm, | ||
real | eps, | ||
bool | xtrans, | ||
phmatrix | xp | ||
) |
Evaluate a triangular system, or
.
alower | Set if ![]() |
aunit | Set if ![]() |
atrans | Set if ![]() ![]() |
a | Matrix containing the lower or upper triangular part of ![]() |
tm | Truncation mode. |
eps | Truncation accuracy. |
xtrans | set if ![]() ![]() |
xp | Input matrix ![]() a->rc . |
void triangularmul_hmatrix_amatrix | ( | bool | alower, |
bool | aunit, | ||
bool | atrans, | ||
pchmatrix | a, | ||
bool | xtrans, | ||
pamatrix | xp | ||
) |
Evaluate a triangular system or
.
alower | Set if ![]() |
aunit | Set if ![]() |
atrans | Set if ![]() ![]() |
a | Matrix containing the lower or upper triangular part of ![]() |
xtrans | set if ![]() ![]() |
xp | Input matrix ![]() a->rc . |
void triangularmul_hmatrix_avector | ( | bool | alower, |
bool | aunit, | ||
bool | atrans, | ||
pchmatrix | a, | ||
pavector | xp | ||
) |
Evaluate a triangular system or
.
alower | Set if ![]() |
aunit | Set if ![]() |
atrans | Set if ![]() ![]() |
a | Matrix containing the lower or upper triangular part of ![]() |
xp | Input vector ![]() a->rc . |
void triangularsolve_hmatrix_avector | ( | bool | alower, |
bool | aunit, | ||
bool | atrans, | ||
pchmatrix | a, | ||
pavector | x | ||
) |
Solve a triangular system or
.
alower | Set if ![]() |
aunit | Set if ![]() |
atrans | Set if ![]() ![]() |
a | Matrix containing the lower or upper triangular part of ![]() |
x | Input vector ![]() |
void trunc_rkmatrix | ( | pctruncmode | tm, |
real | eps, | ||
prkmatrix | r | ||
) |
Truncate an rkmatrix, .
tm | Truncation mode. |
eps | Truncation accuracy ![]() |
r | Source matrix, will be overwritten by truncated matrix. |
void uppersolve2_hmatrix_hmatrix | ( | bool | aunit, |
bool | atrans, | ||
pchmatrix | a, | ||
pctruncmode | tm, | ||
real | eps, | ||
bool | xtrans, | ||
phmatrix | x | ||
) |
Solve an upper triangular system using accumulators, or
.
aunit | Set if ![]() |
atrans | Set if ![]() ![]() |
a | Matrix containing the lower triangular part of ![]() |
tm | Truncation mode. |
eps | Truncation accuracy. |
xtrans | set if ![]() ![]() |
x | Input matrix ![]() a->rc . |
Solve an upper triangular system using accumulators, or
.
aunit | Set if ![]() |
atrans | Set if ![]() ![]() |
a | Matrix containing the upper triangular part of ![]() |
xtrans | set if ![]() ![]() |
xa | Accumulator representation of ![]() |