Clark Kent LU Decomposition (`sksparse.klu`)¶

The sksparse.klu module provides an interface to the SuiteSparse KLU package, which computes basic linear algebra operations for sparse, square matrices.

The main function of this module is to compute the LU factorization of a sparse matrix \(A\) with a fill-reducing permutation \(P\) and \(Q\), and row-scaling \(R\), and offset \(F\) such that:

\[LU + F = RPAQ.\]

These factors can then be used to solve linear systems of the form \(Ax = b\).

The klu module exposes most of the capabilities of the KLU package including:

Computation of the LU factors with fill-reducing permutation for both real and complex sparse matrices \(A\), in any format supported by scipy.sparse.
An interface for using this decomposition to solve problems of the form \(Ax = b\).
The ability to perform the fill-reduction analysis once, and then re-use it to efficiently decompose many matrices with the same pattern of non-zero entries.

Quickstart¶

If \(A\) is a sparse matrix then the following code computes the LU decomposition of \(A\):

from sksparse.klu import klu_factor
L, U, p, q, r, F, _ = klu_factor(A)
# L @ U + F == r[:, np.newaxis] * A[p[:, np.newaxis], q]

See the example below for a demonstration of the effect of the fill-reducing permutation.

Once the factorization has been computed, it can be used to solve a square linear system:

from sksparse.klu import klu_factor
A = ...  # a sparse, square matrix
b = ...  # right-hand side
x = klu_solve(A, b)

Example¶

To see how to use the LU factorization, we can load a sparse matrix from the SuiteSparse Matrix Collection and compute its ordering.

# Part of the scikit-sparse project.
# Copyright (C) 2025 Bernard Roesler. All rights reserved.
# See pyproject.toml for full author list and LICENSE.txt for license details.
# SPDX-License-Identifier: BSD-2-Clause
#
# =============================================================================
#     File: klu_example.py
#  Created: 2025-11-17 10:38
# =============================================================================

"""An example of using KLU to compute the klu_factor factorization of a
sparse matrix.
"""

from pathlib import Path

import matplotlib.pyplot as plt
import numpy as np
from numpy.testing import assert_allclose
from scipy.io import mmread

from sksparse.klu import klu_factor
from sksparse.umfpack import umf_factor

# Load the west0479 matrix (downloaded from the SuiteSparse Matrix Collection:
# <https://sparse.tamu.edu/HB/west0479>)
filepath = Path("data") / "west0479.mtx"
A = mmread(filepath, spmatrix=False)  # read the matrix
A = A.tocsc()

# Compute the LU factorization
# NOTE that KLU does not support a "natural" ordering, so use UMFPACK here
L, U, p, q, rscale = umf_factor(A, ordering_method="none")
Lp, Up, pp, qp, rscalep, Fp, _ = klu_factor(A, ordering="AMD")

# Make sure the factorization is correct
PRAQ = (rscale[:, np.newaxis] * A).tocsc()[p[:, np.newaxis], q]
PRAQp = rscalep[:, np.newaxis] * A[pp[:, np.newaxis], qp]
assert_allclose((L @ U).toarray(), PRAQ.toarray(), atol=1e-9)
assert_allclose((Lp @ Up + Fp).toarray(), PRAQp.toarray(), atol=1e-9)

# Plot the original and permuted matrices
plt.rcParams.update({"font.size": 10})

fig, axs = plt.subplots(num=1, nrows=2, ncols=2, clear=True)
fig.set_size_inches((6, 6), forward=True)
fig.set_constrained_layout(True)
fig.suptitle("KLU Example: LU Factors of west0479")

ax = axs[0, 0]
ax.spy(A, markersize=1)
ax.set_title(r"Original Matrix $A$")

ax = axs[0, 1]
ax.spy(PRAQ, markersize=1)
ax.set_title(r"Permuted Matrix $PAQ$")

ax = axs[1, 0]
ax.spy(L + U, markersize=1)
ax.set_title(r"Original LU Factors $L + U$")
ax.set_xlabel(f"{L.nnz + U.nnz:,} non-zeros")

ax = axs[1, 1]
ax.spy(Lp + Up, markersize=1)
ax.set_title(r"Permuted LU Factors $L_p + U_p$")
ax.set_xlabel(f"{Lp.nnz + Up.nnz:,} non-zeros")

for ax in axs.flat:
    ax.set_aspect("equal")
    ax.set_xticks([])
    ax.set_yticks([])

plt.show()
fig.savefig("klu_example.svg")

This figure shows the effect of AMD ordering that reduces the fill-in of the LU factorization of a sparse matrix.

KLU Example with AMD Ordering — The number of non-zeros in the LU factorization of the original matrix (left) and the permuted matrix (right) using AMD ordering.¶

Function Interface¶

For users who want to directly compute the factorization without needing to manipulate the KLUFactor object, the klu module provides the klu_factor() function that performs both the symbolic analysis and the numeric factorization in one step. It returns the factorization object as a single output, or the individual factors and permutations can be unpacked as multiple outputs.

To directly solve a linear system without explicitly computing the factorization, the klu_solve() function can be used. This function performs the factorization internally and returns the solution to the linear system.

Object Interface¶

For more advanced usage, users can instantiate the KLUFactor class. This class can be instantiated directly using its constructor, or more conveniently using the klu_factor() function.

When instantiated directly, the constructor performs a symbolic analysis of the matrix, but does not compute the numeric factorization. The numeric factorization is then performed by calling the KLUFactor.factorize() method.

The klu_factor() function performs both the symbolic analysis and the numeric factorization in one step, and returns an instance of the KLUFactor class. The resulting KLUFactor object can then be used to solve linear systems using its KLUFactor.solve() method. The KLUFactor.factorize() method can be called again to factor a new matrix with the same sparsity pattern.

Customization and Statistics¶

The klu_factor() accepts several optional parameters that control the behavior of the factorization, such as the fill-reducing ordering method to use. The KLUControl object can be instantiated directly to customize the factorization options before passing it to the KLUFactor constructor or the klu_factor() function. Any options not explicitly set by the user will be given the KLU default values.

After a factorization is computed, statistics about the factorization can be accessed via the KLUFactor.info attribute, which is an instance of the KLUInfo class. This object contains various performance metrics and statistics about the factorization process.

Exceptions and Warnings¶

Warnings issued by KLU are converted into Python warnings of type KLUWarning. The module will also issue a SparseEfficiencyWarning if the input matrix is not a csc_array (note that the csc_matrix class is not supported, as it will be deprecated and is not recommended for use in new code).

Errors detected by KLU or by our wrapper code are converted into exceptions of type KLUError or an appropriate subclass. See the Warnings and Exceptions for details.

Clark Kent LU Decomposition (sksparse.klu)¶