HLIBpro  3.0
Real Valued Double Precision

Admissibility condition

hpro_ admcond_t hpro_d_ admcond_alg (const hpro_ adm_t crit, const double eta, const hpro_d_ matrix_t S, const hpro_ permutation_t row_perm, const hpro_ permutation_t col_perm, int *info)
 

Cluster Trees

Root of a cluster tree with stored index permutations.

hpro_ clustertree_t hpro_d_ clt_build_bsp_nd (const hpro_ coord_t coord, const hpro_d_ matrix_t S, const hpro_ bsp_t bsptype, const unsigned int nmin, int *info)
 
hpro_ clustertree_t hpro_d_ clt_build_alg (const hpro_d_ matrix_t S, const hpro_ alg_t algtype, const unsigned int nmin, int *info)
 
hpro_ clustertree_t hpro_d_ clt_build_alg_nd (const hpro_d_ matrix_t S, const hpro_ alg_t algtype, const unsigned int nmin, int *info)
 
hpro_ clustertree_t hpro_d_ clt_build_alg_part (const hpro_d_ matrix_t S, const unsigned int *partition, const hpro_ alg_t algtype, const unsigned int nmin, int *info)
 

Vectors

size_t hpro_d_ vector_size (const hpro_d_ vector_t v, int *info)
 
hpro_d_ vector_t hpro_d_ vector_copy (const hpro_d_ vector_t v, int *info)
 
size_t hpro_d_ vector_bytesize (const hpro_d_ vector_t v, int *info)
 
void hpro_d_ vector_free (hpro_d_ vector_t v, int *info)
 
double hpro_d_ vector_entry_get (const hpro_d_ vector_t x, const size_t i, int *info)
 
void hpro_d_ vector_entry_set (const hpro_d_ vector_t x, const size_t i, const double f, int *info)
 
hpro_d_ vector_t hpro_d_ vector_build (const size_t size, int *info)
 
void hpro_d_ vector_import (hpro_d_ vector_t v, const double *arr, int *info)
 
void hpro_d_ vector_export (const hpro_d_ vector_t v, double *arr, int *info)
 
void hpro_d_ vector_fill (hpro_d_ vector_t x, const double f, int *info)
 
void hpro_d_ vector_fill_rand (hpro_d_ vector_t x, int *info)
 
void hpro_d_ vector_assign (hpro_d_ vector_t y, const hpro_d_ vector_t x, int *info)
 
void hpro_d_ vector_scale (hpro_d_ vector_t x, const double f, int *info)
 
void hpro_d_ vector_axpy (const double alpha, const hpro_d_ vector_t x, hpro_d_ vector_t y, int *info)
 
double hpro_d_ vector_dot (const hpro_d_ vector_t x, const hpro_d_ vector_t y, int *info)
 
double hpro_d_ vector_norm2 (const hpro_d_ vector_t x, int *info)
 
double hpro_d_ vector_norm_inf (const hpro_d_ vector_t x, int *info)
 
void hpro_d_ vector_permute (hpro_d_ vector_t v, const hpro_ permutation_t perm, int *info)
 

Matrices

void hpro_d_ linearoperator_free (hpro_d_ linearoperator_t A, int *info)
 
hpro_d_ vector_t hpro_d_ linearoperator_range_vector (const hpro_d_ linearoperator_t A, int *info)
 
hpro_d_ vector_t hpro_d_ linearoperator_domain_vector (const hpro_d_ linearoperator_t A, int *info)
 
hpro_d_ linearoperator_t hpro_d_ perm_linearoperator (const hpro_ permutation_t P, const hpro_d_ linearoperator_t A, const hpro_ permutation_t R, int *info)
 
void hpro_d_ linearoperator_apply (const hpro_d_ linearoperator_t A, const hpro_d_ vector_t x, hpro_d_ vector_t y, const hpro_ matop_t matop, int *info)
 
void hpro_d_ linearoperator_apply_add (const double alpha, const hpro_d_ linearoperator_t A, const hpro_d_ vector_t x, hpro_d_ vector_t y, const hpro_ matop_t matop, int *info)
 
size_t hpro_d_ matrix_rows (const hpro_d_ matrix_t A, int *info)
 
size_t hpro_d_ matrix_cols (const hpro_d_ matrix_t A, int *info)
 
hpro_d_ matrix_t hpro_d_ matrix_copy (const hpro_d_ matrix_t A, int *info)
 
hpro_d_ matrix_t hpro_d_ matrix_copy_acc (const hpro_d_ matrix_t A, const hpro_ acc_t acc, const int coarsen, int *info)
 
void hpro_d_ matrix_copyto (const hpro_d_ matrix_t A, hpro_d_ matrix_t B, int *info)
 
void hpro_d_ matrix_copyto_acc (const hpro_d_ matrix_t A, hpro_d_ matrix_t B, const hpro_ acc_t acc, const int coarsen, int *info)
 
hpro_d_ matrix_t hpro_d_ matrix_copy_blockdiag (const hpro_d_ matrix_t A, const unsigned int lvl, const size_t blocksize, int *info)
 
hpro_d_ matrix_t hpro_d_ matrix_copy_blockdiag_acc (const hpro_d_ matrix_t A, const hpro_ acc_t acc, const unsigned int lvl, const size_t blocksize, int *info)
 
hpro_d_ matrix_t hpro_d_ matrix_copy_nearfield (const hpro_d_ matrix_t A, const int without_farfield, int *info)
 
void hpro_d_ matrix_restrict_nearfield (hpro_d_ matrix_t A, const int delete_farfield, int *info)
 
void hpro_d_ matrix_nearfield_to_crs (const hpro_d_ matrix_t A, size_t *nnz, int **rowptr, int **colind, double **coeffs, int *info)
 
size_t hpro_d_ matrix_bytesize (const hpro_d_ matrix_t A, int *info)
 
void hpro_d_ matrix_free (hpro_d_ matrix_t A, int *info)
 
int hpro_d_ matrix_is_sym (const hpro_d_ matrix_t A, int *info)
 
int hpro_d_ matrix_is_herm (const hpro_d_ matrix_t A, int *info)
 
int hpro_d_ matrix_is_complex (const hpro_d_ matrix_t A, int *info)
 
double hpro_d_ matrix_entry_get (const hpro_d_ matrix_t A, const size_t i, const size_t j, int *info)
 
hpro_d_ vector_t hpro_d_ matrix_row_vector (const hpro_d_ matrix_t A, int *info)
 
hpro_d_ vector_t hpro_d_ matrix_col_vector (const hpro_d_ matrix_t A, int *info)
 
hpro_d_ matrix_t hpro_d_ matrix_import_crs (const size_t rows, const size_t cols, const size_t nnz, const int *rowptr, const int *colind, const double *coeffs, const int sym, int *info)
 
hpro_d_ matrix_t hpro_d_ matrix_import_ccs (const size_t rows, const size_t cols, const size_t nnz, const int *colptr, const int *rowind, const double *coeffs, const int sym, int *info)
 
hpro_d_ matrix_t hpro_d_ matrix_import_dense (const size_t rows, const size_t cols, const double *D, const int sym, int *info)
 
