Source-linked AI summary

CAST: Canonical Approximate Schur Tree for Approximate Cholesky on Graphs

Meher Chaitanya, Cameron Musco, Aristides Gionis

arXiv:2609.09255v1stat.MLcs.LGmath.NA

TL;DR

Repeated graph workloads need reusable preconditioners, but exact elimination creates dense Schur cliques. CAST replaces each clique with a leverage-score weighted random spanning tree, while CAST-ρ splits terminals to reduce local sampling variability. The guarantees are local, and experiments indicate CAST-1 is generally faster while CAST-2 helps when added edge contributions are inexpensive.

  • Problem

    Repeated Laplacian or SDDM solves require sparse approximate factors because eliminating a pivot creates a dense Schur clique.

  • Method

    CAST samples a weighted random spanning tree of each Schur clique with inverse-marginal reweighting, while CAST-ρ samples after splitting each neighbor into ρ equal-share copies.

  • Results

    CAST-ρ remains unbiased and connected with exact O(ρd)-time sampling and a 1/ρ local second-moment bound; experiments identify CAST-1 as the faster default and CAST-2 as beneficial in favorable degree regimes.

  • Takeaways & Limitations

    Leverage-score marginals provide the unique local minimax one-tree estimator, while increasing ρ trades lower certified sampling variability against construction work and downstream fill.

  • Takeaways & Limitations

    The guarantees certify a single local Schur-clique replacement and do not track error accumulation or provide a worst-case spectral guarantee for the complete factorization.

Abstract

from arXiv · show

Graph-data workloads such as diffusion estimation, ranking, semi-supervised learning, and network optimization often solve many Laplacian or symmetric diagonally dominant M-matrix (SDDM) systems with the same coefficient matrix. Approximate Cholesky preconditioners eliminate vertices one at a time and store the resulting sparse approximate factorization, the \emph{factor}, whose construction cost is amortized across these solves. But eliminating a vertex, the \emph{pivot}, creates a dense Schur-complement clique among its $d$ active neighbors. We introduce CAST (Canonical Approximate Schur Tree), which replaces this clique with a weighted random spanning tree sampled directly from it. Every realization is connected and contains exactly d-1 edges, while reweighting each selected edge by the reciprocal of its tree-inclusion probability makes the update unbiased. The distribution is independent of the ordering of the pivot neighbors, and we prove that its leverage-score marginals minimize the largest normalized reweighted-edge contribution among unbiased inverse-marginal one-tree estimators. We also introduce CAST-$ρ$, which replaces each pivot neighbor with $ρ$ copies, each carrying a $1/ρ$ share of that neighbor's incident weight, samples a weighted random spanning tree on the expanded clique, and contracts the copies back to the original neighborhood. The resulting update remains unbiased and connected, can be sampled exactly in $O(ρd)$ time, and satisfies a $1/ρ$ bound on the second moment of the normalized local Schur error. Increasing $ρ$ therefore reduces certified local sampling variability, but may increase construction cost and downstream fill. Empirically, we observe that CAST-1 is the faster default, whereas CAST-2 is preferable when its additional edge contributions remain inexpensive.

1 Introduction

Many graph workloads repeatedly solve systems with a fixed operator, making reusable preconditioner construction central. CAST replaces dense Schur cliques with canonical randomized tree updates, while CAST-ρ adds splitting to reduce local sampling variability subject to construction and fill trade-offs.

  • Motivation: Repeated graph-data workloads reuse a fixed Laplacian or SDDM operator across many right-hand sides, so factor construction can be amortized across solves.Examples include ranking, diffusion estimation, label propagation, spectral methods, and network optimization.
  • Problem setting: Exact elimination creates a Θ(d2)-edge Schur clique, whereas approximate Cholesky replaces it with a sparse surrogate containing O(d) edges.Connected sparse updates preserve the clique’s Laplacian nullspace structure while avoiding quadratic fill.
  • CAST: CAST samples a weighted random spanning tree of each Schur clique and uses inverse-marginal reweighting, yielding an order-independent connected update.Its edge-inclusion probabilities equal clique leverage scores, which uniquely minimize the largest normalized reweighted-edge contribution among unbiased one-tree estimators.
  • CAST-ρ: CAST-ρ splits each pivot neighbor into ρ equal-share copies, samples a tree on the expanded clique, and contracts the copies while remaining unbiased and connected.The local second-moment bound scales as 1/ρ, but splitting can increase construction work and downstream fill.
  • Scope: The guarantees are local to a single Schur-clique replacement and do not provide a worst-case spectral guarantee for the complete approximate factorization.The full-factor behavior remains outside the stated guarantee.

