H2Lib  3.0
Functions
h2arith.h File Reference
#include "clusteroperator.h"
#include "h2matrix.h"
#include "h2compression.h"
#include "h2update.h"

Go to the source code of this file.

Functions

ph2matrix build_from_block_lower_h2matrix (pcblock b, pclusterbasis rb, pclusterbasis cb)
 builds a lower triangular H2-matrix based on block tree b More...
 
ph2matrix build_from_block_upper_h2matrix (pcblock b, pclusterbasis rb, pclusterbasis cb)
 builds a upper triangular H2-matrix based on block tree b More...
 
void tolower_h2matrix (ph2matrix a)
 Clears all upper diagonal blocks of a h2matrix. More...
 
void add_hmatrix_h2matrix (pchmatrix a, ph2matrix B, pclusteroperator rwf, pclusteroperator cwf, ptruncmode tm, real tol)
 $ B \gets B + a $ More...
 
pamatrix convert_h2matrix_amatrix (bool atrans, pch2matrix a)
 converts an h2matrix to a new amatrix More...
 
void fastaddmul_clusterbasis_amatrix (pcclusterbasis cb, pamatrix Yt, pamatrix Y)
 computes $ Y = cb\rightarrow V \cdot Yt_{|cb\rightarrow k \times Yt\rightarrow cols} + Y $
 
prkmatrix convert_uniform_rkmatrix (bool atrans, pcuniform a)
 converts an uniform matrix $ a $ to a new rkmatrix
uses fastaddmul_clusterbasis_amatrix More...
 
phmatrix convert_h2matrix_hmatrix (pch2matrix a)
 converts an h2matrix to a new hmatrix
More...
 
prkmatrix mul_h2matrix_rkmatrix (pch2matrix a, bool btrans, pch2matrix B, real tol)
 computes $ a \cdot B $ and converts it to a new rkmatrix More...
 
void addmul_h2matrix (field alpha, pch2matrix a, bool btrans, pch2matrix B, ph2matrix C, pclusteroperator rwf, pclusteroperator cwf, ptruncmode tm, real tol)
 $ C \gets C + alpha * a * B $ More...
 
ph2matrix mul_h2matrix (field alpha, ph2matrix a, ph2matrix B, ph2matrix h2, ptruncmode tm, real tol)
 computes $ a \cdot B $ and converts it to a new h2matrix More...
 
void invert_h2matrix (ph2matrix a, pclusteroperator arwf, pclusteroperator acwf, ph2matrix B, pclusteroperator brwf, pclusteroperator bcwf, ptruncmode tm, real tol)
 $ B \gets a^{-1} $ More...
 
void lowersolve_h2matrix_avector (bool unit, bool atrans, pch2matrix a, pavector x)
 computes $ x \gets T^{-1} \cdot x $, using forward substitution More...
 
void uppersolve_h2matrix_avector (bool unit, bool atrans, pch2matrix a, pavector x)
 computes $ x \gets T^{-1} \cdot x $, using backward substitution More...
 
void lowersolve_h2matrix_amatrix (bool unit, bool atrans, pch2matrix a, bool xtrans, pamatrix X)
 computes $ X \gets T^{-1} \cdot X $ or $ X^* \gets T^{-1} \cdot X^* $, using forward substitution More...
 
void uppersolve_h2matrix_amatrix (bool unit, bool atrans, pch2matrix a, bool xtrans, pamatrix X)
 computes $ X \gets T^{-1} \cdot X $ or $ X^* \gets T^{-1} \cdot X^* $, using backward substitution More...
 
void lowersolve_amatrix_h2matrix (bool aunit, bool atrans, pcamatrix a, bool xytrans, pch2matrix X, ph2matrix Y, pclusteroperator rwf, pclusteroperator cwf, ptruncmode tm, real tol)
 computes $ Y \gets T^{-1} \cdot X $ or $ Y^* \gets T^{-1} \cdot X^* $, using forward substitution More...
 
void uppersolve_amatrix_h2matrix (bool aunit, bool atrans, pcamatrix a, bool ytrans, ph2matrix Y, pclusteroperator rwf, pclusteroperator cwf, ptruncmode tm, real tol)
 computes $ Y \gets T^{-1} \cdot Y $ or $ Y^* \gets T^{-1} \cdot Y^* $, using backward substitution More...
 
void lowersolve_h2matrix (bool aunit, bool atrans, pch2matrix a, bool xytrans, ph2matrix X, pclusteroperator xrwf, pclusteroperator xcwf, ph2matrix Y, pclusteroperator yrwf, pclusteroperator ycwf, ptruncmode tm, real tol)
 computes $ Y \gets T^{-1} \cdot X $ or $ Y^T \gets T^{-1} \cdot X^T $, using forward substitution More...
 
void uppersolve_h2matrix (bool aunit, bool atrans, pch2matrix a, bool ytrans, ph2matrix Y, pclusteroperator yrwf, pclusteroperator ycwf, ptruncmode tm, real tol)
 computes $ Y \gets T^{-1} \cdot Y $ or $ Y^T \gets T^{-1} \cdot Y^T $, using backward substitution More...
 
void lrdecomp_h2matrix (ph2matrix a, pclusteroperator arwf, pclusteroperator acwf, ph2matrix L, pclusteroperator lrwf, pclusteroperator lcwf, ph2matrix R, pclusteroperator rrwf, pclusteroperator rcwf, ptruncmode tm, real tol)
 computes the LR-decomposition $ L \cdot R \gets a $ of a, L has unit diagonal More...
 
void lrsolve_h2matrix_avector (pch2matrix L, pch2matrix R, pavector x)
 $ x \gets R^{-1} \cdot L^{-1} \cdot x $ More...
 
void init_cholesky_h2matrix (ph2matrix a, pclusteroperator*prwf, pclusteroperator*pcwf, ptruncmode tm)
 initialises a, rwf and cwf for choldecomp_h2matrix More...
 
void choldecomp_h2matrix (ph2matrix a, pclusteroperator arwf, pclusteroperator acwf, ph2matrix L, pclusteroperator lrwf, pclusteroperator lcwf, ptruncmode tm, real tol)
 computes the Cholesky decomposition $ L \cdot L^* \gets$ of a and stores the result in L and R More...
 
void cholsolve_h2matrix_avector (pch2matrix a, pavector x)
 $ x \gets L^{-*} \cdot L^{-1} \cdot x $ More...
 

Detailed Description

Author
Knut Reimer
Date
20010-2013

This is the $\mathcal{H}^2$-matrix-arithmetics module of the H2Lib package. It uses the $\mathcal{H}^2$-matrix-structure of h2matrix.h.