klu_factor

sksparse.klu.klu_factor(A, *, control=None, **kwargs)[source]

Compute the LU factorization of a sparse matrix using KLU.

This is a convenience function that creates a KLUFactor object, computes the numeric factorization, and returns the resulting object.

Parameters:

  • A ((M, N) numpy.ndarray or sparse array) – The input matrix to factorize.

  • control (KLUControl, optional) – An optional KLUControl object to set the factorization parameters. If not provided, default parameters are used.

  • **kwargs – Additional keyword arguments passed to the KLUControl constructor.

Returns:

KLUFactor – The LU factorization of the input matrix.

Raises:

KLUSingularMatrixWarning – If the matrix is exactly singular.

See also

KLUFactor, klu_solve

Added in version 0.5.0.

Examples

See: Davis, Timothy A. (2006). Direct Methods for Sparse Linear Systems, p 74 (Figure 5.1)

>>> import numpy as np
>>> from scipy import sparse
>>> from sksparse.klu import klu_factor
>>> N = 8
>>> rows = np.array(
...    [0, 1, 2, 3, 4, 5, 6, 3, 6, 1, 6, 0, 2, 5, 7, 4, 7, 0, 1, 3, 7, 5, 6],
...    dtype=np.int32,
...)
>>> cols = np.array(
...    [0, 1, 2, 3, 4, 5, 6, 0, 1, 2, 2, 3, 3, 4, 4, 5, 5, 6, 6, 6, 6, 7, 7],
...    dtype=np.int32,
...)
>>> vals = np.ones(len(rows), dtype=np.float64)
>>> vals[:7] = np.arange(1, 8, dtype=np.float64)  # make diagonal entries non-unit
>>> A = sparse.csc_array((vals, (rows, cols)), shape=(N, N))
>>> A
<Compressed Sparse Column sparse array of dtype 'float64'
        with 23 stored elements and shape (8, 8)>
>>> # Compute the LU factorization
>>> f = klu_factor(A)
>>> f
<KLUFactor numeric factor of dtype 'float64' with 'int32' indices:
    L: (8, 8) with 16 stored elements
    U: (8, 8) with 17 stored elements>
>>> L, U, p, q, r, F, _ = f  # unpack the factorization
>>> LUF = (L @ U + F).toarray()
>>> RPAQ = (r[:, np.newaxis] * A[p[:, np.newaxis], q]).toarray()
>>> np.allclose(LUF, RPAQ)
True
>>> # Solve a linear system
>>> expect_x = np.arange(N, dtype=np.float64)
>>> b = A @ expect_x
>>> x = f.solve(b)
>>> np.allclose(x, expect_x)
True