2 Preliminaries and Problem Setting

The paper studies reusable approximate Cholesky preconditioners for Laplacian and SDDM systems, focusing on sparse unbiased replacements for dense Schur cliques. Its local design objective is to preserve connectivity while controlling normalized sampling contributions.

  • System model: An SDDM matrix can be written as A = LG + Γ, combining a graph Laplacian with a positive-semidefinite diagonal term.Positive diagonal support on each connected component makes the matrix positive definite.
  • Factor reuse: Laplacian and SDDM systems arise in repeated solves with fixed operators and changing right-hand sides, allowing one approximate factor to be reused.The factor is constructed once and applied as a preconditioner across q solves, with construction amortized as q increases.
  • Schur replacement: Eliminating a pivot of degree d creates a Schur clique with Θ(d2) edges, so approximate methods seek sparse random updates satisfying E[bKv] = Kv.This replacement avoids filling the matrix while preserving unbiasedness.
  • Local objective: For a spanning-tree distribution, inverse-marginal reweighting assigns selected edge e weight we/pe, and local risk is governed by the largest normalized contribution τe/pe.The leverage score τe is the edge weight multiplied by its effective resistance in the local clique.
  • Design criterion: CAST chooses weighted spanning-tree marginals pe = τe, and the transfer-current theorem identifies these marginals as the local minimax choice.The design depends on the weighted Schur clique rather than the processing order of neighboring edges.

3 CAST: Canonical Approximate Schur Trees

CAST constructs sparse Schur updates by sampling weighted spanning trees, optionally after splitting each terminal into equal-share copies. Product-form conductances enable exact O(ρd)-time sampling without materializing the expanded clique.

  • Expanded construction: CAST-ρ replaces each terminal with ρ copies carrying conductance ai/ρ, samples a tree on the expanded product-form clique, and contracts copies back to terminals.Edges within one copy block become self-loops and disappear; cross-block edges become terminal contributions.
  • Contracted update: The auxiliary tree has exactly ρd −1 edges, so contraction yields at most ρd −1 nonzero terminal-edge contributions while retaining connected support.For ρ > 1, the contracted update need not itself be a terminal tree.
  • CAST-1: For ρ = 1, CAST samples a weighted random spanning tree of the original Schur clique and assigns each selected edge its inverse-marginal conductance.The resulting edge marginals equal leverage scores, and the reweighted edge matches the series conductance of the eliminated two-hop path.
  • Exact sampling: Weighted Prüfer codes sample the expanded tree exactly in O(ρd) time without materializing the dense clique.Terminal probabilities are sampled according to ai/a, followed by uniform selection among each terminal’s copies; decoding and contraction are linear in the expanded size.
  • Preconditioner: Applying the local update at every pivot produces CAST-CHOL while maintaining a residual graph and recording each pivot star in the approximate factor.Residual degrees depend on fill created by earlier updates, motivating empirical study of ρ’s downstream effects.

4 Local Theory of the CAST Update

The CAST-ρ update is unbiased and connected for every splitting factor, while CAST-1 uniquely minimizes the largest normalized edge contribution and larger ρ reduces certified local sampling variability.

  • CAST samples and contracts one weighted spanning tree on ρd auxiliary copies in O(ρd) time without materializing the dense clique.
  • For general ρ, splitting bounds each normalized auxiliary-edge contribution and the local second moment at order 1/ρ, while the guarantee remains local.Each further doubling of ρ halves the certified bound; the analysis does not track error accumulation, construction cost, or downstream fill.
  • CAST-ρ preserves the exact Schur update in expectation for every integer ρ ≥1.
  • Every CAST-ρ realization has connected terminal support, with nullspace span{1_d} and rank d − 1.
  • For ρ = 1, leverage-score marginals uniquely minimize the maximum normalized reweighted-edge contribution among unbiased inverse-marginal one-tree estimators.
  • CAST-1 attains local risk 1, with clique-edge inclusion probability τ_ij = (a_i + a_j)/a and reweighted conductance a_i a_j/(a_i + a_j).

