Restarted Halpern PDHG Method

• The paper introduces rHPDHG, a matrix-free, primal-dual method that achieves accelerated two-stage linear convergence for large-scale linear programming and saddle-point problems.
• The method leverages Halpern fixed-point iterations with an outer restarting schedule and reflection steps to reduce dependence on pessimistic global condition constants.
• GPU-aware implementations and adaptive parameter tuning enhance its practical performance, ensuring robust convergence in both feasible and infeasible LP settings.

The Restarted Halpern Primal-Dual Hybrid Gradient (rHPDHG) Method is a matrix-free, primal–dual first-order optimization algorithm designed to solve large-scale linear programming (LP) and related saddle-point problems. By incorporating Halpern-type fixed-point iterations, an outer restarting schedule, and—in some variants—a reflection step, rHPDHG achieves accelerated two-stage linear convergence, improved practical performance, and efficient parallelization. This article delineates the algorithmic principles, theoretical foundations, condition measures, extensions, implementation details, and empirical performance of rHPDHG and related advances.

1. Algorithmic Foundation and Structure

The rHPDHG method addresses LPs of the form: $\min_{x \geq 0,~Ax = b} c^\top x,$ using a primal–dual splitting approach. The core innovation over classical PDHG algorithms is the adoption of Halpern fixed-point anchoring. Specifically, within each inner epoch (indexed by $k$), the update is: $$z^{k+1} = \frac{k+1}{k+2} T(z^k) + \frac{1}{k+2} z^0,$$ where $T$ denotes the standard PDHG operator and %%%%2%%%% is the anchor point (typically the initial point or the output of a prior epoch). After a prescribed number of inner iterations, a restart is triggered by re-anchoring to the current best iterate. This process is repeated across outer epochs, yielding a doubly-nested iterative scheme.

A further enhancement, termed reflected rHPDHG (or $\mathrm{r^2HPDHG}$), applies Halpern iteration to the reflection of the PDHG operator, i.e., using $R = 2T - I$ in place of $T$, thereby taking more aggressive steps.

The restart criterion can be fixed frequency or adaptive, based on monitoring the fixed-point residual: $r(z) = \|z - T(z)\|_P,$ where $P$ is a parameter matrix scaling the norm. Adaptive rules include restarting upon sufficient decrease of the residual, lack of local progress, or reaching a maximum inner-epoch length.

2. Theoretical Convergence and Two-Stage Behavior

A rigorous analysis of rHPDHG establishes global accelerated linear convergence for feasible and bounded LPs. In particular, for a restart period $k_*\geq 2e/\alpha\eta$ (with $\alpha\eta$ a sharpness parameter), the following holds for the restart iterates: $$\mathrm{dist}(z^{n+1,0}, \mathcal{Z}^*) \leq \left(\frac{1}{e}\right)^{n+1} \mathrm{dist}(z^{0,0}, \mathcal{Z}^*),$$ leading to an iteration complexity $O((1/\alpha\eta)\log(1/\varepsilon))$ for target accuracy $\varepsilon$ (Lu et al., 23 Jul 2024).

The convergence admits a refined two-stage characterization:

• Stage I: Active Set Identification—The algorithm identifies the correct support structure by partitioning variable indices into nonbasic ($N$), nondegenerate basic ($B_1$), and degenerate ($B_2$). The duration depends on the nondegeneracy metric

$$\delta = \min \left\{ \min_{i\in N} \frac{c_i+A_i^\top y^*}{\|A\|_2},~ \min_{i \in B_1} x^*_i \right\},$$

resulting in a complexity $O(\|A\|_2 / (\alpha_{L_1} \delta))$, with $\alpha_{L_1}$ a local sharpness constant (Lu et al., 23 Jul 2024).

• Stage II: Local Linear Convergence—Once identification is achieved, the algorithm enjoys accelerated local linear convergence with a rate dictated by another local sharpness $\alpha_{L_2}$, and complexity $O((\|A\|_2/\alpha_{L_2})\log(1/\varepsilon))$.

This refined analysis replaces prior global Hoffman-type constants—often exponentially pessimistic—with practical, locally sharp rates (Lu et al., 2023).

3. Condition Measures: Sharpness, Limiting Error Ratio, and Geometric Quantities

The complexity of rHPDHG is fundamentally governed by intrinsic condition measures:

• LP sharpness ($\mu_p$ and $\mu_d$) quantifies the sensitivity of the objective to deviation from the optimum, capturing how rapidly the objective increases away from the optimal set (Xiong et al., 2023, Xiong et al., 2023).
• Limiting Error Ratio measures local geometric feasibility properties, defined as

$$\theta_p(x) = \frac{\mathrm{Dist}(x, \mathcal{P})}{\mathrm{Dist}(x, \mathbb{R}^n_+)}$$

