---
title: 'PDCS: Primal-Dual Conic Solver'
url: https://www.emergentmind.com/topics/primal-dual-conic-programming-solver-pdcs
type: topic
---

# PDCS: Primal-Dual Conic Solver

Primal-Dual Conic Programming Solver (PDCS) is a first-order, matrix-free solver for large-scale conic optimization problems. It was introduced as a primal-dual large-scale conic programming solver with GPU enhancements and is built around the restarted primal-dual hybrid gradient method, augmented by adaptive reflected Halpern restarts, adaptive step-size selection, adaptive weight adjustment, diagonal rescaling, and bijection-based projection onto rescaled cones [2505.00311]. A subsequent technical note documents the implementation and use of PDCS, including its CPU and GPU backends, solver options, and interfaces with JuMP and CVXPY [2603.15504].

## 1. Mathematical formulation and conic scope

PDCS solves conic programs of the form
\[
\min \ \langle c,x\rangle \quad \text{s.t. } Gx-h \in K,\; l \le x \le u,\; x \in K_p,
\]
where \(G \in \mathbb{R}^{m\times n}\), \(c \in \mathbb{R}^n\), \(h \in \mathbb{R}^m\), \(l,u\) are componentwise lower and upper bounds on a subset of variables, \(K\) is a primal cone, and \(K_p\) is the cone on the primal variable block [2603.15504]. In the 2025 presentation, the same solver family is written in a closely related form,
\[
\min_{x=(x_1,x_2)} \langle c,x\rangle \quad \text{s.t.}\quad Gx-h\in \mathcal K_d^*,\; l\le x_1\le u,\; x_2\in \mathcal K_p,
\]
with the feasible region defined by a product of cones and bound constraints [2505.00311].

The cone \(K\) is assumed to be a Cartesian product of basic cone blocks. The technical note states that PDCS supports blocks including the zero cone \(\{0\}\), the nonnegative orthant \(\mathbb{R}_+^d\), the second-order cone, the exponential cone \(K_{\exp}\), the rotated second-order cone, and the dual exponential cone [2603.15504]. The main paper emphasizes the currently supported problem classes: linear programs, second-order cone programs, convex quadratic programs via SOCP reformulations, and exponential cone programs [2505.00311].

This formulation places PDCS in the general conic-programming setting rather than in the narrower LP-only or QP-only regime. The main paper explicitly treats LPs, SOCPs, convex QPs, and exponential-cone models as target classes, which covers applications such as Fisher market equilibrium problems, Lasso regression, and multi-period portfolio optimization [2505.00311].

## 2. Saddle structure and the PDHG kernel

PDCS is organized around the primal-dual saddle formulation
\[
\min_{x\in X}\max_{y\in Y} L(x,y) := \langle c,x\rangle - \langle y, Gx-h\rangle,
\]
with
\[
X := [l,u]\times K_p,\qquad Y:=K_d,
\]
so that a saddle point \((x,y)\) corresponds to an optimal primal solution \(x\) and an optimal dual solution \((y,c-G^\top y)\) [2603.15504]. The main paper presents the same convex-concave structure as
\[
\min_{x\in \mathcal X}\max_{y\in \mathcal Y}\; \mathcal L(x,y) :=\langle c,x\rangle - \langle y,Gx-h\rangle,
\]
with \(\mathcal X=[l,u]\times \mathcal K_p\) and \(\mathcal Y=\mathcal K_d\) [2505.00311].

The core iteration is a PDHG step. With primal step size \(\tau\) and dual step size \(\sigma\), PDCS computes
\[
x^+ = \operatorname{Proj}_{[l,u]\times K_p}\left(x - \tau(c - G^\top y)\right),
\]
\[
y^+ = \operatorname{Proj}_{K_d}\left(y + \sigma\left(h - G(2x^+ - x)\right)\right).
\]
This is the extrapolated PDHG or Chambolle-Pock update used as the solver’s base mechanism [2603.15504]. The main paper gives the same update in the notation
\[
\hat{x}=\operatorname{Proj}_{[l,u]\times\mathcal K_p}\!\bigl(x-\tau(c-G^\top y)\bigr), \qquad
\hat{y}=\operatorname{Proj}_{\mathcal K_d}\!\bigl(y+\sigma(h-G(2\hat{x}-x))\bigr),
\]
and identifies sparse matrix-vector multiplication with \(G\) and \(G^\top\) and projection onto cones as the essential per-iteration primitives [2505.00311].

PDCS monitors primal feasibility, dual feasibility, and primal-dual gap through absolute and relative error quantities. The technical note defines absolute errors for the primal constraint residual, the dual residual, and the primal-dual gap, and then relative errors normalized by \(1+\|h\|\), \(1+\|c\|\), and \(1+|\langle c,x\rangle|+|\langle y,h\rangle + l^\top\lambda_1-u^\top\lambda_1|\), respectively [2603.15504]. The main paper summarizes the same monitoring strategy through
\[
\text{maxErr}(z)=\max\{err_{\mathrm{rel,p}},err_{\mathrm{rel,d}},err_{\mathrm{rel,gap}}\},
\]
which functions as the solver’s aggregate optimality indicator [2505.00311].