5 Experimental Evaluation

The evaluation compares CAST and prior approximate-Cholesky variants under a common elimination and PCG framework across benchmark collections. CAST-1 is consistently faster in the base comparison, whereas CAST-2 helps in some low-degree or clique-dominated settings but can lose when high-degree contributions create downstream fill.

  • Experimental setup: The evaluation reuses one factor for q = 250 right-hand sides under a common elimination and PCG framework.
  • Experimental setup: The corpus contains 201 SDDM and Laplacian systems, with 190 systems having paired comparative measurements.
  • Cross-collection comparison: Selecting the faster CAST configuration per system yields a 1.9× arithmetic-mean speedup over the correspondingly faster AC configuration across 190 paired systems.
  • Base granularity: CAST-1 improves on AC across every collection, including a mean 1.18× speedup on 128 Chimera-IPM systems and 1.17× on 17 SuiteSparse matrices.
  • Doubled granularity: CAST-2 is 4.43× faster than AC2 on Spielman-IPM systems, but AC2 is faster on 101 of 128 Chimera systems when high-degree pivots dominate.
  • Choosing ρ: The authors conclude that ρ is a structural choice: splitting helps when extra contributions remain confined, but costs when they propagate as downstream fill.

6 Related Work

Prior approximate-Cholesky methods replace dense elimination fill with sparse randomized updates, often preserving connected support. CAST retains this principle while using an order-independent weighted spanning-tree distribution that supports exact local sampling without materializing dense cliques.

  • CAST relative to AC(k): CAST differs from sequential edge-pairing methods by sampling directly from the weighted random spanning-tree distribution, making its marginals order-independent and minimax-optimal.Its edge inclusion probabilities equal leverage scores.
  • Exact local sampling: Weighted Prüfer codes enable exact sampling on the star-induced clique in O(ρd) time without forming either the original or expanded dense clique.The expanded clique has endpoint-product conductances, allowing hierarchical symbol sampling and linear-time decoding and contraction.
  • Empirical comparison: CAST-1 reduces total reuse workload by 14–15% relative to AC, while CAST-2 reaches 10−8 residual tolerance in one to two PCG iterations and averages 4.43× faster than AC2 on Spielman-IPM systems.The CAST-2 advantage is reported when additional fill remains controlled.

B.3 Proof of Theorem 4.2

The proof establishes that the contracted CAST-ρ update is unbiased and connected. It does so by matching expected cross-terminal edge coefficients to the exact Schur clique and by contracting a connected auxiliary spanning tree.

  • Unbiasedness: Unbiasedness follows by summing each auxiliary edge contribution weighted by its tree-inclusion probability.Within-block edges become self-loops, while cross-block contributions remain in the terminal graph.
  • Unbiasedness: For distinct terminals, aggregated expected cross-block conductances equal the exact Schur-clique edge coefficient.The identity uses the copy conductances and the relation that each terminal’s copy conductances sum to its incident pivot weight.
  • Unbiasedness: Therefore, the expected contracted Laplacian equals the exact Schur complement: E[bKv] = Kv.Both matrices are Laplacians on the terminal neighborhood and are determined by their edge coefficients.
  • Connectivity: Contracting the connected auxiliary spanning tree preserves connectivity among terminals, while within-block edges become discardable self-loops.Every cross-block edge has positive conductance, so the quotient support retains the required terminal connections.
  • Connectivity: For ρ = 1 the update is itself a spanning tree, whereas for ρ > 1 its contracted support remains connected but need not be a tree.Copy contraction can create parallel terminal edges that are aggregated after contraction.

B.5 Proof of Corollary 4.4