hpro_d_ matrix_t hpro_d_ matrix_build_sparse (const hpro_ blockclustertree_t bct, const hpro_d_ matrix_t S, const hpro_ acc_t acc, int *info)
 
hpro_d_ matrix_t hpro_d_ matrix_build_dense (const hpro_ blockclustertree_t bct, const hpro_d_ matrix_t D, const hpro_ lrapx_t lrapx, const hpro_ acc_t acc, int *info)
 
hpro_d_ matrix_t hpro_d_ matrix_build_coeff (const hpro_ blockclustertree_t bct, const hpro_d_ coeff_t f, void *arg, const hpro_ lrapx_t lrapx, const hpro_ acc_t acc, const int sym, int *info)
 
hpro_d_ matrix_t hpro_d_ matrix_build_identity (const hpro_ blockclustertree_t bct, int *info)
 
hpro_d_ matrix_t hpro_d_ matrix_assemble_block (const size_t brows, const size_t bcols, hpro_d_ matrix_t *submat, const int copy, int *info)
 
hpro_ permutation_t hpro_d_ matrix_row_perm_i2e (const hpro_d_ matrix_t A, int *info)
 
hpro_ permutation_t hpro_d_ matrix_col_perm_i2e (const hpro_d_ matrix_t A, int *info)
 
hpro_ permutation_t hpro_d_ matrix_row_perm_e2i (const hpro_d_ matrix_t A, int *info)
 
hpro_ permutation_t hpro_d_ matrix_col_perm_e2i (const hpro_d_ matrix_t A, int *info)
 
void hpro_d_ matrix_permute (hpro_d_ matrix_t A, const hpro_ permutation_t row_perm, const hpro_ permutation_t col_perm, int *info)
 
hpro_d_ linearoperator_t hpro_d_ perm_matrix (const hpro_ permutation_t P, const hpro_d_ matrix_t A, const hpro_ permutation_t R, int *info)
 
hpro_d_ linearoperator_t hpro_d_ matrix_product2 (const hpro_ matop_t matop1, const hpro_d_ linearoperator_t A1, const hpro_ matop_t matop2, const hpro_d_ linearoperator_t A2, int *info)
 
hpro_d_ linearoperator_t hpro_d_ matrix_product3 (const hpro_ matop_t matop1, const hpro_d_ linearoperator_t A1, const hpro_ matop_t matop2, const hpro_d_ linearoperator_t A2, const hpro_ matop_t matop3, const hpro_d_ linearoperator_t A3, int *info)
 
hpro_d_ linearoperator_t hpro_d_ matrix_sum2 (const double alpha1, const hpro_ matop_t matop1, const hpro_d_ linearoperator_t A1, const double alpha2, const hpro_ matop_t matop2, const hpro_d_ linearoperator_t A2, int *info)
 
hpro_d_ linearoperator_t hpro_d_ matrix_sum3 (const double alpha1, const hpro_ matop_t matop1, const hpro_d_ linearoperator_t A1, const double alpha2, const hpro_ matop_t matop2, const hpro_d_ linearoperator_t A2, const double alpha3, const hpro_ matop_t matop3, const hpro_d_ linearoperator_t A3, int *info)
 
hpro_d_ matrix_t hpro_d_ matrix_to_dense (const hpro_d_ matrix_t A, int *info)
 
hpro_d_ matrix_t hpro_d_ matrix_to_rank (const hpro_d_ matrix_t A, const hpro_ acc_t acc, int *info)
 
hpro_d_ matrix_t hpro_d_ matrix_approx_rank_aca (const hpro_d_ matrix_t A, const hpro_ acc_t acc, int *info)
 
void hpro_d_ matrix_print_ps (const hpro_d_ matrix_t A, const char *filename, const int options, int *info)
 

Algebra

void hpro_d_ matrix_mulvec (const double alpha, const hpro_d_ matrix_t A, const hpro_d_ vector_t x, const double beta, hpro_d_ vector_t y, const hpro_ matop_t matop, int *info)
 
void hpro_d_ matrix_transpose (const hpro_d_ matrix_t A, int *info)
 
void hpro_d_ matrix_scale (const double f, const hpro_d_ matrix_t A, int *info)
 
void hpro_d_ matrix_add (const double alpha, const hpro_d_ matrix_t A, const double beta, hpro_d_ matrix_t B, const hpro_ acc_t acc, int *info)
 
void hpro_d_ matrix_mul (const double alpha, const hpro_ matop_t matop_A, const hpro_d_ matrix_t A, const hpro_ matop_t matop_B, const hpro_d_ matrix_t B, const double beta, hpro_d_ matrix_t C, const hpro_ acc_t acc, int *info)
 
hpro_d_ matrix_t hpro_d_ matrix_inv (hpro_d_ matrix_t A, const hpro_ acc_t acc, int *info)
 
hpro_d_ vector_t hpro_d_ matrix_inv_diag (hpro_d_ matrix_t A, const hpro_ acc_t acc, int *info)
 
hpro_d_ linearoperator_t hpro_d_ matrix_factorise (hpro_d_ matrix_t A, const hpro_ acc_t acc, int *info)
 
hpro_d_ linearoperator_t hpro_d_ matrix_factorise_inv (hpro_d_ matrix_t A, const hpro_ acc_t acc, int *info)
 
hpro_d_ linearoperator_t hpro_d_ matrix_lu (hpro_d_ matrix_t A, const hpro_ acc_t acc, int *info)
 
hpro_d_ linearoperator_t hpro_d_ matrix_ldl (hpro_d_ matrix_t A, const hpro_ acc_t acc, int *info)
 
hpro_d_ linearoperator_t hpro_d_ matrix_ll (hpro_d_ matrix_t A, const hpro_ acc_t acc, int *info)
 
hpro_d_ linearoperator_t hpro_d_ matrix_chol (hpro_d_ matrix_t A, const hpro_ acc_t acc, int *info)
 
hpro_d_ linearoperator_t hpro_d_ matrix_waz (hpro_d_ matrix_t A, const hpro_ acc_t acc, int *info)
 
hpro_d_ linearoperator_t hpro_d_ matrix_lu_inv (hpro_d_ matrix_t A, const hpro_ acc_t acc, int *info)
 
hpro_d_ linearoperator_t hpro_d_ matrix_ldl_inv (hpro_d_ matrix_t A, const hpro_ acc_t acc, int *info)
 
hpro_d_ linearoperator_t hpro_d_ matrix_ll_inv (hpro_d_ matrix_t A, const hpro_ acc_t acc, int *info)
 
hpro_d_ linearoperator_t hpro_d_ matrix_chol_inv (hpro_d_ matrix_t A, const hpro_ acc_t acc, int *info)
 
void hpro_d_ matrix_solve_lower_left (const hpro_d_ matrix_t L, const hpro_ matop_t op_L, hpro_d_ matrix_t A, const hpro_ acc_t acc, const int pointwise, const int unit_diag, int *info)
 
void hpro_d_ matrix_solve_lower_right (hpro_d_ matrix_t A, const hpro_d_ matrix_t L, const hpro_ matop_t op_L, const hpro_ acc_t acc, const int pointwise, const int unit_diag, int *info)
 
