Arrow–Hurwicz Iteration

• Arrow–Hurwicz iteration is a primal-dual descent-ascent method for solving convex–concave saddle-point problems by updating primal variables through descent and dual variables through ascent.
• The method operates in both continuous forms—yielding differential equations—and discrete forms using projections, proximal updates, and semi-implicit schemes for enhanced stability.
• Extensions such as damping, regularization, and coupled corrections address oscillatory issues, making it crucial for applications in distributed optimization, Navier–Stokes solvers, and sphere packing.

Arrow–Hurwicz iteration denotes a class of primal-dual descent-ascent methods for saddle-point problems and constrained optimization. In its canonical form, it evolves a primal variable by descent in the Lagrangian and a dual variable by ascent in the constraint residual; in continuous time this gives a differential system, while in discrete time it gives stationary splitting, projection, proximal, or semi-implicit iterations. Recent work treats Arrow–Hurwicz both as a classical method in its own right and as a structural template underlying distributed least-squares solvers, incompressible Navier–Stokes iterations, sphere-packing algorithms, and modern primal-dual schemes such as PDHG (Niederländer, 2022, Liu et al., 2017, Geredeli et al., 2022, Li et al., 2024).

1. Canonical saddle-point structure

The basic setting is a convex-concave Lagrangian. In the Hilbert-space formulation for linearly constrained convex minimization,

$\inf\{\,f(x)\mid Ax-b=0\,\},$

the associated Lagrangian is

$L(x,\lambda)=f(x)+\langle \lambda,Ax-b\rangle_Y,$

and saddle points satisfy

$$\nabla f(\bar x)+A^*\bar\lambda=0,\qquad A\bar x-b=0.$$

The continuous-time Arrow–Hurwicz system is then

$$\dot x+\nabla f(x)+A^*\lambda=0,\qquad \dot\lambda+b-Ax=0,$$

that is, primal steepest descent coupled to dual ascent in the feasibility defect (Niederländer, 2022).

A more general projected form arises for convex-concave saddle problems

$$\inf_{u\in C}\sup_{p\in K}L(u,p),\qquad L(u,p)=\langle Au,p\rangle+F(u)-G(p),$$

with closed convex sets $C\subset V$ and $K\subset W$. The discrete Arrow–Hurwicz update is a projected ascent in $p$ and a projected descent in $u$, with an extrapolated primal variable in the generalized scheme studied for nonsmooth saddle systems (Phan et al., 2017). In the linear special case

$$L(u,p)=\langle Au,p\rangle+\langle f,u\rangle-\langle g,p\rangle,$$

the scheme reduces to a Chambolle–Pock-type update, convergent under

$L(x,\lambda)=f(x)+\langle \lambda,Ax-b\rangle_Y,$0

(Phan et al., 2017).

The same template appears in network optimization. For a distributed least-squares problem with global system

$L(x,\lambda)=f(x)+\langle \lambda,Ax-b\rangle_Y,$1

equation $L(x,\lambda)=f(x)+\langle \lambda,Ax-b\rangle_Y,$2 is assigned to node $L(x,\lambda)=f(x)+\langle \lambda,Ax-b\rangle_Y,$3, each node stores $L(x,\lambda)=f(x)+\langle \lambda,Ax-b\rangle_Y,$4, and consensus is enforced by the graph Laplacian $L(x,\lambda)=f(x)+\langle \lambda,Ax-b\rangle_Y,$5. The distributed problem becomes

$L(x,\lambda)=f(x)+\langle \lambda,Ax-b\rangle_Y,$6

with

$L(x,\lambda)=f(x)+\langle \lambda,Ax-b\rangle_Y,$7

The resulting continuous-time Arrow–Hurwicz–Uzawa-type flow is

$L(x,\lambda)=f(x)+\langle \lambda,Ax-b\rangle_Y,$8

so the primal variable is the stack of node states and the dual variable enforces consensus (Liu et al., 2017).

2. Continuous-time dynamics and asymptotic theory

For the classical differential system, the governing operator

$L(x,\lambda)=f(x)+\langle \lambda,Ax-b\rangle_Y,$9

