Differentiation Through Black-Box Quadratic Programming Solvers
Authors: Connor W. Magoon
Fengyu Yang
Noam Aigerman
Shahar Z. Kovalsky
Venue: NeurIPS 2025 poster ()
Keywords: optimization
differentiable optimization
quadratic programming
Primary area: optimization
TL;DR: We introduce dQP, a solver-agnostic, modular framework for differentiable quadratic programming (QP), enabling plug-and-play differentiation with any solver and demonstrating advantages in large-scale sparse problems.

Abstract: Differentiable optimization has attracted significant research interest, particularly for quadratic programming (QP). Existing approaches for differentiating the solution of a QP with respect to its defining parameters often rely on specific integrated solvers. This integration limits their applicability, including their use in neural network architectures and bi-level optimization tasks, restricting users to a narrow selection of solver choices. To address this limitation, we introduce **dQP**, a modular and solver-agnostic framework for plug-and-play differentiation of virtually any QP solver. A key insight we leverage to achieve modularity is that, once the active set of inequality constraints is known, both the solution and its derivative can be expressed using simplified linear systems that share the same matrix. This formulation fully decouples the computation of the QP solution from its differentiation. Building on this result, we provide a minimal-overhead, open-source implementation (<https://github.com/cwmagoon/dQP>) that seamlessly integrates with over 15 state-of-the-art solvers. Comprehensive benchmark experiments demonstrate dQP’s robustness and scalability, particularly highlighting its advantages in large-scale sparse problems.

Differentiation Through Black-Box Quadratic Programming Solvers
Connor W. Magoon"| Fengyu Yang" Noam Aigerman', Shahar Z. Kovalsky' Department of Mathematics, University of North Carolina at Chapel Hill,
2 Université de Montréal and Mila Equal contribution
Abstract
Differentiable optimization has attracted significant research interest, particularly for quadratic programming (QP). Existing approaches for differentiating the solu- tion of QP with respect its defining parameters often rely on specific integrated solvers. This integration limits their applicability, including their use in neural network architectures and bi-level optimization tasks, restricting users selection of sol choices. To address this limitation, we introduce dQP, a modular and solver-agnostic framework for plug-and-play differentiation of virtually any QP solver. A key insight we leverage achieve modularity is that, once the active set of inequality constraints is known, both the solution and its derivative can be expressed using simplified linear systems that share the same matrix. This formu- lation fully decouples the computation of the QP solution from its differentiation. Building this result, we provide tation (https //github. com/ cwmagoon/ (QP) that seamlessly integrates with over 15 State-of-the-ee solvers. Comprehensive benchmark experiments demon- strate dQP's robustness and scalability, particularly highlighting its advantages large-scale sparse problems.
narrow
-overal-overhedd, open-source implemen-
1 Introduction
Computational methods often rely on solving optimization problems, i.e., finding optimum of an objective function subject to constraints. This has led to the development of numerous high- performance open-source and commercial solvers, particularly for constrained convex optimization [14 27). Recently, there has been growing interest in differentiable optimization, which enables gradient-based learning through optimization layers that show promise across applications such as
image classification [8]. optimal transport [97 (88, zero-sum games (75). tessellation [35), control
9 40), decision- makingg [110), robotics (62). biology (116), and NLP [III].
Differentiable optimization focuses on computing the derivatives of the optimal solution to an optimization problem with respect to parameters defining the problem's 's objective and constraints.
Rooted classical sensitivity analysis and parametric programming 999 100 49 73 26 93),
differentiable optimization has gained momentum through the idea that optimization problems encoding prior domain knowledge can be embedded as parameterized "layers" within neural networks training such layers requires computing the derivative (gradient) of the optimization problem's solution during backpropagation [
This paper focuses on differentiating solutions to Quadratic Programs (QPs) core class of convex problems involving the minimization of quadratic objective subject linear inequality constraints.
Despite significant research efforts, existing approaches for differentiating QPs fail to fully leverage state-of-the -ortt--rtt--ttttt solvers. General- peneral-purpose purpose methods that support broad class of optimization problems
39th Conference on Neural Information Processing Systems (NcurIPS 2025).

---

(a)
(b)
10a
10
OptNet (qpth)
OPI ayer (ProxQP)
SCOPTH (OSQP))
10a
RPQP (OSQP)
Layes (SCS)
Duality Gap
Time (8)
Durs (Gurobi)
10
10a
10
10a
10
10
10
10
10
Problem Size
Problem Size
Figure 1: Comparison of differentiable QP methods for projection onto the probability simplex, evaluated by (a) total solve and differentiation time and (b) solution aceuracy (duality gap). For moderately sized problems, our approach, using the Gurobi QP solver, outperforms existing methods on both metrics. As problem size increases, our method remains efficient, while others become intractable.
often include QPs, but they sacrifice efficiency and robustmess for generality, underperforming compared to QP-specific methods. On the other hand, recent QP-specific approaches are tightly coupled with particular or proprietary solvers (e.g., for batching small QPs), limiting their flexibility 15,17
This tight integration restricts broader applicability. Solving QPs efficiently and reliably often requires advanced solvers like Gurobi [58] and MOSEK [10], developed over years of academic and commercial effort. These solvers handle scales and complexities that non-industrial implementations cannot. However, no single solver is optimal for all problems, making solver selection critical. To our knowledge, no existing approach provides robust, efficient, and fully flexible framework and implementation for differentiating QPs across solvers.
We introduce dQP, modular, colver-agnostic framework that transforms any QP solver into differentiable layer. We strengthen modularity by introducing an explicit characterization QP's primal-dual solution in terms of its active constraint set. Our key insight is that once the active inequalities are known, both the solution and its derivative can be computed from simplified linear systems that share the from differentiation but also enables differentiation of solvers that output only primal solutions, by completing the corresponding dual variables negligible computational cost.
a
matrix. This formulation not only fully decouples optimization
implement dQP as an open-source, hinimal--vveheall layer on top of qpsolvers [3), enabling seamless integration with over 15 commercial and free solvers, and straightforward extension to additional ones. Our modular approach allows users to select the best solver for their task while pro- viding differentiabilltyy. We evaluate dQP benchmark of over 2,000 diverse QPs and demonstrate superior performance compared to existing methods, particularly in large-scale structured sparse problems, demonstrated in Figure T
Our contributions:
1. We design and implement modular differentiable layer compatible with any QP solver. Our open-source code is available at /https://githul: com/ cwmagoon/ dQP
We demonstrate cate-of-thhe-at performance in solving and differentiating large-scale, sparse QPs.
We introduce an explicit, simplified characterization QP's solution and gradient in terms of its active constraints, which underpins our modular framework.
2 Related Works
Differentiable Optimization. OptNet 5
differentiates QPs through their optimality conditions, and focuses on small dense problems for GPU batching. They solve the full Jacobian system efficiently by reusing the factorization employed in their custom interior-point method. However, noted in

---

[7]] this comes at the cost of ill ill-conditioning due to symmetrization. More recent differentiable QP methods include Dit- and SCQPTH [009] 30 29 which use first-order ADMM and approximately differentiate the fixed point map, and QPLayer 17 which accommodates infeasibility with extended conservative Jacobians. Similarly to OptNet, several of these are tightly integrated with custom algorithms, often to enable access to internal computations required for differentiation. Alt-Diff is coupled with a custom ADMM method, SCQPTH reimplements OSQP, and QPLayer is built on ProxQP. Several works have noticed the importance of the active constraint set in differentiating constrained optimization problems S 55 89], as well as in the context of quadratic programming
86 83 A common observation is that the algebraic system obtained through implicit differentiation can be simplified by removing rows corresponding to inactive constraints. Some works have additionally observed that backpropagation itself can be cast as distinct consttii-- QP parameterized by the incoming gradients g [6]. However, existing approaches have not used these observations create differentiable layer that supports arbitrary black-box QP solvers, missing the opportunity to fully decouple optimization and differentiation. Other classes of optimization problems, such convex cone programs 2 and mixed-integer programs [89], have also been differentiated.
19.
Some frameworks 21,24, 92 of optimization problems, but their support of QP has significant limitations. CVXPYLayers (I] reformulates the QP into a cone program to use diffep internally (2] and, as a result, does not support specialized QP solvers, relying exclusively on the cone solvers SCS, ECOS, and Clarabel. The framework Theseus [92] directly handles only unconstrained problems and similarly lacks support for QP-specific solvers. TorchOpt (96) and Optax (39] support unconstrained optimizers used in meta-learning, but require user-defined optimality conditions and user-supplied solvers to handle problems such as QPs. JAXopt (24] includes an implicit differentiation wrapper for CVXPY, which requires symbolic compilation of the QP and does not support sparse matrices; for sparse QPs, JAXopt includes a differentiable re-implementation of OSQP (similar to SCQPTH). Altogether, these drawbacks lead to subpar performance which is reported for some generic methods in previous work on differentiable QPs (e.g., (17). For completeness, Table in Appendix A summarizes 24 relevant methods, including those discussed above.
39
96 22, 90: [022 provide differentiable interfaces to broader classes
Implicit Layers. Optimization layers are a class of implicit layers that use implicit differentiation to compute gradients of solution mappings that lack flosed-fomm expressions [43). Another class of implicit layer are deep equilibrium models, which are defined by fixed-point mappings and can be
interpreted as infinitely deep networks [14 66 44 59
113
15 Similar techniques extend beyond
algebraic equations to neural ODEs, where the adjoint state method from parametric PDE control is applied [76] [15
18 34). Implicit differentiation also plays key role in bi-level programming 3 and meta-learning, where it enables optimization of the outer learning loop [50] 61 95 [01. Alternative methods bypass implicit differentiation using approximations for example,
36
72 [66
63
by applying automatic differentiation to iterative algorithms via loop unrolling [19 20 8T 105),
or by differentiating single iteration or using Jacobian-fiee techniques for fixed- point mappings
53 52 25
Sensitivity Analysis and Parametric Programming. There is extensive mathematical theory on the local behavior of optimization problems under perturbations [99 100], particularly regarding solution sensitivity and stability 47 49 48
26 While the implicit function theorem is central to this analysis, applying it to optimization problems requires an intermediate step: reformulating the problem through its optimality conditions, which demands several regularity assumptions. Sensitivity theory supports applications such as multi-parametric programming [93, including model predictive control, where problems must often be solved repeatedly for many parameter values, increasing computational costs. To address this, [2 observed that QPs admit closed-form solutions when the active set is known in advance, enabling offline precomputation. The parameter space can be partitioned into regions with fixed active sets [[06), where the solution is stable under perturbations.
73
This idea continues to inform modern methods 45
[2]
3 Approach
We begin by formulating the problem and establishing the basic theory QP differentiation. Then, we connect this foundation to our central theoretical observation and conclude with our straightforward algorithm derived from it. We note that various subparts of our discussion have been used to develop methods for differentiating QPs, see Section However, no single work has fully developed the
3

---

explicit theory we present or leveraged it to practically implement fully modular differentiable QP layer.
3.1 Problem Setup: Differentiating QPs
We consider a quadratic program in standard form,
z (0) arg min
2² p²"P(0)= (())²
(1)
A(0)= 1(0) s.t.
C(0) & d(0),
where P E R ,9 E RRx,,,bb R2, C R=* and d E m are smoothly parameterized by some 0 E RR, We assume that the QP is feasible and strictly convex (i.e., P > 0). To simplify notation, in the following we omit 0.
This work focuses on computing (g)" = 36 (0), the derivative of the optimal point of the QP ① with respect to the parameters 9. Intuitively, this derivative quantifies the change in the optimal point of the QP in response to change in its parameters f. Our goal is to compute 0g) efficiently and independently of the method used to approximate the optimal point 2 2²()).
An important use case, highlighted recent work on differentiable optimization, involves incorporat- ing a QP of the form ① as the E-th layer of neural network. In this setting, the input to the QP layer, Te, serves as the parameter vector 0 that defines the QP's objective and constraints, and the output the optimal solution ø4+1 2* (xz). Generally, training such a network requires backpropagating gradients through each layer, which involves computing the Jacobian 8x. For QP layer, this
80441
Jacobian corresponds to (0)) the derivative of the optimal point with respect to the problem parameters. This derivative is also crucial in descent-baas based methods for solving bi-level optimization problems [36]
3.2 KKT Conditions and Sensitivity Analysis
Our goal is to differentiate QPs using only the solution returned by black-box numerical solver. To do so, we first establish the necessary theoretical foundations. As is standard in optimization, the required derivatives
obtained via sensitivity analysis of the KKT conditions. This section distills key ideas from optimization theory, sensitivity and parametric analysis, and differentiable programming, framing them in the context of QPs to support the development of dQP.
Optimality Conditions. The first-order Karush- Kuhn-Tucker (KKT) conditions 66 provide useful algebraic characterization of the optimal points of constrained optimization problems. For the QP ( the KKT conditions take the form,
14
Pz + + A + CT 1t¹ =0 Ar* - 0, Cz" - -d 0, 4¹¹ 2 0
D( (C2 - d) 0.
where D(a") = (p"), and the additional variables X E IRP and n² E IR" are the optimal dual variables of the linear equalities and inequalities, respectively. The primal-dual solution of the QP (I)
is defined by (²) = (=" (0), A (0), μd² (0)). Under strict convexity and feasibility, the QP (T has unique solution (#)) and the KKT conditions (2) are necessary and sufficient for its optimality.
Active Set and Complementary Slackness. The last equation in (2). the nonlinear complementary slackness condition, plays central role in work. Intuitively, it encodes the two possible states of each inequality constraint in (T). (Oz* - d)j I O. Fither (i) the constraint is active, meaning it holds with equality (Cz" d), 0, in which case n² 2 0; (ii) it is inactive, satisfied with strict inequality, which case 14 = 0. Notably, a constraint inactive, the same optimal solution would be obtained even if that constraint were removed. We denote the active by
(() = [j (C(0)=*( (0) - -d(0); 0}.
Derivatives via Sensitivity Analysis. differentiate QPs, we apply the Basic Sensitivity Theorem (Theorem 2.1 in (46)), which underpins differentiation of the KKT conditions with respect to 0.

---

Differentiability at requires the additional assumption of strict complementary slackness, ruling out the degenerate case where both (C2* - -d, and μ³, 0, thus ensuring that the active set remains unchanged under small perturbations of 0. Under this assumption, the primal-dual point (0) = (=" (0), X (0), μr" (0)) is differentiable in neighborhood of optimal for the QP [ uniquely satisfies the KKT conditions ( and maintains strict complementary slackness. Crucially, the active set (() remains fixed in this neighborhood.
With the active set stable, the equality conditions in locally characterize ( Implicit differentia- tion of these yields the Jacobians of the solution Ogp in terms of the linear system,
P
A-
C
00221 og7
TopP** doq + 8 + OOC² n²) 0g Oob D(n")(0002 0yd)
0
D(
D(C²*
d)
Bop
Under the assumptions of the Basic Sensitivity Theorem, the linear system 3 is invertible. It
degenerates exactly in the presence of weakly active constraints nę = (Cz* - d), = O, for which the QP is non-differentiable (see, e.g., (5)) For any inactive constraint j & J, the dual variable 23 vanishes, and thus the corresponding rows and columns of (3 can be removed, yiclding the simplified reduced form,
P
AT
C
Ogz
OgP2* + doq + 00 + 0,CT ogAz* 0gb o,C ayd,
A
0
0
(4)
C₁
0
0
8epu
where ,(,, CJ and d, denote restriction to rows corresponding to active inequality constraints j E J.
3.3 Extracting Derivatives from a QP Solver's Solution
With the above theory, derive main theoretical results and introduce dQP, a straightforward algorithm for efficient and robust differentiation of any black-box QP solver.
Our approach stems from two straightforward yet powerful insights: (1) given the primal solution of a QP, the active set can be easily identified; (2) once the active set known, both the primal- dual optimal point and its derivatives can be derived explicitly in closed-form. Furthermore, these quantities can computed efficiently via single matrix factorization of reduced-dimension symmetric system.
These observations lead to simple algorithm: first, solve the optimization problem using any QP solver; then, iden- tify the active set from the solution and solve linear system to compute the derivatives. Consequently, we can define "backward pass" for any layer that uses a QP solver, allowing for the seamless integration of any solver best suited to the problem, thus leveraging years of search and development invested in state-of-the-art QP solvers.
Figure 2: Illustration of active set dif- ferentiation. Left: a QP is shown by its quadratic level sets and polyhedral fea- sible set; the solution lies on a facet of the boundary; perturbations of the con- straints lead to perturbations in the so- lution. Right: the perturbation of the solution remains the same when inactive constraints are eliminated.
Locally Equivalent Equality-Constraiie QP. Consider the QP (T] and its optimal point (0). along with the set (0) of active constraints, see Section 3.2 We define the reduced cqustity-onsttrained quadratic program, obtained by removing inactive inequalities and converting active inequality constraints into equality constraints,
z" (0) = arg min
T T +
(5)
A(0) C(0))(() = d(0))(() (())
b(0)
s.t.
Under the assumptions of Section 3.2 this simpler QP is, in fact, locally equivalent to the QP (1. illustrated in Figure Moreover, it provides an explicit expression for both the primal-dual optimal point and its derivatives:
5

---

Theorem 3.1. The QP (5) is locally equivalent to the reduced equality-constrainee QP
and its
solution (0) (()) (8)) admits the explicit form
P
A
A [
A
(6)
CJ
0
dj
Furthermore, the optimal point can be explicitly differentiated to obtain
Ogz
P
AT
CJ
-ooq Ogb [ad
0g OpA OBCr
Ag
CgC
2
Oa
0
0
A
(7)
-
0op
CJ
0 0
0
0
proof of this Theorem, based on the Basic Sensitivity Theorem 46
provided in Appendix B
along with calculation of the derivatives using differential matrix calculus 991
78 We note that this
result closely related to analyses studied in multi cuulti-parametiic programming 21 93 106
12
Notably, for quadratic programming, the Basic Sensitivity Theorem allows us bypass implicit differentiation techniques [70]. We emphasize that this does not imply that the general solution or its active set admit a closed-form expression. Rather, while we perform explicit differentiation, the implicit function theorem remains key in establishing the local equivalence between the original and reduced problems - the resulting derivatives are the same, though the derivations differ. The derivatives 7 are obtained via ordinary (explicit) differentiation of the closed-fon solution to the reduced OP (6), whereas those in (4) arise from implicit differentiation of the full nonlinear KKT conditions 2 followed by pruning inactive constraints. This distinction, formalized in Theorem 3.11 underscores a critical computational insight: once a black-box solver returns the primal solution, the active set can be identified and the derivatives computed directly via (7) Moreover, if the solver returns only the primal solution, the corresponding dual variables can be recovered using (6). Since both (6) and (7) rely on the same KKT matrix K,,, single matrix factorization (e.g.. via SuperL.U (74)) suffices for both steps, adding negligible overhead. These insights culminate in the core algorithm of dQP, summarized in Algorithm !
Explicit Differentiation.
Numerical Computation. Our approach enables compact and efficient gradient computation. The linear system in (7). used compute both derivatives and dual variables, symmetric and of reduced size. In contrast, implicit differentiation of the full KKT conditions (I) yields signif- icantly larger, asymmetric system B). Beyond simplifying the derivative computation, our ap- proach enables the use of fast, specialized linear solvers that exploit the reduced systems sym- metric indefinite KKT matrix structure (e.g., g., using an LDL factorization as in QDLDL [08 37)
Empirically, we observe that the reduced linear system is often significantly better condi- tioned than its full counterpart. Figure 3 illus- trates this with an example of a QP governed
105
full
reduced
Condition
by two parameters 0 (0,, (0),0)) of (106), cal-
03
culated using DAQP [33. Eliminating inactive constraints improves conditioning significantly.
Figure 3 also highlights key challenge: com- puting derivatives near non-differentiable singu- larities where the active set changes and some constraints become weakly active. In such re- gions, implicit differentiation becomes severely ill-conditioned. Our method is also affected, primarily through difficulty in determining the correct active set an approximate solution. Several strategies have been proposed for improving active set identification (3 84 28). Our implementation includes an optional refinement heuristic (Appendix D), though in all our experiments (Section 4), simple hard thresholding of the primal
01
01
Figure 3: (a) Active-set parameter space, coloring regions where the active set is constant. (b) Condi- tion number of the full and reduced linear systems along a line in parameter space.
residual rj = (C2² d), 2 -E proved robust.
Implementation. Our open-source implementation is available at https://gi thub. com/ cwmagoon/ dQp We implement dQP (Algorithm() [ as a fully differentiable PyTorch module [88]. providing intuitive interface for integrating differentiable QPs into machine learning and bi-level

---

Algorithm - dQP: Differentiation through Black-box Quadratic Programming Solvers
Input: p.g, P,q.A,b, A, C, d, and tolerance CJ
Output: 2", _==,"" jr" and 0øz", OgX", Ōoje"
1: Solve QP ① with any solver for the primal solution 2" (and X", , if available)
2: Compute the active set by thresholding:
J = (j (Cz"
- -d); 2 -EJ
P A²
C3
3: Factorize the reduced KKT system matrix:
K,
0
O 0
4: Compute X", u=" (if not obtained in step (1)):
X
KJ
di
MJ
0a="
-99 agb Ogd,
8g P
8g
8,0
5: Compute the derivatives:
0.
KJ²
8g A
0
X
86p
CbCs
0
Table 1: Performance of differentiable QP methods for 129 sparse problems from the Maros-Meszaros (MM) dataset 779..
Solver
Full
Subset Problezs Solved Metlods
Probs
Arg Irdd [=xx]
Arg Bldd [cxs]
Ang Total [mx]]
/og BeeS/Tbral [94
Accuracy [duality sap]
e Pribs Sulved
Ang Fird
Arg Bidd [(s]]
Aug Toxal
Avg Boud/toel
Accuracy
Salved
[es]
[ea]]
(%)
[duslity sap)
IQP ((Phradtmma)))
996
1467
57%
1.75
OPLayer (PruOPP)
15889
19721
2.21
2838
433
2910
1:27
OpeNeer Ciget)
39329
619
9:36
(Dsgge))
10
optimization workflows. The implementation supports both dense and sparse problems end-to-end, with appropriate QP and linear solver support. As a PyTorch module, ( must be rendered backward pass; this is described in Appendix C To ensure modularity, the forward pass supports any QP solver interfaced via the open-source gpsolvers library (32), which provides lightweight access to over 15 commercial and free open-source solvers and easily supports the integration of new ones. We likewise offer flexibility in the choice of linear solver for differentiation, including support for large-scale sparse solvers like Pardiso [103] and symmetric indefinite solvers like QDLDL [088] 37). For users unsure of which solver to use, dQP includes a profiling tool to help identify the best-performing QP solver for given problem. Additional implementation details such as constraint normalization, handling non non-differentiabiiity. and warm-starting options for bi-level optimization are discussed in Appendix D
4 Experimental Results
We have extensively tested dQP to ensure its robustness and evaluate its performance against existing methods for differentiable quadratie programming. Notably, we demonstrate dQP's strengths in handling large-scale structured and sparse problems, thus complementing custom differentiable GPU- batched solvers such as OptNet [5), which are optimized for solving many small, dense problems simultaneously. Given this focus, and considering the limited availability of state-of-the-e GPU- batchable QP solvers, we conduct our experiments on CPUs, similar to prior works [17] 29 109 I. Our evaluation includes large benchmark of over 2,000 challenging sparse and dense QPs taken from public and randomly generated problem datasets, designed to test dQPs robustness and performance. Additional experimental details provided in Appendix
Modularity and Performance. We tested dQP on the QP Benchmark suite 31, focusing on 129 large-scale sparse problems from the standard Maros-Meszaros (MM) dataset (79]. which are widely used as stress tests for QP solvers. We compared dQP against other differentiable QP methods available in the PyTorch framework: OptNet (5). QPLayer [17]. SCQPTH (29). and CVXPYLayers Ii using the authors' open-source implementations, each paired with its respective solver. We did not include other general frameworks that either lack direct support for generic QPs in PyTorch, showed subpar performance our preliminary tests (consistent with findings reported in prior work (17), or rely on forward passes built diffcp/CVXPPLLayerr (see Table 4.. For dQP's forward pass, we selected the top performfgnnng QP solver for each problem, as identified by QP Benchmark [3I] to demonstrate the importance of enabling differentiation for the solver that is best suited to each task.
The scatter plot in Figure 4 shows total runtime (forward + backward), duality gap (accuracy), and problem dimension (illustrated by point size) for each differentiable solver and problem. Aggregate
7

---

10a
Duality Gap
OpeNee (ipeh)
30%
QPLayer (Prox)P)
60%
10
e SCCOPT (OSQe)
43%
CVXyYLaees 14%
Ours (QPRanchark)) 100%
10.
10a
10
10' 10'
Time (s)
Figure 4: Accuracy versus total forward/backward solve for the Maros-Meszarge dataset 79]. Each point represents solved problem; point size illustrates dimension; problems solved solely by dQP are circled. The legend shows percentages of success rates; the solvers dQP used and their counts are PIQP 100, Gurobi 9, ProxQP 9, Clarabel 7, OSQP 2, MOSEK 1, QPALM 1.
performance across the full dataset and the subset of problems solvable by all methods is summarized T CVXPYLayers failed on all but small-scale problems and is therefore excluded from the reported statistics. The MM dataset proved especially challenging, with OptNet and SCQPTH contrast, dQP solved all MM problems and was the only method succeed on 38 of them (circled in the figure). Moreover, dQP achieved the best total forward/backward runtime and accuracy in 80% and 78% of all problems, respectively. performed particularly well on large-scale instances (dimension over 1000), being fastest and most accurate in 97% and 93% of such cases. Additional experiments on 450 random dense QPs and 625 sparse QPs (dimensions 10 to 10+) are provided in Appendix E.1.3 showing similarly strong performance.
solving fewer than 50% of the problems.
Scalability. We evaluated dQP on large-scale sparse problems that are highly structured. a regime where state-of-the state-of-the-art QP solvers have significant advantage over less optimized solvers. Specifically, we tested dQP and other available differentiable QP solvers on two prototypical projection layers expressed as constrained QPs:
(z arg min ||. 21? s.t. 05
(8)
and
Tnn) arg min > ||, =, s.t. || - 2j+1 1l S 1.
Results for P1 are shown in Figure ! demonstrating dQP's scalability compared to OptNet and QPLayer. Other methods fail on all but small problems (see Appendix [E.I.]]. In dimensions greater than 2000, dQP outperforms competing methods by 2-3 orders of magnitude in both speed and accuracy. Competing methods are limited to dense calculations and fail in dimensions beyond 10a. We note that P1 is the projection onto the probability simplex, also known as SparseMAX, for which more efficient, non-Qu non-QQ-baseed methods exist [80). Results for P2, representing projection onto "chains" with bounded links, exhibit similar scalability and are detailed in Appendix E.1.2
Table 2: Performance analysis for projection onto the probability simplex, formulated in Addi- tional details are provided in Appendix F.I.T and Table
Problem
450
4600
1000)
aQP (Gurcbii)
Time [mns]
Accuracy
1.07 x 10
1.63
8.88 10-10
1.92
2,28 x 10
3.13
1.47 x 10
5.00
2.72 x 10
16.46
9.55 x 10
49.00
6.67 10
476 64
OptNer (cpch)
Accuracy
4.04 10 2.92
10-#
x 10
10-a 70856.c3
Time [ms]
5.19
37.06
182.99
$514.8
QPLayer (ProxQP)
Time [ms]]
Accurcy
9.53 x 16
0.29
3.65 x 1
1661
4.16 x 10
77.56
2.19 x 10
751.14
1.16 10
79514.19
1.94 10-3
946174.88
BPQP (OSQP)
Time [ms]
Accurccy
2.04
0.32
5.64
0.76
9.66
6.81
3.04 x
21.08
1.76 10
358.78
10 1407.05
CVXPYLayers (SCS)
Accunccy
1.31 10
9.47 x 10
1.51 10-3
2.70 x 10-
Time [ms]
1.97
13.22
156.37
662.60
SCOPTH (0SQP*)
Accuracy
Time [ms]
6.19 x -
14.75
8

---

Table 3: Performance analysis for projection onto chains, formulated in 9 Additional details are provided in Appendix E.112 and Table
Metric
Problem Size
200
500
1000
2000
4000
10000
100000
(Gurobi)
Time [mx]
Accuracy
2.73 10-
6.15
2.02 10
3.79 x 10-"
26.19
9.16 x 10-6
47.94
2.64 x 10"5
88.35
4.29 10
224.89
2.81 10
2432.64
OptNet (upsh)
Time (ms)
Accuracy
6.97 10-
26.38
1.75 10
169.64
9.22 x 10
901.25
2.43 x 10-7
5481.56
2.00 10-7
34799.98
1.98 10-
571710.06
QPlayer (ProxQP)
Accuracy Time [ms]
10-5 8:04
78 10
81.58
577.11
10-
2.97 10-
3996 92
30748.07
10-4
471649.91
10-1
SCQPTH (OSQP")
Accunccy
(ms]
1.67 10-5
13,20
2 83 10"
67.55
4.76 10-5
407.13
6.64 x 10-5
2755.28
7,80 10-5
10028.11 15
1.21 10
195462.18
CVXPYLayers (SCS)
Time [ms)
Accuracy
129.26
10-1
2. 18 10-1
683.34
2.97 x 10-1
2048.34
Bi-Level Geometry Optimization. We further demonstrate the scalability of our approach on bi-level optimization problem introduced in 669
M arg min 1
(10)
M
s.t. (M") (M)) arg arg (tr MD) s.t. Bu u, CMu = 0}
In this problem E RRnx2 represents the n vertex coordinates of triangular mesh, M E R parameterized Laplacian, and B, u and C' encode boundary conditions. (²MM) is the dual variable corresponding to linear inequalities of the lower-level problem. The results described in [69] imply that if μ ((M) vanishes then v² (M*). at the optimal Laplacian M represents a straight-edge intersection-fre -free drawing of the mesh [112] The solution to (IO) is visualized in Figure S (]a) for an example large-scale ant mesh. Additionally, Figure (b shows that dQP scales more favorably with problem size compared to OptNet, QPLayer and SCQPTH: in particular, only dQP supports problems with over 10d vertices.
(b)
(opth) syer
SCOPTH (OSQ!")
(PIQP)
Mcan Time
(Rack)
10a
10a
10a Vertices
10
Figure 5: (a) Planar embedding of a large-scale ant mesh (15k vertices). Zoom-ins show: (top) a non-injective harmonic map with edge overlaps highlighted in red; (bottom) dQP's injective solution to the bi-level problem (IOp. (b) Scalability: solver runtimes as mesh size increases for a synthetic problem.
Conclusion
We introduce the dQP framework for differentiating QPs by leveraging the local equivalence between QP and simpler equality-constrained problem. dQP provides a straightforward differentiable interface for any QP solver, yielding an efficient QP layer that can be readily integrated into neural architectures, among other applications. We note that our current method does not yet support full parallelization or GPU acceleration, as state-of oftte-of-the-ant the-art sparse and scalable QP solvers with GPU support are still lacking. We recognize these as important challenges and limitations to address in future work. We see the present work as an important step toward developing similar solver-agnostic differentiable layers for other popular optimization problems (e.g., semidefinite programming), which plan explore in future research.

---

Acknowledgments
The authors gratefully acknowledge partial support from NSF grant DMS-2152289, FRG: Collab- orative Research: Mathematical and Statistical Analysis of Compressible Data on Compressive Networks; the NSERC Discovery Grant 2GGIN-2202-04600 Practical Neural Geometry Processing; and the FRQNT Etablissement de la relève professorale grant 365040, Calcul rapide léger des déformations à I' l'ddde de réseaux neuronaux.
References
[1] Agrawal, A., Amos, B., Barratt, S., Boyd. S., Diamond, S., and Kolter, J. Z. Differentiable convex optimization layers. In Advances Neural Information Processing Systems, 2019.
Agrawal, A., Barratt, S., Boyd, S., Busseti, E., and Moursi, W. M. Differentiating through cone program. Journal of Applied and Numerical Optimization, 1(2):107-115, 2019.
[3] Alesiani, F. Implicit bilevel optimization: Differentiating through bilevel optimization pro- gramming. In Proceedings of the AAAI Conference on Artificial Intelligence, volume 37, pp. 14683-14699 2023.
[4] Amos, B. Differentiable optimization- modeling for machine learning. Ph. Thesis, 2019.
[5] Amos, B. and Kolter, J. Optnet:
Differentiable optimization as a layer in neural networks.
In International conference on machine learning, pp. 136-145. PMLR, 2017.
[6] Amos, B. and Kolter, J. Z. cvtpy httpe://gghabbcccmm////sss/////p master qpth/ http///pthh//oeee///// solvers/ cvxpy. py 2017.
[7] Amos, B. and Kolter, J. Z. gpth-sparse. hasterr/ptth/aollvrr///dd/p
.yy 2017.
[8] Amos, B., Xu, I.,, and Kolter, J. 7.. Input convex neural networks. In Precup, D. and Teh, Y. W. (eds.), Proceedings of the 34th International Conference on Machine Learning, volume 70 of Proceedings of Machine Learning Research, pp. 146-155. PMLR, 06-11 Aug 2017.
[9] Amos, B., Jimenez, I., Sacks, J., Boots, B., and Kolter, J. Z. Differentiable mpe for end-to-end planning and control. In Bengio, S., Wallach, H., Larochelle, H.. Grauman, K., Cesa-Bianchi, N., and Garnett, R. (eds.), Advances in Neural Information Processing Systems, volume 31. Curran Associates, Inc., 2018.
[10] Andersen, E. D. and Andersen, K. D. The mosek interior point optimizer for linear program- ming: an implementation of the homogeneous algorithm. In High performance optimization, pp. 197-232. Springer, 2000.
[11] Andrychowicz, M., Denil, M., Gomez, S., Hoffman, M. W., Pfau, D., Schaul, T., Shillingford, B., and De Freitas, N. Learning to learn by gradient descent by gradient descent. Advances in neural information processing systems, 29, 2016.
[12] Arstrröm D. and Axehill, D. A high-performant multi solver, 2024.
multi- parametric quadratic programming
[13] Arnström, D., Bemporad, A., and Axchill, D. A dual active-set solver for embedded quadratic programming using recursive LDL" updates. IEEE Transactions on Automatic Control, 67(8): 4362-4369, 2022.
[14] Bai, S., Kolter, J. Z., and Koltun, V. Deep equilibrium models. In Wallach, H., Larochelle, H., Beygelzimer, A., d'Alché-Buc, F., Fox, E., and Garnett, (eds.), Advances in Neural Information Processing Systems, volume 32. Curran Associates, Inc., 2019.
[15] Bai, S., Koltun, V., and Kolter, J. Z. Multiscale deep equilibrium models. In Larochelle, H., Ranzato, M., Hadsell, R., Balcan, M., and Lin, H. (eds.), Advances Neural Information Processing Systems, volume 33, pp. 5238-5250 Curran Associates, Inc., 2020.
10

---

[16] Bambade, A., Schramm, F., El Kazdadi, S., Caron, S., Taylor, A. B., and Carpentier, J. Proxqp: an efficient and versatile quadratic programming solver for real-time roboties applications and beyond. 2023.
[17] Bambade, A., Schramm, F., Taylor, A. B., and Carpentier, J. Leveraging augmented--agrangian techniques for differentiating over infeasible quadratic programs in machine learning. In The Twelfth International Conference on Learning Representations, ICLR 2024, Vienna, Austria, May 7-11, 2024. OpenReview.net, 2024.
[18] Beatson, A., Ash, J. T., Roeder, G., Xue, T., and Adams, R. P. Learning composable energy surrogates for pde order reduction. In Proceedings of the 34th International Conference Neural Information Processing Systems, NIPS '20, Red Hook, NY, USA, 2020. Curran Associates Inc. ISBN 9781713829546.
[19] Belanger, D. and McCallum, A. Structured prediction energy networks. In International Conference on Machine Learning, pp. 983-992. PMLR, 2016.
[20] Belanger, D., Yang, B., and McCallum, A. End-to-end learning for structured prediction energy networks. In International Conference Machine Learning, pp. 429-439. PMLR, 2017.
[21] Bemporad, A., Morari, M.. Dua, V., and Pistikopoulos. E. N. The explicit linear quadratic regulator for constrained systems. Automatica, 38(1):3-20, 2002. ISSN 0005-1098.
[22] Besançon, M., Dias Garcia, J., Legat, B., and Sharma, A. Flexible differenntabbee optimization via model transformations. INFORMS Journal Computingg, 36(2):456 2024.
[23] Blondel, M. and Roulet, V. The elements of differentiable programming. arXiv preprint arXiv:2403. 14606, 2024.
[24] Blondel, M., Berthet, Q., Cuturi, M., Frostig, R., Hoyer, S., Llinares-ópezz F., Pedregosa, F., and Vert, J.-P. Efficient and modular implicit differentiation. Advances in neural information processing systems, 35:5230-5244 2022.
[25] Bolte, J., Pauwels, E., and Vaiter, S. One-step differentiation of iterative algorithms. In Teirty-seventh Conference on Neural Information Processing Systems, 2023.
[26] Bonnans, J. F. and Shapiro, A. Perturbation analysis of optimization problems. Springer Science & Business Media, 2013.
[27] Boyd, S. and Vandenberghe, L.. Convex optimization. Cambridge university press, 2004.
[28) Burke, J. V. and Moré, J. On the identification of active constraints. SIAM Journal Numerical Analysis, 25 1988. ISSN 00361429.
[29] Butler, A. Seqpth: an efficient differentiable splitting method for convex quadratic program- ming. 08 2023.
[30] Butler, A. and Kwon, R. H. Efficient differentiable quadratic programming layers: an admm approach. Computational Optimization and Applications, 84(2):449 476, 2023.
[31] Caron, S., Zaki, A., Otta, P., Amström, D., Carpentier, J., Yang, F., and Leziart, P.-
gpbenchmark: Benchmark for quadratic programming solvers available in Python. https://github. htthsshhhh..cmmcoeeeee hasn/qqpolverse wwrpep/wllelesooooooom cctn://cppo...sssppppperecrrrrrr 2024.
[32] Caron, S. et al. QPSOLVERS: Quadratic Programming Solvers in Python. https://githe
com/ http///llerr//psolllece/ March 2024.
[33] Cartis, C. and Yan, Y. Active-set prediction for interior point methods using controlled perturbations. Computational Optimization and Applications, 63(3))633 684, Apr 2016. ISSN 1573-2894.

---

[34] Chen, R. T. Q.. Rubanova, Y., Bettencourt, J., and Duvenaud, D. K. Neural ordinary differential equations. In Bengio, S., Wallach, H., Larochelle, H.. Grauman, K., Cesa-Bianchi, N., and Garnett, R. (eds.), Advances in Neural Information Processing Systems, volume 31. Curran Associates, Inc., 2018.
[35] Chen, R. T. Q. Amos, B., and Nickel, M. Semi-discrete normalizing flows through differ- entiable voronoi tessellation. In ICLR Workshop on Deep Generative Models for Highly Structured Data, 2022.
[36] Colson, B., Marcotte, P., and Savard, G. An overview of bilevel optimization. Annals of operations research, 153:235-256, 2007.
[37] Davis, T. A. Algorithm 849: concise sparse cholesky factorization package. ACM Transac- tions Mathematical Software (TOMS), (144)5587-591 2005.
[38] de Avila Belbute-Peres, F., Smith, K., Allen, K., Tenenbaum, J., and Kolter, J. Z. End-to-end differentiable physics for learning and control. In Bengio, S., Wallach, H., Larochelle, H., Grauman, K., Cesa-Bianchi, N., and Garnett, R. (eds.), Advances in Neural Information Processing Systems, volume 31. Curran Associates, Inc., 2018.
[39] DeepMind, Babuschkin, I., Baumli, K., Bell, A., Bhupatiraju, S., Bruce, J., Buchlovsky, P., Budden, D., Cai, T., Clark, A., Danihelka, I., Dedieu, A., Fantacci, C., Godwin, J., Jones, C., Hemsley, R., Hennigan, T., Hessel, M., Hou, S., Kapturowski, S., Keck, T., Kemaev, I., King. M., Kunesch, M., Martens, L., Merzic, H., Mikulik, V., Norman, T., Papamakarios, G., Quan, J., Ring, R., Ruiz, F., Sanchez, A., Sartran, L., Schneider, R., Sezener, E., Spencer, S., Srinivasan, S., Stanojević, M., Stokowiee, W., Wang, L.,, Zhou, G., and Viola, F. The DeepMind JAX Ecosystem, 2020. URL /
com/ woon://oglledeeee
[40] Ding. S., Wang, J., Du, Y., and Shi, Y. Reduced policy optimization for continuous control with hard constraints. In Proceedings of the 37th International Conference on Neural Information Processing Systems, NIPS '23, Red Hook, NY, USA, 2024. Curran Associates Inc.
[41] Du, X., Aigerman, N., Zhou, Q., Kovalsky, S. Z., Yan, Y., Kaufman, D. M., and Ju, T. Lifting simplices to find injectivity. ACM Transactions on Graphics, 39(4). 2020.
[42] Dunning, I., Huchette, J., and Lubin, M. Jump: A modeling language for mathematical optimization. SIAM Review, 59(2):955-320, 2017. doi: 10.1117/15M1020555
[43] Duvenaud, D., Kolter, J. Z., and Johnson, M. Deep implicit layers autorial-neury odes, deep equilibirum models, and beyond. Neural Information Processing Systems Tutorial, 2020.
[44] El Ghaoui, L., Gu, F., Travacca, B., Askari. A., and Tsai, A. Implicit deep learning. SIAM Journal on Mathematics of Data Science, 3(3):930-955) 958, 2021.
[45] Ferreau, H. J., Kirches, C., Potschka, A., Bock, H. G., and Diehl, M. qpoases: A parametric active-set algorithm for quadratic programming. Mathematical Programming Computation, 327-363, 2014.
[46] Fiacco, A. Sensitivity analysis for nonlinear programming using penalty methods. Mathe- mutical Programming, 10(1):287-311, Dec 1976. ISSN 1436-4646.
[47] Fiacco, A. V. Introduction to sensitivity and stability analysis non linear programming. New York: Academic Press,, 1983.
[48] Fiacco, A. and McCormick, G. P. Nonlinear programming: sequential unconstrained minimization techniques. John Wiley & Sons, New York, NY, USA, 1968. Reprinted by SIAM Publications in 1990.
[49] Fiacco, Anthony Vand Ishizuka, Y. Sensitivity and stability analysis for nonlinear program- ming. Annals of Operations Research, 27(1):215-235, 1990.
[50] Finn, C. B. Learning to learn with gradients. University of California, Berkeley, 2018.
12

---

[51] Frison, G. and Diehl, M. Hpipm: high-performnn performance quadratic programming framework for model predictive control. IFAC-PapersOnLine, 53(2):6563 6569, 2020. ISSN 2405-8963. st IFAC World Congress.
[52] Fung, S. W., Heaton, H.. Li. Q.. McKenzie, D., Osher, S. and Yin, W. Jfb: Jacobian backpropagation for implicit networks. In AAAI Conference on Artificial Intelligence, 2021.
[53] Geng, Z., Zhang, X.-Y., Bai, S., Wang, Y., and Lin, Z. On training implicit models. In Beygelzimer, A., Dauphin, Y., Liang, P., and Vaughan, J. W. (eds.), Advances in Neural Information Processing Systems, 2021.
[54] Goldfarb, D. and Idnani, A. A numerically stable dual method for solving strictly convex quadratic programs. Mathematical Programming, 27(1):1-33, Sep 1983. ISSN 1436-4646.
[55] Goulart, P. J. and Chen, Y. Clarabel: An interior-poimt solver for conic programs with quadratic objectives, 2024.
[56] Gould, S., Fernando, B., Cherian, A., Anderson, P., Cruz, R. S., and Guo, E. On differentiating parameterized argmin and argmax problems with application to bi-level optimization. arXiv 2016.
preprint 160.v:1607.0544
[57] Gould, S., Hartley, R., and Campbell, D. Deep declarative networks. IEEE Trunsactions Puttern Analysis and Machine Intelligence, (00)3398-40004 2022.
[58] Gurobi Optimization, LLC. Gurobi Optimizer Reference Manual, 2024.
[59] Gurumurthy, S., Bai, S., Manchester, Z., and Kolter, inference and input optimization in equilibrium networks. Advances in Neural Information Processing Systems, 34:16818 16832, 2021.
[60] Hermans, B., Themelis, A., and Patrinos, P. QPALM: A Newton-type Proximal Augmented Lagrangian Method for Quadratic Programs. In 58th IEEE Conference on Decision and Control, Dec. 2019.
[61] Hochreiter, S., Younger, A. S., and Conwell, R. Learning learn using gradient descent. In Artificial Neural Networks-ICANN 2001: International Conference Vienna, Austria, August 21-25, 2001 Proceedings 11, pp. 87-94. Springer, 2001.
[62] Holmes, C., Dümbgen, F., and Barfoot, T. D. Sdprlayers: Certifiable backpropagation through polynomial optimization problems in robotics, 2024.
[63] Hospedales, T., Antoniou, A., Micaelli, P., and Storkey, A. Meta-learning in neural networks: A survey. IEEE transactions on pattern analysis and machine intelligence, 44(9):5149-510 2021.
[64] Huangfu, Q. and Hall, J. Parallelizing the dual revised simplex method. Mathematical Programming Computation, 10(1):119-144 Mar 2018. ISSN 1867-2957.
[65] Karush, W. Minima of functions of several variables with inequalities as side constraints. Master's thesis, Department of Mathematics, University of Chicago, 1939.
[66] Kawaguchi, K. On the theory of implicit deep learning: Global convergence with implicit layers. In International Conference on Learning Representations, 2021.
[67] Kingma, D. P. and Ba, J. Adam: A method for stochastic optimization, 2017.
[68] Kovalsky, S. Z.,, Glasner, D., and Basri, R. A global approach for solving edge-matching puzzles. SIAM Journal on Imaging Sciences, 8(2):16-933,, 2015.
[69] Kovalsky, S. Z., Aigerman, N., Daubechies, I., Kazhdan, M., Lu, J., and Steinerberger, S. Non-convex planar harmonic maps. arXiv preprint arXiv:2001.0133 2020.
[70] Krantz, S. G. and Parks, H. R. The Implicit Function Theorem. Birkhäuser Boston, MA, 2012.
13

---

[71] Kuhn, H. W. and Tucker, A. W. Nonlinear programming, In Neyman, J. (ed.), Proceedings of the Second Berkeley Symposium on Mathematical Statistics and Probability, pp. 481-492. University of California Press, 1951.
[72] Kunisch, K. and Pock, T. A bilevel optimization approach for parameter learning in variational models. SIAM Journal on Imaging Sciences, 6(2):938-99) 9838,983,2222 2013.
[73] Lee, G., Tam, N., and Nguyen, D. Y. Quadratic programming and affine variational inequalities, qualitative study. Springer New York, NY, 01 2010.
[74] Li, X. An overview of superlu: Algorithms, implementation, and user interface. ACM Trans. Math. Softw., 31(3):302-325 September 2005. ISSN 0098-3500.
[75] Ling, C. K., Fang, F., and Kolter, J. 7. What game are we playing? end-to-end learning in normal and extensive form games. arXiv preprint urXiv:1805.022
2018.
[76] Lions, J.-L. Optimal control of systems governed by partial differential equations. Springer- Verlag, Berlin, 1971. ISBN 3340051155. Translation of Contrôle optimal de systèmes gouvernés par des équations aux dérivées partielles.
[77] Liu, Z., Liu, L., Wang, X., and Zhao, P. LFWLayer: Differentiable frank-wolfe optimization layer. In The Second Tiny Papers Track at ICLR 2024, 2024.
[78] Magnus, J. R. and Neudecker, H. Matrix Differential Calculus. New York, 1988.
[79] Maros, I. and Mészáros, C. repository of convex quadratic programming problems. Opti- mization methods and software, 11(1-4):671-681) 1999.
[80] Martins, A. and Astudillo, R. From softmax to sparsemax: A sparse model of attention and multi-label classification. In Balcan, M. F. and Weinberger, K. Q. (eds.), Proceedings The 33rd International Conference on Machine Learning, volume 48 of Proceedings of Machine Learning Research, pp. 1614-1623, New York, New York, USA, 20-22 Jun 2016. PMLR.
[81] Metz, L., Maheswaranathan, N., Nixon, J., Freeman, D., and Sohl- Shhk-Dickstei,, J. Understanding and correcting pathologies in the training of learned optimizers. In International Conference on Machine Learning, pp. 4556-4565 PMLR, 2019.
[82] Narciso, D. A., Pappas, I., Martins, F,, and Pistikopoulos, E. N. A new solution strategy for multiparametric quadratic programming. Computers Chemical Engineering, 164:107882, 2022. ISSN 1998-1354,
[83] Niculae, V., Martins, A., Blondel, M., and Cardie, C. SparseMAP: Differentiable sparse structured inference. In Dy, J. and Krause, A. (eds.), Proceedings of the 35th International Conference on Machine Learning, volume 80 of Proceedings of Machine Learning Research, pp. 3799-3808. PMLR, 10-15 Jul 2018
[84] Oberlin, and Wright, S. J. Active set identification in nonlinear programming. SIAM Journal on Optimization, 17(2):577-29, 2006. Copyright- Copyright] 2006 Society for Industrial and Applied Mathematics; Last updated- 2023- 12-04.
[85] O' O'Donoguue, B., Chu, E., Parikh, N., and Boyd, S. Conic optimization via operator splitting and homogeneous self-dual embedding. Journal of Optimization Theory and Applications, 14463):0422-068, 2016.
[86] Pan, J., Ye, Z.,, Yang, X., Yang, X., Liu, W., Wang, L.,, and Bian, J. BPQP: A differentiable convex optimization framework for efficient end-to-end learning. In The Thirty-eighth Annual Conference on Neural Information Processing Systems, 2024. URL Pet//forumiid-VWKYaay/7v
[87] Pandala, A. G., Ding. Y., and Park, H.-W. gpswift: A real-time sparse quadratic program solver for robotic applications. IEEE Robotics and Automation Letters, 4(4):3355 2019.
[88] Paszke, A. et al. PyTorch: an imperative style, high-performance deep learning library. Curran Associates Inc., Red Hook, NY, USA, 2019.
14

---

[89] Paulus, A., Rolinek, M., Musil, V., Amos, B., and Martius, G. Comboptnet: Fit the right np-hard problem by learning integer programming constraints. In Meila, M. and Zhang, (eds.). Proceedings of the 38th International Conference on Machine Learning, volume 139 of Proceedings of Machine Learning Research, pp. 8443-844 8453. PMLR, 18-24 Jul 2021.
[90] Paulus, A., Martius, G., and Musil, LPGD: general framework for backpropagation through embedded optimization layers. In Salakhutdinov, R., Kolter, Z., Heller, K., Weller, A., Oliver, N., Scarlett, J.. and Berkenkamp, (eds.), Proceedings of the 41st International Conference on Machine Learning, volume 235 of Proceedings of Machine Learning Research, pp. 19989-40011 PMIR, 21-27 Jul 2024.
[91] Petersen, K. B. and Pedersen, M. The matrix cookbook. Technical University of Denmark, ((15):510, 2008.
[92] Pineda, L., Fan, T., Monge, M., Venkataraman, S., Sodhi, P., Chen, R. T. Q., Ortiz, J., DeTonee D., Wang, A., Anderson, S., Dong, J., Amos, B., and Mukadam, M. Theseus: A library for differentiable nonlinear optimization. In Koyejo, S., Mohamed, S., Agarwal, A., Belgrave, D., Cho, K., and Oh, (eds.), Advances in Neural Information Processing Systems, volume 35, pp. 18801-3888. 8. Curran Associates, Inc., 2022.
[93] Pistikopoulos, E. N., Diangelakis, N. A., and Oberdieck, R. Multi-parametric optimization and control. John Wiley & Sons, 2020.
[94] Pogančić, M. V., Paulus, A., Musil, V., Martius, G., and Rolinek, M. Differentiation of blackbox combinatorial solvers. In International Conference on Learning Representations, 2020.
[95] Rajeswaran, A., Finn, C., Kakade, S. M., and Levine, S. Meta-leaming with implicit gradients. Advances in neural information processing systems, 32, 2019.
[96] Ren*, J., Feng", X., Liu*, B., Pan*, X., Fu, Y., Mai, L., and Yang. Y. Torchopt: An efficient library for differentiable optimization. Journal of Machine Learning Research, 24(367):1-14 2023.
[97] Rezende, D. and Racanière, S. Implicit riemannian concave potential maps, 2021.
[98] Richter-Powell, J., Lorraine, J., and Amos, B. Input convex gradient networks. arXiv preprint arXiv:2111.12187, 2021.
[99] Rockafellar, R. T. Convex Analysis. Princeton University Press, Princeton, 1970. ISBN 9781400873173.
[100] Rockafellar, T. and Wets, Variational Analysis, volume 317. Springer Berlin, Heidelberg. 01 1998. ISBN 978-3-540-62772-2
[101] Sambharya, R., Hall, G., Amos, B., and Stellato, B. Learning to warm-start fixed-point optimization algorithms. Journal of Machine Learning Research, 25(166): 2024.
[102] Schaller, M. and Boyd, S. Code generation for solving and differentiating through convex optimization problems, 2025.
[103] Schenk, O. and Gärtner, K. Solving unsymmetric sparse systems of linear equations with pardiso. Future Generation Computer Systems, (0(3):475-4877 2004.
[104] Schwan, R., Jiang, Y., Kuhn, D., and Jones, C. N. PIQP: A proximal interior-point quadratic programming solver. In 2023 62nd IEEE Conference on Decision and Control (CDC), pp. 1088-1095 2023.
[105] Scieur, D., Gidel, G.. Bertrand, Q., and Pedregosa, F. The curse of unrolling: Rate of diffcrentiating through optimization. In Oh, A. H., Agarwal, A., Belgrave, D., and Cho, K. (eds.), Advances in Neural Information Processing Systems, 2022.
[106] Spjøtvold, J., Kerrigan, E. C., Jones, C. N., Tøndel, P., and Johansen, T. A. On the facet-to- facet property of solutions to convex parametric quadratic programs. Automatica, 42(12): 2209-2214, 14, 2006. ISSN 0005-1098.
15

---

[107] Stellato, B. and O' Donoghue, B. Osqpth: A interface to qdldl for parametric qp solvers. https://githnm com/ ooqp osqpth 2023.
[108] Stellato, B., Banjac, G., Goulart, P., Bemporad, A., and Boyd, S. OSQP: an operator splitting solver for quadratic programs, Mathematical Programming Computation, 12(4):637-672, 2020.
[109] Sun, H., Shi, Y., Wang, J., Tuan, H. D., Poor, H. V., and Tao, D. Alternating differentiation for optimization layers. arXiv preprint arXiv:2210.01800 2022.
[110] Tan, Y., Terekhov, D., and Delong, A. Learning linear programs from optimal decisions. In Larochelle, H., Ranzato, M., Hadsell, R., Balcan, M.. and Lin, H. (eds.), Advances in Neural Information Processing Systems, volume 33, pp. 19738-19749) Curran Associates, Inc., 2020.
[111] Thayaparan, M., Valentino, M., Ferreira, D., Rozanova, J., and Freitas, A. Diff-explainer: Differentiable convex optimization for explainable multi-hop inference. Transactions of the Association for Computational Linguistics, 10:1 103-1119 2022.
[112] Tutte, W. T. How to draw graph. Proceedings of the London Mathematical Society, 3(1): 743-767, 1963.
[113] Winston, E. and Kolter, J. Z. Monotone operator equilibrium networks. In Larochelle, H., Ranzato, M., Hadsell, R., Balcan, M., and Lin, H. (eds.), Advances in Neural Information Processing Systems, volume 33, pp. 10718 10718-10728. Curran Associates, Inc., 2020.
[114] Wright, S. J. Numerical optimization, 2006.
[115] Xue, T., Beatson, A., Adriaenssens, S., and Adams, Amortized finite element analysis for fast PDE-constrained optimization. In III, H. D. and Singh, A. (eds.), Proceedings of the 37th International Conference on Machine Learning, volume 119 of Proceedings of Machine Learning Research, pp. 10638-10647. PMLR, 13-18 Jul 2020.
[116] Zhang, H., Nguyen, D. H., and Tsuda, K. Differentiable optimization layers enhance gnn-based mitosis detection. Scientific Reports, 13(1):14306, 2023.
16

---

A Existing Methods for Differentiable Optimization
Table provides an overview of existing differentiable optimization methods for quadratic program- ming (QP), more general conic programs (CP) and other optimization problems.
Table 4: Comparison of differentiable optimization libraries and layers.
Name
Programs
Solvers
Peatures
Limitations
OptNet 1
spth
GPU batchable, Re-uses factorization forward backward solve.
Solver specilic. Spunity supported.
QPLayer
ProcQP
Suppurts infeasible QPs. Tighily cou- plas to Prodge.
specifiic Spunity supported.
SCQPTH 20
QP
OSQP
Deffrrntiitts ADMM updates. Sup- ports infeusibilit detection and auto- matic scaling.
Solver specificc Spansty not suppertd.
Diflep H
BCOS, Clurabel
Supports spurse programs (CP).
Supports specific cunic solvers.
CVXPYLayers
dificp
problem fomuuation
Inherits difiep limitations.
choice: Supports sparsiy.
BPQP X7
QRSOCP
osge
Computes derivatives by solving a sec- and equaltyy conetty-conttiinned QP.
Published implesenniaion is incom- plete.
Alt 00
Custom ADMM
Difieremiates ADMM updates.
Sup-
Repuriedly slow Dk]
solutions.
Difiereniates updates.
constraints only.
LPGD X0
diflip
Modifies difTep hundle degenerate derivatives.
Inherits diffep limitations.
CVXPYgan 1
QP. SOCP
DSQU
Accelarates CVXPYlayers cooppi- Salver specific. lation problem Commulation.
3AXogt-geneal
User supplied
User supplied
Defferemtiaes afhitrary implicit fune-
User manually provide KKT.
JAXogt OSQP
QP
Custom GPU OSQP
GPU baschable.
Solver specific.
CVXPP-CVXPPQP
QP
CVXPY
Flexible problem fonmulation and sulver choice.
Inherits CVXPY limiitsions, Spunity not supparted.
SDPRLayen
X
CP
CVXPY
Flexible problem fommulation and sulver choice.
Inheriss CVXPY limitations.
Torchoye
User supplied
User supplied
Difeerentiates arbittryy implicit fune- tiuss.
User must monuaily prncide
Optox 30
Specialived
Specializod
Supports specialized prbllemshatser. projections.
DdftUpt.j
Various
JuMP supported solvers
42
Supparts variety problems.
JaMP provides support only for CDSMO. OSQP, lpopt
Thescus 22
NLS
CHOLMOD, cudaL.U, BaSpoCho
Supports nos-lincar squares (NLS). GPU batchuble. Supports spinity.
hard constraints
OptNee-Spunec 2
QP
apth
Published implementation is incom- plcte.
OptNet CVXPY 0
CVXPY
Published implementation plete.
OSQPT: QP
OSQP
Published implementation incom- plate.
DFWLayer 17
QP
Frank-Walic
Automatiaally difierentiates am unnolled Frank Walle solver.
Solver specific.
ComboppNNe x9
Gurobi
Differestiation combinatorial salvers
Backprop
MIP.
Dijkstra
Differentiaioon
solvers
Remarks. This table includes all relevant differentiable optimization methods known to us at the time of writing. In some cases, code was available online without an associated publication. Among all listed methods that directly support generic QP without modification, only dQP, QPLayer, SDPRLayers, and DiffOpt.jl return both the dual optimal point of QP and its derivatives. Notably, dQP is the sole QP method that supports a wide variety of solvers with nanimal--vvrread through gpsolvers, in contrast to a significant number of methods that build on top of CVXPY/CVXPYLLayee which requires a compilation step (e.g., SDPRLayers, JAXopt-CXXPPy diffep which has limited solvers compared to ordinary, non-differentiable CVXPY (e.g., LPGD). Several methods are available only outside the Py Torch framework (e.g., JAXopt, Optax, DiffOpttjl).
OptNet-CVXPP) and/or
17

---

B Proof of Theorem 1
In this section provide proof of Theorem which restate below:
Theorem B.1. The QP (5) is locally equivalent to the reduced equality-constrained QP and its solution (0) (2" (0), (0), (0) admits the explicit form
P
AT
CJ 0
A
A
(11)
CJ
d
Furthermore, the optimal point can be explicitly differentiated to obtain
0g=
P
A
OgP
0gATT
OgC
Oog Oob Oads
Og
A
AgA
x
(12)
Oop C
0
OoC,
0
Proof. We begin by establishing that the QP ( and the equality-constrer reduced QP
(5) are
equivalent. For any 0 satisfying the assumptions of the theorem, the QP (I) has unique solution characterized by the KKT system
P(0)=* g(0) A(0) (0)) C(0)"
(0))= -0)- b(0) = 0
c(0)="(0) d(0) S 0
(13)
(0(0) 2 0
D(" (() d(0)) 0.
Complementarity implies that active constraints j E J(0) have n (0)) > and therefore must be satisfied with equality (C(0)=* d(0); O, while inactive constraints j E J(0) have 4² (0); and thus can be eliminated, without altering the solution. Therefore, the unique solution (#) = (=" (0), X (0), μ² (0)) of (13) is also the unique solution of the reduced system
P(O)=* - (0) A(O) (0) ((O)Tro)) (0) = 0
A(0)= (0) b(0) = 0
(14)
(()) (0) d(0))(p) O.
which are exactly the KKT conditions of the -quality-constrained reduced QP (5). Uniqueness of
solution then implies that [1] and (5) are pointwise equivalent at O. Beyond this, since P. 4, A, C, Basic Sensitivity Theorem [46) asserts that the primal-dual solution ç*(8 (0) for ( is differentiable function of 0 in neighborhood of 0 defined implicitly through the KKT's equality conditions, and that the active set (0) is fixed in this neighborhood.
are smoothly parameterized by
Thus, ( and (5) are locally equivalent in neighborhood of the parameter O.
144 implies that the reduced primal-dual solution (G) (0), X- (0), (J))) satisfies K (0) (f), where
P(O)
A(O)T
C
(())
K)(0)
A(O)
m)(0)
biO) dd
(15)
C(0))(0)
Under the assumptions of the theorem, the reduced KKT matrix K2( is invertible and
Ç = Kj Us,
(16)
yielding 66. Moreover, since 7(0) is constant in local neighborhood, the Basic Sensitivity Theorem establishes that (0) differentiable. Using the derivative of the matrix inverse 78 explicitly differentiate [6) to obtain
91]. we
O K (0,K) los KK (06u)) -K (O'(00)))C} x)'(0o)))
(17)
yielding
18

---

C Backpropagation
Like other differentiable QP layers implemented within automatic differentiation frameworks such as PyTorch [88) we do not directly return the derivative ōoς*. This is in part because dQP directly receives the QP parameters P. q. A, b, C, d and not 8, and in backpropogation we are not concerned with 0. This is rather accounted for in the next step outside dQP, usually by automatic differentiation. Secondly, backpropagation requires that we compute a so-called Jacobian-vector product which are products of the Jacobians with an "incoming" gradient of a quantity or loss & that depends S*, This requires less computation and does not require the formation of 3-tensor. Since Ç) = K V. a formal matrix-vecto multiplication, the secobian-vector product well-known,
Vuze= = ) Vezl,
(18)
and
vK,t=-Vizz v CJ"
(19)
with respect to Kj, u, respectively. Although backpropagation introduces transposition, the re-use of factorization from solving for the active duals is unaffected. This follows from the symmetry of the reduced KKT matrix which simplifies (18) into V, Kj' Vezl. Next, we extract the gradients with respect to the parameters by the chain rule. This amounts to tracking their position in the blocks and accounting for symmetry constraints, is helpful to write (dg, d,, dup) = -V so that we express
1 (dz T + d
vpl
Uge = d2
2
Vxl - dx='T + Xº
U&l = -dx
(20)
Vc, - dp,z² =T #āde
Va vaye=-daap
similar to OptNet [5]. We note that the gradient with respect to P is constrained to lie within the subspace of symmetric matrices. Similarly, if the matrices P, A, are sparse, then we project the
gradient to lie within the non-zero entries, which can be implemented efficiently in Equations 20
Although the above argument is for scalar loss É, the same approach is naturally adapted Ç is mapped to vector in the immediate next layer.
D Implementation Details
Choosing a solver. Since our work enables users to choose any QP solver as the front-end for their differentiable QP applications, include a simple diagnostic tool for quantitatively measuring solver performance. We present an example result in Figure 6 for the cross geometry experiment in Figure 9 finding PIQP, OSQP, and QPALM M to be the most efficient. For this reason, choose PIQP in the geometry experiments. We also include tools for checking the solution and gradient accuracy.
Tolerances. In addition to the active set tolerance eJ. QP solvers often support additional user- provided tolerances. These include the primal residual which measures violations of feasibility, the dual residual which measures violations of stationary, and for some solvers also the duality gap, which provides a direct handle on solution accuracy. We inherit the structure of qpsolvers for setting custom tolerances different QP solvers, though set heuristic default which is sufficient for many of the experiments in this work.
Convexity and Feasibility, Two key assumptions of our method are strict convexity and feasibility. However, these are often violated in practice. We include optional checks that P is symmetric positive definite. On the other hand, we do not perform any special handling for infcasibility a limitation of our method compared to, for example, QPLayer [17.
For non-differentiable problems, we solve for the derivatives in the least- squares sense, plugging the system into qpsolvers which can handle squsrss or a standard least- squares solver. We attempt to anticipate weakly active constraints which cause non non-differentiabiliy by measuring the norms the primal residual and the dual. Additionally, the reduced KKT is non-invertible if the active dual solution is not unique and so we check necessary condition: the total number of active constraints plus the number of equalities must be less than the dimension. these checks are passed, we attempt the standard lincar solve and pass to least-squares if it fails.
Non-differentiable Points.
19

---

10
abb
10
10
10-a
Time
10a
10ah
Clarabel
HIGHS
ProxQP
CVXOPT
QP Solver
Figure 6: Evaluating the best QP solver for the cross geometry problem (Figure 9 using our diagnostic tool. The solution tolerance regimes are varied between cabs 10 10
Normalization. Some problems have large variations in scale between different rows within the constraints. This influences the primal residual and thus the active which is determined by comparing with an absolute threshold tolerance. To address this issue for these problems, we include an optional differentiable normalization step the constraints before Algorithm Under this choice, the resulting relative primal residual becomes the scale-invariant distance to the constraint.
[ is carried out.
Equality Constraints. While we include equality constraints in our general formulation, they are not required.
Warm-Start. Since qpsolvers supports warm-starting, we inherit it as an option and store data in the PyTorch module from previous outer iterations, which can be used as initialization. This is useful for bi-level optimization problems where the input 0 changes little between outer iterations. We did use this feature in our bi-level optimization experiment.
Fixed Parameters. For fixed parameters, do not compute the corresponding derivative. This saves the cost of unwrapping the linear solve as in (20) and saves the memory to form the loss gradients, which are matrices for P,A,C. C.
Active Set Refinement. Inaccuracy in a solution may lead to instability the active set near weakly active constraints, degrading the gradient quality. To show this, we repeat the experiment in Figure 3 which has simple polyhedral active set parameter space. One setup where instability appears is illustrated in Figure 7 where
absolute solver tolerance eabs 10 and much tighter active tolerance e = 10 Qualitatively, the active set at each solution is severely degraded, even for points away from the boundaries where the set changes. We provide optional heuristic algorithm
(a)
1.5
(b)
1.5
%
0,
-1.5
-1.5
-1.5
Figure 7: The set-up in figure 3 with looser solver tolerance éabs 10 active tolerance e, 10 and solver PIQP. (a) The computed active set is degraded due to the inaccurate solution. (b) Our heuristic active set refinement algorithm recovers the ground truth active sets.
to address this, which recovers the desired set this problem. First, we order the constraints by increasing residual and select initial active set from the tolerance ey. Then, progressively add
20

---

constraints by checking if the residual of the system 1ostemoffor Çj decreases, and greedily accepting until adding constraints no longer improves the residual. At each step, we keep the primal solution from the forward fixed, and solve for the new active dual variables. While this algorithm works well on simple examples, more sophisticated and efficient techniques may desired for harder problems. We did not use this refinement algorithm in any of our experiments.
QP Solvers. Throughout this work, we use number of QP solvers available in qpsolvers including Clarabel [55), DAQP [3]. Gurobi [58], HiGHS (64] HPIPM [51] MOSEK [00, OSQP [08]. PIQP
104. ProxQP [6], QPALM [60] qpSWIFT
87 quadprog
54 and SCS 85
E
Experimental Details
For completeness and reproducibilityy we include additional details on experiments. We run all experiments and methods CPU, including methods that support GPU such as OptNet.
E.1 Performance Evaluation
All experiments in this section were run on a Macbook Air with Apple M2 chips, 8 cores, and 16GB RAM.
In our QP benchmark experiments, we evaluate the solution accuracy using the primal residual rp (the maximum error equality and inequality constraints), dual residual ra (the maximum error on the dual feasibility condition), and duality gap ra (the difference between primal and dual optimal values).
Tp max (I -l |C (Cz-d[+) d)+)
Td P2z + q A CT 1110c
Ty = P=" P= + q" + 6rt x d'ul
Throughout our experiments, we present results for the duality gap to indicate the solution accuracy since, for a strongly convex QP, zero duality gap rg 0 is a necessary and sufficient condition for optimality.
For the forward, we set the absolute residual tolerance to tolerance ç 10 We run each problem separately with batch size 1.
10 and the active constraint
In our benchmark, we regard a problem as successfully solved it meets the following criteria:
1. The solve time is less than practical 800s time limit.
2. The primal residual, dual residual, and duality gap are less than 1.0. This a coarse check, less stringent than the imposed tolerances.
3. The differentiation is executed, and does not lead to a fatal error (e.g. due to non-invertibility of a linear system).
Experimental results averaged over 5 independent samples.
Since SCQPTH does not support equality constraints, we convert them into a corresponding set of inequality constraints.
E.1.1 Projection onto the probability simplex
For projection onto the probability simplex, formulated in Pi, equation we set x E R with I; (O, 1) drawn randomly from standard normal distribution. The dataset, with 500 problems, has dimensions n E (10, 20, 50, 100, 220, 450, 1000, 2100, 4600, 10000, 100000). For n S 4600, each dimension contains 50 problems and 25 problems for n > 4600. Gurobi is chosen as dQP's forward solver. Figure 1 shows the median performance within the 1/4 and 3/4 quantiles for each dimension. SCQPTH failed for all problems with n > 50. The statisties Table that dQP outperforms competing methods for differentiable QP in both forward and backward times.
21

---

Table 5: Time and accuracy performance statistics for projection onto the probability simplex.
Solver
Merric
Problem Size
20
100
450
10001
4600
10000
100000
Accumeey
1:07 10 "
8.88 10
2.26 10
1.47 10
2.72 x 10
9.55 x 10
6.67 x: 10
(Chundes)
Forward [ms]
Buckward [ms]
0.24
1.38
1:68
2.66
4.37
0.60
15.83
2.58
42:21
420.91
Total [ms]]
1.92
3.13
5.06
18:46
Accuracy
4.04 x10
4.24 x 10-
1.64 10
2.07 10
3.55 10
6.08 x 10
Folled
OpcNet (gpth
Forward [ms]]
Bockwaad [ms]
2.72
0.20
4.72
0.46
33.46
3.99
165.50
17.48
7788.73
720.43
65976.15
4958, 74
Total [ms]
2.92
5.19
37.66
182.99
8514.65
70836.43
Accarsyy
0.53 10
3.65 10-3
10
10
1.16 10-3
10-5
Failed
QPLayer (ProxQP)
Forwand [ms]]
Backward [ms]
0.14
0.14
0.37
66.73
10.85
657.88
91.72
71774.25
7591.93
77831.58
Total [ass]
0.29
1.61
77.56
751.14
79314.49
946174.68
Accunsc
2.04 10
5.54 x 100
9.66 10
3,04 x: 10
1.76 10
4.02 10
Failed
(OSOP)
Forward [ms]
Bockward [ms]
0.11
0.21
0.40
0.34
3.88
1.90
14.45
7,10
226.27
192.73
712.20
Totall [ms]
21.88
358.78
Accuracy
1.31 10
9.47 x 10-3
1.51 10-1
2.70 10
Fuiled
Failed
Failed
CVXPYLayens (SCS)
Forwand [ms]
Backward [ms]
0.96
1.30
10.51
2.92
131.60
23,83
$37.26
123.25
Total [mis]
1.97
13.42
156.57
662.60
Accurney
10
Failed
Failed
Failed
Failed
Pailed
Failed
SCOPTH (OSQP))
Backvard [ms]]
[ms]
14.35
0.40
Total (ams]
14.75
Table 6: Time and accuracy performance statistics for projection onto chains.
Solver
Mouric
Problem Size
200
200
1000
2000
4000
10000
100000
Accemacy
10
2.02 10
3.79 10
2.64 10-
10
2,81 10
JQP (Gurobi)
Backward [ms]
[ms]
0.49
5.66
12.04
0.98
24.41
1.74
44.79
3.18
82.57
5.81
209.79
14.69
2263.54
172.80
Total [mns]
6.15
12.99
26.19
47.94
86. 35
224.89
2452.64
Accuracy
97 10
1.75 10-7
9.22. x 10
2.43 10
2.60 10
1,98 10
Fuiled
OptNect (upth)
Forward [ms]
156.49
845.24
32505
530702.00
Rackward [ms]
1.98
13.06
61.41
365.02
2206.20
35438.33
Total [mas]
25.38
169.64
907.25
5491.56
34709.98
571710.06
Accercy
8.46 10
9.78 x 10
1.82 x 10
2.97 10
6.55 x 10
1,03 x 10
Failed
QPLayer (ProsQP)
Farward [ms]
6.60
09:50
505.05
3684.47
20221 57
414295,22
[ms]
512.96
3533
57219.08
Total [ms)
8:04
81.93
377.11
2096.92
30748.07
471649.91
Accuracy
1.07 10
2.83 x 10-5
4.70 x 10-5
6.61 x 10
7,80 x 10-5
1.21 x 10
Failed
(OSQP")
Forward [ms]
Backoard [ms]
10.02
3.17
39.49
25,
236.01
170.88
1617.88
1126.46
8058 89
837
129085.97
65507 05
[ms]
67.55
407.13
2755.28
16628. 15
195:62.18
Accemcy
1,69 10
2.18 x 10
2.97 x 10-1
Falled
Pailed
Falled
Failed
CyxPYLayers (SCS)
Backward [ms]]
Forward [ms]
36.27
94.70
$10.04
167.21
1644.47
404.12 12
Total [mms]
683.34
2048.34
E.1.2 Projection onto chain
For projection onto chains with links of length bounded by in oc-norm, as formulated in Pe. equation ( the input point cloud 21,..., E Rd set with x; N/ 100Ia), where the number 100. By varying the dimension of the vector, d, we generated 300 problem dataset with dimensions n E (200,500, 2000,2000,4000 1,000,000,,0000,,,11 contains 50 problems; for > 4000 each dimension contains 25 problems. Gurobi is chosen as dQP's forward solver. Figure 8 and Table [Jble/sddeooostraat the solvers have performance similar to the projection onto the probability simplex shown in Figure[] [ in terms of efficiency. Additionally, dQP successfully solves large-scale problems that other solvers fail to solve.
of points
10000, 100000). For n S 4000, each dimension
E.1.3 Random sparse/dense problems
We generated two datasets of random QPs: sparse and dense.
For the random sparse dataset, P L, where L is the standard Laplacian matrix of k-nearest graph (k = 3). Entries of C and A are filled by AN(0,) 1) random numbers with density of 5 x 10-4; entirely zero rows are avoided. The vectors d and are generated so that the constraints are feasible d = C1 + Il and b- A1. The dataset contains 625 problems with dimensions spanning E {100, 220, 450, 1000, 2100, 4600, 10000), where m = n and p = n/2. For n S 4600, each dimension contains 100 problems; for n 4600 each dimension contains 25 problems. The KKT
22

---

(a)
(b)
10a
OptNer
QPLayer (Prex(P)
10a
SCOPTH (OSQP)
10a
CaypyLayers (SCS)
Duality Gap
Time (8)
Ours (Carubi)
104
10"
10
10a
10a
10'
10
10
10
10
Problem Size
Problem Size
Figure Time and accuracy performance for projection onto chains.
systems in these problems tend to be ill-conditioned. We use Gurobi as dQP's forward solver and employ least-squares solver for the backward. In our experiments, shown in Table 7 OptNet and CVXPYLayers fail on all problems and SCQPTH is substantially slower, failing for n 2 4600. dQP demonstrates superior accuracy and efficiency
QPLayer.
For the random dense dataset, P QTQ + 100+1, where Q E RO with Qy ~ U(0, 1). C' Rmx" with Cij (000, 1), and A R with Ag ~ U(0, 1). The constraint vectors are chosen in the same way as the sparse case. The dataset contains 450 problems with dimensions spanning
(10, 20, 50, 100, 220, 450, 1000, 2100, 4600). dimension contains 50 problems. We use DAQP and ProxQP as dQP's forward solvers. In the regime of these dense problems, qpth and ProxQP the forward solvers used by OptNet and QPLayer, respectively-are highly competitive with other available solvers, both free and commercial. As such, do expect dQP to yield substantial gains. Indeed, Table T shows that our method is comparable to OptNet and QPLayer in both runtime and accuracy. For smaller dimensions (n < 1000), DAQP offers higher accuracy, while for larger problems, ProxQP proves more efficient.
n, and p n/2. Each
Table 7: Time and accuracy performance statistics random sparse problems.
Solver
Metrc
Problem Size
220
450
1000
2100
10000
Accuracy
4.46
9.23 10
1.34 10
6.89 10
1.34 10
3.16 x 10
3.43
dQP (Gurobi)
Furward [ms]
Backward [ms]
2.57
1.79
2.86
3.44
5.53
4.73
11.07
9.70
24.03
60.65
2446.70
309.10
143209 89
9364.61
Total [ms]
4,37
6,33
10.28
20.72
90.01
2760.07
151471.27
OptNet (qpti)
Failed
Acceracy
6.46 10-4
1.25 x 10-5
1.69 x 10-6
3.04 10
6.12 10-5
1.77 x 10-3
7.82 10-3
OPLayer (ProxQP)
Forward [ms]
1.04
5:47
$1.11
235.00
2268 24
23507 22
199009 91
Rackward [ms]
1.17
[ms]]
27133.19
Accency
2.21 10
3.96 x 10
1.64 x 10
2,86 10
L.32 10
Failed
Failed
SCQPTH (OSQP")
Backward [ms]
Farward (ms]
4.48
1.37
8.36
4.92
25.000
27.64
4153.39
177.8
700227 33
1558 04
[mn)]
81211.50
CVXPYLayess (SCS)
Pailed
Failed
E.2 Bi-Level Geometry Optimization
The geometry experiments cores.
Intel(R) Core(TIM) 17-8850H CPU @ 2.60GHz with 6
We include supplementary example in F which illustrates the boundary conditions (inequality constraints in (10) of (68) using blue arrows (corresponding to applying the Laplacian to output vertices) and cones at points where the shape is locally non-convex. If the arrows fall outside the cones as shown in (b), then the map is non-invertible; instead they lie strictly inside the cones then the map is invertible shown in (c). The process is also flexible: by adding regularizing term M M bl that measures the distance to the combinatorial Laplacian M, shown in (d),
All TM IMallF
23

---

Table 8: Time and accuracy performance statistics on random dense problems.
Solver
Metric
Problem Size
100
2100
4600
Accuracy
1.59
1.20 10
2.35 10
4.08 x 10
5.26 x 10
Failed
4QP (daqp)
Forward [ms]
0.20
1.31
131.35
1115.62
10065,77
Backward [ms]
0.14
0.48
11.22
$6,90
313.90
Total [ms]
144.91
10379.68
Accuracy
4.71 10-"
6.42 10"5
9.11 10
7.26 x 10"1
1.3 x 10
4.25 10-1
Forward [ms]
2.54
2583.82
26408.12
JQP (prexqp)
Backward [ms]
0.17
1.85
13.53
70.25
385.04
3569.77
Total [ms]
0.46
4.32
73.68
455.22
2938.93
29771.33
Accurccy
6.80 x 10-*
2.51 10
3.80 10
3.51 x 10-7
2 2.80 x: 10-6
3.34 x 10-5
Forward [ms]
7,09
78.56
463.04
3176.59
29387.34
OptNet (qpth)
Backward [ms]
0.23
0.45
5.55
29.30
185.80
1540.07
Total [ms]
3.22
7.56
84,20
491.57
3362.25
30931,00
Accuracy
3.08 10
6.88 10
3.98 10
1.31 x 10
1,35 x 10
1 1.48 10
Forward [ms]
0.14
0.99
43.11
407.77
3973.89
43740.91
QPLayer (ProxQP)
Backward [mx]
0.15
0.34
9.67
74.24
601.25
5781.58
Total [ms]
0.29
1.35
52.99
482.17
4575.13
49558.44
Accuracy
3.48 x 10
4.62 10
4.32 x 10
8.54 x 10
1.83 x 10
2.26 x 10-3
SCOPTH (OSQP")
Forward [ms]
10.01
26.72
120.12
664.74
68022.36
384565.40 40
Backward [ms]
0.47
1.28
26.80
184.22
1733.72
15203.05
Total [ms]
10.50
27.90
147.25
850.37
8550.04
399699.87
Accuracy
5.40 x
Failed
Failed
Failed
Failed
Failed
CVXPYLayers (SCS)
Forward [ms]
3.15
Backward [ms]
1.08
Total [ms]
4.22
we can enhance map quality in the sense that the triangles change their shape less with respect to the input mesh (a).
The upper-level loss for the bi-level cross experiment is shown in Figure [0 where the unregularized loss is driven the desired tolerance, accompanied by sudden changes in the active set as the dual variables are driven to zero. We terminate the optimization at convergence, once all constraints are inactive to guarantee injective map. For the regularized optimization (Figure TO())), choose the regularization hyperparameter to be À = 10 after sample testing. Despite the regularization, the dual loss can still be driven 0 to reach injective map.
To optimize over Laplacians M, we directly parameterize the space Laplacians; we impose that the diagonals are the absolute row sums during optimization and that the off-diagonals are negative. We also constrain M have the same sparsity pattern as the combinatorial Laplacian. Lastly, since the Laplacian M is the quadratic term in I the resulting QP is not strictly convex because the Laplacian does not have strictly positive eigenvalues. To address this, we perturb M by a small scaling of the identity 10 10
Throughout the geometry experiments, we use the same solution tolerance Cats 10-5 and active tolerance 10-4 with the forward solver PIQP determined by our tool for choosing the best solver described Appendix[D [ For the upper-level optimization, we use the Adam optimizer with We initialize the bi-level optimization with M..
learning rate 10-2 [67]
We report only forward (QP) time for OptNet, QPLayer and SCQPTH in Figure 51b) because OptNet and SCQPTH do not output, nor differentiate the duals, and while QPLayer does, it suffers poor scaling from dense operations. The backward timing that we report for dQP excludes any contribution coming from the set-up of the parameterized Laplacian and directly report the time to solve the reduced KKT and extract the gradients with respect to M.
aa
Figure 9: Computing invertible map into non-convex plus sign. A naive mapping from (a) the square into (b) the plus without cone constraints has inversions (red). Including cone constraints (blue) and performing the bi-level optimization leads (c) an invertible map which can be regularized to (d) an enhanced quality map.
24

---

(a)
(b)
10
10
10
10
Loss
10
Loss
10
10a
10a
10a
10a
104
104
0
20
40
60
80
100
120
0
50
100
150 200 250 300 350 Iteration
Iteration
Figure 10: The evolution of the loss for the mappings the square into the cross. Iterations for which the active set changes denoted with vertical red lines. (a) Without regularization, the loss is driven monotonically to the tolerance. (b) With a competing regularizing loss term (dashed), convergence to the tolerance is slowed but not prevented.
The (Figure and (Figure ((a) meshes and boundary constraints are obtained from the datasets in 41. We create the synthetic mesh refinement example in Figure 5 b) by perturbing the corner of square mesh.
25

---

NeurIPS Paper Checklist
1. Claims
Question: Do the main claims made in the abstract and introduction accurately reflect the paper's contributions and scope?
Answer: [Yes]
Justification: See Section 4
Guidelines:
The answer NA means that the abstract and introduction do not include the claims made in the paper.
The abstract and/or introduction should clearly state the claims made, including the contributions made in the paper and important assumptions and limitations. A No or
NA answer to this question will not perceived well by the reviewers.
The claims made should match theoretical and experimental results, and reflect how much the results
be expected to generalize other settings.
is fine to include aspirational goals as motivation as long is clear that these goals are not attained by the paper.
2. Limitations
Question: Does the paper discuss the limitations of the work performed by the authors?
Answer: [Yes]
Justification:
We discuss limitations our work in Section
Guidelines:
The answer NA means that the paper has no limitation while the answer No means that the paper has limitations, but those not discussed in the paper.
The authors are encouraged to create separate "Limitations" section in their paper.
The paper should point out any strong assumptions and how robust the results are to violations of these assumptions (e.g., independence assumptions, noiseless settings, model well-specification. asymptotic approximations only holding locally). The authors should reflect how these assumptions might be violated in practice and what the implications would be.
The authors should reflect on the scope of the claims made, e.g., if the approach was only tested on few datasets or with few runs. In general, empirical results often depend on implicit assumptions, which should be articulated.
The authors should reflect the factors that influence the performance of the approach. For example, facial recognition algorithm may perform poorly when image resolution low images are taken in low lighting. Or speech-to-text system might be used reliably to provide closed captions for online lectures because it fails to handle technical jargon.
The authors should discuss the computational efficiency of the proposed algorithms and how they scale with dataset size.
If applicable, the authors should discuss possible limitations of their approach address problems of privacy and faimess.
While the authors might fear that complete honesty about limitations might be used by reviewers as grounds for rejection, a worse outcome might be that reviewers discover limitations that aren't 't acknowledged in the paper. The authors should use their best judgment and recognize that individual actions in favor of transparency play an impor- tant role in developing norms that preserve the integrity of the community. Reviewers will be specifically instructed to not penalize honesty concerning limitations.
3. Theory assumptions and proofs
Question:
For each theoretical result, does the paper provide the full set of assumptions and complete (and correct) proof?
Answer: [Yes]
26

---

The proof of our primary theoretical result, Theorem 3.1 is provided in Appendix [
Justification:
Guidelines:
The answer NA means that the paper does not include theoretical results.
All the theorems, formulas, and proofs in the paper should be numbered and cross- referenced.
All assumptions should be clearly stated referenced in the statement of any theorems.
The proofs can either appear in the main paper or the supplemental material, but if they appear in the supplemental material, the authors are encouraged to provide short proof sketch provide intuition.
Inversely, any informal proof provided in the core of the paper should be complemented by formal proof's provided in appendix supplemental material.
Theorems and Lemmas that the proof relies upon should be properly referenced.
4. Experimental result reproducibility
Does the paper fully disclose all the information needed to reproduce the main ex- perimental results the paper to the extent that affects the main claims and/or conclusions of the paper (regardless of whether the code and data are provided or not)?
Question:
Answer: [Yes]
Justification: Algorithm and implementation details are provided in Section 3 and exper- imental details are given in Appendix Our open-source implementation will be pub- licly released. An anonymized version is available for review at:
com/
Guidelines:
The answer NA means that the paper does not include experiments. If the paper includes experiments, a No answer to this question will not perceived well by the reviewers: Making the paper reproducible is important, regardless of whether the code and data are provided or not.
If the contribution is a dataset and/or model, the authors should describe the steps taken make their results reproducible verifiable.
Depending on the contribution, reproducibility can be accomplished in various ways. For example, if the contribution is a novel architecture, describing the architecture fully might suffice, or if the contribution a specific model and empirical evaluation, may be necessary to either make it possible for others to replicate the model with the same dataset, or provide access to the model. In general. releasing code and data often one good way to accomplish this, but reproducibility can also be provided via detailed instructions for how to replicate the results, access to a hosted model (e.g., in the case of large language model), releasing of a model checkpoint, or other means that are appropriate to the research performed.
While NeurlPS does not require releasing code, the conference does require all submis- sions to provide some reasonable avenue for reproducibility, which may depend on the nature of the contribution. For example
(a) If the contribution is primarily a new algorithm, the paper should make it clear how to reproduce that algorithm.
(b) the contribution is primarily a new model architecture, the paper should describe the architecture clearly and fully.
(c) If the contribution is a new model (e.g., a large language model), then there should either be way to access this model for reproducing the results or way to reproduce the model (e.g., with an open- open-source dataset or instructions for how to construct the dataset).
(d) We recognize that reproducibility may be tricky in some cases, in which authors are welcome to describe the particular way they provide for reproducibility. In the case of closed-source models, may be that access to the model is limited in some way (e.g., to registered users), but it should be possible for other researchers have some path reproducing or verifying the results.
27

---

5. Open access to data and code
Question: Does the paper provide open access to the data and code, with sufficient instruc- tions faithfully reproduce the main experimental results, described in supplemental material?
Answer: [Yes]
Justification: As note in Sections and complete open-source implementation will be publicly released. An anonymized version is available for review at: https
//github. wtw.chu..ooooceec dQP- andn/dP
Guidelines:
The answer NA means that paper does not include experiments requiring code.
Please see the NeurlPS code and data submission guidelines (https://nips. cc/
for more details.
While we encourage the release of code and data, we understand that this might be possible, so "No" is an acceptable answer. Papers cannot be rejected simply for not including code, unless this is central to the contribution (e.g., for new open-source benchmark).
The instructions should contain the exact command and environment needed to run to reproduce the results. See the NeurIPS code and data submission guidelines (https
//nips.cc/ httips..ccppbllcgguu.eee
for more details.
The authors should provide instructions data access and preparation, including how to access the raw data, preprocessed data, intermediate data, and generated data, etc.
The authors should provide scripts to reproduce all experimental results for the new proposed method and baselines. only subset of experiments are reproducible, they should state which ones are omitted from the script and why.
At submission time, to preserve anonymity, the authors should release anonymized versions (if applicable).
Providing as much information as possible in supplemental material (appended to the paper) recommended, but including URLs to data and code is permitted.
6. Experimental setting/details
Question: Does the paper specify all the training and test details (e.g., data splits, hyper- parameters, how they were chosen, type of optimizer, etc.) necessary to understand the results?
Answer: [Yes]
Justification: Experimental details are provided in Section and Appendix E Our open- source implementation reproduces the experiments described in the manuscript.
Guidelines:
The answer NA means that the paper does not include experiments.
The experimental setting should be presented in the of the paper to level of detail that is necessary to appreciate the results and make sense of them.
The full details can be provided either with the code, in appendix, or as supplemental material.
7. Experiment statistical significance
Question: Does the paper report error bars suitably and correctly defined or other appropriate information about the statistical significance of the experiments?
Answer: [Yes]
Justification: Illustrations of deviations in aggregated experimental results are included in Figures
3 and g
Guidelines:
The answer NA means that the paper does not include experiments.
The authors should answer "Yes" if the results are accompanied by emoor bars, confi- dence intervals, or statistical significance tests, least for the experiments that support the main claims of the paper.
28

---

The factors of variability that the error bars are capturing should be clearly stated (for example, train/test split, initialization, random drawing of some parameter, or overall run with given experimental conditions).
The method for calculating the error bars should be explained (closed form formula, call to a library function, bootstrap, etc.)
The assumptions made should be given (e.g.. Normally distributed errors).
It should be clear whether the error bar is the standard deviation the standard error of the mean.
It is OK to report 1-sigma error bars, but one should state it. The authors should preferably report a 2-sigma error bar than state that they have 96% CI, if the hypothesis of Normality of errors is not verified.
For asymmetric distributions, the authors should be careful not show in tables figures symmetric error bars that would yield results that are out range (e.g. negative error rates).
If error bars are reported in tables or plots, The authors should explain the text how they were calculated and reference the corresponding figures or tables in the text.
8. Experiments compute resources
For each experiment, does the paper provide sufficient information on the com- puter resources (type of compute workers, memory, time of execution) needed reproduce the experiments?
Question:
Answer: [Yes]
Justification:
Implementation details are provided in Appendix D
Guidelines:
The answer NA means that the paper does not include experiments.
The paper should indicate the type of compute workers CPU or GPU, internal cluster, or cloud provider, including relevant memory and storage.
The paper should provide the amount compute required for each of the individual experimental runs as well estimate the total compute.
The paper should disclose whether the full research project required more compute than the experiments reported in the paper (e.g., preliminary or failed experiments that didn't make into the paper).
9. Code ethics
Question: Does the research conducted the paper conform, in every respect, with the NeurlPS Code of Ethics httptt/ttttt.ea.........//aalaat
Answer: [Yes]
Justification:
This paper conforms with the NeurIPS Code of Ethics.
Guidelines:
The answer NA means that the authors have not reviewed the NeurlPS Code of Ethics.
If the authors answer No, they should explain the special circumstances that require deviation from the Code of Ethics.
The authors should make sure to preserve anonymity (e.g., there a special consid- cration due to laws or regulations their jurisdiction).
10. Broader impacts
Question: Does the paper discuss both potential positive societal impacts and negative societal impacts of the work performed?
Answer: [NA]
Justification: As a methodological paper, our work has no direct societal impact. Down- stream effects depend how the method is applied in specific domains.
Guidelines:
The answer NA means that there is no societal impact of the work performed.

---

If the authors answer NA or No, they should explain why their work has no societal impact or why the paper does not address societal impact.
Examples of negative societal impacts include potential malicious or unintended uses (e.g., disinformation, generating fake profiles, surveillance), fairness considerations (e.g., deployment of technologies that could make decisions that unfairly impact specific groups), privacy considerations, and security considerations.
The conference expects that many papers will be foundational research and not tied to particular applications, let alone deployments. However, if there is direct path to any negative applications, the authors should point it out. For example, it is legitimate to point out that improvement in the quality of generative models could be used to generate deepfakes for disinformation. On the other hand, it is not needed to point out that generic algorithm for optimizing neural networks could enable people to train models that generate Deepfakes faster.
The authors should consider possible harms that could arise when the technology is being used as intended and functioning correctly, harms that could arise when the technology is being used as intended but gives incorrect results, and harms following from (intentional or unintentional) misuse of the technology.
If there are negative societal impacts, the authors could also discuss possible mitigation strategies (e.g., gated release of models, providing defenses in addition attacks, mechanisms for monitoring misuse, mechanisms to monitor how system learns from feedback over time, improving the efficiency and accessibility of ML).
11. Safeguards
Does the paper describe safeguards that have been put in place for responsible release of data models that have high risk for misuse (e.g., pretrained language models, image generators, or scraped datasets)?
Question:
Answer: [NA]
Justification: [NA]
Guidelines:
The answer NA means that the paper poses no such risks.
Released models that have high risk for misuse or dual-use should be released with necessary safeguards to allow for controlled use of the model, for example by requiring that users adhere to usage guidelines or restrictions to access the model or implementing safety filters.
Datasets that have been scraped from the Internet could pose safety risks. The authors should describe how they avoided releasing unsafe images.
We recognize that providing effective safeguards challenging, and many papers do not require this, but we encourage authors to take this into account and make best faith effort.
12. Licenses for existing assets
Are the creators or original owners of assets (e.g.. code, data, models), used in the paper, properly credited and the license and terms of use explicitly mentioned and properly respected?
Question:
Answer: [Yes]
Justification:
All related work is cited and acknowledged.
Guidelines:
The answer NA means that the paper does not use existing assets.
The authors should cite the original paper that produced the code package dataset.
The authors should state which version of the asset used and, possible, include URL
The name of the license (e.g., CC-BY 4.0) should be included for each asset.
For scraped data from a particular source (e.g., website), the copyright and terms of service of that source should be provided.
30

---

If assets are released, the license, copyright information, and terms of use in the package should be provided. For popular datasets, paperswi thcode. has curated licenses for some datasets. Their licensing guide can help determine the license of dataset.
batasets.papersstt
For existing datasets that are re-packaged, both the original license and the license of the derived asset (if it has changed) should be provided. If this information is available online, the authors encouraged to reach out the asset's s creators.
13. New assets
Question: Are assets introduced in the paper well documented and is the documentation provided alongside the assets?
Answer: [Yes]
Justification: Our open-source implementation will be publicly released. An anonymized version is available for review at: http:///gtttuuu..oodope/t--oo wtt.ccmmm///pp.
Guidelines:
The answer NA means that the paper does not release assets.
Researchers should communicate the details of the dataset/code/model as part of their submissions via structured templates. This includes details about training, license, limitations, etc.
The paper should discuss whether and how consent was obtained from people whose asset is used.
At submission time, remember anonymize your assets (if applicable). You can either create an anonymized URL or include an anonymized zip file.
14. Crowdsoureing and research with human subjects
Question: For crowdsourcing experiments and research with human subjects, does the paper include the full text of instructions given to participants and screenshots, if applicable, as well details about compensation (if any)?
Answer: [NA]
Justification: [NA]
Guidelines:
The answer NA means that the paper does not involve crowdsourcing nor research with human subjects.
Including this information in the supplemental material is fine, but if the main contribu- tion of the paper involves human subjects, then much detail as possible should be included the main paper.
According to the NeurlPS Code of Ethics, workers involved in data collection, curation, other labor should be paid at least the minimum wage in the country of the data collector.
15. Institutional review board (IRB) approvals or equivalent for research with human subjects
Question: Does the paper describe potential risks incurred by study participants, whether such risks were disclosed the subjects, and whether Institutional Review Board (IRB) approvals (or an equivalent approval/review based on the requirements of your country or institution) were obtained?
Answer: [NA]
Justification: [NA]
Guidelines:
The answer NA means that the paper does involve crowdsourcing nor research with human subjects.
Depending the country in which research is conducted, IRB approval (or equivalent) may be required for any human subjects research. If you obtained IRB approval, you should clearly state this in the paper.
31

---

We recognize that the procedures for this may vary significantly between institutions and locations, and we expect authors adhere the NeurIPS Code of Ethics and the guidelines for their institution.
For initial submissions, do not include any information that would break anonymity (if applicable), such as the institution conducting the review.
16. Declaration of LLM usage
Does the paper describe the usage of LLMs if it is an important, original, or non-standard component of the core methods in this research? Note that if the LLM is used only for writing, editing, formatting purposes and does not impact the core methodology. scientific rigorousness, or originality of the research, declaration is not required.
Question:
Answer: [NA]
Justification: [NA]
Guidelines:
The answer NA means that the method development in this research does not involve LLMs as any important, original, or non-standard components.
Please refer to our LLM policy (https: /ttpps////eurpss.... /cc/////eeeeeeecs///22///// for what should or should be described.
32


Reviews:
[{"review_id": "XCEsMY4BZo", "reviewer": "Reviewer_13br", "summary": "This paper provides a modular approach, called dQP, for differentiable optimization through QP. Unlike existing approaches, dQP can be equipped with any solver. Once the solution is identified, the QP system is reduced to equality constraint problems by restricting to the active constraints. Therefore, the differentiation operation is reduced to solving a simple linear system. Experiments show dQP is very promising, particularly for large-scale sparse systems.", "questions": "Can your framework be extended to handle linear programs or, more generally, convex programs in which the quadratic objective is replaced by a linear or other convex function? If not, could you clarify the specific technical limitations that prevent such an extension?", "limitations": "Yes, the limitation has been properly addressed.", "rating": 5, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 3, "clarity": 3, "significance": 3, "originality": 2, "strengths_and_weaknesses": "# Strengths\nThe manuscript is well structured and easy to follow, with each theoretical idea accompanied by intuitive explanations and supporting figures.  The experiments convincingly demonstrate both performance and scalability; in particular, dQP is the only method in the comparison that successfully solves very large‐scale QPs while maintaining competitive accuracy.\n\n\n# Weaknesses \n\nA potential limitation is the migration to GPU implementation. The current implementation is CPU-only. Porting the method to GPUs will almost certainly require iterative solvers for the reduced KKT system, as I highly doubt the speedup of direct linsysm solver in GPU; yet such iterative solvers often trade speed for numerical precision. Recent first-order, GPU-friendly QP solvers (e.g., PDQP (https://github.com/jinwen-yang/PDQP.jl).) suggest one possible path forward, but PDQP is written in Julia and cannot be plugged into PyTorch without significant engineering. A native GPU backend may therefore need to be developed from scratch.\n\nAnother remaining issue is how to handle degeneracy. The paper introduces an active-set refinement heuristic that appears effective in practice, but a deeper analysis of how solver accuracy influences active-set identification—especially in degenerate or near-degenerate regimes—would strengthen the contribution. While full theoretical treatment may be beyond the current scope, even an empirical study of pathological cases would add confidence in the method’s robustness.", "was_revised": true, "final_justification": "The rebuttal have addressed the proposed issue and i am generally ok. I will keep my score", "rebuttal": "We thank the reviewer for their thorough summary, review, and questions. \n\nWe agree that integration with GPU-based sparse iterative solvers for QP and linear systems will be highly desirable once such tools become widely available, as we noted in our conclusion. In the meantime, our work contributes a scalable CPU-based solution that enables the solution of large-scale problems beyond the reach of existing approaches.\n\nSome problems inherently exhibit sensitivity or ill-conditioning in the active set, leading to near-degeneracy in the KKT conditions. Nevertheless, as noted in lines 243–245, we used a fixed threshold to determine the active set across all our experiments. Our implementation provides supplementary tools to help users select solvers, tune parameters, and optionally apply an active-set refinement heuristic (Appendix D). The latter was not used in our experiments. \n\n> Can your framework be extended to handle linear programs or, more generally, convex programs in which the quadratic objective is replaced by a linear or other convex function?\n\nA1: Linear programs are challenging because their solutions are generically at extremal points, leading to complexities in differentiation that we, and most other differentiable QP solvers, do not directly address. We are currently working on extensions of the framework to more general classes of convex optimization, like described by the reviewer."}, {"review_id": "ROBxHgSqIu", "reviewer": "Reviewer_gwAp", "summary": "The paper investigates how to compute gradients through black-box Quadratic Programming (QP) solvers to efficiently incorporate QP layers into neural networks or bilevel optimization tasks. It introduces dQP, a modular framework that enables automatic differentiation of any QP solver, transforming it into a differentiable layer. The core theoretical insight is based on a key observation: given a known active constraint set, the gradient of the QP solution can be explicitly expressed as a simplified linear system, significantly reducing computational complexity.", "questions": "* How is the duality gap defined in the context of dQP? Why is it used as a convergence metric? While the paper uses the duality gap to report solution quality, a precise formal definition (in this context) and explanation of why it reliably measures convergence across solver backends would be helpful. Does it account for both primal and dual feasibility in the reduced KKT system?\n\n* While dQP is clearly more flexible, I am still unclear on the primary motivation from a user perspective. Why would users choose dQP over existing approaches like OptNet, which already perform quite well? In other words, what is the key intuition or motivation for dQP?", "limitations": "N/A", "rating": 5, "confidence": 5, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 3, "clarity": 4, "significance": 2, "originality": 2, "strengths_and_weaknesses": "# Strengths\n* The paper provides an extensive and rigorous investigation of differentiable quadratic programming (QP), situating its contributions within a rich body of related work. The authors evaluate their method across a wide variety of benchmarks, including dense and large-scale sparse QPs, and demonstrate dQP’s superior scalability and robustness.\n* A major strength is the explicit formulation of the differentiation strategy based on the active set. The authors clearly connect sensitivity analysis, the KKT conditions, and the derivation of a reduced KKT system. This makes the framework both conceptually clean and practically efficient.\n*  The proposed dQP framework is notable for supporting a wide range of existing QP solvers via a modular design. It decouples differentiation from the solver itself, enabling plug-and-play usage with over 15 solvers, including commercial and open-source ones, which significantly broadens the applicability of differentiable optimization.\n\n\n# Weakness\n* The paper reports timing comparisons in milliseconds. However, the differences—especially in the context of total runtime—may amount to only a few seconds even for large-scale problems. It would help if the authors could clarify the scenarios or applications (e.g., bi-level learning or large-scale simulation) where these runtime savings are practically consequential.\n\n* The paper demonstrates excellent performance on a wide benchmark, but there is little discussion of potential failure modes (e.g., sensitivity to incorrect active set identification, degeneracy, or ill-conditioning near singularities). While Appendix D mentions heuristics, it would strengthen the paper to include a quantitative evaluation of robustness in such edge cases.\n\nOverall, I think this is a well-written paper currently.", "was_revised": true, "final_justification": "The author's response well addressed my concern, so I will maintain my original score and vote for acceptance.", "rebuttal": "We thank the reviewer for their thorough summary and review. We address their questions below.\n\n> How is the duality gap defined in the context of dQP? Why is it used as a convergence metric? While the paper uses the duality gap to report solution quality, a precise formal definition (in this context) and explanation of why it reliably measures convergence across solver backends would be helpful. Does it account for both primal and dual feasibility in the reduced KKT system?\n\nA1: We define the duality gap and residuals for QP in Appendix E.1 Lines 742-748. We use the duality gap measure for accuracy as it readily provides an upper bound on the objective error (with respect to the true unknown optimal objective value) [27]. The duality gap is routinely used in benchmarks (e.g., QPBenchmark) and in the stopping criteria of many solvers. We separately measure primal and dual feasibility, and while they are not reported, they are used along with the duality gap, as described in Appendix E.1 Lines 751-756,  to determine “solver success.” We will clarify these points in the main text.\n\n[27] Boyd and Vandenberghe, Convex Optimization, Pages 241-242\n\n> While dQP is clearly more flexible, I am still unclear on the primary motivation from a user perspective. Why would users choose dQP over existing approaches like OptNet, which already perform quite well? In other words, what is the key intuition or motivation for dQP?\n\nA2: Users should choose dQP for large-scale, sparse problems, where its ability to leverage a wide range of commercial solvers (as noted by the reviewer) can yield orders-of-magnitude speedups over other differentiable QP methods (Figure 1). The advantage of this flexibility is also evident in the MM experiment (Figure 4). While OptNet performs well on some problems, it fails on 70% of the MM instances (Figure 4), underperforms on the rest (Table 1), and scales poorly for sparse problems (Figures 1, 5; Tables 2-3, 5-7). Other differentiable solvers face similar limitations.\n\nThis also connects to the reviewer’s point about runtime differences: though small in isolation, these differences accumulate and become crucial (or even prohibitive) in bilevel optimization and learning applications. We thank the reviewer for pointing this out and will clarify it in the revision."}, {"review_id": "S3SNAcZTwK", "reviewer": "Reviewer_EkQp", "summary": "This work studies how to compute the gradients of QP solutions with respect to QP problem parameters, making QP solve a differentiable layer that can be potentially integrated into a larger computational graph. The problem is a classic, extensively studied topic in sensitivity analysis. The proposed method obtains QP solutions from an external QP solver and analyzes the KKT conditions that connect the QP solution to the problem parameters. One of the key ingredients in the technical method is the reduced dimensionality analyzed from the active set of QP constraints. The paper provides an open-source implementation and evaluates its performance against several differentiable QP implementations in PyTorch.", "questions": "Perhaps my only major concern with this paper is the novelty of the technical method. I will first list my understanding of the contribution #1 (lines 57-58) in the paper, which seems to focus on two key claims:\n1. The paper derives the gradients of a QP solver from its primal solution without dependency on the exact QP solver itself. This leads to the solver-agnostic feature of dQP, which can be used with any existing QP solvers.\n2. The paper uses the set of active constraints at the solution to simplify the dimension of the system solved for gradient computation, which seems to be rooted in the standard implicit function theorems/adjoint methods applied to the KKT system of the QP problem.\n\nFor #1: Decoupling backward differentiation from its forward solver is not a novel thing in sensitivity analysis, e.g., previous works have explained it when differentiating Navier-Stokes flows [Stuck 2012]. This fact is a natural result of the implicit function theorem and the classic adjoint method (Bradley et al. [2004]).\n\nFor #2: Pruning inactive constraints is a standard practice in differentiable programming when a solver contains constraints. A concrete example is differentiable simulation with a Columb contact and friction model [Li et al. 2022], where inactive constraints from the complementarity conditions in the Columb model are pruned in gradient computation.\n\nTherefore, while I agree that dQP may be novel for differentiating QP problems, I feel it is also an expected and predictable result from existing theories, as similar patterns from common practice in sensitivity analysis/differentiable programming can be found in existing works that differentiate other numerical solvers. This shakes my confidence in the significance of the proposed method. My questions to the authors are:\n1. Could you summarize what new knowledge this paper brings to differentiable programming/sensitivity analysis?\n2. Is there anything QP-specific that can be highlighted or exploited in developing your technical method?\n\nReferences\n- Andrew Bradley, PDE-constrained optimization and the adjoint method, 2004.\n- Arthur Stuck, Adjoint Navier-Stokes methods for hydrodynamic shape optimisation, 2012.\n- Yifei Li, et al. DiffCloth: Differentiable Cloth Simulation with Dry Frictional Contact, 2022.", "limitations": "I appreciate that the paper is upfront about the lack of GPU support. I am actually a bit surprised to learn that “state-of-the-art sparse and scalable QP solvers with GPU support are still lacking.” I thought QP algorithms are largely based on solving linear systems of equations, and large-scale linear solvers have GPU support these days. I am taking this sentence at face value and won't let this claim affect my evaluation of this work negatively, but I am curious to learn more about this claim if the authors are willing to elaborate on it.", "rating": 4, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 3, "clarity": 3, "significance": 2, "originality": 2, "strengths_and_weaknesses": "Overall, I like the direction this paper explores, and I want to commend the open-source implementation.\n\n**Strengths**\n\n1. I appreciate that the paper offers an open-source implementation of the proposed algorithm. I think the engineering effort and the modular design look valuable to me.\n\n2. The experiments on many baselines and benchmark problems seem comprehensive. I think these experiments have covered some widely accepted and fairly strong numerical solvers, e.g., Gurobi, Pardiso, and cvx.\n\n**Weaknesses**\n\nThe novelty and significance of the contribution “QPs can be explicitly differentiated using only the primal solution via a reduced dimension locally equivalent linear system determined by the active set” is unclear to me. I will leave my specific questions in the “Questions” section below.", "was_revised": true, "final_justification": "The rebuttal's answers to my questions are more or less expected. Technical innovation is not the strength of this work. Instead, it finds and fills a gap in differentiating QP with a well-integrated system. This was roughly what I thought when I was writing my original review, so the rebuttal does not significantly affect my rating either. I am comfortable with maintaining my borderline accept rating.", "rebuttal": "We thank the reviewer for their thorough summary and review. We address their questions below.\n\n> Decoupling backward differentiation from its forward solver is not a novel thing in sensitivity analysis … Could you summarize what new knowledge this paper brings to differentiable programming/sensitivity analysis? Is there anything QP-specific that can be highlighted or exploited in developing your technical method?\n\nA1: We agree with the reviewer that we did not introduce the general idea of separating backward differentiation from the forward solver. However, in the context of QPs, prior methods often tightly integrate differentiation with specific solvers, e.g., to enable information sharing between the forward and backward. As we acknowledge in Section 3 and the related work section, components of our approach (e.g., pruning inactive constraints) have appeared before, but no prior work has synthesized these ideas into a fully modular differentiable QP layer. Our contribution lies in combining these elements to provide a clear, flexible framework, and in making the key observation that, for QPs, a single symmetric solve over the active set can yield both dual variables if needed and derivatives. Practically, we also demonstrate that for certain problems, especially large-scale, sparse ones, the benefits of using efficient forward solvers far outweigh any gains from solver-specific integration, supporting the case for modularity over tight coupling.\n\n> I appreciate that the paper is upfront about the lack of GPU support. I am actually a bit surprised to learn that “state-of-the-art sparse and scalable QP solvers with GPU support are still lacking.” I thought QP algorithms are largely based on solving linear systems of equations, and large-scale linear solvers have GPU support these days. I am taking this sentence at face value and won't let this claim affect my evaluation of this work negatively, but I am curious to learn more about this claim if the authors are willing to elaborate on it.\n\n\nA2: We are happy to elaborate, though we are not experts on this topic. Until recently, Gurobi maintained that the potential benefits of GPU-based solvers were minimal, guiding our own view. A key challenge is that large-scale problems of interest are typically sparse, and efficiently exploiting sparsity on GPUs remains a difficult and developing area. Recently, Gurobi and NVIDIA began collaborating on a GPU-based solver, “cuOpt.” Gurobi documents the latest and previous updates on their efforts to use GPUs in optimization in their forum post, “Does Gurobi support GPUs?”"}, {"review_id": "XPkEsc8Yhy", "reviewer": "Reviewer_dAwU", "summary": "This manuscript focuses on differentiation of quadratic programming problems. They note that there are many QP solvers and many previous attempts try to force the QP into some specific solver technique and differentiate that approach. But, this can make it so you can't always use the best QP solver. Thus instead of forcing everything to use a potentially not so great solver, the authors derive an adjoint rule so that a black-box QP solver can be used and a single linear solve (via LDT factorization) can reverse for the derivative. They then go on to show that if you use the right QP solver for each problem its faster than using a not so good one.", "questions": "1. You support 15 commercial solvers and 1 open source one with your new package: in Figure 4 when you say \"ours\", which QP solver is that using? It's ambiguous as to what exactly was used to solve there since that refers to your package and differentiation rule, but since there are 16 solvers it's not clear which one, or if that makes a difference.\n2. Why not benchmark with the 16 different solvers together to better explain in which regime the given methods do better/worse?\n3. Are there any domains where the other differentiation techniques do better? Obviously not in the large scaling case, but maybe there are cases that could be highlighted where not treating it as a black-box actually enables some optimization?\n4. The method in the paper assumes zero residual of the forward solve and uses the solution value for the derivative. If there is an error in the forward solve, how does that error propagate into the derivative? Does it scale like the square root of the condition number of the Jacobian matrix being factorized? Does using the LDL / Cholesky could be improving the numerical accuracy by decreasing the growth here vs LU/QR?", "limitations": "Yes", "rating": 4, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 3, "clarity": 3, "significance": 3, "originality": 3, "strengths_and_weaknesses": "The strength and the weakness is kind of the same thing here... it's pretty obviously the way to do this? I can easily see the mindset of the authors here. You scour the literature and notice that things like JaxOpt is differentiating the QP by transforming it into a convex optimization and using CVXPY, which yes you can always solve a QP by using a convex optimization solver because QPs are a subset but no one does that because specializing on the quadratic properties makes the solve much faster and more robust. And Alt-Diff's custom ADMM method is clearly way slower than actually good QP solvers like commercial solvers (e.g. Gurobi). So clearly the right thing to do is to just define the derivative of a QP solve and use a good solver and it will just be a whole lot faster than these bad QP solvers. So they did that and got the expected results. QED.\n\nIs that novel enough? I don't really know. What's more surprising really is that the other authors were allowed to get away without this as the standard baseline test. The ADMM method and other papers, they should have been asked \"test your method against defining the implicit differentiation rule over Gurobi and show the performance difference\", and when it failed to be anywhere close to efficient that should've been the end of the story. In most other areas of differentiable programming, this has always been the baseline. You'd have to test differentiating through an ODE solver, nonlinear solver, etc. vs deploying a black-box rule that covers all solvers, and these black-box rules are implemented over lots of solvers in different domains already. I double checked the literature here and yes, it seems in QP solvers people just didn't do this, which is really odd, so the authors are doing the right thing that should've been done the first time like in all of the other domains. Is that novel enough for a NeurIPS publication? It is clearly the right thing to do and the best method, but most practitioners would have guessed that before reading it and its standard exercise to derive and define that adjoint rule, so there's nothing that is really shocking here other than the fact that this wasn't already done on QPs. Thus I am on the fence as to novelty, it's a strength that it's the clear best method to do but it's a weakness that it's not unobvious.\n\nThe fact that the adjoint is a linear solve is not surprising: the implicit differentiation theorem proves that all such cases of solvers which iterate to convergence has an adjoint rule that is given as a linear system. And the derivation of this linear system from the KKT equations is straightforward and standard. That's all run of the mill graduate-level homework these days. But, by subsetting it to the active set, noticing that gives a symmetric solve and thus can use an LDL factorization, that is big insight that makes this a lot more efficient over the naive application of the implicit differentiation rule. Thus I would point to this as the novelty of the paper, not what is mentioned in the title or abstract about supporting any black-box QP solver and that being the right way to do this, but this trick within the adjoint definition is a nice little trick.", "was_revised": true, "final_justification": "I think it's pretty clear from the author responses that all of the good properties from the method are likely to come from just using a better QP solver and nothing really to do with their algorithm. Which is fine, but it should be better clarified in the paper. That makes it borderline, but since it's the right algorithm to do and QP solvers aren't doing this, it is at least acceptable.", "rebuttal": "We thank the reviewer for their thorough summary and review. We address their questions below.\n> You support 15 commercial solvers and 1 open source one with your new package: in Figure 4 when you say \"ours\", which QP solver is that using? It's ambiguous as to what exactly was used to solve there since that refers to your package and differentiation rule, but since there are 16 solvers it's not clear which one, or if that makes a difference.\n\nA1: In Figure 4, “Ours (QPBenchmark)” means that for each problem we select the QP solver  that performed the best in the open-source benchmark, QPBenchmark. We noted this in  L277-279 and specified the solver counts in the caption.\n\nWe use the open-source QPSolvers as a unified interface, allowing dQP to support all solvers it provides, currently including three commercial solvers (Gurobi, MOSEK, NPPro) and several open-source ones. QPSolvers is regularly updated with new solvers, which dQP inherits automatically. \n\nWe will clarify these points in the text.\n\n> Why not benchmark with the 16 different solvers together to better explain in which regime the given methods do better/worse?\n\nA2: We found that including all ~2000 dQP variants with every supported solver in Figure 4 would have shifted the focus toward QP solver comparisons (covered extensively by QPBenchmark) and obscured the modularity benefits. We will revise the caption of Figure 4 and related text to better clarify.\n\n> Are there any domains where the other differentiation techniques do better? Obviously not in the large scaling case, but maybe there are cases that could be highlighted where not treating it as a black-box actually enables some optimization?\n\nA3: Yes, as we mention in L263-265. For example, OptNet’s reuses the internal factorizations in their custom GPU-batchable forward solver in their backward, which is not possible with black-box solvers. As a result, they can solve very small dense QPs in batches without a backward differentiation bottleneck, but cannot scale up to large, sparse problems. \n\n> The method in the paper assumes zero residual of the forward solve and uses the solution value for the derivative. If there is an error in the forward solve, how does that error propagate into the derivative?\n\nA4: The analysis of error propagation is complicated by the combinatorial quality of the active set that is crucial to dQP’s algorithm. If the problem is well-conditioned and the solvers used are sufficiently accurate so that the active set identified is locally stable, the error propagation to the derivatives boils down to the sensitivity of linear systems to perturbations, characterized by the condition number like the reviewer notes [24].\n\n[24] Blondel et. al., Efficient and Modular Implicit Differentiation, Theorem 1"}]

Author rebuttal:

Meta-review: This paper introduces dQP, a modular and solver-agnostic framework for plug-and-play differentiation of virtually any QP solver. The authors' open-source implementation integrates more than 15 state-of-the-art solvers. Extensive experiments comparing against other strong differentiable solvers demonstrate the accuracy and robustness of dQP on small- and large-scale settings, and in both sparse and dense regimes.

In some sense, the submission does not present anything conceptually new. The principles of designing QP layers via implicit differentiation through the KKT conditions, of active set differentiation, of using Jacobian-vector products, and of decoupling the forward pass from the backward pass are not new. (While the authors call their differentiation technique "explicit differentiation," this seems to be at odds with how the differentiable layers literature uses this term - this is still differentiation through the KKT conditions at the solution point, albeit using an equivalent problem, rather than differentiating through the unrolled iterations of the QP solver.) In addition, the solver presented is CPU-only and is not evaluated/usable in batched settings.

Despite this, the reviewers felt the paper makes a strong contribution through its software package. In particular, it is beneficial to allow the use of many state-of-the-art solvers, allowing the right solver for the problem to be chosen and for in-built tricks such as warm-starting to be leveraged. The package is indeed benchmarked thoroughly in non-batched settings, demonstrating good performance. The active-set refinement strategy is important in ill-conditioned settings, and other differentiable solvers do not necessarily provide this. There is support for sparse solvers, allowing for associated speedups, whereas the GPU-based solvers in the literature only support dense representations.

Overall, the strengths of the practical contribution (software package and extensive benchmarking) seem to outweigh the limited conceptual novelty, providing a case for acceptance. However, the authors should consider significantly updating the narrative of the paper to make the true contributions clearer (notably, to indicate that the main contribution and source of benefit is the software package that enables better/diverse solvers to be used in the forward pass, rather than anything about the backward pass algorithm).
Decision: Accept (poster)