File numerics/src/tools/CSparseMatrix.h

Go to the source code of this file

Structure definition and functions related to sparse matrix storage in Numerics.

Defines

CS_INT
csn
css
NCOMPLEX
NSM_NCOL_CSR(mat)
NSM_NROW_CSR(mat)
NSM_UNKNOWN_ERR(func, orig)

Typedefs

typedef struct cs_di_numeric cs_din
typedef struct cs_di_symbolic cs_dis
typedef struct cs_dl_numeric cs_dln
typedef struct cs_dl_symbolic cs_dls
typedef struct cs_di_sparse CSparseMatrix

Functions

int CSparseMatrix_aaxpby(const double alpha, const CSparseMatrix *A, const double *x, const double beta, double *y)

Matrix vector multiplication : y = alpha*A*x+beta*y.

Return

0 if A x or y is NULL else 1

Parameters
  • [in] alpha: matrix coefficient

  • [in] A: the sparse matrix

  • [in] x: pointer on a dense vector of size A->n

  • [in] beta: vector coefficient

  • [inout] y: pointer on a dense vector of size A->n

CSparseMatrix *CSparseMatrix_alloc_for_copy(const CSparseMatrix *const m)

Allocate a CSparse matrix for future copy (as in NSM_copy)

Return

an newly allocated matrix

Parameters
  • m: the matrix used as model

int CSparseMatrix_check_csc(CSparseMatrix *T)

Check if the given triplet matrix is properly constructed (col and row indices are correct)

Return

0 if the matrix is fine, 1 otherwise

Parameters
  • T: the sparse matrix to check

int CSparseMatrix_check_triplet(CSparseMatrix *T)

Check if the given triplet matrix is properly constructed (col and row indices are correct)

Return

0 if the matrix is fine, 1 otherwise

Parameters
  • T: the sparse matrix to check

int CSparsematrix_chol_factorization(CS_INT order, const CSparseMatrix *A, CSparseMatrix_factors *cs_chol_A)

compute a Cholesky factorization of A and store it in a workspace

Return

1 if the factorization was successful, 1 otherwise

Parameters
  • order: control if ordering is used

  • A: the sparse matrix

  • cs_chol_A: the parameter structure that eventually holds the factors

CS_INT CSparseMatrix_chol_solve(CSparseMatrix_factors *cs_chol_A, double *x, double *b)

reuse a Cholesky factorization (stored in the cs_chols_A) to solve a linear system Ax = b

Return

0 if failed, 1 otherwise

Parameters
  • cs_chol_A: contains the Cholesky factors of A, permutation information

  • x: workspace

  • [inout] b: on input RHS of the linear system; on output the solution

void CSparseMatrix_copy(const CSparseMatrix *const A, CSparseMatrix *B)

Copy a CSparseMatrix inside another CSparseMatrix.

Reallocations are performed if B cannot hold a copy of A

Parameters
  • [in] A: a CSparseMatrix

  • [inout] B: a CSparseMatrix

void CSparseMatrix_free_lu_factors(CSparseMatrix_factors *cs_lu_A)

Free a workspace related to a LU factorization.

Parameters
  • cs_lu_A: the structure to free

int CSparsematrix_lu_factorization(CS_INT order, const CSparseMatrix *A, double tol, CSparseMatrix_factors *cs_lu_A)

compute a LU factorization of A and store it in a workspace

Return

1 if the factorization was successful, 1 otherwise

Parameters
  • order: control if ordering is used

  • A: the sparse matrix

  • tol: the tolerance

  • cs_lu_A: the parameter structure that eventually holds the factors

int CSparseMatrix_max_abs_by_columns(const CSparseMatrix *A, double *max)
int CSparseMatrix_max_by_columns(const CSparseMatrix *A, double *max)
CSparseMatrix *CSparseMatrix_new_from_file(FILE *file)
int CSparseMatrix_print_in_file(const CSparseMatrix *A, int brief, FILE *file)

print a matrix to a text file

Parameters
  • A: matrix to print

  • brief: if positive, print only a portion of the matrix

  • file: file descriptor

int CSparseMatrix_scal(const double alpha, const CSparseMatrix *A)

Multiply a matrix with a double alpha*A > A.

Parameters
  • alpha: the coefficient

  • A: the matrix

CS_INT CSparseMatrix_solve(CSparseMatrix_factors *cs_lu_A, double *x, double *b)

reuse a LU factorization (stored in the cs_lu_A) to solve a linear system Ax = b

Return

0 if failed, 1 otherwise

Parameters
  • cs_lu_A: contains the LU factors of A, permutation information

  • x: workspace

  • [inout] b: on input RHS of the linear system; on output the solution

CSparseMatrix *CSparseMatrix_spfree_on_stack(CSparseMatrix *A)

Free space allocated for a SparseMatrix.

note : cs_spfree also free the cs_struct this fails when the struct is allocated on the stack.

Return

NULL on success

Parameters
  • A: the sparse matrix

CS_INT CSparseMatrix_symmetric_zentry(CSparseMatrix *T, CS_INT i, CS_INT j, double x)

Add an entry to a symmetric triplet matrix only if the absolute value is greater than DBL_EPSILON.

Return

integer value : 1 if the absolute value is less than DBL_EPSILON, otherwise the return value of cs_entry.

Parameters
  • T: the sparse matrix

  • i: row index

  • j: column index

  • x: the value

CS_INT CSparseMatrix_to_dense(const CSparseMatrix *const A, double *B)
CS_INT CSparseMatrix_zentry(CSparseMatrix *T, CS_INT i, CS_INT j, double x)

Add an entry to a triplet matrix only if the absolute value is greater than DBL_EPSILON.

Return

integer value : 1 if the absolute value is less than DBL_EPSILON, otherwise the return value of cs_entry.

Parameters
  • T: the sparse matrix

  • i: row index

  • j: column index

  • x: the value

struct CSparseMatrix_factors
#include <CSparseMatrix.h>

Information used and produced by CSparse for an LU factorization.

Public Members

CS_INT n

size of linear system

csn *N

numerics factorization

css *S

symbolic analysis