is monotone and continuous, hence maximally monotone, and the Cauchy problem is well posed for every initial datum (Niederländer, 2022). A nonexpansiveness property holds: for two solutions $$\nabla f(\bar x)+A^*\bar\lambda=0,\qquad A\bar x-b=0.$$0 and $$\nabla f(\bar x)+A^*\bar\lambda=0,\qquad A\bar x-b=0.$$1,

$$\nabla f(\bar x)+A^*\bar\lambda=0,\qquad A\bar x-b=0.$$2

The velocity norm is non-increasing, and if the saddle set is nonempty then trajectories are bounded (Niederländer, 2022).

A central identity is the “Lagrangian identity”

$$\nabla f(\bar x)+A^*\bar\lambda=0,\qquad A\bar x-b=0.$$3

More generally, for any saddle point $$\nabla f(\bar x)+A^*\bar\lambda=0,\qquad A\bar x-b=0.$$4,

$$\nabla f(\bar x)+A^*\bar\lambda=0,\qquad A\bar x-b=0.$$5

This yields boundedness, integrability of the primal-dual gap, and existence of limits of distances to saddle points (Niederländer, 2022).

In the merely convex case, the Cesàro averages

$$\nabla f(\bar x)+A^*\bar\lambda=0,\qquad A\bar x-b=0.$$6

satisfy

$$\nabla f(\bar x)+A^*\bar\lambda=0,\qquad A\bar x-b=0.$$7

and converge weakly to a saddle point. Under strong convexity, the actual trajectory satisfies

$$\nabla f(\bar x)+A^*\bar\lambda=0,\qquad A\bar x-b=0.$$8

with weak convergence of $$\nabla f(\bar x)+A^*\bar\lambda=0,\qquad A\bar x-b=0.$$9 to a saddle point (Niederländer, 2022). If $$\dot x+\nabla f(x)+A^*\lambda=0,\qquad \dot\lambda+b-Ax=0,$$0 or $$\dot x+\nabla f(x)+A^*\lambda=0,\qquad \dot\lambda+b-Ax=0,$$1 is bounded below, this strengthens to strong convergence and improved decay such as

$$\dot x+\nabla f(x)+A^*\lambda=0,\qquad \dot\lambda+b-Ax=0,$$2

When $$\dot x+\nabla f(x)+A^*\lambda=0,\qquad \dot\lambda+b-Ax=0,$$3 is $$\dot x+\nabla f(x)+A^*\lambda=0,\qquad \dot\lambda+b-Ax=0,$$4, $$\dot x+\nabla f(x)+A^*\lambda=0,\qquad \dot\lambda+b-Ax=0,$$5 is bounded, and either $$\dot x+\nabla f(x)+A^*\lambda=0,\qquad \dot\lambda+b-Ax=0,$$6 or $$\dot x+\nabla f(x)+A^*\lambda=0,\qquad \dot\lambda+b-Ax=0,$$7 is bounded below, the paper derives exponential-type decay, including the underdamped, critically damped, and overdamped regimes familiar from damped harmonic oscillators (Niederländer, 2022).

A regularized variant augments the flow by vanishing Tikhonov terms: $$\dot x+\nabla f(x)+A^*\lambda=0,\qquad \dot\lambda+b-Ax=0,$$8 Under $$\dot x+\nabla f(x)+A^*\lambda=0,\qquad \dot\lambda+b-Ax=0,$$9 assumptions on $$\inf_{u\in C}\sup_{p\in K}L(u,p),\qquad L(u,p)=\langle Au,p\rangle+F(u)-G(p),$$0, the differential inequalities

$$\inf_{u\in C}\sup_{p\in K}L(u,p),\qquad L(u,p)=\langle Au,p\rangle+F(u)-G(p),$$1

and nonemptiness of the saddle set, solutions converge strongly to the minimum-norm saddle point

$$\inf_{u\in C}\sup_{p\in K}L(u,p),\qquad L(u,p)=\langle Au,p\rangle+F(u)-G(p),$$2

If, in addition,

