umf_solve

sksparse.umfpack.umf_solve(A, b, *, trans='N', rhs_batch_size=100, control=None, **kwargs)[source]

Solve a linear system using UMFPACK.

This is a convenience function that creates a UMFFactor object, computes the numeric factorization, and solves the linear system.

Parameters:
  • A ((N, N) numpy.ndarray or sparse array) – The input matrix.

  • b ((N,) or (N, K) numpy.ndarray or sparse array) – The right-hand side vector or matrix.

  • trans (str, optional) – The type of system to solve. Possible values are:

    • N: solve \(A x = b\) (default)

    • T: solve \(A^{\top} x = b\)

    • H: solve \(A^{H} x = b\)

    Note

    If \(A\) is real, then T and H are equivalent.

  • rhs_batch_size (int, optional) – If b is a 2D sparse array, this parameter controls the number of columns to be solved simultaneously. A larger number will increase memory consumption by converting more columns at a time to dense arrays, but may improve runtime.

  • control (UMFControl, optional) – The control parameters to use for the factorization. If not provided, default parameters are used.

  • kwargs (keyword arguments, optional) – Additional keyword arguments to pass to UMFControl if control is not provided.

Returns:

x ((N,) or (N, K) numpy.ndarray or sparse array) – The solution vector or matrix of the same type and shape as the input right-hand side b.

Raises:

UMFPACKSingularMatrixWarning – If the matrix is detected to be singular to working precision. In that case, the solution will have infinite or NaN values, but other entries may still be valid.

Added in version 0.5.0.

Notes

The underlying UMFPACK solver can only handle 1D dense array inputs. If the RHS b is a 2D array, this method will solve each column independently. If b is dense, there is a slight performance gain (~5% in time) by passing it as a Fortran-contiguous array (e.g. by using numpy.asfortranarray()), since the columns are then stored contiguously in memory. Otherwise, each column will be copied to a temporary Fortran-contiguous buffer before solving.

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.umfpack import umf_solve
>>> 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)>
>>> # Solve a linear system
>>> expect_x = np.arange(N, dtype=np.float64)
>>> b = A @ expect_x
>>> x = umf_solve(A, b)
>>> np.allclose(x, expect_x)
True