- The paper introduces a scalable framework integrating Krylov subspace methods with distributed strategies (Block Jacobi and BBD) to overcome sparse linear solve bottlenecks in interior point methods.
- It employs adaptive preconditioning—using sparse filtering, diagonal correction, and preconditioner reuse—to enhance convergence and reduce factorization costs.
- Numerical results demonstrate significant speedups over PETSc and PARDISO, validating its robust performance on both well-conditioned and ill-conditioned large-scale problems.
SDSL-Solver: Scalable Distributed Sparse Linear Solvers for Large-Scale Interior Point Methods
Motivation and Bottlenecks in Large-Scale Interior Point Methods
Interior Point Methods (IPMs) are central in solving large-scale convex optimization problems, including linear and quadratic programming. At each iteration, the computational bottleneck lies in the sparse linear system solve, which constitutes more than 70% of the iteration time, especially as problem scale reaches millions of variables. Direct solvers (LU/Cholesky) incur prohibitive fill-in and poor parallel scalability, rendering them infeasible for industrial-scale applications. Krylov subspace iterative methods provide scalable alternatives, but their efficacy is dependent on high-quality preconditioners, which are challenging to construct for ill-conditioned systems. The need for distributed, scalable solvers optimized for such demands sets the foundation for the SDSL-Solver framework.

Figure 1: Workflow of an interior point method, with sparse linear system solve dominating computational cost.
SDSL-Solver Framework Architecture
SDSL-Solver integrates Krylov subspace iterative solvers with two complementary distributed parallel strategies—Block Jacobi and Bordered Block Diagonal (BBD)—along with numerics-based sparse filtering, diagonal correction for preconditioner improvement, and intelligent preconditioner reuse. The framework dynamically selects the appropriate distributed method based on runtime diagnostics: Block Jacobi is suited for well-conditioned, diagonally dominant systems; BBD, leveraging Schur complement techniques, is employed for ill-conditioned or globally coupled matrices. Adaptive preconditioning is constructed through sparse filtering (dropping negligible off-diagonal entries), MC64-based diagonal scaling/maximization, and diagonal correction. Preconditioner reuse amortizes the expensive construction across iterations. All methods are orchestrated within a hybrid MPI+OpenMP computing model.

Figure 2: SDSL-Solver architecture: dynamic selection between Block Jacobi and BBD distributed strategies, with adaptive preconditioning and reuse mechanisms.
Methodological Innovations
Distributed Parallelization Methods
Block Jacobi partitions the matrix row-wise, with independent diagonal blocks preconditioned via ILU/IC factorization—enabling communication-free sparse triangular solves. It suits matrices where diagonal dominance is preserved and minimal inter-block coupling exists.
BBD leverages nested dissection for block structure, preserves inter-block coupling via border/interface systems, and coordinates Schur complement solves across processes. This is essential for general or ill-conditioned matrices where local preconditioners are insufficient for convergence.
Adaptive Preconditioning Techniques
Sparse Filtering Algorithm constructs high-quality sparser preconditioners by dropping off-diagonal elements below a threshold Ï„ relative to diagonal entries, significantly reducing factorization time while retaining spectral properties necessary for rapid Krylov convergence.
Diagonal Correction introduces a positive perturbation δ to each diagonal entry after MC64 scaling, improving the preconditioner's condition number, especially in late IPM stages, thereby ensuring convergence under extreme spectral disparity.
Preconditioner Reuse leverages structural stability between IPMs iterations, reusing symbolic/numerical (and, where feasible, both) factorizations to further reduce computational overhead, and—in some cases—improving solver robustness via implicit regularization.
Numerical Results and Scalability
Comprehensive benchmarks with matrix dimensions ranging from 104 up to 5×106 were executed on multi-node X86 and Kunpeng clusters. SDSL-Solver achieved the following strong numerical results:
- Block Jacobi: On well-conditioned benchmarks, SDSL-Solver surpassed PETSc by a geometric mean speedup of 6.23×, and a geometric mean speedup of 97.54× over single-node PARDISO.
- BBD: On ill-conditioned or globally coupled problems, speedups of 7.77× (PETSc) and 5.85× (PARDISO) were attained.
- In numerous benchmarks (e.g., PageRank-type matrices), direct solvers failed due to memory or integer indexing limits, while SDSL-Solver demonstrated robust scalability and convergence.
- The filtering approach enabled up to 90% reduction in factorization cost, maintaining IPM convergence for aggressive thresholds (τ=10−3).
- Diagonal correction enabled SDSL-Solver to reach optimality on 15 benchmarks for which direct solvers (including MKL PARDISO) failed.
- Preconditioner reuse delivered speedups up to δ0 across IPM iterations and reduced the aggregate IPM step count by implicit regularization effects.
- End-to-end wall-clock time for large-scale optimization was reduced by up to δ1 over direct solvers.

Figure 3: Strong scaling of SDSL-Solver on X86 cluster, demonstrating linear speedup trend up to four nodes.
Practical and Theoretical Implications
The results demonstrate the necessity of architecture-aware distributed solvers synergized with adaptive preconditioning for unlocking the scalability of IPMs in industrial-size optimization. SDSL-Solver's modular approach allows flexible adaptation to matrix structure, delivers robust performance even as spectral properties deteriorate across IPM iterations, and enables high-performance on commodity multi-node architectures. Beyond LP/QP, the methodology is extensible to other convex programming frameworks that rely on iterative sparse linear solves.
The theoretical implications suggest that further advancements could be realized in (i) dynamic thresholding for filtering via spectral analysis, (ii) recursive interface decomposition for BBD to mitigate interface solve bottlenecks, and (iii) automated switching strategies driven by convergence statistics rather than fixed heuristics.
Future Directions
- Accelerator Support: Implementation on GPU/NPU and integration with high-bandwidth memory architectures will further reduce solution time for very large problems.
- Adaptive Threshold Selection: Automated filtering threshold choice, possibly using eigenvalue distribution estimates, can reduce manual tuning and improve solver robustness.
- Extension to Nonlinear Domains: Support for second-order cone and semidefinite programming will extend SDSL-Solver's applicability across a wider optimization landscape.
Conclusion
SDSL-Solver provides a scalable, distributed, and numerically robust framework for sparse linear solves in IPMs. Its integration of Krylov subspace methods, adaptive preconditioning, and dynamic distributed strategies delivers strong speedups and superior stability compared to state-of-the-art direct and iterative solvers. These techniques collectively enable industrial-scale optimization, with promising directions for further performance improvements and broader applicability in future AI-driven optimization research.
(2604.23979)