void hpro_d_ matrix_solve_upper_right (hpro_d_ matrix_t A, const hpro_d_ matrix_t U, const hpro_ acc_t acc, const int pointwise, const int unit_diag, int *info)
 
void hpro_d_ matrix_solve_diag_left (const hpro_d_ matrix_t D, const hpro_ matop_t op_D, hpro_d_ matrix_t A, const hpro_ acc_t acc, const int pointwise, int *info)
 
void hpro_d_ matrix_solve_diag_right (hpro_d_ matrix_t A, const hpro_d_ matrix_t D, const hpro_ matop_t op_D, const hpro_ acc_t acc, const int pointwise, int *info)
 
void hpro_d_ matrix_add_identity (hpro_d_ matrix_t A, const double lambda, int *info)
 
double hpro_d_ matrix_norm_frobenius (const hpro_d_ matrix_t A, int *info)
 
double hpro_d_ matrix_norm_frobenius_diff (const hpro_d_ matrix_t A, const hpro_d_ matrix_t B, int *info)
 
double hpro_d_ matrix_norm_spectral (const hpro_d_ matrix_t A, int *info)
 
double hpro_d_ matrix_norm_spectral_inv (const hpro_d_ matrix_t A, int *info)
 
double hpro_d_ matrix_norm_spectral_diff (const hpro_d_ matrix_t A, const hpro_d_ matrix_t B, int *info)
 
double hpro_d_ matrix_norm_inv_approx (const hpro_d_ matrix_t A, const hpro_d_ matrix_t B, int *info)
 
double hpro_d_ linearoperator_norm_spectral (const hpro_d_ linearoperator_t A, int *info)
 
double hpro_d_ linearoperator_norm_spectral_inv (const hpro_d_ linearoperator_t A, int *info)
 
double hpro_d_ linearoperator_norm_spectral_diff (const hpro_d_ linearoperator_t A, const hpro_d_ linearoperator_t B, int *info)
 
double hpro_d_ linearoperator_norm_inv_approx (const hpro_d_ linearoperator_t A, const hpro_d_ linearoperator_t B, int *info)
 

Solver

void hpro_d_ solver_solve (const hpro_ solver_t solver, const hpro_d_ linearoperator_t A, hpro_d_ vector_t x, const hpro_d_ vector_t b, const hpro_d_ linearoperator_t W, hpro_ solve_info_t *solve_info, int *info)
 
void hpro_ ds_solver_solve (const hpro_ solver_t solver, const hpro_d_ linearoperator_t A, hpro_d_ vector_t x, const hpro_d_ vector_t b, const hpro_s_ linearoperator_t W, hpro_ solve_info_t *solve_info, int *info)
 

Input/Output

hpro_d_ matrix_t hpro_d_ hformat_load_matrix (const char *filename, int *info)
 
hpro_d_ linearoperator_t hpro_d_ hformat_load_linearoperator (const char *filename, int *info)
 
hpro_d_ vector_t hpro_d_ hformat_load_vector (const char *filename, int *info)
 
void hpro_d_ hformat_save_linearoperator (const hpro_d_ linearoperator_t A, const char *filename, int *info)
 
void hpro_d_ hformat_save_matrix (const hpro_d_ matrix_t A, const char *filename, int *info)
 
void hpro_d_ hformat_save_vector (const hpro_d_ vector_t v, const char *filename, int *info)
 
hpro_d_ matrix_t hpro_d_ samg_load_matrix (const char *filename, int *info)
 
void hpro_d_ samg_save_matrix (const hpro_d_ matrix_t A, const char *filename, int *info)
 
hpro_d_ vector_t hpro_d_ samg_load_vector (const char *filename, int *info)
 
void hpro_d_ samg_save_vector (const hpro_d_ vector_t x, const char *filename, int *info)
 
hpro_d_ matrix_t hpro_d_ matlab_load_matrix (const char *filename, const char *matname, int *info)
 
void hpro_d_ matlab_save_matrix (const hpro_d_ matrix_t M, const char *filename, const char *matname, int *info)
 
hpro_d_ vector_t hpro_d_ matlab_load_vector (const char *filename, const char *vecname, int *info)
 
void hpro_d_ matlab_save_vector (const hpro_d_ vector_t v, const char *filename, const char *vecname, int *info)
 
hpro_d_ matrix_t hpro_d_ pltmg_load_matrix (const char *filename, int *info)
 
hpro_d_ matrix_t hpro_d_ hb_load_matrix (const char *filename, int *info)
 
void hpro_d_ hb_save_matrix (const hpro_d_ matrix_t M, const char *filename, int *info)
 
hpro_d_ matrix_t hpro_d_ mm_load_matrix (const char *filename, int *info)
 
hpro_d_ matrix_t hpro_d_ load_matrix (const char *filename, int *info)
 
hpro_d_ vector_t hpro_d_ load_vector (const char *filename, int *info)
 

Detailed Description

This modules defines double precision real valued functions for the usage of 𝖧𝖫𝖨𝖡𝗉𝗋𝗈 in C programs.

Function Documentation

◆ admcond_alg()

hpro_ admcond_t hpro_d_ admcond_alg ( const hpro_ adm_t  crit,
const double  eta,
const hpro_d_ matrix_t  S,
const hpro_ permutation_t  row_perm,
const hpro_ permutation_t  col_perm,
int *  info 
)

create admissibility condition based on algebraic connectivity in S

  • row_perm and col_perm define the (optional) mapping between internal to external ordering for row and column index sets, respectively (see hlib_ct_perm_i2e); set to NULL if not present

◆ clt_build_alg()

hpro_ clustertree_t hpro_d_ clt_build_alg ( const hpro_d_ matrix_t  S,
const hpro_ alg_t  algtype,
const unsigned int  nmin,
int *  info 
)

build clustertree via algebraic clustering based on given sparse matrix S and partitioning type algtype.

◆ clt_build_alg_nd()

hpro_ clustertree_t hpro_d_ clt_build_alg_nd ( const hpro_d_ matrix_t  S,
const hpro_ alg_t  algtype,
const unsigned int  nmin,
int *  info 
)

build clustertree via algebraic clustering with nested dissection based on given sparse matrix S and partitioning type algtype.

◆ clt_build_alg_part()

hpro_ clustertree_t hpro_d_ clt_build_alg_part ( const hpro_d_ matrix_t  S,
const unsigned int *  partition,
const hpro_ alg_t  algtype,
const unsigned int  nmin,
int *  info 
)

build clustertree via algebraic clustering based on matrix S and use a predefined partition for first step

  • partition must have size "rows(S)" and contain ids ∈ [0,k-1], where k is the total number of partitions
  • index "i" will go into group partition[i]

◆ clt_build_bsp_nd()

hpro_ clustertree_t hpro_d_ clt_build_bsp_nd ( const hpro_ coord_t  coord,
const hpro_d_ matrix_t  S,
const hpro_ bsp_t  bsptype,
const unsigned int  nmin,
int *  info 
)