$$\inf_{u\in C}\sup_{p\in K}L(u,p),\qquad L(u,p)=\langle Au,p\rangle+F(u)-G(p),$$3

then the zero is unique and stronger distance estimates follow (Battahi et al., 2024). For $$\inf_{u\in C}\sup_{p\in K}L(u,p),\qquad L(u,p)=\langle Au,p\rangle+F(u)-G(p),$$4, $$\inf_{u\in C}\sup_{p\in K}L(u,p),\qquad L(u,p)=\langle Au,p\rangle+F(u)-G(p),$$5, the framework yields explicit decay bounds, and for $$\inf_{u\in C}\sup_{p\in K}L(u,p),\qquad L(u,p)=\langle Au,p\rangle+F(u)-G(p),$$6 the key residual and velocity estimates simplify to $$\inf_{u\in C}\sup_{p\in K}L(u,p),\qquad L(u,p)=\langle Au,p\rangle+F(u)-G(p),$$7 rates for the quantities controlled by the theorem (Battahi et al., 2024).

3. Discrete iterations, projections, and step-size restrictions

The explicit projected scheme for

$$\inf_{u\in C}\sup_{p\in K}L(u,p),\qquad L(u,p)=\langle Au,p\rangle+F(u)-G(p),$$8

updates $$\inf_{u\in C}\sup_{p\in K}L(u,p),\qquad L(u,p)=\langle Au,p\rangle+F(u)-G(p),$$9 by

$C\subset V$0

and $C\subset V$1 by

$C\subset V$2

followed by primal extrapolation. Convergence is proved under

$C\subset V$3

and

$C\subset V$4

These conditions explicitly couple the Lipschitz moduli of $C\subset V$5 and $C\subset V$6, the operator norm $C\subset V$7, and the product $C\subset V$8 (Phan et al., 2017).

A semi-implicit variant replaces the descent direction by the filtered term

$C\subset V$9

so that the primal update becomes

$K\subset W$0

The corresponding convergence condition is

$K\subset W$1

notably without any dependence on $K\subset W$2. The paper identifies this removal of $K\subset W$3 from the admissible parameter range as the key robustness advantage of the semi-implicit scheme (Phan et al., 2017).

For distributed least squares over a fixed connected graph, Euler discretization of the Arrow–Hurwicz–Uzawa flow yields

$K\subset W$4

$K\subset W$5

In shifted stacked form,

$K\subset W$6

There exists a critical step size

$K\subset W$7

such that $K\subset W$8 gives exponential convergence of all node states to the least-squares solution, whereas $K\subset W$9 yields divergence for some initial conditions. If the graph/equation spanning condition fails, divergence can occur for any $p$0 (Liu et al., 2017).

4. Oscillation, nonconvergence, and stabilization

A recurrent theme in the literature is that the classical Arrow–Hurwicz mechanism can be only marginally stable. In the distributed least-squares flow, the shifted dynamics are linear,

$p$1

with Lyapunov function

$p$2

The flow converges to consensus least squares if and only if for every Laplacian eigenvector $p$3,

$p$4

If this spanning condition fails, the system matrix acquires a purely imaginary eigenvalue and oscillatory solutions occur (Liu et al., 2017).

A similar phenomenon appears in the proximal Arrow–Hurwicz discretization of

$p$5

where the primal step is followed by a dual step using $p$6. For the example

$p$7

the iteration becomes

$p$8

which is the symplectic-Euler discretization of

$p$9

The associated Hamiltonian

$u$0

is conserved by the continuous flow, so trajectories lie on closed level sets rather than converging to the saddle. The paper uses this to explain the periodic behavior and six-point closed loop observed for proximal Arrow–Hurwicz (Li et al., 2024).

In sphere packing, classical Arrow–Hurwicz is likewise reported to be too oscillatory. For the smooth two-sphere one-dimensional constraint

$u$1

the linearized algorithm has characteristic polynomial

$u$2

so the eigenvalues are purely imaginary and the equilibrium is a center. The damped multi-step variant instead yields a cubic polynomial

$u$3

with sufficient stability condition

$u$4