## 3. Adaptive acceleration, restart, and primal-dual balancing

The distinctive features of PDCS lie in the mechanisms wrapped around the base PDHG iteration. The technical note lists five enhancements: adaptive step-size selection, adaptive reflected Halpern iteration, adaptive restart, primal weight updates, and diagonal preconditioning or rescaling [2603.15504]. The main paper groups the same ideas under adaptive step-size selection, adaptive reflected Halpern restarts, adaptive restart, and adaptive primal weight adjustment [2505.00311].

For step-size selection, PDCS parameterizes the primal and dual steps by a primal weight \(\omega\) and a scalar \(\eta\),
\[
\tau = \frac{\eta}{\omega}, \qquad \sigma = \eta \omega,
\]
and adapts \(\eta\) by a line-search heuristic inspired by PDLP [2603.15504]. This removes the need to rely exclusively on conservative worst-case choices such as \(\tau\sigma \le 1/\|G\|^2\), while retaining the PDHG structure [2603.15504].

After the raw PDHG update, PDCS applies a reflected Halpern iteration. The main paper writes the update as
\[
z^{t,k+1} = \frac{k+1}{k+2}\Bigl((1+\beta^{t,k})\hat z^{t,k}-\beta^{t,k} z^{t,k}\Bigr) +\frac{1}{k+2}z^{t,0},
\]
anchoring the current inner-loop dynamics at the restart point \(z^{t,0}\) while also extrapolating through the reflection term [2505.00311]. The reflection coefficient is selected adaptively as
\[
\operatorname{Proj}_{[0,1]}\!\left(-0.1\log_{10}(\text{maxErr}) + 0.2\right),
\]
so that the method is more conservative when residuals are large and more aggressive as residuals shrink [2505.00311].