build cluster tree using binary space partitioning and nested dissection based on given sparse matrix S

  • for period see HPRO_D_PF(clt_build_bsp)
Parameters
coordcoordinates of indices
Ssparse mat. for connectivity
bsptypepartitioning type
nminminimal cluster size

◆ ds_solver_solve()

void hpro_ ds_solver_solve ( const hpro_ solver_t  solver,
const hpro_d_ linearoperator_t  A,
hpro_d_ vector_t  x,
const hpro_d_ vector_t  b,
const hpro_s_ linearoperator_t  W,
hpro_ solve_info_t *  solve_info,
int *  info 
)

mixed precision solving: A is double precision but W is single precision

◆ hb_load_matrix()

hpro_d_ matrix_t hpro_d_ hb_load_matrix ( const char *  filename,
int *  info 
)

read matrix from file filename

◆ hb_save_matrix()

void hpro_d_ hb_save_matrix ( const hpro_d_ matrix_t  M,
const char *  filename,
int *  info 
)

save matrix M to file filename

◆ hformat_load_linearoperator()

hpro_d_ linearoperator_t hpro_d_ hformat_load_linearoperator ( const char *  filename,
int *  info 
)

read matrix from file filename

◆ hformat_load_matrix()

hpro_d_ matrix_t hpro_d_ hformat_load_matrix ( const char *  filename,
int *  info 
)

read matrix from file filename

◆ hformat_load_vector()

hpro_d_ vector_t hpro_d_ hformat_load_vector ( const char *  filename,
int *  info 
)

read vector from file filename

◆ hformat_save_linearoperator()

void hpro_d_ hformat_save_linearoperator ( const hpro_d_ linearoperator_t  A,
const char *  filename,
int *  info 
)

save linear operator A to file filename

◆ hformat_save_matrix()

void hpro_d_ hformat_save_matrix ( const hpro_d_ matrix_t  A,
const char *  filename,
int *  info 
)

save matrix A to file filename

◆ hformat_save_vector()

void hpro_d_ hformat_save_vector ( const hpro_d_ vector_t  v,
const char *  filename,
int *  info 
)

save vector v to file filename

◆ linearoperator_apply()

void hpro_d_ linearoperator_apply ( const hpro_d_ linearoperator_t  A,
const hpro_d_ vector_t  x,
hpro_d_ vector_t  y,
const hpro_ matop_t  matop,
int *  info 
)

computes \( y := op(A) \cdot x \)

◆ linearoperator_apply_add()

void hpro_d_ linearoperator_apply_add ( const double  alpha,
const hpro_d_ linearoperator_t  A,
const hpro_d_ vector_t  x,
hpro_d_ vector_t  y,
const hpro_ matop_t  matop,
int *  info 
)

computes \( y := y + \alpha op(A) \cdot x \)

◆ linearoperator_domain_vector()

hpro_d_ vector_t hpro_d_ linearoperator_domain_vector ( const hpro_d_ linearoperator_t  A,
int *  info 
)

create vectors from the domain of the linear operator A

◆ linearoperator_free()

void hpro_d_ linearoperator_free ( hpro_d_ linearoperator_t  A,
int *  info 
)

free resources coupled with linear operator A

◆ linearoperator_norm_inv_approx()

double hpro_d_ linearoperator_norm_inv_approx ( const hpro_d_ linearoperator_t  A,
const hpro_d_ linearoperator_t  B,
int *  info 
)

compute inversion error |I-BA|_2, where B is supposed to be A^-1

◆ linearoperator_norm_spectral()

double hpro_d_ linearoperator_norm_spectral ( const hpro_d_ linearoperator_t  A,
int *  info 
)

compute spectral norm of matrix

◆ linearoperator_norm_spectral_diff()

double hpro_d_ linearoperator_norm_spectral_diff ( const hpro_d_ linearoperator_t  A,
const hpro_d_ linearoperator_t  B,
int *  info 
)

compute relative spectral norm of A-B, e.g. |A-B|_2/|A|_2

◆ linearoperator_norm_spectral_inv()

double hpro_d_ linearoperator_norm_spectral_inv ( const hpro_d_ linearoperator_t  A,
int *  info 
)

compute spectral norm of inverse matrix

◆ linearoperator_range_vector()

hpro_d_ vector_t hpro_d_ linearoperator_range_vector ( const hpro_d_ linearoperator_t  A,
int *  info 
)

create vectors from the range of the linear operator A

◆ load_matrix()

hpro_d_ matrix_t hpro_d_ load_matrix ( const char *  filename,
int *  info 
)

read matrix from file filename

◆ load_vector()

hpro_d_ vector_t hpro_d_ load_vector ( const char *  filename,
int *  info 
)

read vector from file filename

◆ matlab_load_matrix()

hpro_d_ matrix_t hpro_d_ matlab_load_matrix ( const char *  filename,
const char *  matname,
int *  info 
)

read matrix named matname in Matlab V7 format from file filename; if matname == NULL or matname == "", the first matrix in filename will be read

◆ matlab_load_vector()

hpro_d_ vector_t hpro_d_ matlab_load_vector ( const char *  filename,
const char *  vecname,
int *  info 
)

read vector named vname from file filename in Matlab V7 format

◆ matlab_save_matrix()

void hpro_d_ matlab_save_matrix ( const hpro_d_ matrix_t  M,
const char *  filename,
const char *  matname,
int *  info 
)

save matrix M named matname in Matlab V7 format to file filename

◆ matlab_save_vector()

void hpro_d_ matlab_save_vector ( const hpro_d_ vector_t  v,
const char *  filename,
const char *  vecname,
int *  info 
)

save vector named vname to file filename in Matlab V7 format

◆ matrix_add()

void hpro_d_ matrix_add ( const double  alpha,
const hpro_d_ matrix_t  A,
const double  beta,
hpro_d_ matrix_t  B,
const hpro_ acc_t  acc,
int *  info 
)

compute B := alpha A + beta B with block-wise precision acc

◆ matrix_add_identity()

void hpro_d_ matrix_add_identity ( hpro_d_ matrix_t  A,
const double  lambda,
int *  info 
)

compute A := A + lambda I

◆ matrix_approx_rank_aca()

hpro_d_ matrix_t hpro_d_ matrix_approx_rank_aca ( const hpro_d_ matrix_t  A,
const hpro_ acc_t  acc,
int *  info 
)

convert matrix to lowrank format using approximation by ACA with accuracy defined by acc

  • currently only dense matrices supported

◆ matrix_assemble_block()

hpro_d_ matrix_t hpro_d_ matrix_assemble_block ( const size_t  brows,
const size_t  bcols,
hpro_d_ matrix_t *  submat,
const int  copy,
int *  info 
)

assemble brows × bcols block matrix out of sub blocks submat

  • sub matrices have to be stored column wise in submat, e.g.,

    \begin{eqnarray*} A_{00} & A_{01} \\ A_{10} & A_{11} \end{eqnarray*}

    is stored { A_00, A_10, A_01, A_11 },
  • if copy is non-zero, sub matrices are copied, otherwise used directly in new block matrix,
  • index sets of sub matrices are adjusted to be consistent in new block matrix,
  • if sub matrices are H-matrices, permutations from internal to external ordering (and vice versa) are also adjusted.
