Python reference

See the Python guide for installation, examples and backend support.

apxchol.factorize(A, *, backend='auto', num_threads=0, seed=None, partitioner=None, keep_factor=True, factor_storage_dtype=None, **advanced)

Build the reusable approximate-Cholesky factor of A.

Parameters:
  • A (scipy.sparse matrix or array) –

    The square, real symmetric matrix in A @ x = b, typically a graph Laplacian or SDDM matrix. Nonempty rows need positive diagonals. Supply both triangles; sparse inputs are converted to sorted, duplicate-free CSC without modifying the caller’s matrix.

    Use laplacian() to convert adjacency input first. A Laplacian RHS must sum to zero within each connected component.

    Positive off-diagonal pairs are lumped onto the diagonal for preconditioning, subject to lump_max_mass (see Solver.lumped). PCG and residual checks still use the original A.

  • backend ({"auto", "cpu", "gpu"}) – Auto may prefer CPU for small systems when it supports the requested storage; otherwise it checks GPU compatibility. Read the resolved choice from Solver.backend. An explicit GPU request fails when the installed build, hardware or options cannot support it.

  • num_threads (int) – CPU thread ceiling for setup and every solve/apply. Positive values override the calling OpenMP task’s setting; zero inherits it. Small workloads may use fewer threads. Each operation restores the setting.

  • seed (int or None) – RNG seed for the randomized clique sampling (factorization is deterministic per seed at one thread only). None uses the C++ default.

  • partitioner (str or None) – Independent-set selector: “block_greedy”, “priority_greedy”, or “baumann_kyng”. None uses the C++ default. Unknown names raise at factorization time.

  • keep_factor (bool) – Retain an extra host factor copy for save(), chol(), L and D export; enabled by default. False saves this memory while still allowing repeated solves, P, factor_nnz and fill_ratio.

  • factor_storage_dtype (numpy dtype, dtype name, or None) –

    np.float16 or np.float32 selects triangular-solve storage; dtype names also work. None uses the backend default: FP32 on CPU/Metal, FP16 on CUDA.

    FP16 stores scaled off-diagonals with FP32 diagonals and scales. It does not change the exported factor’s FP32 type or CPU arithmetic’s FP64 type. Unsupported FP16 requests fail; FP64 storage is unavailable.

  • **advanced – Passed straight through to the core options: degree_quantile, degree_multiplier, degree_tiebreak, exact_clique_max_degree, residual_peel (“natural” | “min_degree” | “bk_serial”), cpu_schedule (“auto” | “levels”), drop_tolerance (in [0, 1]; 0 disables dropping), regularization (nonnegative diagonal fraction; 0 disables), lump_max_mass (fraction in [0, 1]; 0 rejects positive off-diagonals), and verbose (default False).

Return type:

Solver

Notes

Unknown option keys raise ValueError.

degree_multiplier only takes effect when degree_quantile=0. Omitting the quantile uses the backend’s default. Metal requires a quantile in (0, 1) and includes all vertices tied at its degree cutoff.

apxchol.solve(A, b, rtol=None, maxiter=None, *, x0=None, out=None, **opts)

One-shot convenience: factorize A, then solve A x = b.

x0 and out have the same meaning as in Solver.solve(). The factor is discarded afterwards, so this passes keep_factor=False (no factor export) unless overridden.

Return type:

SolveResult

apxchol.laplacian(A)

Convert a weighted adjacency matrix to its graph Laplacian.

Forms L = D - W using absolute edge weights, with self-loops removed. D contains the weighted vertex degrees; isolated vertices have zero rows. Pass this result to factorize() or solve(), which expect the matrix in L @ x = b rather than adjacency input.

Parameters:

A (scipy.sparse matrix or array) – Square, real, symmetric weighted adjacency matrix, in any sparse format. Explicit zeros and self-loops are dropped; duplicate entries are summed. A is not modified.

Returns:

The float64 Laplacian, ready to pass to factorize().

Return type:

scipy.sparse.csc_matrix

Examples

>>> A = scipy.io.mmread("com-Amazon.mtx").tocsc()   # adjacency
>>> res = apxchol.solve(apxchol.laplacian(A), b)
class apxchol.Solver(A, **opts)

