CSparseMatrix (functions)


siconos.numerics.CSparseMatrix_aaxpy(double alpha, CSparseMatrix *A, array_like (np.float64, 1D)x, double beta, array_like (np.float64, 1D)y) → int[source]

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

Parameters:
  • alpha – matrix coefficient
  • A – the sparse matrix
  • x – pointer on a dense vector of size A->n
  • beta – vector coefficient
  • y – pointer on a dense vector of size A->n
Returns:

0 if A x or y is NULL else 1


siconos.numerics.CSparseMatrix_alloc_for_copy(CSparseMatrix * m) → CSparseMatrix *[source]

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

Parameters:m – the matrix used as model
Returns:an newly allocated matrix

siconos.numerics.CSparseMatrix_check_csc(CSparseMatrix *T) → int[source]

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

Parameters:T – the sparse matrix to check
Returns:0 if the matrix is fine, 1 otherwise

siconos.numerics.CSparseMatrix_check_triplet(CSparseMatrix *T) → int[source]

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

Parameters:T – the sparse matrix to check
Returns:0 if the matrix is fine, 1 otherwise

siconos.numerics.CSparseMatrix_copy(CSparseMatrix * A, CSparseMatrix *B) → None[source]

Copy a CSparseMatrix inside another CSparseMatrix.

Reallocations are performed if B cannot hold a copy of A

Parameters:
  • A – a CSparseMatrix
  • B – a CSparseMatrix

siconos.numerics.CSparseMatrix_free_lu_factors(CSparseMatrix_lu_factors *cs_lu_A) → None[source]

Free a workspace related to a LU factorization.

Parameters:cs_lu_A – the structure to free

siconos.numerics.CSparseMatrix_new_from_file(FILE *file) → CSparseMatrix *[source]

siconos.numerics.CSparseMatrix_print_in_file(CSparseMatrix *A, int brief, FILE *file) → int[source]

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

siconos.numerics.CSparseMatrix_solve(CSparseMatrix_lu_factors *cs_lu_A, array_like (np.float64, 1D)x, array_like (np.float64, 1D)b) → CS_INT[source]

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

Parameters:
  • cs_lu_A – contains the LU factors of A, permutation information
  • x – workspace
  • b – on input RHS of the linear system; on output the solution
Returns:

0 if failed, 1 otherwise


siconos.numerics.CSparseMatrix_spfree_on_stack(CSparseMatrix *A) → CSparseMatrix *[source]

Free space allocated for a SparseMatrix.

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

Parameters:A – the sparse matrix
Returns:NULL on success

siconos.numerics.CSparseMatrix_to_dense(CSparseMatrix * A, array_like (np.float64, 1D)B) → CS_INT[source]

siconos.numerics.CSparseMatrix_zentry(CSparseMatrix *T, CS_INT i, CS_INT j, double x) → CS_INT[source]

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

Parameters:
  • T – the sparse matrix
  • i – row index
  • j – column index
  • x – the value
Returns:

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


siconos.numerics.CSparsematrix_lu_factorization(CS_INT order, CSparseMatrix *A, double tol, CSparseMatrix_lu_factors *cs_lu_A) → int[source]

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

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
Returns:

1 if the factorization was successful, 1 otherwise