Parameters
browsnumber of block rows
bcolsnumber of block columns
submatarray of sub matrices
copycopy matrices if true

◆ matrix_build_coeff()

hpro_d_ matrix_t hpro_d_ matrix_build_coeff ( const hpro_ blockclustertree_t  bct,
const hpro_d_ coeff_t  f,
void *  arg,
const hpro_ lrapx_t  lrapx,
const hpro_ acc_t  acc,
const int  sym,
int *  info 
)

build H-matrix over block cluster tree bct of precision acc out of dense matrix given by a matrix coefficient function f by using a low-rank approximation for admissible blocks

Parameters
bctblock ct. defining indexsets
fcoeff. function defining mat.
argadd. argument to coeff. fn.
lrapxlow-rank approximation method
accaccuracy of the approximation
symif non-zero, matrix is sym.

◆ matrix_build_dense()

hpro_d_ matrix_t hpro_d_ matrix_build_dense ( const hpro_ blockclustertree_t  bct,
const hpro_d_ matrix_t  D,
const hpro_ lrapx_t  lrapx,
const hpro_ acc_t  acc,
int *  info 
)

build H-matrix over block cluster tree bct with contents defined by given dense matrix D with low-rank approximation method lrapx and block-wise accuracy acc

Parameters
bctblock ct. defining indexsets
Ddense mat. to be converted
lrapxlow-rank approximation method
accaccuracy of the approximation

◆ matrix_build_identity()

hpro_d_ matrix_t hpro_d_ matrix_build_identity ( const hpro_ blockclustertree_t  bct,
int *  info 
)

build H-matrix over block cluster tree bct representing identity

◆ matrix_build_sparse()

hpro_d_ matrix_t hpro_d_ matrix_build_sparse ( const hpro_ blockclustertree_t  bct,
const hpro_d_ matrix_t  S,
const hpro_ acc_t  acc,
int *  info 
)

build H-matrix over block cluster tree bct with contents defined by given sparse matrix S

Parameters
bctblock ct. defining indexsets
Ssparse mat. to be converted
accaccuracy of the approximation

◆ matrix_bytesize()

size_t hpro_d_ matrix_bytesize ( const hpro_d_ matrix_t  A,
int *  info 
)

return size of matrix in bytes

◆ matrix_col_vector()

hpro_d_ vector_t hpro_d_ matrix_col_vector ( const hpro_d_ matrix_t  A,
int *  info 
)

create vectors matching the column indexsets of matrix A

◆ matrix_cols()

size_t hpro_d_ matrix_cols ( const hpro_d_ matrix_t  A,
int *  info 
)

return number of columns of matrix

◆ matrix_copy()

hpro_d_ matrix_t hpro_d_ matrix_copy ( const hpro_d_ matrix_t  A,
int *  info 
)

return copy of matrix A

◆ matrix_copy_acc()

hpro_d_ matrix_t hpro_d_ matrix_copy_acc ( const hpro_d_ matrix_t  A,
const hpro_ acc_t  acc,
const int  coarsen,
int *  info 
)

return copy of matrix A with accuracy acc; if coarsen != 0, copy is coarsend

◆ matrix_copy_blockdiag()

hpro_d_ matrix_t hpro_d_ matrix_copy_blockdiag ( const hpro_d_ matrix_t  A,
const unsigned int  lvl,
const size_t  blocksize,
int *  info 
)

copy blockdiagonal part of first lvl levels or blocks of size blocksize of matrix

◆ matrix_copy_blockdiag_acc()

hpro_d_ matrix_t hpro_d_ matrix_copy_blockdiag_acc ( const hpro_d_ matrix_t  A,
const hpro_ acc_t  acc,
const unsigned int  lvl,
const size_t  blocksize,
int *  info 
)

copy blockdiagonal part of first lvl levels or blocks of size blocksize of matrix with given accuracy acc

◆ matrix_copy_nearfield()

hpro_d_ matrix_t hpro_d_ matrix_copy_nearfield ( const hpro_d_ matrix_t  A,
const int  without_farfield,
int *  info 
)

copy nearfield part of matrix

  • if without_farfield is non-zero, farfield matrices (e.g., low-rank) will be skipped, otherwise set to rank 0

◆ matrix_copyto()

void hpro_d_ matrix_copyto ( const hpro_d_ matrix_t  A,
hpro_d_ matrix_t  B,
int *  info 
)

copy matrix A to B

◆ matrix_copyto_acc()

void hpro_d_ matrix_copyto_acc ( const hpro_d_ matrix_t  A,
hpro_d_ matrix_t  B,
const hpro_ acc_t  acc,
const int  coarsen,
int *  info 
)

copy matrix A to B with accuracy acc; if coarsen != 0, B is coarsend

◆ matrix_entry_get()

double hpro_d_ matrix_entry_get ( const hpro_d_ matrix_t  A,
const size_t  i,
const size_t  j,
int *  info 
)

get single entry at position (i,j) in matrix A

◆ matrix_factorise()

hpro_d_ linearoperator_t hpro_d_ matrix_factorise ( hpro_d_ matrix_t  A,
const hpro_ acc_t  acc,
int *  info 
)

factorise matrix A up to block-wise precision acc

  • depending on form of A, e.g. if unsymmetric, symmetric or hermitian, an appropriate factorisation method is chosen
  • A will be overwritten by the factors
  • the return value is a matrix which can be used to evaluate the factorisation, e.g. for matrix-vector mult. (this is not possible with A after factorisation as it holds just the data, albeit, it is still needed)

◆ matrix_factorise_inv()

hpro_d_ linearoperator_t hpro_d_ matrix_factorise_inv ( hpro_d_ matrix_t  A,
const hpro_ acc_t  acc,
int *  info 
)

same as hlib_matrix_factorise but return matrix object for evaluation of inverse of A.

◆ matrix_free()

void hpro_d_ matrix_free ( hpro_d_ matrix_t  A,
int *  info 
)

free resources coupled with matrix A

◆ matrix_import_ccs()

hpro_d_ matrix_t hpro_d_ matrix_import_ccs ( const size_t  rows,
const size_t  cols,
const size_t  nnz,
const int *  colptr,
const int *  rowind,
const double *  coeffs,
const int  sym,
int *  info 
)

import given real sparse matrix in CCS format into internal sparse matrix format

  • the data is copied into internal representation
Parameters
rowsnumber of rows and
colscolumns of the matrix
nnznumber of non-zero coeff.
colptrarray of column pointers
rowindarray of row indices
coeffsarray of matrix coefficients
symif non-zero, matrix is sym.

◆ matrix_import_crs()

hpro_d_ matrix_t hpro_d_ matrix_import_crs ( const size_t  rows,
const size_t  cols,
const size_t  nnz,
const int *  rowptr,
const int *  colind,
const double *  coeffs,
const int  sym,
int *  info 
)

import given real sparse matrix in CRS format into internal sparse matrix format

  • the data is copied into internal representation