Reusable approximate-Cholesky solver: factor built once, solve many b.

Options (all keyword-only) are documented on factorize().

classmethod from_csc(data, indices, indptr, **opts)

Build from NumPy-compatible CSC arrays, without requiring SciPy.

The operator is square with len(indptr) - 1 rows. Indices are zero-based integers. Unsorted entries and duplicates are accepted; the native importer sorts and sums them without changing the input. Options are the same as factorize().

save(path)

Save the operator and retained factor to a portable binary file.

Requires keep_factor=True. Existing destinations are rejected. The file contains no RHS, solution or device-specific workspace.

classmethod load(path, **opts)

Restore a solver from an archive, without requiring SciPy.

Options use the constructor’s defaults. The stored factor is adopted, not factorized again; backend, storage and preparation settings apply. backend="auto" uses CPU for adoption; "gpu" requests upload. Tolerance, iteration budget and warm starts remain solve arguments.

property shape

Shape of the square operator, as (n, n).

property backend: str

"cpu" or "gpu".

Type:

Backend selected at construction

property sddm: bool

Whether the solver uses SDDM rather than Laplacian handling.

property factor_storage_dtype: dtype

Resolved triangular-solve storage dtype; CPU arithmetic remains float64.

property lumped: int

Number of positive off-diagonal entries moved during preconditioning.

Zero for a Laplacian or SDDM matrix. Each positive pair is moved to the diagonal of the preconditioner and counts as two entries. PCG and the reported residual still use the original, unchanged matrix.

solve(b, rtol=None, maxiter=None, *, x0=None, out=None)

Solve A @ x = b with PCG, reusing this solver’s factor.

b has shape (n,) or (n, nrhs). Matrix inputs use the native batched solver and return one convergence result per column. Omitted rtol and maxiter use the compiled native defaults.

x0 is an initial solution estimate with the RHS shape; an already-converged guess takes zero iterations. out must have that shape and be writable, C- or F-contiguous float64. It is returned as SolveResult.x.

The final residual is always checked against the original matrix and RHS. Inspect SolveResult.converged for success.

Return type:

SolveResult

apply(r)

Apply the preconditioner to a vector or RHS matrix, preserving shape.

Return type:

ndarray

aslinearoperator()

Return a SciPy LinearOperator for another iterative solver.

That solver controls convergence. With a GPU-backed factor, each application exchanges vectors with the host; solve() keeps the GPU iteration on the device.

property P: ndarray

P[original_vertex] = position in the order.

Always available (the permutation survives even keep_factor=False). The returned array is read-only; .copy() it if you need to modify.

Type:

Elimination permutation

chol()

Return a copy of the FP32 lower-triangular factor in permuted order.

The result is a scipy.sparse.csc_matrix with A[p][:, p] ≈ G @ G.T, where p = np.argsort(self.P). Export requires keep_factor=True (the default), regardless of backend or internal factor storage dtype.

For a Laplacian with k connected components, the factor represents the rank-n-k subspace. Each component’s last eliminated column carries a placeholder diagonal.

property L

Unit-lower-triangular factor in permuted space (G with unit diagonal).

A[p][:, p] ≈ L @ diags(D) @ L.T with p = np.argsort(self.P); the same Laplacian caveat as chol() applies. A fresh matrix is returned on every access.

property D: ndarray

Diagonal of the L·D·L^T form, i.e. diag(chol())**2 (fresh copy).

property factor_nnz: int

nnz of the lower-triangular factor G (diagonal included).

property fill_ratio: float

(2 * factor_nnz - n) / nnz(A).

Symmetric convention: the numerator counts the factor G’s pattern reflected to a full symmetric one (the diagonal once), matching nnz(A) for a full symmetric A.

Type:

Fill relative to the input

class apxchol.SolveResult(x, iters, residual, converged)

Solution and convergence information for the original system.

residual is checked against the original matrix and RHS. For a matrix RHS, iters, residual and converged contain one entry per column; for a vector RHS they are scalars.

Parameters:
  • x (ndarray)

  • iters (int | ndarray)

  • residual (float | ndarray)

  • converged (bool | ndarray)

x: ndarray
iters: int | ndarray
residual: float | ndarray
converged: bool | ndarray