For the nonsmooth constraint $u$5, the undamped equilibrium is already asymptotically stable in the two-sphere one-dimensional case, with fastest convergence at $u$6, while the damped version is stable under

$u$7

(Degond et al., 2016).

Across these settings, damping, implicit filtering, and vanishing regularization recur as remedies for purely imaginary spectrum, periodic orbits, and slow decay (Degond et al., 2016, Phan et al., 2017, Battahi et al., 2024).

5. Major applications and specialized formulations

Arrow–Hurwicz iteration is best understood as a reusable saddle-point mechanism whose concrete form depends on the objective, the constraint operator, and the discretization.

Setting	Arrow–Hurwicz formulation	Reported result
Linearly constrained convex minimization	Continuous primal-dual flow $u$8, $u$9	$$L(u,p)=\langle Au,p\rangle+\langle f,u\rangle-\langle g,p\rangle,$$0 ergodic gap decay; stronger rates under curvature (Niederländer, 2022)
Network least squares	AHU flow with consensus constraint $$L(u,p)=\langle Au,p\rangle+\langle f,u\rangle-\langle g,p\rangle,$$1	Necessary and sufficient graph-dependent convergence criterion (Liu et al., 2017)
Steady Navier–Stokes FEM	Decoupled velocity solve plus pressure correction	Grad-div and Anderson acceleration improve convergence (Geredeli et al., 2022)
Divergence-free VEM Navier–Stokes	A-H iteration on nonlinear saddle system from VEM discretization	Geometric convergence with mesh-independent contraction factor (Du et al., 16 Jul 2025)
Sphere packing	Damped, second-order, multi-step Arrow–Hurwicz	Better stability than classical AHA in tested cases (Degond et al., 2016)

In incompressible steady Navier–Stokes discretizations, the attraction of Arrow–Hurwicz is the decoupling of velocity and pressure. A typical finite-element form computes $$L(u,p)=\langle Au,p\rangle+\langle f,u\rangle-\langle g,p\rangle,$$2 from a convection-diffusion-type problem,

$$L(u,p)=\langle Au,p\rangle+\langle f,u\rangle-\langle g,p\rangle,$$3

followed by a pressure correction

$$L(u,p)=\langle Au,p\rangle+\langle f,u\rangle-\langle g,p\rangle,$$4

Grad-div stabilization adds

$$L(u,p)=\langle Au,p\rangle+\langle f,u\rangle-\langle g,p\rangle,$$5

to the momentum equation, improving incompressibility and nonlinear convergence. Under $$L(u,p)=\langle Au,p\rangle+\langle f,u\rangle-\langle g,p\rangle,$$6 and

$$L(u,p)=\langle Au,p\rangle+\langle f,u\rangle-\langle g,p\rangle,$$7

the method becomes equivalent to the iterated Picard penalty method, with a linear convergence rate bounded by the known IPP rate under

$$L(u,p)=\langle Au,p\rangle+\langle f,u\rangle-\langle g,p\rangle,$$8

(Geredeli et al., 2022).

The numerical consequences are substantial. For the $$L(u,p)=\langle Au,p\rangle+\langle f,u\rangle-\langle g,p\rangle,$$9 lid-driven cavity, earlier Arrow–Hurwicz results reported around $L(x,\lambda)=f(x)+\langle \lambda,Ax-b\rangle_Y,$00 iterations with Taylor–Hood elements, whereas the stabilized Scott–Vogelius version converged in about $L(x,\lambda)=f(x)+\langle \lambda,Ax-b\rangle_Y,$01 iterations under good parameters. The same study reports successful computations at $L(x,\lambda)=f(x)+\langle \lambda,Ax-b\rangle_Y,$02 and $L(x,\lambda)=f(x)+\langle \lambda,Ax-b\rangle_Y,$03, and finds that the combination of grad-div stabilization and Anderson acceleration is the best-performing variant (Geredeli et al., 2022).

For divergence-free mixed virtual element discretizations of steady incompressible Navier–Stokes, the Arrow–Hurwicz iteration is analyzed directly at the algebraic level. Starting from a discrete Stokes solve, the update reads