Parameters
rowsnumber of rows and
colscolumns of the matrix
nnznumber of non-zero coeff.
rowptrarray of row pointers
colindarray of column indices
coeffsarray of matrix coefficients
symif non-zero, matrix is sym.

◆ matrix_import_dense()

hpro_d_ matrix_t hpro_d_ matrix_import_dense ( const size_t  rows,
const size_t  cols,
const double *  D,
const int  sym,
int *  info 
)

import given real dense matrix in column major format into internal dense matrix format

  • the data is copied into internal representation
Parameters
rowsnumber of rows and
colscolumns of the matrix
Daccuracy of the approximation
symif non-zero, matrix is sym.

◆ matrix_inv()

hpro_d_ matrix_t hpro_d_ matrix_inv ( hpro_d_ matrix_t  A,
const hpro_ acc_t  acc,
int *  info 
)

compute inverse of A with block-wise precision acc; A will be overwritten with A^-1 upon exit

  • the return value is either A, or a new representation for the inverse of A (A is still needed then)

◆ matrix_inv_diag()

hpro_d_ vector_t hpro_d_ matrix_inv_diag ( hpro_d_ matrix_t  A,
const hpro_ acc_t  acc,
int *  info 
)

compute diagonal of inverse of A with local accuracy acc

  • A will be overwritten during computation

◆ matrix_is_complex()

int hpro_d_ matrix_is_complex ( const hpro_d_ matrix_t  A,
int *  info 
)

return 1 if matrix A is complex valued and 0 otherwise

◆ matrix_is_herm()

int hpro_d_ matrix_is_herm ( const hpro_d_ matrix_t  A,
int *  info 
)

return 1 if matrix A is hermitian, e.g. A = A^H

◆ matrix_is_sym()

int hpro_d_ matrix_is_sym ( const hpro_d_ matrix_t  A,
int *  info 
)

return 1 if matrix A is symmetric, e.g. A = A^T

◆ matrix_ldl()

hpro_d_ linearoperator_t hpro_d_ matrix_ldl ( hpro_d_ matrix_t  A,
const hpro_ acc_t  acc,
int *  info 
)

LDL factorisation of matrix A up to block-wise precision acc

  • A will be overwritten by L*D*L^T, which only allows matrix vector mult.
  • in case of hermitian matrices, LDL^H is computed
  • for return value see "hlib_matrix_lu"

◆ matrix_ldl_inv()

hpro_d_ linearoperator_t hpro_d_ matrix_ldl_inv ( hpro_d_ matrix_t  A,
const hpro_ acc_t  acc,
int *  info 
)

LDL factorisation of matrix A up to block-wise precision acc

  • A will be overwritten by L*D*L^T, which only allows matrix vector mult.
  • in case of hermitian matrices, LDL^H is computed
  • for return value see "hlib_matrix_lu_inv"

◆ matrix_ll()

hpro_d_ linearoperator_t hpro_d_ matrix_ll ( hpro_d_ matrix_t  A,
const hpro_ acc_t  acc,
int *  info 
)

Cholesky factorisation of matrix A up to block-wise precision acc

  • A will be overwritten by L*L^T, which only allows matrix vector mult.
  • in case of hermitian matrices, L*L^H is computed
  • for return value see "hlib_matrix_lu"

◆ matrix_ll_inv()

hpro_d_ linearoperator_t hpro_d_ matrix_ll_inv ( hpro_d_ matrix_t  A,
const hpro_ acc_t  acc,
int *  info 
)

Cholesky factorisation of matrix A up to block-wise precision acc

  • A will be overwritten by L*L^T, which only allows matrix vector mult.
  • in case of hermitian matrices, L*L^H is computed
  • for return value see "hlib_matrix_lu_inv"

◆ matrix_lu()

hpro_d_ linearoperator_t hpro_d_ matrix_lu ( hpro_d_ matrix_t  A,
const hpro_ acc_t  acc,
int *  info 
)

LU factorise matrix A up to block-wise precision acc

  • A will be overwritten by L*U, which only allows matrix vector mult.
  • the return value is a matrix which can be used to evaluate LU, e.g. for matrix-vector multiplication (this is not possible with A after factorisation as it holds just the data, albeit, it is still needed)

◆ matrix_lu_inv()

hpro_d_ linearoperator_t hpro_d_ matrix_lu_inv ( hpro_d_ matrix_t  A,
const hpro_ acc_t  acc,
int *  info 
)

LU factorise matrix A up to block-wise precision acc

  • A will be overwritten by L*U, which only allows matrix vector mult.
  • the return value is a matrix which can be used to evaluate (LU)^-1, e.g. for matrix-vector multiplication (this is not possible with A after factorisation as it holds just the data, albeit, it is still needed)

◆ matrix_mul()

void hpro_d_ matrix_mul ( const double  alpha,
const hpro_ matop_t  matop_A,
const hpro_d_ matrix_t  A,
const hpro_ matop_t  matop_B,
const hpro_d_ matrix_t  B,
const double  beta,
hpro_d_ matrix_t  C,
const hpro_ acc_t  acc,
int *  info 
)

compute C := alpha op(A) op(B) + beta C with block-wise precision acc

◆ matrix_mulvec()

void hpro_d_ matrix_mulvec ( const double  alpha,
const hpro_d_ matrix_t  A,
const hpro_d_ vector_t  x,
const double  beta,
hpro_d_ vector_t  y,
const hpro_ matop_t  matop,
int *  info 
)

compute y := alpha op(A) x + beta y, where op(A) is either A or A^T

Parameters
xsource vector of dimension columns(op(A))
ydest. vector of dimension rows(op(A))
matopmatrix operation, e.g. transpose

◆ matrix_nearfield_to_crs()

void hpro_d_ matrix_nearfield_to_crs ( const hpro_d_ matrix_t  A,
size_t *  nnz,
int **  rowptr,
int **  colind,
double **  coeffs,
int *  info 
)

returns nearfield part of given matrix in CRS format

  • CRS is stored in nnz, rowptr, colind and coeffs
Parameters
nnznumber of non-zero coeff.
rowptrarray of row pointers
colindarray of column indices
coeffsarray of matrix coefficients

◆ matrix_norm_frobenius()

double hpro_d_ matrix_norm_frobenius ( const hpro_d_ matrix_t  A,
int *  info 
)

compute Frobenius norm of matrix

◆ matrix_norm_frobenius_diff()

double hpro_d_ matrix_norm_frobenius_diff ( const hpro_d_ matrix_t  A,
const hpro_d_ matrix_t  B,
int *  info 
)

compute Frobenius norm of A-B

◆ matrix_norm_inv_approx()

double hpro_d_ matrix_norm_inv_approx ( const hpro_d_ matrix_t  A,
const hpro_d_ matrix_t  B,
int *  info 
)

compute inversion error |I-BA|_2, where B is supposed to be A^-1

◆ matrix_norm_spectral()

double hpro_d_ matrix_norm_spectral ( const hpro_d_ matrix_t  A,
int *  info 
)

compute spectral norm of matrix

◆ matrix_norm_spectral_diff()

double hpro_d_ matrix_norm_spectral_diff ( const hpro_d_ matrix_t  A,
const hpro_d_ matrix_t  B,
int *  info 
)

