On the Infimal Sub-differential Size of Primal-Dual Hybrid Gradient Method and Beyond

Linear programming (LP) [5, 16] is a fundamental problem in mathematical optimization and operations research. The traditional LP solvers are essentially based on either simplex method or barrier method. However, it is highly challenging to further scale up these two methods. The computational bottleneck of both methods is solving linear equations, which does not scale well on modern computing resources, such as GPUs and distributed systems. Recently, a numerical study [2] demonstrates that an enhanced version of Primal-Dual Hybrid Gradient (PDHG) method, called PDLP¹¹1The solver is open-sourced at Google OR-Tools. , can reliably solve LP problems to high precision. As a first-order method, the computational bottleneck of PDHG is matrix-vector multiplication, which can efficiently leverage GPUs and distributed system. It sheds light on solving huge-scale LP where the data need to be stored in a distributed system.

This work is motivated by this recent line of research on using PDHG for solving large-scale LP. The convergence guarantee of PDHG, for generic convex-concave problems or for LP, has been extensively studied in the literature [66, 10, 11, 1, 4]. Most of the previous convergence analyses on PDHG are based on either the primal-dual gap or the distance to the optimal solution setting. However, it can be highly non-trivial to compute these metrics while running the algorithms. For primal-dual gap, the iterated solutions when solving LP are often infeasible to the primal and/or the dual problem until the algorithm identifies an optimal solution; thus, the primal-dual gap at any iterates may likely be infinity. For the distance to the optimal solution set, it is not computable unless the optimal solution set is known. Furthermore, even if these values can be computed, they often oscillate dramatically in practice, which makes it hard to evaluate the solving progress. Motivated by these issues, the goal of this paper is to answer the following question:

Is there an evaluatable and monotonically decaying progress metric for analyzing PDHG?

We provide an affirmative answer to the above question by proposing a new progress metric for analyzing PDHG, which we dub infimal sub-differential size (IDS). IDS essentially measures the distance between $0$ and the sub-differential set of the objective. It is a natural extension of the gradient norm to potentially non-smooth problems. Compared with other progress measurements for PDHG, such as primal-dual gap and distance to optimality, IDS always has a finite value and can be computable directly only using the information of the current solution without the need of knowing the optimal solution set. More importantly, IDS decays monotonically during the solving process, making it a natural progress measurement for PDHG. The design of IDS also take into consideration the geometry of PDHG by using a natural norm of PDHG.

Furthermore, most previous PDHG analyzes [10, 11] focus on the ergodic rate, that is, the performance of the average PDHG iterates. On the other hand, it is frequently observed that the last iteration of PDHG has comparable and sometimes superior performance compared to the average iteration in practice [2]. Our analysis on IDS provides an explanation of such behavior. We show that IDS of the last iterate of PDHG has $\mathcal{O}(1/k)$ convergence rate for convex-concave problems. As a direct byproduct, the primal dual gap at the last iterate has $\mathcal{O}(1/\sqrt{k})$ convergence rate, which is inferior to the $\mathcal{O}(1/k)$ convergence of the average iterate [11]. This explains why the average iterate can be a good option. Furthermore, we show that if the problem further has metric sub-regularity, a condition that is satisfied by many optimization problems, IDS and primal dual gap both enjoy linear convergence, whereas average iterates can often have sublinear convergence. This, accompanied with other results under metric sub-regularity in literature, explains the numerical observation that the last iteration may have superior behaviors.

More formally, we consider a generic minimax problem with a bilinear interaction term:

where $f:\mathbb{R}^{n}\rightarrow(-\infty,+\infty]$, $g:\mathbb{R}^{m}\rightarrow(-\infty,+\infty]$ are simple proper lower semi-continuous convex functions and $A\in\mathbb{R}^{m\times n}$ is a matrix. Note that functions $f(x)$ and $g(y)$ can encode the constraints $x\in\mathcal{X}$ and $y\in\mathcal{Y}$ by using an indicator function. Throughout the paper, we assume that there is a finite optimal solution to the minimax problem (1).The primal-dual form of linear programming can be formulated as an instance of (1):

Beyond LP, problems of form (1) are ubiquitous in statistics and machine learning. For example,

• •

  In image processing, total variation (TV) denoising model [7] plays a central role in variational methods for imaging and it can be reformulated as $\min_{x\geq 0}\max_{y:\|y\|_{\infty}\leq 1}\mathcal{L}(x,y)=\frac{\lambda}{2}\|x-f\|_{2}^{2}+y^{T}\nabla x$, where $\nabla$ is the differential matrix.