The corollary proves a lower bound on the worst normalized reweighted-edge contribution and identifies equality uniquely with leverage-score marginals. The weighted random spanning-tree distribution attains this optimum.

  • Structural properties: Connected support guarantees rank d − 1 and nullspace span{1_d} for the sampled Laplacian.The result follows from the connected-support characterization of graph-Laplacian nullspaces.
  • Structural properties: Every spanning-tree realization contains exactly d − 1 edges, fixing the expected total number of selected edges.The proof uses the identity summing edge indicators over each realization.
  • Minimax guarantee: The maximum normalized contribution satisfies R_K(D) ≥ 1 for every unbiased inverse-marginal one-tree estimator.The bound follows because the relevant quantity is a positive weighted average of the ratios τ_e/p_e.
  • Minimax guarantee: Equality holds only when p_e = τ_e for every edge, and leverage-score marginals achieve equality under weighted random spanning-tree sampling.The transfer-current theorem supplies p_e = τ_e for that distribution.

B.7 Proof of Lemma 4.6

The lemma analyzes dependence among sampled tree edges through the transfer-current matrix and bounds the resulting covariance structure. This supports controlling the second moment of the local Schur error.

  • Local error structure: Within-block auxiliary edges contribute zero after contraction, while cross-block edges produce rank-one positive-semidefinite local error matrices.For a cross-block edge, the error matrix is formed from the terminal difference vector and its conductance.
  • Local error structure: The local edge-error contribution is bounded by 1/ρ times the pivot projection, and its squared contribution inherits the corresponding scalar factor.The bound uses copy conductances ai/ρ and aj/ρ together with the pivot effective resistance.
  • Edge dependence: Tree-edge indicators form a determinantal point process whose kernel is the symmetric positive-semidefinite transfer-current matrix.Its diagonal entries equal the edge-inclusion probabilities p_e.
  • Covariance bound: The covariance matrix of edge indicators satisfies C ⪯ diag(p).The proof uses determinantal pair correlations and the Schur product theorem.
  • Second-moment control: The second-moment argument combines unbiasedness with the covariance bound and positive-semidefinite Gram structure.The proof compares the covariance matrix with diag(p) through a Frobenius inner product.

C.1 SuiteSparse benchmark

On 17 nontrivial SuiteSparse matrices, CAST-1 generally outperformed AC, while CAST-2 showed more mixed gains against AC2 and benefited from factor reuse. Aggregate comparisons exclude trivial or non-finishing cases and reflect both construction and repeated-solve costs.

  • Excluded cases: bcsstm25 was excluded because the 10−8 tolerance was below attainable double-precision accuracy for some right-hand sides, while all methods converged in one iteration at 10−5.This identifies the failure as a matrix–tolerance issue rather than a solver-specific failure.
  • Base-granularity comparison: CAST-1 was faster than AC on 14 of 17 matrices, with an average speedup of 1.17×.The largest gains were 1.55× on nos7, 1.54× on shallow_water2, and 1.48× on nos6.
  • Base-granularity comparison: CAST-1’s aggregate improvement arose primarily from lower solve costs, although it was slower on apache1, jnlbrng1, and bcsstm24.On apache1, CAST-1 required 54 iterations per solve versus 28 for AC.
  • Doubled-granularity comparison: CAST-2 was faster than AC2 on 11 of 17 matrices, with an average speedup of 1.069×.Its strongest reported gains included 1.42× on shallow_water2, 1.27–1.31× on the ecology matrices, and 1.22× on apache1.
  • Cost decomposition: CAST-2 construction cost 2.09 µs/nnz versus 0.63 µs/nnz for AC2 on one matrix, making CAST-1 preferable there.On ecology1, CAST-2 nevertheless achieved a 1.31× total-time speedup through fewer iterations and lower per-iteration application cost.
  • Factor reuse: Factor reuse changes the cost comparison: AC2 was 51% more expensive than AC at q = 1, while aggregate results at q = 250 modestly favored CAST-2.The additional construction cost of doubled-granularity variants is amortized as the number of right-hand sides increases.

C.2 IPM Sequences on Chimera Graphs