compute relative spectral norm of A-B, e.g. |A-B|_2/|A|_2

◆ matrix_norm_spectral_inv()

double hpro_d_ matrix_norm_spectral_inv ( const hpro_d_ matrix_t  A,
int *  info 
)

compute spectral norm of inverse matrix

◆ matrix_permute()

void hpro_d_ matrix_permute ( hpro_d_ matrix_t  A,
const hpro_ permutation_t  row_perm,
const hpro_ permutation_t  col_perm,
int *  info 
)

reorder matrix A entries with given row and column permutations

  • only supported for dense matrices
  • NULL permutation will be handled as identity

◆ matrix_print_ps()

void hpro_d_ matrix_print_ps ( const hpro_d_ matrix_t  A,
const char *  filename,
const int  options,
int *  info 
)

print matrix A in postscript format to file filename options (may be ORed) :

  • HPRO_MATIO_SVD : print singular value decomposition in each block
  • HPRO_MATIO_ENTRY : print each entry of matrix
  • HPRO_MATIO_PATTERN : print sparsity pattern (non-zero entries)
Parameters
optionsoptions defining output

◆ matrix_product2()

hpro_d_ linearoperator_t hpro_d_ matrix_product2 ( const hpro_ matop_t  matop1,
const hpro_d_ linearoperator_t  A1,
const hpro_ matop_t  matop2,
const hpro_d_ linearoperator_t  A2,
int *  info 
)

represents product of matrices \(A_1 \cdot A_2\)

  • all matrix objects are only referenced, not copied

◆ matrix_product3()

hpro_d_ linearoperator_t hpro_d_ matrix_product3 ( const hpro_ matop_t  matop1,
const hpro_d_ linearoperator_t  A1,
const hpro_ matop_t  matop2,
const hpro_d_ linearoperator_t  A2,
const hpro_ matop_t  matop3,
const hpro_d_ linearoperator_t  A3,
int *  info 
)

represents product of matrices \(A_1 \cdot A_2 \cdot A3\)

  • all matrix objects are only referenced, not copied

◆ matrix_restrict_nearfield()

void hpro_d_ matrix_restrict_nearfield ( hpro_d_ matrix_t  A,
const int  delete_farfield,
int *  info 
)

restrict given matrix to nearfield part

  • if delete_farfield is non-zero, farfield matrices (e.g., low-rank) will be removed, otherwise set to rank 0

◆ matrix_row_perm_e2i()

hpro_ permutation_t hpro_d_ matrix_row_perm_e2i ( const hpro_d_ matrix_t  A,
int *  info 
)

return mapping from external (user) to internal (in H-matrix) ordering for rows and columns.

  • returned object is reference to internal permutation: do NOT free
  • if matrix has no permutation, the return value is NULL

◆ matrix_row_perm_i2e()

hpro_ permutation_t hpro_d_ matrix_row_perm_i2e ( const hpro_d_ matrix_t  A,
int *  info 
)

return mapping from internal (in H-matrix) to external (user) ordering for rows and columns.

  • returned object is reference to internal permutation: do NOT free
  • if matrix has no permutation, the return value is NULL

◆ matrix_row_vector()

hpro_d_ vector_t hpro_d_ matrix_row_vector ( const hpro_d_ matrix_t  A,
int *  info 
)

create vectors matching the row indexsets of matrix A

◆ matrix_rows()

size_t hpro_d_ matrix_rows ( const hpro_d_ matrix_t  A,
int *  info 
)

return number of rows of matrix

◆ matrix_scale()

void hpro_d_ matrix_scale ( const double  f,
const hpro_d_ matrix_t  A,
int *  info 
)

scale matrix A by factor f

◆ matrix_solve_diag_left()

void hpro_d_ matrix_solve_diag_left ( const hpro_d_ matrix_t  D,
const hpro_ matop_t  op_D,
hpro_d_ matrix_t  A,
const hpro_ acc_t  acc,
const int  pointwise,
int *  info 
)

solve op(D) B = A with diagonal D and known A

  • B will overwrite A
  • if pointwise is true, evaluation is performed pointwise, otherwise blockwise, e.g. diagonal leaf matrices are considered full and not diagonal

◆ matrix_solve_diag_right()

void hpro_d_ matrix_solve_diag_right ( hpro_d_ matrix_t  A,
const hpro_d_ matrix_t  D,
const hpro_ matop_t  op_D,
const hpro_ acc_t  acc,
const int  pointwise,
int *  info 
)

solve B op(D) = A with diagonal D and known A

  • B will overwrite A
  • if pointwise is true, evaluation is performed pointwise, otherwise blockwise, e.g. diagonal leaf matrices are considered full and not diagonal

◆ matrix_solve_lower_left()

void hpro_d_ matrix_solve_lower_left ( const hpro_d_ matrix_t  L,
const hpro_ matop_t  op_L,
hpro_d_ matrix_t  A,
const hpro_ acc_t  acc,
const int  pointwise,
const int  unit_diag,
int *  info 
)

solve L B = A with lower triangular L and known A

  • B will overwrite A
  • if pointwise is true, evaluation is performed pointwise, otherwise blockwise, e.g. diagonal leaf matrices are considered full and not triangular
  • if unit_diag is true, the pointwise/blockwise diagonal is considered to be the identity

◆ matrix_solve_lower_right()

void hpro_d_ matrix_solve_lower_right ( hpro_d_ matrix_t  A,
const hpro_d_ matrix_t  L,
const hpro_ matop_t  op_L,
const hpro_ acc_t  acc,
const int  pointwise,
const int  unit_diag,
int *  info 
)

solve B op(L) = A with lower triangular L and known A

  • B will overwrite A
  • if pointwise is true, evaluation is performed pointwise, otherwise blockwise, e.g. diagonal leaf matrices are considered full and not triangular
  • if unit_diag is true, the pointwise/blockwise diagonal is considered to be the identity

◆ matrix_solve_upper_right()

void hpro_d_ matrix_solve_upper_right ( hpro_d_ matrix_t  A,
const hpro_d_ matrix_t  U,
const hpro_ acc_t  acc,
const int  pointwise,
const int  unit_diag,
int *  info 
)

solve B U = A with upper triangular U and known A

  • B will overwrite A
  • if pointwise is true, evaluation is performed pointwise, otherwise blockwise, e.g. diagonal leaf matrices are considered full and not triangular
  • if unit_diag is true, the pointwise/blockwise diagonal is considered to be the identity

◆ matrix_sum2()

hpro_d_ linearoperator_t hpro_d_ matrix_sum2 ( const double  alpha1,
const hpro_ matop_t  matop1,
const hpro_d_ linearoperator_t  A1,
const double  alpha2,
const hpro_ matop_t  matop2,
const hpro_d_ linearoperator_t  A2,
int *  info 
)

represents product of matrices \(A_1 + A_2\)

  • all matrix objects are only referenced, not copied

◆ matrix_sum3()