Adaptive restart is based on the normalized duality gap
\[
\rho(r;z) := \frac{1}{r}\sup_{z'\in B(r;z)} \left[L(x',y)-L(x,y')\right],
\]
or, in practice, a diagonal-surrogate approximation \(\rho_N\) computed efficiently by bisection [2603.15504]. Restart is triggered by sufficient decay, necessary decay combined with lack of local progress, or a hard cap on the inner-loop length [2603.15504]. For cones such as exponential and second-order cones, the note states that normalized duality-gap computation can become numerically delicate; in that case PDCS switches to a fallback criterion based on a weighted KKT error,
\[
\mathrm{KKT}_\omega^2 = \omega^2 \mathrm{err}_{\mathrm{abs},p}^2 + \omega^2 \mathrm{err}_{\mathrm{abs},d}^2 + \mathrm{err}_{\mathrm{abs},gap}^2
\]
in the technical-note notation [2603.15504].

The primal weight \(\omega\) is updated only at restarts. It balances the relative scale of primal and dual motion, and the technical note reports instability when \(\omega>10^5\) or \(\omega<10^{-5}\), in which case the solver resets it to its initial value [2603.15504]. In the main paper, the initialization rule is
\[
\omega = \frac{\|c\|_2}{\|h\|_2}
\]
when both norms exceed \(10^{-10}\), and \(\omega=1\) otherwise [2505.00311].

## 4. Rescaled cone projections, software realization, and GPU implementation

A major implementation layer in PDCS is diagonal rescaling. The technical note writes the rescaled matrix as
\[
G = D_1 G D_2,
\]
with positive diagonal matrices \(D_1,D_2\), and transforms the conic program accordingly [2603.15504]. The main paper gives the corresponding transformed data
\[
\hat G = D_1 G D_2,\qquad \hat c = cD_2,\qquad \hat h=D_1 h,\qquad \hat l=D_2^{-1}l,\qquad \hat u=D_2^{-1}u,
\]
together with transformed cones \(\hat{\mathcal K}_{d^*}=D_1\mathcal K_{d^*}\) and \(\hat{\mathcal K}_p=D_2^{-1}\mathcal K_p\) [2505.00311]. The stated purpose is conditioning improvement without loss of sparsity.

Rescaling complicates projection because transformed cones are no longer always standard. PDCS therefore uses bijection-based root-finding methods to compute projections onto rescaled cones [2505.00311]. For the rescaled second-order cone, projection is reduced to a scalar root-finding problem. For the rescaled exponential cone, the main paper uses Moreau decomposition and parameterizes the projection by a scalar \(\rho\), again leading to a bijection search [2505.00311]. The technical note states that, for general cones, projections onto transformed cones can still be reduced to efficient root-finding procedures and identifies this point as especially important for second-order and exponential cones [2603.15504].

The GPU implementation, cuPDCS, is designed to keep as much computation on the GPU as possible and to reduce CPU-GPU transfers [2505.00311]. The main paper describes three levels of parallelization for cone projection: grid-wise, block-wise, and thread-wise. It further states that zero, nonnegative, and exponential cone projections are cheap and work well with thread-wise parallelism, whereas SOCP projection may use thread-wise, block-wise, or grid-wise strategies depending on workload size [2505.00311]. This architectural specialization is directly tied to the solver’s stated emphasis on sparse matrix-vector multiplication as the core computational operation.

PDCS is publicly available at `https://github.com/ZikaiXiong/PDCS` under the Apache License 2.0 [2603.15504]. The technical note documents three interfaces—JuMP, CVXPY, and a direct Julia API—and recommends the GPU backend whenever possible [2603.15504]. It also states that the CVXPY interface currently supports the GPU backend only [2603.15504].

## 5. Empirical behavior and application domains

The main paper reports experiments on CBLIB benchmark problems, Fisher market equilibrium problems, Lasso problems, and multi-period portfolio optimization problems [2505.00311]. Its headline empirical claim is that cuPDCS is generally more efficient than state-of-the-art commercial solvers and other first-order methods on large-scale conic program applications, and that it exhibits better scalability, efficiency, and robustness than other first-order methods on CBLIB, with the advantages more pronounced in large-scale, lower-accuracy settings [2505.00311].

The performance profile is not uniform across scales. On the mostly small-scale CBLIB benchmark, barrier or interior-point methods are usually much faster, but cuPDCS is reported to be more stable than other first-order methods and often competitive or better on medium and large instances [2505.00311]. In Fisher market equilibrium, an exponential-cone application, cuPDCS beats SCS(GPU), is often competitive with or better than COPT and MOSEK on large instances, becomes particularly strong when the number of nonzeros grows beyond roughly \(2\times 10^7\) to \(4\times 10^7\), and is especially good at low accuracy [2505.00311].

For Lasso reformulated as SOCP, the main paper states that cuPDCS is the best performer across all tested solvers, including commercial and GPU barrier solvers, and remains effective up to \(m=n=7.5\times 10^5\) [2505.00311]. It also notes that runtime increases only modestly when tightening tolerance from \(10^{-3}\) to \(10^{-6}\), and that nearly linear convergence is observed empirically on this class [2505.00311]. In multi-period portfolio optimization, cuPDCS consistently outperforms SCS(indirect) and ABIP, is more scalable than interior-point methods on large instances, and continues to run effectively when MOSEK* becomes memory-limited for large horizons [2505.00311].

The technical note frames these results in operational terms. It identifies PDCS as especially attractive for large-scale conic optimization problems with millions of variables and constraints, where interior-point methods are typically too memory-intensive, and recommends keeping preconditioning, adaptive step-size selection, and adaptive restart enabled [2603.15504]. A plausible implication is that PDCS is intended primarily for throughput-oriented conic workloads rather than for the fastest possible high-accuracy termination on small to medium instances.

## 6. Position within the conic-solver landscape

PDCS belongs to the first-order primal-dual branch of conic optimization. Its immediate algorithmic ancestry is PDHG, but the technical note explicitly positions it in the same family as PDLP, cuPDLP, cuPDLP-C, HPR-LP, and rHPDHG, while extending the restarted-PDHG idea from LPs to more general conic settings [2603.15504]. This distinguishes it from homogeneous self-dual embedding solvers such as the ADMM-based method of O’Donoghue, Chu, Parikh, and Boyd, which solves a cone-subspace feasibility problem and is designed to return either primal-dual solutions or certificates of infeasibility or unboundedness [1312.3039].

Within first-order conic methods, PDCS is related in spirit to projection-based primal-dual schemes such as PIPG and xPIPG. PIPG is a first-order primal-dual method for conic optimization that guarantees \(O(1/k)\) convergence of both primal-dual gap and constraint violation in the general convex case, improves the primal-dual gap to \(O(1/k^2)\) under strong convexity, and improves the constraint violation to \(O(1/k^3)\) [2108.10260]. xPIPG extends this line by adding extrapolation and automatic infeasibility detection, with iterates that either satisfy primal-dual optimality conditions asymptotically or generate proofs of primal or dual infeasibility [2203.04188]. PDCS differs by emphasizing restarted PDHG, adaptive reflection, diagonal rescaling, and GPU-oriented cone-projection kernels rather than the specific update and certificate machinery of PIPG-type methods.

The principal contrast is with interior-point methods. Recent work on asymmetric long-step primal-dual IPMs with dual centering targets conic optimization problems for which the dual formulation is simpler than the primal one and retains best-known polynomial-time complexity bounds while using only Cholesky factorization [2503.10155]. Separate 2025 results sharpen the complexity analysis of symmetric primal-dual interior-point algorithms by proving exact and asymptotically tight bounds on a key scaling measure, which in turn implies the standard \(O(\vartheta^{1/2}\ln(1/\epsilon))\) iteration complexity for broad classes of primal-dual IPMs [2509.10263]. Against that background, PDCS occupies the matrix-free, sparse-linear-algebra end of the solver spectrum. This suggests a division of labor: interior-point solvers remain the reference approach for rapidly reaching high accuracy, whereas PDCS is designed for very large sparse conic programs in which sparse matrix-vector multiplication, cone projection, restart heuristics, and GPU execution dominate the practical cost model.

Source: https://www.emergentmind.com/topics/primal-dual-conic-programming-solver-pdcs