$L(x,\lambda)=f(x)+\langle \lambda,Ax-b\rangle_Y,$04

$L(x,\lambda)=f(x)+\langle \lambda,Ax-b\rangle_Y,$05

Under the discrete smallness condition $L(x,\lambda)=f(x)+\langle \lambda,Ax-b\rangle_Y,$06 and

$L(x,\lambda)=f(x)+\langle \lambda,Ax-b\rangle_Y,$07

the paper proves

$L(x,\lambda)=f(x)+\langle \lambda,Ax-b\rangle_Y,$08

with $L(x,\lambda)=f(x)+\langle \lambda,Ax-b\rangle_Y,$09 independent of the mesh size (Du et al., 16 Jul 2025).

6. Relation to PDHG, Chambolle–Pock, and modern primal-dual methods

Contemporary work often interprets Arrow–Hurwicz as the baseline member of a broader family of primal-dual algorithms. In the generalized Lasso setting,

$L(x,\lambda)=f(x)+\langle \lambda,Ax-b\rangle_Y,$10

is rewritten as

$L(x,\lambda)=f(x)+\langle \lambda,Ax-b\rangle_Y,$11

The proximal Arrow–Hurwicz algorithm performs

$L(x,\lambda)=f(x)+\langle \lambda,Ax-b\rangle_Y,$12

$L(x,\lambda)=f(x)+\langle \lambda,Ax-b\rangle_Y,$13

PDHG keeps the same primal update but uses the extrapolated point

$L(x,\lambda)=f(x)+\langle \lambda,Ax-b\rangle_Y,$14

in the dual step. Expanded in increment form,

$L(x,\lambda)=f(x)+\langle \lambda,Ax-b\rangle_Y,$15

$L(x,\lambda)=f(x)+\langle \lambda,Ax-b\rangle_Y,$16

The terms $L(x,\lambda)=f(x)+\langle \lambda,Ax-b\rangle_Y,$17 and $L(x,\lambda)=f(x)+\langle \lambda,Ax-b\rangle_Y,$18 are the coupled $L(x,\lambda)=f(x)+\langle \lambda,Ax-b\rangle_Y,$19-correction and $L(x,\lambda)=f(x)+\langle \lambda,Ax-b\rangle_Y,$20-correction emphasized by the paper (Li et al., 2024).

A high-resolution ODE model for PDHG is

$L(x,\lambda)=f(x)+\langle \lambda,Ax-b\rangle_Y,$21

with second-derivative corrections retained in the full derivation. Discrete PDHG is identified as the implicit Euler discretization of this system, whereas proximal Arrow–Hurwicz corresponds to a mixed explicit-implicit Euler discretization of the lower-resolution system obtained by dropping the $L(x,\lambda)=f(x)+\langle \lambda,Ax-b\rangle_Y,$22 coupling terms (Li et al., 2024). This gives a precise dynamical explanation for a widely observed algorithmic contrast: Arrow–Hurwicz-type discretizations may preserve periodic structure, while PDHG’s coupled correction acts as a small but essential perturbation that yields convergence.

The resulting rates are correspondingly different. For PDHG, the averaged iterates satisfy an ergodic $L(x,\lambda)=f(x)+\langle \lambda,Ax-b\rangle_Y,$23 bound, and when one component is strongly convex the averaged primal iterates satisfy

$L(x,\lambda)=f(x)+\langle \lambda,Ax-b\rangle_Y,$24

in the one-step-size setting described by the paper (Li et al., 2024). The discrete numerical-error Lyapunov quantity is also monotone decreasing, which the paper treats not as a nuisance term but as a dissipative mechanism.

This suggests that Arrow–Hurwicz iteration is less a single fixed algorithm than a foundational primal-dual architecture. Chambolle–Pock extrapolation, PDHG’s coupled correction, semi-implicit filtering, grad-div stabilization, Anderson acceleration, damping, and Tikhonov regularization all modify the same descent-ascent core in order to recover stronger contraction, suppress oscillation, or improve robustness (Phan et al., 2017, Geredeli et al., 2022, Li et al., 2024, Battahi et al., 2024).