py-mkl-pardiso

Python pybind11 wrapper for the Intel oneMKL PARDISO sparse direct solver.

pymklpardiso exposes the real-valued subset of Intel's PARDISO sparse direct solver to Python via pybind11. It works with SciPy sparse matrices in CSR format and NumPy arrays.

Supported platforms: Linux (x86_64), Windows (AMD64).

Installation

pip install py-mkl-pardiso

Or install from source (requires MKL):

git clone https://github.com/bodono/py-mkl-pardiso.git
cd py-mkl-pardiso
pip install -e ".[test]"

Building from source requires:

• Intel oneMKL (set MKLROOT if not auto-detected)
• A C++17 compiler
• Python >= 3.10
• pybind11 >= 2.12

Quick start

import numpy as np
import scipy.sparse as sp
from pymklpardiso import PardisoSolver, MTYPE_REAL_SYM_POSDEF

# Build a symmetric positive-definite matrix (upper triangle, CSR)
A_full = np.array([
    [4.0, 1.0],
    [1.0, 3.0],
])
A_upper = sp.csr_matrix(np.triu(A_full))
A_upper.sort_indices()

# Create solver — analyzes, factors, and is ready to solve
solver = PardisoSolver(A_upper, MTYPE_REAL_SYM_POSDEF)

b = np.array([1.0, 2.0])
x = solver.solve(b)
print(x)  # [0.09090909 0.63636364]

Refactoring workflow

When the sparsity pattern stays the same but values change (e.g., in an iterative algorithm), use refactor() to skip symbolic analysis:

solver = PardisoSolver(A_upper, MTYPE_REAL_SYM_POSDEF)
for new_values in value_generator:
    solver.refactor(new_values)
    x = solver.solve(b)

API reference

PardisoSolver(A, mtype, iparms=None, msglvl=0)

Create a PARDISO solver instance. The constructor extracts the CSR sparsity pattern from A, applies any iparms overrides, and runs symbolic analysis

• numeric factorization so the solver is ready to call solve().

Parameter	Type	Default	Description
A	sparse CSR	(required)	Square sparse matrix (any object with indptr, indices, data, shape). For symmetric types, pass only the upper triangle in CSR format.
mtype	int	(required)	Matrix type (see constants below).
iparms	dict	None	Optional {index: value} iparm overrides.
msglvl	int	0	Message level (0 = silent, 1 = print statistics).

Matrix type constants

Constant	Value	Description
MTYPE_REAL_STRUCT_SYM	1	Real structurally symmetric
MTYPE_REAL_SYM_POSDEF	2	Real symmetric positive definite
MTYPE_REAL_SYM_INDEF	-2	Real symmetric indefinite
MTYPE_REAL_NONSYM	11	Real nonsymmetric

Core methods

solver.solve(b) Solve Ax = b. Accepts 1D (n,) or 2D (n, nrhs) arrays. Returns the solution as a new NumPy array (Fortran-contiguous for 2D).

solver.solve_into(b, x) Solve Ax = b writing into pre-allocated x. For 2D arrays, both b and x must be Fortran-contiguous.

solver.refactor(values) Re-factorize with new nonzero values (phase 22 only). Does not re-run symbolic analysis. Raises if symbolic analysis is invalid; use factor() to re-analyze from scratch. values must match the stored sparsity pattern exactly.

solver.factor(values) Re-analyze and re-factorize with new values (phases 11 + 22). This always runs fresh symbolic analysis. Use this for error recovery or when iparm changes require fresh symbolic analysis. values must match the stored sparsity pattern exactly.

Analysis invalidation: symbolic analysis (phase 11) becomes invalid when an iparm value changes via set_iparm() or set_iparm_all(), or after release() / reset(). Changing numeric values alone does not invalidate analysis. After a successful factor(), subsequent refactor() calls continue to work, including with value-dependent analysis settings such as iparm[10] = 1. When analysis is invalid, refactor() raises; call factor() to re-analyze and recover.

Other methods

Method	Description
solver.release()	Free PARDISO internal memory.
solver.n	Matrix dimension (property).
solver.nnz	Number of nonzeros (property).
solver.mtype	Matrix type (property).
solver.set_perm(perm)	Set fill-reducing permutation.
solver.clear_perm()	Clear permutation.
solver.has_perm()	Whether a permutation is set.
solver.set_iparm(idx, value)	Set a single iparm entry.
solver.get_iparm()	Get all 64 iparm values.
solver.get_iparm_value(idx)	Get a single iparm value.
solver.set_iparm_all(iparm)	Set all 64 iparm values.
solver.set_msglvl(msglvl)	Change message level.
solver.run_phase(phase)	Run an arbitrary PARDISO phase.
solver.run_phase_into(phase, b, x)	Run a phase with RHS/output arrays.

iparm notes

• iparm[0] is locked to 1 (user-supplied parameters).
• iparm[34] is locked to 1 (zero-based indexing).
• See the MKL PARDISO iparm documentation for all parameters.

License

MIT