• •

  In statistics, LASSO [57, 12] is a popular regression model, whose primal-dual formulation is $\min_{x}\max_{y}\mathcal{L}(x,y)=\lambda\|x\|_{1}+y^{T}Ax-(\frac{1}{2}\|y\|_{2}^{2}+b^{T}y)$.

• •

  In machine learning, support vector machine (SVM) [15, 6] is a classic supervised learning method for classification and can be formulated as $\min_{x}\max_{-1/n\leq y_{i}/b_{i}\leq 0}\mathcal{L}(x,y)=\frac{\lambda}{2}\|x\|_{2}^{2}+y^{T}Ax-\frac{1}{n}\sum_{i=1}^{n}\frac{ny_{i}}{b_{i}}$.

Primal-dual hybrid gradient method (PDHG, a.k.a. Chambolle and Pock algorithm [10]) described in Algorithm 3 is one of the most popular algorithms for solving the structured minimax problem (1). It has been extensively used in the field of image processing: total-variation image denoising [7], multimodal medical imaging [38], computation of nonlinear eigenfunctions [28], and many others. More recently, PDHG also serves as the base algorithm for a new large-scale LP solver PDLP [2, 3, 4].

The contribution of the paper can be summarized as follows:

• •

  We propose to identify the “right” progress metric when studying an optimization algorithm. We propose a new progress metric IDS for analyzing PDHG, and show that it monotonically decays (see Section 2).

• •

  We show that IDS converges to $0$ with $\mathcal{O}(1/k)$ convergence rate when solving a convex-concave minimax problem (1) (Theorem 1). As a by-product, this shows $\mathcal{O}(1/\sqrt{k})$ non-ergodic convergent rate of the duality gap for PDHG when solving convex-concave minimax problems.

• •

  We show that IDS converges to $0$ with a linear convergence rate if the problem satisfies metric sub-regularity, which is satisfied by many applications.

• •

  We extend the above results of PDHG to other classic algorithms for minimax problems, including proximal point method, alternating direction method of multipliers and linearized alternating direction method of multipliers.

• •

  The proofs of the above results are surprisingly simple, which is in contrast to the recent literature on the last iterate convergence for cousin algorithms.

In this section, we introduce a new progress measurement, which we dub infimal sub-differential size (IDS), for PDHG. IDS is a natural extension of the squared norm of gradient for smooth optimization problems to non-smooth optimization problems. Compared to other progress measurements for PDHG, such as primal-dual gap and distance to optimality, IDS always has a finite value and is computable directly without the need to know the optimal solution set. More importantly, unlike the other metric that may oscillate over time, the IDS monotonically decays along the iteration of the PDHG, further suggesting that the IDS is a natural progress measurement for PDHG. Here is the formal definition of IDS:

Note that the first-order optimality condition of (1) is $0\in\mathcal{F}(z)$. In other words, if the original minimax problem (1) is feasible and bounded, the minimization problem

has optimal value $0$, and furthermore (3) and (1) share the same optimal solution set, thus IDS is a valid progress measurement.

When measuring the distance between $0$ and the set $\mathcal{F}(z)$, a crucial part is the choice of the norm. The simplest norm which is also used in practice is the $\ell_{2}$ norm. However, this may not be the natural norm for all algorithms. To understand the intuition, Figure 1 plots the trajectory of PDHG and proximal point method (PPM) to solve a simple unconstrained bilinear problem

We can observe that the iterations of two algorithms exhibit spiral patterns toward the unique saddle point $(0,0)$. Nevertheless, in contrast to the circular pattern of PPM, the spiral of PDHG spins in a tilting elliptical manner. The fundamental reason for such an elliptical structure is the asymmetricity of the primal and dual steps in PDHG. The elliptical structures of the iterates can be captured by the choice of norms in the analysis. This is consistent with the fact that while the analysis of PPM [55] utilizes $\ell_{2}$ norm, the analysis of PDHG [11, 4, 34] often utilizes a different norm defined by $P_{s}=\begin{pmatrix}\frac{1}{s}I_{n}&-A^{T}\\ -A&\frac{1}{s}I_{m}\end{pmatrix}$. The $P_{s}$ has already been used in analyzing PDHG, for example, [34] studies the contraction of PDHG under $P_{s}$ norm and motivate variants of PDHG from this perspective. We here propose to use $P_{s}^{-1}$ in the sub-differentiable (dual) space.