hpro_d_ linearoperator_t hpro_d_ matrix_sum3 ( const double  alpha1,
const hpro_ matop_t  matop1,
const hpro_d_ linearoperator_t  A1,
const double  alpha2,
const hpro_ matop_t  matop2,
const hpro_d_ linearoperator_t  A2,
const double  alpha3,
const hpro_ matop_t  matop3,
const hpro_d_ linearoperator_t  A3,
int *  info 
)

represents product of matrices \(A_1 + A_2 + A3\)

  • all matrix objects are only referenced, not copied

◆ matrix_to_dense()

hpro_d_ matrix_t hpro_d_ matrix_to_dense ( const hpro_d_ matrix_t  A,
int *  info 
)

convert matrix to dense format

◆ matrix_to_rank()

hpro_d_ matrix_t hpro_d_ matrix_to_rank ( const hpro_d_ matrix_t  A,
const hpro_ acc_t  acc,
int *  info 
)

convert matrix to lowrank format using best approximation (SVD) with truncation defined by acc

◆ matrix_transpose()

void hpro_d_ matrix_transpose ( const hpro_d_ matrix_t  A,
int *  info 
)

transpose matrix A

◆ matrix_waz()

hpro_d_ linearoperator_t hpro_d_ matrix_waz ( hpro_d_ matrix_t  A,
const hpro_ acc_t  acc,
int *  info 
)

Factorise matrix A into \(WAZ=I\), e.g., factorization of inverse of A

  • A will be overwritten by \(W^{-1} Z^{-1}\).
  • the return value is a linear operator which can be used to evaluate \(A^{-1}\), e.g. for matrix-vector multiplication (this is not possible with A after factorisation as it holds just the data, albeit, it is still needed)

◆ mm_load_matrix()

hpro_d_ matrix_t hpro_d_ mm_load_matrix ( const char *  filename,
int *  info 
)

read matrix from file filename

◆ perm_linearoperator()

hpro_d_ linearoperator_t hpro_d_ perm_linearoperator ( const hpro_ permutation_t  P,
const hpro_d_ linearoperator_t  A,
const hpro_ permutation_t  R,
int *  info 
)

represents permuted operator \(P \cdot A \cdot R\)

◆ perm_matrix()

hpro_d_ linearoperator_t hpro_d_ perm_matrix ( const hpro_ permutation_t  P,
const hpro_d_ matrix_t  A,
const hpro_ permutation_t  R,
int *  info 
)

represents permuted matrix \(P \cdot A \cdot R\)

  • all objects (P, A and R) are only referenced, not copied

◆ pltmg_load_matrix()

hpro_d_ matrix_t hpro_d_ pltmg_load_matrix ( const char *  filename,
int *  info 
)

read matrix from file filename

◆ samg_load_matrix()

hpro_d_ matrix_t hpro_d_ samg_load_matrix ( const char *  filename,
int *  info 
)

read matrix A from file "\a filename" in SAMG format

  • expecting format definition in basename(filename)>.frm

◆ samg_load_vector()

hpro_d_ vector_t hpro_d_ samg_load_vector ( const char *  filename,
int *  info 
)

read vector from file filename in SAMG format

  • expecting format definition in basename(filename)>.frm

◆ samg_save_matrix()

void hpro_d_ samg_save_matrix ( const hpro_d_ matrix_t  A,
const char *  filename,
int *  info 
)

save matrix A to file "\a filename" in SAMG format

  • also writing format definition to basename(filename)>.frm

◆ samg_save_vector()

void hpro_d_ samg_save_vector ( const hpro_d_ vector_t  x,
const char *  filename,
int *  info 
)

save vector to file filename in SAMG format

◆ solver_solve()

void hpro_d_ solver_solve ( const hpro_ solver_t  solver,
const hpro_d_ linearoperator_t  A,
hpro_d_ vector_t  x,
const hpro_d_ vector_t  b,
const hpro_d_ linearoperator_t  W,
hpro_ solve_info_t *  solve_info,
int *  info 
)

solve A x = b using solver with optional preconditioning, e.g. W A x = W b if W != NULL; information about the solving process is stored in solve_info

◆ vector_assign()

void hpro_d_ vector_assign ( hpro_d_ vector_t  y,
const hpro_d_ vector_t  x,
int *  info 
)

copy content of vector x to vector y x and y have to be of equal type

◆ vector_axpy()

void hpro_d_ vector_axpy ( const double  alpha,
const hpro_d_ vector_t  x,
hpro_d_ vector_t  y,
int *  info 
)

compute y := y + a x

◆ vector_build()

hpro_d_ vector_t hpro_d_ vector_build ( const size_t  size,
int *  info 
)

create scalar vectors of size size initialised with 0

◆ vector_bytesize()

size_t hpro_d_ vector_bytesize ( const hpro_d_ vector_t  v,
int *  info 
)

return size of matrix in bytes

◆ vector_copy()

hpro_d_ vector_t hpro_d_ vector_copy ( const hpro_d_ vector_t  v,
int *  info 
)

copy vectormatrix

◆ vector_dot()

double hpro_d_ vector_dot ( const hpro_d_ vector_t  x,
const hpro_d_ vector_t  y,
int *  info 
)

compute dot product <x,y> = x^H * y

◆ vector_entry_get()

double hpro_d_ vector_entry_get ( const hpro_d_ vector_t  x,
const size_t  i,
int *  info 
)

get single entry at position i in vector x

◆ vector_entry_set()

void hpro_d_ vector_entry_set ( const hpro_d_ vector_t  x,
const size_t  i,
const double  f,
int *  info 
)

set single entry at position i in vector x

◆ vector_export()

void hpro_d_ vector_export ( const hpro_d_ vector_t  v,
double *  arr,
int *  info 
)

copy data from vector v into given C arrays

  • arr has to be of size hlib_vector_size( v )

◆ vector_fill()

void hpro_d_ vector_fill ( hpro_d_ vector_t  x,
const double  f,
int *  info 
)

fill vector x with constant value f

◆ vector_fill_rand()

void hpro_d_ vector_fill_rand ( hpro_d_ vector_t  x,
int *  info 
)

fill vector x with random values

◆ vector_free()

void hpro_d_ vector_free ( hpro_d_ vector_t  v,
int *  info 
)

free resources coupled with vector v

◆ vector_import()

void hpro_d_ vector_import ( hpro_d_ vector_t  v,
const double *  arr,
int *  info 
)

copy data from C array arr into vector v

  • arr has to be of size hlib_vector_size( v )

◆ vector_norm2()

double hpro_d_ vector_norm2 ( const hpro_d_ vector_t  x,
int *  info 
)

compute euklidean norm of x

◆ vector_norm_inf()

double hpro_d_ vector_norm_inf ( const hpro_d_ vector_t  x,
int *  info 
)

compute infinity norm of x

◆ vector_permute()

void hpro_d_ vector_permute ( hpro_d_ vector_t  v,
const hpro_ permutation_t  perm,
int *  info 
)

reorder vector entries with given permutation

◆ vector_scale()

void hpro_d_ vector_scale ( hpro_d_ vector_t  x,
const double  f,
int *  info 
)

scale vector x by f

◆ vector_size()

size_t hpro_d_ vector_size ( const hpro_d_ vector_t  v,
int *  info 
)

return size of vectors