On Chimera-IPM sequences, CAST-1 consistently reduces total workload relative to AC, while CAST-2’s benefit depends on pivot-degree structure and reuse. The experiments measure repeated Newton-step solves, including construction and verified PCG solves.

  • Benchmark and protocol: 128 Chimera-IPM systems use five graph instances, with 23–28 Newton-step systems per instance and five duality-gap targets.The systems are weighted graph Laplacians with n = 100,000 vertices, and performance is aggregated over Newton steps and targets within each graph instance.
  • CAST-1 results: At q = 250, CAST-1 is faster than AC on all 128 systems, with per-instance speedups from 1.13× to 1.24×.Normalized solve cost is lower for CAST-1 on all five instances; construction costs are equal or nearly equal on i1 and i2 and higher on i3–i5.
  • CAST-2 results: CAST-2 is fastest on i3, achieving 1.16× speedup over AC2 and winning on all 27 systems.i3 has the smallest mean and maximum pivot degrees among the highlighted instances.
  • CAST-2 trade-offs: CAST-2 construction costs rise with high-degree pivots, reaching approximately 10× CAST-1 on i4 and i5 and 13.1× on i1.A degree-d pivot uses 2d − 1 auxiliary tree edges before aggregation, compared with d − 1 for CAST-1, increasing downstream terminal-edge contributions.
  • Choosing the splitting factor: At q = 250, CAST-2 outperforms CAST-1 only on i3; on instances with dmax ≥92, it is 1.3–4.6× slower.The results support an empirical degree-aware choice: ρ = 2 when reuse amortizes construction, and ρ = 1 when high-degree pivots make fill dominant.

C.3 IPM Sequences on Spielman Graphs

On Spielman IPM sequences, uniformly small elimination degrees make CAST-2 inexpensive and highly effective, while extreme-weight near-tree structure makes CAST-1 less consistent. The preferred splitting factor depends on the trade-off between convergence and downstream construction cost.

  • 40 systems span 3.4 × 10^5 to 2.1 × 10^7 vertices, with up to 6.4 × 10^7 nonzeros.
  • Every pivot has degree at most five, and the mean pivot degree is 2.00 at every scale.Minimum-degree elimination avoids large pivots despite high-degree input hubs.
  • CAST-2 converges in one to two PCG iterations on every recorded solve, with worst verified relative residual 5.5 × 10^-9.At q = 250, mean speedups over AC2 range from 4.04–4.72×.
  • CAST-1 and CAST-2 produce factors containing approximately 0.66 times as many nonzeros as the input, with similar construction costs.CAST-2's improved convergence does not increase measured factor size on this family.
  • CAST-1 is less consistent on extreme-weight near-tree systems, while CAST-2 remains inexpensive when degree d ≤ 5.At k = 100, one unlucky CAST-1 draw causes AC to win the aggregate despite CAST-1 being faster on nine of ten systems.
  • A uniformly light degree distribution favors ρ = 2 when factors are reused, whereas high-degree pivots favor ρ = 1 because fill and construction costs dominate.This is presented as an empirical selection heuristic, not an evaluated adaptive algorithm.

C.4 Sachdeva-star stress tests

Sachdeva-star instances stress approximate Cholesky with unweighted, clique-dominated structure whose difficulty grows with instance size. CAST-2 is substantially more stable and faster than the base-granularity methods, although it does not reduce these systems to one or two iterations.

  • Sachdeva-star graphs are unweighted and structurally challenging because clique-block eliminations create dense local neighborhoods.Their difficulty contrasts with the extreme-weight, near-tree Spielman family.
  • AC's mean iteration count rises from 50.8 at k = 50 to 1,120.8 at k = 600, while CAST-1 rises from 37.9 to 352.9.
  • Across the same range, AC2 rises from 27.6 to 65.0 iterations per solve, while CAST-2 rises from 25.2 to 44.7.The doubled-granularity variants are substantially more stable.
  • CAST-2 is faster than AC2 on every instance, with speedup increasing from 1.23× at k = 50 to 2.46× at k = 600.At the largest instance, CAST-2's total-time advantage over AC is approximately 41.4×.
  • CAST-2 is already the fastest method at q = 1, so its advantage does not depend on amortizing construction cost.
  • CAST-2's iteration growth remains mild because its extra contributions fall within already dense clique blocks, unlike sparse neighborhoods with heavy degree tails.The observed improvement is consistent with reduced local sampling variability, but the local theory does not imply a global iteration bound.
Loading 2609.09255v1…