where $\mathcal{P}$ is the feasible set. The limiting value as $x\to x^*$ informs algorithm complexity near the optimum.

• Geometric Condition Numbers such as $\Phi$ in (Xiong, 5 Oct 2024), defined explicitly from the optimal basis and solution, allow for accessible, closed-form iteration bounds.

A central insight is the reciprocal relationship between sharpness and the iteration bound: as the LP becomes "sharper" (larger $\mu_p$), the algorithm converges more rapidly. The two-stage complexity splits accordingly: an initial phase determined by sharpness and stability under data perturbation, followed by a faster, locally conditioned phase controlled by the minimal nonzero solution component and the conditioning of the optimal basis matrix.

4. Handling Infeasibility and Certificates

rHPDHG, like its base PDHG, is effective for both feasible and infeasible LPs. For infeasible cases, the method converges to the infimal displacement vector $v$ (the minimum-norm element of $\mathrm{cl}\,\mathrm{Range}(T-I)$). The updated Halpern and reflection strategies yield accelerated linear convergence in extracting infeasibility certificates, even in the absence of strong regularity conditions (Lu et al., 23 Jul 2024, Jiang et al., 2023). This accelerates and stabilizes the identification of invalid LP instances.

5. Extensions: Reflective and Parameter-Free Variants

The reflected rHPDHG extension ($\mathrm{r^2HPDHG}$) applies Halpern averaging with the reflected operator $R=2T-I$, reducing the required inner-epoch length by half and improving per-epoch contraction rates: $$z^{k+1} = \frac{k+1}{k+2} (2T(z^k) - z^k) + \frac{1}{k+2} z^0$$ All convergence guarantees of rHPDHG extend, with constants improved by a factor of 2 (Lu et al., 23 Jul 2024, Lu et al., 18 Jul 2025).

Related methods, such as adaptive PDHG with automatic stepsize tuning (Goldstein et al., 2013) and relaxed, parameter-free variants with line search (McManus et al., 21 Mar 2025), provide alternate avenues for removing user input dependence on problem constants, using residual balancing or operator-based line searches.

6. GPU-Aware Implementation and Empirical Evaluation

The matrix-free, matrix–vector-multiplication structure of rHPDHG makes it ideal for GPU acceleration. Reference implementations (e.g., HPDLP in Julia/CUDA (Lu et al., 23 Jul 2024), cuPDLP+ (Lu et al., 18 Jul 2025)) exploit parallelism and tailored restart/anchor management. In large-scale tests on MIPLIB LP relaxations:

• rHPDHG and $\mathrm{r^2HPDHG}$ realize per-instance speedups of $1.27\times$–$1.33\times$ over previous PDHG-based solvers at moderate accuracy, and $2\times$–$4\times$ at high accuracy or with presolve (Lu et al., 18 Jul 2025).
• The PID-controlled primal weight update in cuPDLP+ dynamically balances primal and dual progress by updating the weight $\omega$ using a logarithmic-scale PID controller:

$$\log (\omega^{n+1}) = \log (\omega^n) - (K_P e^n + K_I \sum_{i=1}^n e^i + K_D (e^n - e^{n-1})),$$

where $e^n$ is the log-scale ratio of primal and dual residuals (Lu et al., 18 Jul 2025).

• Constant stepsize rules (e.g., $\eta=0.99/\|A\|_2$) remove the need for per-iteration spectral estimation.
• Empirical performance matches or exceeds that of existing GPU-based first-order solvers; e.g., cuPDLP+, incorporating r^2HPDHG, achieves substantial gains in solve time and robustness.

7. Practical Guidelines and Perspective

For effective application of rHPDHG:

• Exploit row preconditioning and tuned stepsize ratios to moderate the iteration bound’s dependence on condition numbers and sharpness (Xiong et al., 2023).
• Use adaptive or residual-based restart criteria to maintain fast local convergence after identification of the optimal structure (Lu et al., 23 Jul 2024, Lu et al., 18 Jul 2025).
• Apply reflection (r^2HPDHG) to further speed up contraction, especially in ill-conditioned regimes.
• Warm-start or combine with learning-based approaches (e.g., PDHG-Net) for very large instances or in hybrid architectures (Li et al., 4 Jun 2024).
• For challenging test cases, manipulate the disparity ratio of the optimal variables to control the iteration count in Stage I (Xiong, 1 Jan 2025).

In summary, the restarted Halpern PDHG method unifies fixed-point averaged iteration, restart-based acceleration, local geometric conditioning, and scalable, parallel-first-order optimization into a framework with theoretically optimal and empirically robust convergence for large-scale LP and convex saddle-point problems. Recent advances in solution certification, parameter-free variants, step-size scheduling, and GPU-tailored implementations have made rHPDHG and its extensions competitive across diverse domains and datasets.