In previous works, the convergence guarantee of PDHG is often stated with respect to the duality gap [10, 11], $\max_{\hat{x}\in\mathcal{X},\hat{y}\in\mathcal{Y}}\mathcal{L}(x,\hat{y})-\mathcal{L}(\hat{x},y)$, or the distance to the optimal solution set [4, 26], $dist(z,\mathcal{Z}^{*})$. Unfortunately, these values are barely computable or informative because the iterates are usually infeasible until the end of the algorithms (thus, the duality gap at the iterates is often infinity) and the optimal solution set is unknown to the users beforehand. In practice, people instead use residual as a performance metric because of its computational efficiency. For example, the KKT residual, i.e., the violation of primal feasibility, dual feasibility, and complementary slackness, is widely used in the LP solver as the performance metric. Indeed, the KKT residual of LP is upper bounded by IDS upto a constant. More formally, consider the primal-dual form of standard LP

and the corresponding KKT residual for $z=(x,y)$, $x\geq 0$ is given by $\left\|\begin{pmatrix}Ax-b\\ [A^{T}y-c]^{+}\\ [c^{T}x-b^{T}y]^{+}\end{pmatrix}\right\|_{2}$.
Then it holds for any $z=(x,y)$, $x\geq 0$ that

where $R$ is an upper bound for the norm of iterates $z_{k}$ and $\rho_{r}(z)$ is the normalized duality gap defined in [4, Equation (4a)]. The first inequality follows from [4, Proposition 5] and the second inequality uses [4, Lemma 4]. Therefore, IDS under $l_{2}$ norm is an upper bound of KKT residual for LP. Since all norms in finite-dimensional Hilbert space are equivalent (upto a constant), this showcases that the IDS under $P_{s}$ norm also provides an upper bound on KKT residual.

Since $\mathcal{F}(z)$ is in the dual space, we utilized the dual norm of $P_{s}$ in the definition of IDS. In fact, for the simple bilinear problem in Figure 1, the use of $P_{s}^{-1}$ norm is the key to guarantee the monotonic decay of IDS for PDHG iterates. Indeed, this monotonicity result holds for any convex-concave problems and can be formalized in the below proposition:

We end this section by making two remarks on the proof of Proposition 1:

In this section, we present the sublinear convergence rate of IDS. Our major result is stated in Theorem 1. As a direct consequence, the result implies sublinear rate of non-ergodic iterates on duality gap.

As a direct consequence of Theorem 1, we can obtain $\mathcal{O}(\frac{1}{\sqrt{k}})$ convergence rate on the primal-dual gap for the last iteration of PDHG. The last iteration of PDHG has a slower convergence rate compared to the average iteration, which has $\mathcal{O}(\frac{1}{k})$ rate on the primal-dual gap [11]. This is consistent with the recent discovery that the last iteration has a slower convergence than the average iteration in a smooth convex-concave saddle point for the extragradient method [30], Douglas-Rachford splitting [19, 18], and primal-dual coordinate descent method [27].

In previous sections, we show that the IDS monotonically decays along iterates of PDHG and has a sublinear convergence rate for a convex-concave minimax problem. In this section, we further show that the IDS has linear convergence if the minimax problem further satisfies a regularity condition. This regularity condition is satisfied by many applications, such as unconstrained bilinear problem, linear programming, strongly-convex-strongly-convex problems, etc.

We begin by introducing the following regularity condition.

Recently, the metric sub-regularity and related sharpness conditions have been used to analyze primal-dual algorithms [26, 48, 1, 4, 60, 62]. In particular, it has been shown that strongly-convex-strongly-concave problems satisfy this condition globally, and piecewise linear quadratic functions satisfy this condition on a compact set [65, 42]. This includes many examples, such as LASSO, support vector machine, linear programming and quadratic programming on a compact region in the primal-dual formulation. Suppose the step-size satisfies $s\leq\frac{1}{2\|A\|}$, one can show that unconstrained bilinear problem and generic linear program satisfy the regularity condition (9). Below we list a few concrete examples that satisfy metric subregularity:

The proofs of these examples follow directly from [4, 26].

The next theorem presents the linear convergence of IDS for PDHG under metric sub-regularity:

In this section, we demonstrate that the derived sublinear and linear rates of the last iteration of PDHG in previous sections are tight, i.e., there exists an instance on which the obtained convergence rate cannot be improved. The tightness herein does not refer to a lower bound result, i.e., this argument does not rule out primal-dual first-order algorithms that can achieve even better performance in the decay of IDS.