Enforcing convex constraints in Graph Neural Networks
Authors: Ahmed Rashwan
Keith Briggs
Chris Budd
Lisa Maria Kreusser
Venue: NeurIPS 2025 poster ()
Keywords: Graph Neural Networks
Constrained optimisation
Dykstra's projection algorithm
Primary area: deep_learning

Abstract: Many machine learning applications require outputs that satisfy complex, dynamic constraints. This task is particularly challenging in Graph Neural Network models due to the variable output sizes of graph-structured data. In this paper, we introduce ProjNet, a Graph Neural Network framework which satisfies input-dependant constraints. ProjNet combines a sparse vector clipping method with the Component-Averaged Dykstra (CAD) algorithm, an iterative scheme for solving the best-approximation problem. We establish a convergence result for CAD and develop a GPU-accelerated implementation capable of handling large-scale inputs efficiently. To enable end-to-end training, we introduce a surrogate gradient for CAD that is both computationally efficient and better suited for optimization than the exact gradient. We validate ProjNet on four classes of constrained optimisation problems: linear programming, two classes of non-convex quadratic programs, and radio transmit power optimization, demonstrating its effectiveness across diverse problem settings.

Enforcing convex constraints in Graph Neural Networks
Ahmed Rashwan University of Bath
Keith Briggs
Chris Budd University of Bath
Lisa Kreusser
BT Research
University of Bath and Monumo
Abstract
Many machine learning applications require outputs that satisfy complex, dy- namic constraints. This task is particularly challenging in Graph Neural Network models due the variable output sizes of graph-structured data. In this paper, we introduce ProjNet, a Graph Neural Network framework which satisfies input- dependant constraints. ProjNet combines a sparse vector clipping method with the cuaponent-Averaged Dykstra (CAD) algorithm, an iterative scheme for solving the best-approximation problem. We establish convergence result for CAD and develop a GPU-accelerated implementation capable of handling large-scale in- puts efficiently. To enable end-to-end training, we introduce a surrogate gradient for CAD that is both computationally efficient and better suited for optimization than the exact gradient. We validate ProjNet on four classes of constrained op- timisation problems: linear programming, two classes of non-convex quadratic programs, and radio transmit power optimization, demonstrating its effectiveness across diverse problem settings.
Introduction
In many machine learning (ML) applications, model outputs must satisfy certain constraints to en- sure feasibility, safety, interpretability [22 40 A common approach to enforcing constraints through activation functions (e.g.. Sigmoid, Softmax) in the output layer [41], which restrict outputs to predefined set. However, while effective for simple linear box constraints, this ap- proach is not applicable for more complicated types of constraints. Many real-world applications require complex, data-dependent constraints. For example, in roboties, constraints can encode fea- sible movements or collision avoidance [21 [4| 63]. Similarly, in energy systems, power generation limits vary with demand fluctuations [48] and in industrial processes, physical constraints must be respected [46]. Such tasks are often modelled as optimisation problems with convex constraints.
Several approaches have been proposed for enforcing given set of convex constraints C' in ML models, each with distinct strengths and drawbacks. Many existing methods provide only approxi- mate feasibility, typically by incorporating the constraints into the training loss using barrier func- tions or duality [22 27 33). Instead, we will only be considering methods which yield strict feasi- bility guarantees. One common strategy involves first selecting arbitrary reference point within C, and then using an ML model to compute an adjustment vector about this reference, ensuring that the final output remains feasible and accomplishes the desired task; we refer to this as the vec- tor clipping method 33 52 59). Typically, the reference point is obtained using a convex solver, which generally more computationally expensive than evaluating the adjustment vector. While this approach works well for fixed constraint sets, becomes significantly more costly when the constraints vary with the model input, as new reference point must be computed for each instance. This limitation makes training vector clipping models with larger batch sizes impractical.
email: ar30090bath. ac.uk
39th Conference on Neural Information Processing Systems (NcurIPS 2025).

---

A related line of research explores the integration of differentiable optimisation layers within ML models 3 12]. These methods leverage convex solvers capable of solving batches of optimisa- tion problems concurrently, with gradients computed analytically. As a result, optimisation layers are more efficient in handling input-dependent constraints and have applications beyond constraint satisfaction. However, key limitation lies in their scalability to large problem instances. The solvers used typically rely on interior-point methods, which require matrix factorizations that are difficult to efficiently implement on GPUs, particularly for sparse matrices (4). This makes it difficult to use this approach with GNNs, since graph-structured data is typically sparse and leads to heterogenous batches, requiring efficient exploitation of sparsity sustain performance.
In addition, simple iterative algorithms have been embedded in ML models for constraint satisfac- tion, sharing similarities with our approach. A prominent example is the Sinkhorn algorithm (19), which has been used to compute stochastic matrices for learning permutations [24 55) and for satis- fying positive-linear constraints [61]. These iterative algorithms are well-suited for GPU implemen- tations. However, their gradients are typically computed by unrolling the iterative procedure through back- pack-propagation, which can be both computationally expensive and memory-intensive. We argue that unrolling is not always required: since many iterative algorithms-such as Sinkhor--roduce solutions to specific optimisation problems, their gradients can often be derived analytically, leading to more efficient methods of gradient computation [44]
Projection methods have long been used for finding feasible points in convex sets [16]. A classic
example is Von Neumann's 's algorithm (6. which uses alternating projections to find feasible point in the intersection of two convex sets. While this algorithm produces feasible points, it generally does not allow us to characterize which point in is returned. In contrast, Dykstra's projection algorithm (34) can find the unique projection of a point onto the intersection of convex sets. For our purposes, being able to determine the exact point produced by the algorithm crucial, as it enables the computation of surrogate gradient. The Component-Averagel Dykstra (CAD) algorithm (14) a variant of Dykstra's method which exploits problem sparsity by using component averaged projections [15] on the constituent sets C, CAD can be efficiently accelerated using GPU scattering operations (35 and works seamlessly on batches of problems. We present a convergence result for the CAD and demonstrate its application to constraint satisfaction in ML models.
Building on these ideas, we introduce ProjNet, a GNN architecture designed for constraint satis- faction. ProjNet leverages the CAD) algorithm to guarantee feasibility, while incorporating sparse variant of the vector clipping method to improve model expressiveness. The architecture is end-end differentiable and fully GPU accelerated. Although this work focuses linear constraints, most of the methods presented naturally extend to broader classes convex constraints. We demonstrate the effectiveness of ProjNet across four classes of constrained optimization problems.
1.1 Contributions
Our key contributions are as follows:
1. We establish theoretical foundations of our work by proving convergence result for the CAD algorithm, explicitly determining its limit, and incorporate these into an efficient GPU imple- mentation, enabling significant speed-ups and scalability for large-scale applications.
2. We introduce computationally efficiently surrogate gradient for the CAD algorithm, tailored for MI. applications.
3. We present a refined variant of the vector clipping method, which we refer to as sparse vector clipping, which leverages problem sparsity and is well-suited for use in GNN models.
4. We propose ProjNet, GNN architecture for constraint satisfaction problems which integrates the CAD algorithm with sparse vector clipping. Unlike many existing approaches, our method naturally supports batched graph inputs.
5. We demonstrate the effectiveness of our GNN-based approach on four constrained optimisa- tion problems: linear programming, two classes of non-convex quadratic programs, and radio transmit power optimisation.

---

1.2 Overview
This paper is structured as follows. Section describes the individual components of our GNN architecture. We study the CAD algorithm in Section providing a convergence result, surrogate gradient, and numerical experiments. In Section 4 the effectiveness of our method is demonstrated on four classes of optimisation problems. Finally, Section 5attionn5conce concludes the paper.
2 ProjNet: A GNN framework for constraint satisfaction
We propose ProjNet, GNN-based framework for solving constraint satisfaction problems. After introducing some preliminaries and notation in Section 2.1 we construct a constrained input graph ge in Section 2.2 which embeds set of linear output constraints C' into given input graph G. Pro- jNet is then defined as a model GNN, that maps any constrained input graph to a feasible solution, Le., GNNo(9c) C' for all Gc. The architecture enforces constraint satisfaction through two core
components: a projection layer and sparse vector clipping layers, detailed in Sections 2.3 and 2.4
respectively.
2.1 Preliminaries
Convex constraint sets We define a convex, non-empty constraint set the intersection of individual constraint, where each C, = {z E R" 9:((1))
0 k=1 C, C R as
< 0}. Here,
N; c: {l...., n specifies the subset of variables affected by each constraint C. and : RR - R a convex function. Given component j E (1..... n), we write I.j (i: j E N₂) = (1...., m] for the set of constraints affecting. j. and set 15 = |T, We are especially interested in linear constraints where (x) A b, for i & m, and A, denote the rows of some matrix
[A Am] Rmxx and b E RRm. For linear constraints, the constraint set C' {x IR" Ar & b} is polytope. We denote the unique projection of point x onto a convex set C' by (a(x)
a)) = arg min || - =||2.
GNN inputs We assume that all problem instances are given as weighted graphs g = (x, E) where the rows X E B correspond to n nodes G. each with a k-dimensional node feature, and E E R=*n encodes edge weights of all the edges connecting the nodes. The input of the GNN will be a constrained input graph Gc which will be associated with graph G and constraint set C, and will be defined in Section 2.2 We denote the output of the GNN by y I² where each output component y, E I corresponds to the output feature of node i, Our aim is to enforce the constraint y e C where the constraint C is input-dependant.
Batched inputs GNNs can be simultancously applied batch of constrained input graphs as- sociated with batch of input graphs (9I...., Ga) of the form Gj = (X,, E.) for i E {l,.... .......3 8} by interpreting them as one large disconnected graph g (cat(X,,. ... . Xs), diag(E E1..... E3)) Here, cat denotes concatenation and diag computes a block-diagonal matrix. This perspective al- lows GNNs and all other forms of graph computation discussed in this paper to generalise seamlessly to the batched case. Since these batched graphs are larger and more sparse than the individual input graphs Ge, is particularly important to develop methods that exploit their sparsity.
2.2 Constrained input graphs
29 where a linear ob- jective is combined with linear constraints into single graph, we extend this methodology beyond linear objectives to more general tasks represented by some graph G. Specifically, we construct a constrained input graph Go from given input graph G and set of linear constraints Ar This new graph Gc = (X, b, A, E)) is heterogeneous, comprising variable nodes X, constraint nodes b, and two types of edges: E and A. The graph contains n variable nodes with k-dimensional features from the rows of X, and m constraint nodes, each associated with a scalar feature b. Edges in E rep- resent relationships among variable nodes X, while the matrix serves as weighted adjacency matrix capturing connections between variable and constraint nodes. This unified graph represen- tation enables the use of single GNN to handle tasks that incorporate both problem instances and their associated constraints.
Inspired by graph-based formulations in mixed-integer programming [7
3

---

Some of the we methods we present exploit the sparsity of Gc by identifying independent groups of constraints. Specifically, we wish to identify the finest partition P of {L..... mj such that con- straints from different components of P do not share any constrained variables. This partition corre- sponds to the connected components of the bipartite graph (X, b, A) and can be computed efficiently using GPU graph processing software 62] or alternatieely, via simple label propagation scheme outlined in Appendix D.I
2.3 Component-averaged Dykstra (CAD) algorithm
To ensure feasibility, we wish to compute the projection of any vector & E R" onto C. denoted by Pc(x). To achieve this, we use the component-averiee scheme which exploits problem sparsity and is well suited for GPU implementations. CAD belongs to a family of Dykstra-style projection algorithms, and is closely related to the two-set Dykstra algorithm [34] and the Simultaneous Dykstra Algorithm (38). All of these methods operate under common principle: computing the projection Pc onto C' by using the projections Pc, on individual sets C. For m 2, i.e. C' C₁ a C2, the iterates of the two-set Dykstra algorithm are defined as:
Dykstra (CAD) algorithm (14), an iterative
(A) Pc, (z(A) + (*)),
p (k++1)= I(k) +p(k) yit))
(2.1)
x (k+1) = Pcs (v (k) + (k))
(++1) y( + (()) - (b+1).
with initial values x(1)) = x and p()) = (2) O. Assuming is non-empty, the sequence points (K) defined by (2.1 is guaranteed to converge to Pa(x)) [16 34
variation of Dykstra's s algorithm can be extended to m > 2 by using a product space formulation
38 to reduce the constraint setting to the two constraint case (2.1. referred to the Simultane-
ous Dykstra Algorithm which holds for any m E N:
x (k+1)
Σ 1-1
Pc. (k) + p,/h),
(2.2)
(k+1) p
(k) + pi (k)
Pc, (x + (1),,, I6iEm.
pm O. The Simultaneous Dykstra Algorithm is well- suited for GPU implementations as it only relies on simple vector operations and allows for parallel computation of the projections Pc, (x) p (.)) However, the convergence of (2.2) can be slow for sparse problems, as for a given constraint Cj only a small number of components of (A() + p. and Pc, (x) + p; (k)) differ when NN < n, and because this, many components of the iterates Ia only slowly change. To overcome the slow convergence of (2.2) instead of averaging over all constraints as in (2.2). averaging over constraints relevant for a given component has been proposed in [14 15 Then, for any N, the CAD algorithm yields
where (1) x and p
(k+1) 1 Σ
(Pe, (2rl)) + +p!h)))
;
4 iEL;
(2.3)
(k+1) (h) + p Pa. (x. + p (K)) p
1 m.
To the best of our knowledge, no convergence results are available for 2..3
2.3.1 The linear CAD algorithm
Dykstra-type algorithms assume that individual set projections Po, (z) can be directly computed. This is possible for many classes of constraint sets such linear subspaces, half-spaces, and second- order cones [9. Even when Pc: (x) is not known in closed-form, we can use approximate hyperplane projections [1O]. For linear constraints of the form C, = {x R A&x b,). these individual projections are given by
Po, (x) = x + min 0, b, AAT
AT.
(2.4)
||A:||

---

We refer to the CAD algorithm (2.3) with projections defined by (2.4) as the linear CAD algo- rithm [37]. One the key advantages of linear CAD that it naturally lends itself GPU imple- mentations as the sparse multiplications in (2.3) and 2.4) can be efficiently computed using GPU scattering operations [35). This makes the linear CAD algorithm significantly more scalable than traditional optimisers for computing projections onto convex polytopes.
2.4 Sparse vector clipping
Given a feasible point E C, we propose a sparse vector clipping layer that transforms 2 into a leamned output y C by leveraging the sparsity of Gc. We start by learning direction vector G 2) = v E R" first. For each constraint C, we want to determine the maximum scaling factor ac; that + ac,v E C2, i.e. ac; max i 0: ac; t C.). Note that only nodes j E N; are relevant for constraint Ci and as such only (zy) Jen), and (U)))EN, need to be considered in practice when finding ocC,- Computing each factor ac; tractable for many types of convex sets C₁ and for many such scaling factors, closed-form equations are provided in (59].
unconstrained
The overall scaling factor is the minimum over the individual constraints aC ,
The standard vector clipping scheme ensures feasibility by computing out- put y min acte E C, where min(1, acte is the clipped direction. As C' is convex, any point y E C' can be reached from any feasible point & e C' using the above construction. However, when is large, ac can become very small, limiting how far the the output y move from and thus significantly reducing expressiveness.
To address this, we propose to exploit constraint independence. Let P be the partition of the con- straints into independent sets, as defined Section 2.2 For each component p E P, we compute a local scaling factor a₂ miniep OC,, This allows us to compute separate clipped directions for each independently constrained set of variables. More precisely, for each E P, we compute output variables y/ 2j + min(1, ap) Vy for all N, where E p. This definition is valid since, by construction, no two components of P constrain the same variable.
By applying scaling locally within each independent constraint group, this approach leverages the sparsity of C to allow more expressive model outputs. The resulting vector y E remains feasible, and the method is fully differentiable with respect to both the initial point and direction vector v.
To summarise, a sparse vector clipping layer y SVCo(Go;2) procedes follows:
1. Compute unconstrained direction GNN e R
2. For each constraint i, compute individual factors ac,, only depending on (2y)JeN, and (uy)JeN.
For each connected component p P. compute factor Opp mingep OCC,-
4. Return output vector C, where for each component p P and all constraints i p. the variables j E N; are computed as = 2j + min(1, app U/g-
We can stack multiple sparse vector clipping layers to iteratively refine an input 2 E C. For this, let (0) 2 E C denote the input of the first layer and iteratively define (*{++) SVC,, (Go; a: (x()), this leads to a set + feasible points z(U), ,.... 2()) C. Notably, subsequent feasible points (z())) are significantly less computationally expensive than the initial projection 200), which is obtained using the CAD algorithm.
2.5 Network architecture
A common technique for enforcing output constraints in machine learning involves computing an initial (unconstrained) output and then projecting it onto the feasible set C. As shown in (47), this projection-bast bassed strategy is a universal approximator for constrained functions. Our approach, called ProjNet, follows this principle by applying the CAD algorithm solution & C' to the constrained problem. Since z lies on the boundary of C' when w & C. we propose improve model expressiveness by making the interior of C' more accessible using sparse vector clipping layers, introduced in the previous section.
which results in a feasible
Given constraint input graph Go input, ProjNet, visualised in Figure consists of three steps:
5

---

1. GNN: Compute an unconstrained output (9c) E R" using GNN over the con- strained graph Gc-
2. CAD Projection: Compute (o(w) using the GPU-accelerated CAD algorithm.
Sparse Vector Clipping: Compute a sequence of feasible points z(0), .4,)..(1) Cusing / sparse vector clipping layers, starting from 20) with model output z y-
ProjNet fully GPU-accelerated, enabling efficient processing of large-scale inputs Gc. Both the CAD algorithm and sparse vector clipping layers are differentiable with respect to their input, mak- ing ProjNet fully end-to-end differentiable allowing to learn y using standard backpropagation. We will further show in Section 3.3 that it is possible to compute surrogate gradient for projections without resorting unrolling procedures.
C
g
W
9c
(1)
GNN
0
w
(2)
CAD Projection
C x3)
(1)
Vector Clipping
yo 6
y
Figure 1: An illustration of the ProjNet architecture. Shows forward pass on constraint polygon where points are colour coded with the modules used to compute them. Model takes as input graph with n nodes and a set of linear constraints c R", and outputs a feasible point y C.
3 Analysis of the Component-averagee Dykstra (CAD) algorithm
3.1 Convergence of the CAD algorithm
While both the two-set (2.1 and simultaneous 2.2 Dykstra algorithms converge Pa(x), show that the CAD algorithm (2.3) converges to a certain non-orthogonal projection determined by the sparsity structure of C. We show this by using sparse product space formulation, demonstrating that (2.3) can be regarded as a special case of (2.I) under a certain transformation.
Theorem 1. For input x E IR"" the CAD algorithm (2.3) converges to the projection Pt((x) arg min yeo Σj= 1/11 - xT))²
particular, for imput point (fz/ /T)ienn and feasible set a C), the CAD algorithm converges to ((Pe(z)),/ /
((15/VL)
use the CAD algorithm for computing the orthogonal projection Po(x), scale by its inputs by 1/ and rescale the outputs by Vly. We incorporate this in our GPU implementation of CAD outlined in Appendix D.2
We prove Theorem Appendix A.T
Computational efficiency
We compare the linear CAD algorithm's 's runtime against Gurobi [3]] popular commercial solver for computing the projection Pc(xx) on polytopes C. We used the dual barrier method for Gurobi as we found it to be the fastest amongst the available options for the class of problems considered. Our results in Figure demonstrate that the CAD algorithm can solve projection problems up to two orders of magnitude faster than Gurobi when the both the numerical tolerance is relaxed and

---

the initial point is close to C. However, this advantage diminishes significantly when either of these conditions not met. To encourage the unconstrained vector w in ProjNet be closer the feasible set, we add an additional penalty term -| Pc (c) 2 to the loss function during training. The hyperparameter e provides a trade-off between speed and performance: increasing c speeds up the CAD algorithm by moving w closer to C, but may reduced task performance.
10
Gurett (i- 381
CAo (o-38)
CAD [0+3)
30F
305
Variable
Variables
Variaddes
10
10
Figure 2: Comparing runtimes of CAD algorithm and Gurobi for the linear projection problem as a function of output dimension n. å values in the legend are a measure the distance between the initial point and the feasible set. Legend and y-axis are identical for all figures.
3.3 Surrogate gradients
To embed the CAD algorithm into ML models, its derivative required during the backward pass which -according to Theorem[]- is identical to the derivative of Pc. In this subsection, we first present the exact gradient of Pc before introducing surrogate gradient that is both computationally cheaper and more effective for optimization compared to the exact gradient. Surrogate gradients have proven useful in scenarios where the exact gradient either doesn't 't exist or is unsuitable for optimisation, particularly in spiking neural networks [32] 49) and arright-through estimators (7).
For fixed linear constraints C, let H₂ TR be the of all points orthogonal to Pa(a:). Define the tangent cone (53) at Pc(x), denoted by Tμ = IR", as the set of all directions v E IR" that can be approached by moving toward points in C from More formally, is in T₂ if there exist sequences (I:)ieN C and (Ti)ae c R., such that x. x, Ti and ZimE U.
Since Pe the projection onto a polytope, is piecewise-afe affine and therefore differentiable almost everywhere. Its Jacobian eg Pc is given by the linear projection Pr,,H, (26) which can be computed using the CAD algorithm for instance. However, this gradient not appropriate for training MI. models, as it can be expensive to compute and often has low rank, which may lead to poor local optima during optimization.
To remedy, we propose to use a surrogate gradient by projecting onto H2 only, namely
Pc Pf. [
de dr
C, C,
(3.1)
where d, (x - Pclx)l and denotes the identity operator.
Proposition [ Istates some desirable properties this surrogate gradient, proved in Appendix A.2 Proposition 1. Let C = IR" be polytope defined by non-redundant linear inequalities and let E TRn be a point where Pg differentiable. Then the surrogate Jacobian PıI,, defined in 3.1 satisfies the following properties:
Rank guarantee: rank( (ank(PH,)2 PH,)2n while no such bounds exist for (ank(C,,P)).
2. Exactness: PH. PC
x int(C), or Pe(a) int((C.) for all but half-space C,-
Alignment with descent direction: For any v E R", (C(Pc( PH,(v) (v)) (o))1 indi- cating that PH. (v) aligns with the descent direction of the true projected gradient.
4. Local equivalence of gradient steps: There exists 3 0 such that Polx è Pc())) Pelx PHr, (v)) for all e R² with ||l||| < B. This indicates that grudient-descenn steps are identical for surrogate and exact gradients for sufficiently small step-size.
7

---

Table Average objective values and their standard deviation for each application problem. For linear programming, optimality is measured instead. Last row combines programming methods. Data obtained for different training seeds and problem instances.
Type
Method
Linear
Quadratic (ER)
Quadratic (BA)
Transmit power
ProjNet (cs- 0)
0.9973 + 0.0001
6.83 + 1.47
7.73 1.52
3.01 0.23
ProjNet (Ch- = 0.01)
0.9905 + 0.0004
6.78 + 1.31
7.86 + 1.59
2.82 + 0.36
ProjNet (ch- 1)
0.9632 = 0.0015
6.75 1.50
7.67 1.38
1.79 0.31
ML
No SVC (c₁- = 0)
N/A
6.76 1.67
7.72 1.47
2.92 0.26
No CAD
0.2954 = 0.7644
1.93 1.33
1.91 1.16
1.79 0.30
No SVC or CAD
0.3365 + 0.2662
+ +.29
1.43 1.04
1.71 + 0.30
trust-constr
N/A
6.86 + 1.51
7.98
1.68
2.15 + 0.59
Classic
LP DCP /
1.0000 + 0.0000
6.38 - 1.80
7.60 1.70
2.79 + 0.61
Random
0.1585 = 3.6216
0.05 + 1.07
0.00 + 0.80
1.59 + 0.27
Numerical experiments in Appendix T demonstrate that surrogate gradient both faster and more stable than the exact gradient when C' is a polytope. While beyond the scope of this work, we conjecture that the benefits of the surrogate gradient (3.I) extend to general convex constraints C.
Numerical results
To assess the performance of the ProjNet model, test it on four classes of linearly constrained optimisation problems. For comparison, we employ both classical optimization- which detailed for each specific task, well as ML baselines in the form of ablations of the two main components of ProjNet: sparse vector clipping (SVC) and CAD. For reference, we also included method which samples random point in C using the hit-and-run algorithm 43 All methods considered are guaranteed to produce feasible outputs y C.
methods,
The problem instances used in evaluation were randomly generated using different graph distri- butions, as outlined in Appendix as function of problem size.
Figures 3 and we report the runtime for selected methods
Average objective values are shown in Table For linear programming, we present the ratio of achieved-to-oge to-timma objective values, while for the other three applications, we report average ob- jective values. Consequently, the standard deviations for linear programming appear smaller, since they are normalized by each instance's optimal value obtained from an LP solver.
4.1 Linear programming
Linear programs (LPs) are perhaps the simplest non-trivial class of convex problems. They are
well-studied class problems with very broad applications and some highly optimised solvers (28).
Most algorithms for LPs focus on obtaining numerically precise, exact solutions. Such methods are
usually CPU-based they utilise sparse linear algebra computations that not easily accelerated
using GPUs. Being reliant on CPU computation makes such classical methods less scalable to large LPs, for this reason there has been some recent interest in using first-order methods solve LPs (5)
644. but state-of-the-art solvers remain mostly rooted in CPU computation. We consider the problem of finding fast, approximate solutions large-scale LPs of the form max cT
x, s.t. Aa b.
(4.1)
Classical baselines We consider three LP solvers: Gurobi, a widely used commercial solver; PDLP (5). a recent first-order method for large-scale LPs; and PDLP-GPU, its GPU-accelerated implementation. We use the gurobipy interface for Gurobi, OR-Tools [51] for PDLP, and CVXPY 20 for PDLP-GPU, which interfaces with NVIDIA's cuOpt solver.

---

Prodict Sowett
FOLP
Prodiet Iaw
1014
107
187
30
300
10
10
10
Figure 3: Comparing runtime of ProjNet, PDLP, and Gurobi for linear programming. We trained three ProjNet models with different e values shown in legend. Error bars show upper/llwer quar- tiles for each point. Legend and y-axis are identical for all figures.
4.2 Non-convex quadratic programming
We consider quadratic programs of the form
max x' QQ cTx, s.t. Ar < b.
(4.2)
If Q is negative-definite then the problem is convex and can be solved efficiently using interior- point solvers [9]. Alternatively, Dykstra-type algorithms also be directly used solve convex quadratic problems (38), although this is not a popular approach due to its slower convergence. We will instead be considering the general, non-convex case in which is an arbitrary symmetric matrix.
Quadratic problems of the form 4.2 be modelled as weighted graphs g (c. Q). where the linear coefficients e are treated as node features and the matrix Q defines undirected edge weights. To evaluate our approach, we tested two types of graph topologies for G: Erdōs-Rényi (ER) 25 and Barabási-Albers (BA) random graphs.
Classical baselines We tested Difference of convex programming (DCP) method, which yields approximate solutions to (4.2) by solving a sequence of approximating convex quadratic programs [42. and trust-constr, trust region method which yields feasible approximate solutions [3] imple- mented in the minimmie package 660
DCP
DCP comtr Profict to.1)
frait-coratr Trojict 13
trudf-cnattr
Projiett (n.
Veriablka
Vertabbea
Veriebles
Quadratic (ER)
Quadratic (BA)
Transmit power
Figure 4: Comparing runtime of ProjNet, trust-constr, and programming baselines for three classes of optimisation problems. We show two ProjNet models with different values for ey. Precision is set t 10-3 Plots include error bars showing upper/lower quartiles for each point. Y-axis is shared.
4.3 Transmit power optimisation
Given n radio transmitters, let Hl E R2 channel gain matrix and let or > denote background noise. We study the problem of optimising vector of radio transmit powers & Pmax such that the sum channel capacities
Harn G = log (1 Ejmr Hi +0

---

is maximised, subject to individual minimum capacity requirements s; > O:
max
Ci
s.t. 9
C 0 24 & Pm
(4.3)
ZER
This is a linearly constrained problem with non-convex objective which can be shown to be NP- hard in general [18). We focus on feasible approximate solutions to (4.3
Classical baselines We tested trust-constr 13
and Fractional programming (FP) 57 which
uses quadratic transform that yields sequence of convex optimisation problems whose solutions converge to a local minimum to (4.3
4.4 Discussion
As result of its GPU-acceleration, ProjNet demonstrates clear advantage runtime over elassical optimization methods on large-scale problems. When it comes to solution quality, ProjNet is very competitive with classical methods while outperforming them for transmit power optimization.
As anticipated Section 3.2 the hyperparameter Ch. governs the trade-off between accuracy and speed, enabling to tailor ProjNet's performance specific application requirements. Ablation baselines confirm that both the SVC and CAD components offer a meaningful increase in perfor- mance over standard alternatives, with CAD providing the most substantial performance gains.
5 Conclusion
We presented ProjNet, a GNN architecture designed to solve constraint satisfaction problems in- volving large-scale, input-dependent linear constraints. ProjNet leverages the CAD algorithm, for which we present convergence result, a GPU implementation, and an efficient surrogate gradient for training ML models that embed the CAD algorithm.
We evaluated the ProjNet architeturre on four classes of linearly constrained optimisation problems. In comparison to both classical and ML baselines, we demonstrate that ProjNet offers strong and tunable trade-off between computational efficiency and solution quality.
6
Acknowledgments
We gratefully acknowledge the support the EPSRC Programme Grant EP/VO26259/1, "The math- ematics of Deep Learning'.
References
[1] Akshay Agrawal et al. "Differentiable convex optimization layers". In: Advances in Neural Information Processing Systems 32 (2019).
[2] Réka Albert and Albert-Lászl6 Barabási. "Statistical mechanics of complex networks". Rev. Mod. Phys. 74 (1 2002), pp. 47-97.
[3] Brandon Amos and J. Zico Kolter. "Optnet: Differentiable optimization layer in neural networks". In: International Conference on Machine Learning. PMLR. 2017, pp. 136-145.
[4] Hartwig Anzt, Stanimire Tomov, and Jack Dongarra. "On the performance and energy effi- ciency sparse linear algebra on GPUs". In: The International Journal of High Performance Computing Applications 31.5 (2017), pp. 375-390.
[5] David Applegate al. "Practical large-scale linear programming using primal-dual hybrid gradient".. In: Advances in Neural Information Processing Systems 34 (2021), pp. 20243- 20257.
[6] Heinz H. Bauschke and Jonathan M. Borwein. "On the convergence of von Neumann's 's alter- nating projection algorithm for two sets". In: Ser-Valued Analysis (1993), pp. 185-212. [ Yoshua Bengio, Nicholas Léonard, and Aaron Courville. Estimating Propagating Gra- dients Through Stochastic Neurons for Conditional Computation. 2013. aixiv: 1308. 3432
[cs.LG] URL.: https: //ap::////xvvv org/ abs/ 1308, 3432
10

---

[8] Ernesto G Birgin and Marcos Raydan. "Robust stopping criteria for Dykstra's algorithm". In:
SIAM Journal on Scientific Computing 26.4 (2005), pp. 1414-1414.
[9] Stephen P. Boyd and Lieven Vandenberghe. Convex optimization. Cambridge university press, 2004.
[10] Lev M. Bregman et al. "Finding the projection of point onto the intersection of convex sets via projections onto half-spacess". In: Journal of Approximation Theory 124.2 (2003), pp. 194-218.
[1]] Shaked Brody, Uri Alon, and Eran Yahav. "How attentive are graph attention networks?" arXiv preprint arXiv:2105. 14491 (2021).
[12] Mathis Brosowsky al. "Sample-specific output constraints for neural networks". In: Pro- ceedings of the AAAI Conference on Artificial Intelligence. Vol. 35. 8. 2021, pp. 6812-6821. [13] Richard H Byrd, Mary E Hribar, and Jorge Nocedal. "An interior point algorithm for large- scale nonlinear programming". In: SIAM Journal on Optimization 9.4 (1999), pp. 877-900. [14] Yair Censor. "Computational acceleration of projection algorithms for the linear best approx- imation problem". In: Linear Algebra and its Applications 416.1 (2006), pp. 111-123.
[15] Yair Censor, Dan Gordon, and Rachel Gordon. "Component averaging: An efficient iterative parallel algorithm for large and sparse unstructured problems". In: Parallel computing 27.6 (2001), pp. 777-808.
[16] Yair Censor and Maroun Zaknoon. "Algorithms and convergence results of projection meth- ods for inconsistent feasibility problems: review". In: Pure and Applied Functional Analysis
3.4 (2018), pp. 565-586.
[17] Ziang Chen et al. "On Representing Linear Programs by Graph Neural Networks". In: The Eleventh International Conference on Learning Representations. 2023.
[18] Mung Chiang et al. "Power control by geometric programming". In: IEEE Transactions on Wireless Communications 6.7 (2007), pp. (640-26551.
[19] Marco Cuturi. "Sinkhorn distances: Lightspeed computation of optimal transport". In: Ad- vances in neural information processing systems 26 (2013).
[20] Steven Diamond and Stephen Boyd. "CVXPY: A Python-embeddeed modeling language for convex optimization". In: Journal of Machine Learning Research 17.83 (2016), pp. 1-5.
[21] Meng Ding and Guoliang Fan. "Multilayer joint gait-pose manifolds for human gait motion modeling". In: IEFF Transactions on Cybernetics 45.1 (2014), pp. 2413-2424.
[22] Mark Eisen and Alejandro Ribeiro. "Optimal wireless resource allocation with random edge graph neural networks". In: IEEE Transactions on Signal Processing 68 (2020), pp. 2977- 2991.
[23] Joseph M Elble, Nikolaos V Sahinidis, and Panagiotis Vouzis. "GPU computing with Kacz- marz's and other iterative algorithms for linear systems". In: Parallel computing 36.5-6 (2010), pp. 215-231.
[24] Patrick Emami and Sanjay Ranka. "Learning permutations with sinkhorn policy gradient". (2018).
In: arXiv preprint
[25] Paul Erdös and Alfréd Rényi. "On Random Graphs I", In: Publicationes Mathematicae De- brecen 6 (1959), p. 290.
[26] Francisco Facchinei. Finite fiitte-dimessional variational inequalities and complementarity prob- lems. Springer Series in Operations Research and Financial Engineering. Springer, 2003.
[27] Ferdinando Fioretto et al. "Lagrangian duality for constrained deep learning". In: Joint Euro- pean conference on machine learning and knowledge discovery in databases. Springer. 2020, pp. 118-135.
[28] Saul I Gass. Linear programming: methods and applications. Courier Corporation, 2003.
[29] Maxime Gasse et al. "Exact combinatorial optimization with graph convolutional neural net-
works". In: Advances in neural information processing systems 32 (2019).
[30] Paul J. Goulart and Yuwen Chen. Clarabel: An interior-point solver for conic programs with quadratic objectives. 2024. arXiv: 2405. 12762 [math. OC]
URL: https: /amps://arrae aalaai.v.oor......
[31] Gurobi Optimization, LLC. Gurobi Optimizer Reference Manual. 2024. URL: www.
gurobi.com

---

[32] Julia Gygax and Friedemann Zenke. Elucidating the theoretical underpinnings of surrogate
gradient learning in spiking neural networks. 2024, arXiv: 2404. 14964 JCS.NEE URL:
https: /ttps://arxiv. org/ abs/ .200/.14444.4444444
[33] Qingyu Han et al. "A gnn-guided predict-and-search framework for mixed-integer linear pro- gramming". In: arXiv preprint urXiv:2302.05633 (2023).
[34] Shih-Ping Han. "A successive projection method". In: Mathematical Programming 40.1
(1988), pp. 1-14.
[35] Bingsheng He al. "Efficient gather and scatter operations on graphics processors". In: Pro-
ceedings of the 2007 ACM/IEEE Conference on Supercomputing, 2007, pp. 1-12.
[36] Alan J Hoffman. "On approximate solutions of systems of linear inequalities". In: Selected Papers Of Alan J Hoffman: With Commentary. World Scientific, 2003, pp. 174-176.
[37] Alfredo N. Jusem and Alvaro De Pierro. "A simultaneous iterative method for comput- ing projections on polyhedra". In: SIAM Journal on Control and Optimization 25.1 (1987), pp. 231-243.
[38] Alfredo N. Iusem and Alvaro R. De Pierro. "On the convergence of Han's method for con- vex programming with quadratic objective". In: Mathematical Programming 52.1 (1991), pp. 265-284.
[39] Andrei V. Konstantinov and Lev V. Utkin. "A new computationally simple approach for implementing neural networks with output hard constraints". In: Doklady Mathematics. Vol. 108. Suppl 2. Springer. 2023, S233-S241.
[40] James Kotary and Ferdinando Fioretto. Learning Constrained Optimization with Deep Aug- mented Lagrangian Methods. 2024. arXiv: 2403. 03454 [cs LG] URL: https: / //axiv..
org/ abs/ /bs/a//...... 2403. 2403.034544
[41] Yann LeCun, Yoshua Bengio, and Geoffrey Hinton. "Deep learning". In: Nature 521.7553
(2015), pp. 436-444.
[42] Thomas Lipp and Stephen Boyd. "Variations and extension of the convex-concall proce- dure". In: Optimization and Engineering 17 (2016), pp. 263-2877
[43] László Lovász and Santosh Vempala. "Hit-and-run is fast and fun". In: preprint, Microsoft Research (2003).
[44] Giulia Luise et al. "Differential properties of sinkhorn approximation for learning with wasserstein distance". In: Advances in Neural Information Processing Systems 31 (2018).
[45] Martín Abadi al. TensorFlow: Large-Scale Machine Learning on Heterogeneous Systems. Software available from tensorflow.org. 2015. URL: htttss://wwt tensorf org/
[46] Giulio Mattera, Alessandra Caggiano, and I Luigi Nele. "Optimal data-driven control of man- ufacturing processes using reinforcement learning: an application to wire are additive manu- facturing". In: Journal of Intelligent Manufacturing (2024), pp. 1-20.
[47] Youngjae Min, Anoopkumar Sonar, and Navid Azizan. Hard- Constrained Neural Networks with Universal Approximation Guarantees. 2024. arXiv: 2410 10807 [cs. LG] URL:
https: /taps://arxiv... org/ /a////////// aaaa//aaa///0//// 27410.10007 10807
[48] Boran Morvaj, Ralph Evins, and Jan Carmeliet. "Optimization framework for distributed energy systems with integrated electrical grid constraints". In: Applied energy 171 (2016), pp. 296-313.
[49] Emre Neftci, Hesham Mostafa, and Friedemann Zenke "Surrogate gradient learning in spiking neural networks: Bringing the power of muddient-based optimization to spiking neural networks". In: IEEE Signal Processing Magazine 36.6 (2019), pp. 51-63.
[50] Adam Paszke et al. PyTorch: An Imperative Style, Perfo-Performance Deep Learning Library.
2019. arxiv: 1912. 01703 [G]] URL: /tttss://spriiir org/ abs/ 1912. 01703
[51] Laurent Perron and Vincent Furnon. OR-Tools Version 19.11. Google, May 7, 2024. URL:
https://developeee gow.oo..............ooooooooooo. com htmm//ptmmizaaiii.
[52] Chendi Qian and Christopher Morris. "Towards graph neural networks for provably solving convex optimization problems". In: arXiv preprini arXiv:2502. 02400 (2025).
[53] Ralph Tyrrell Rockafellar. "Clarke's 's tangent cones and the boundaries of closed sets in R="
In: Nonlinear Analysis: theory, methods and applications 3 (1979), pp. 1154
[54] Hyunwoo Ryu et al. "Equivariant descriptor fields: Se (3)-equivariant energy-based models for end-to-end visual robotic manipulation learning". In: arXiv preprint orXiv22206.0832/
(2022).
12

---

[55] Rodrigo Santa Cruz et al. "Visual permutation learning". In: JEEE Transactions on Pattern Analysis and Machine Intelligence 41.12 (2018), pp. 3100-3114. Roland Schwan et al. "PIQP: Proximal Interior-Poins Quadratic Programming Solver".
[56]
2023 62nd IEEE Conference Decision and Control (CDC). 2023, pp. 1088-1093. DOE:
10. N0.1109/CDC49755
[57]
Kaiming Shen and Wei Yu. "Fractional programming for communication systems Part I: Power control and beamforming". In: IEEE Transactions on Signal Processing 66.10 (2018), pp. 2616-2630.
[58] B. Stellato et al. "OSQP: operator splitting solver for quadratic programs". In: Mathemal- ical Programming Computation 12.4 (2020), pp. 637-672. DOI: 10.
020
00179-2
URL: https://doi. 00g/10.
s0lss2-2222220 020- 00179-2
[59] Jesus Tordesillas, Jonathan How, and Marco Hutter. "Rayen: Imposition of hard convex constraints on neural networks". In: arXiv preprint TrXiv:2307.08322 (2023).
[60] Pauli Virtanen et al. "SciPy 1.0: Fundamental Algorithms for Scientific Computing in Python". In: Nature Methods 17 (2020), pp- 261-272.
[61] Runzhong Wang et al. "LinSATNet: the positive linear satisfiability neural networks". In: International Conference on Machine Learning. PMLR. 2023, pp. 36605-3666
[62]
Yangzihao Wang al. "Gunrock: A high-performance graph processing library on the GPU". In: Proceedings of the 2t ACM SIGPLAN symposium on principles and practice of parallel programming. 2016, pp. 1-12.
[63] Andrea Maria Zanchettin and Paolo Rocco. "Motion planning for robotic manipulators using robust constrained control". R In: Control Engineering Practice 59 (2017), pp. 127-136.
[64] Mingqiang Zhu and Tony Chan. "An efficient primal-dual hybrid gradient algorithm for total variation image restoration". In: UCLA CAM Report 34.2 (2008).
13

---

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: The effectiveness of the proposed method is demonstrated through combination of numerical and theoretical results.
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 contri- butions made in the paper and important assumptions and limitations. A No or NA answer to this question will not be perceived well by the reviewers.
The claims made should match theoretical and experimental results, and reflect how much the results can be expected to generalise to other settings.
It is fine to include aspirational goals as motivation as long as it 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:
Limitations are shown using experimental results in Sections and 3.2
Guidelines:
The answer NA means that the paper has no limitation while the answer No means that the paper has limitations, but those are not discussed in the paper.
The authors are encouraged to create a separate "Limitations" section their paper.
The paper should point out any strong assumptions and how robust the results are to vi- olations 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 on the factors that influence the performance of the approach. For example, a facial recognition algorithm may perform poorly when image resolution is low images taken in low lighting. Or speech-to-text system might not 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 to address problems of privacy and fairness.
While the authors might fear that complete honesty about limitations might be used by reviewers grounds for rejection, worse outcome might be that reviewers discover lim- itations that aren't I acknowledged in the paper. The authors should use their best judgment and recognize that individual actions in favor of transparency play important role in de- veloping 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 a complete (and correct) proof?
Answer: [es]]
14

---

Justification:
Theoretical claims were presented in the main paper and are proved in Appendix
A.1
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 or 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 provide short proof sketch to provide intuition.
Inversely, any informal proof provided in the core of the paper should be complemented by formal proofs provided in appendix or supplemental material.
Theorems and Lemmas that the proof relies upon should be properly referenced.
Experimental Result Reproducibility
Question: Does the paper fully disclose all the information needed reproduce the main exper- imental results the paper to the extent that it affects the main claims and/or conclusions the paper (regardless of whether the code and data are provided not)?
Answer: [Yes]
Justification:
Implementation described in Appendix
D
Guidelines:
The answer NA means that the paper does not include experiments.
If the paper includes experiments, No answer to this question will not be 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 to make their results reproducible verifiable.
Depending on the contribution, reproducibility can be accomplished in various ways. For example, if the contribution is novel architecture, describing the architecture fully might suffice, if the contribution essary to either make it possible for others to replicate the model with the same dataset, or provide access the model. In general. releasing code and data is 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 a large language model), releasing of a model checkpoint, or other means that are appropriate to the research performed.
specific model and empirical evaluation, may be nec-
While NeurIPS does not require releasing code, the conference does require all submissions to provide some reasonable avenue for reproducibility, which may depend 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) If the contribution is primarily new model architecture, the paper should describe the architecture clearly and fully.
(c) the contribution either be a way to access this model for reproducing the results way reproduce the model (e.g., with an open-source dataset or instructions for how to construct the dataset).
new model (e.g., large language model), then there should
(d) We recognize that reproducibility may be tricky in some cases. in which case authors are welcome describe the particular way they provide for reprodueibility. In the case of closed-source models, it 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 to reproducing verifying the results.
5. Open access to data and code
Question: Does the paper provide open access the data and code, with sufficient instructions to faithfully reproduce the main experimental results, as described in supplemental material?
15

---

Answer:
[No]
Justification: While we are currently not able to provide code, we are hoping to eventually release both code and data.
Guidelines:
The answer NA means that paper does not include experiments requiring code.
Please see the NeurlPS code and data submission guidelines lattps://nine cc/ for more details.
While we encourage the release of code and data, we understand that this might not be pos- sible, 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: Seaiii/guuidee/oeeeeeeiee http///ooSeuuiisseeeeeee sww.uiiiiiiioooccccc.c... for more details.
/
The authors should provide instructions on data access and preparation, including how to access the raw data, preprocessed data, intermediate data, and generated data, etc.
The authors should provide scripts reproduce all experimental results for the new pro- posed method and baselines. If only a subset of experiments are reproducible, they should state which ones are omitted from the script and why.
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 pa- per) is recommended, but including URLs to data and code is permitted.
Experimental Setting/Details
Question: Does the paper specify all the training and test details (e.g., data splits, hyperparame- ters, how they were chosen, type of optimizer, etc.) necessary to understand the results?
Answer: [Yes]
Justification: Discussed in appendix.
Guidelines:
The answer NA means that the paper does not include experiments.
The experimental setting should be presented in the core of the paper to a 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 supplemental material.
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: Our experimental results are equipped with error bars. For instance, the error bars in Figures and and represent upper and lower quartiles.
Guidelines:
The answer NA means that the paper does not include experiments.
The authors should answer "Yes" the results are accompanied by error bars, confidence intervals, or statistical significance tests, at least for the experiments that support the main claims of the paper.
The factors of variability that the error bars are capturing should be clearly stated (for example, train/test split, initialisation, random drawing of some parameter, overall run with given experimental conditions).
The method for calculating the error bars should be explained (closed form formula, call to library function, bootstrap, etc.)
The assumptions made should be given (e.g., Normally distributed errors).
It should be clear whether the error bar the standard deviation or the standard error the mean.
16

---

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 to show in tables figures symmetric error bars that would yield results that are out of range (e.g. negative error rates).
If error bars reported in tables plots, The authors should explain in the text how they were calculated and reference the corresponding figures or tables in the text.
8. Experiments Compute Resources
Question: For each experiment, does the paper provide sufficient information on the computer resources (type of compute workers, memory, time of execution) needed to reproduce the exper- iments?
Answer: [Yes]
Justification: Hardware is mentioned in appendix D The computation involved was light so we did not feel that detailed compute resource discussion was needed.
Guidelines:
The answer NA means that the paper does 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 of compute required for each of the individual exper- imental runs as well as 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 't make it into the paper).
Code Of Ethics
Question:
Does the research conducted in the paper conform, in every respect, with the NeurlPS Code of Ethies httppp//hherrrps..cc/p///////hct
Answer: [Yes]
Justification: We have carefully reviewed the NeurIPS Code Ethics. The research conducted in our paper conforms, in every respect, 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 de- viation from the Code of Ethics.
The authors should make sure to preserve anonymity (e.g., if there a special consideration due to laws regulations in 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: While we are aware of any immediate societal impact of this work. we demon- strated the effectiveness of method that may have some long-term social impacts.
Guidelines:
The answer NA means that there is societal impact of the work performed. If the authors answer NA or No, they should explain why their work has no societal impact 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.. de- ployment of technologies that could make decisions that unfairly impact specific groups), privacy considerations, and security considerations.
17

---

The conference expects that many papers will be foundational research and not tied to particular applications, let alone deployments. However, there is a direct path any negative applications, the authors should point it out. For example, it is legitimate to point out that an 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 (c.g., gated release models, providing defenses in addition to attacks, mecha- nisms for monitoring misuse, mechanisms to monitor how system learns from feedback over time, improving the efficiency and accessibility of ML).
11. Safeguards
Question:
Does the paper describe safeguards that have been put in place for responsible release of data or models that have high risk for misuse (e.g., pretrained language models, image generators, or scraped datasets)?
Answer: [NA]
Justification:
The paper poses no such risks.
Guidelines:
The answer NA means that the paper poses such risks.
Released models that have a high risk for misuse dual-use should released with nec- essary safeguards to allow for controlled use of the model, for example by requiring that users adhere to usage guidelines 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 is 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
Question: Are the creators original owners of assets (e.g., code, data, models), used in the paper, properly credited and are the license and terms of use explicitly mentioned and properly respected?
Answer: [Yes]
Justification: Libraries used were credited in Appendix ??.
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 is used and, if possible, include URL.
The name of the license (e.g., CC-BY 4.0) should be included for each asset.
For scraped data from particular source (e.g., website), the copyright and terms of service of that source should be provided.
If assets are released, the license, copyright information, and terms of use the package should be provided. For popular datasets, paperswi wwwprrswthhhh.....m com/ commdaaasets has curated licenses for some datasets. Their licensing guide can help determine the license a dataset.
For existing datasets that are re-packaged, both the original license and the license of the derived asset (if it has changed) should provided.
If this information is not available online, the authors are encouraged to reach out to the asset's creators.
13. New Assets
18

---

Question: Are new assets introduced in the paper well documented and is the documentation provided alongside the assets?
Answer: [NA] d
Justification: No new assets introduced.
Guidelines:
The answer NA means that the paper does not release new assets.
Researchers should communicate the details of the dataset/code/model as part of their sub- missions 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 to anonymize your assets (if applicable). You either create an anonymized URL or include an anonymized zip file.
14. Crowdsourcing 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, well as details about compensation (if any)?
Answer: [NA]
Justification:
No crowdsourcing nor research with human subjects involved.
Guidelines:
The answer NA means that the paper does not involve crowdsoureing research with human subjects.
Including this information in the supplemental material is fine, but if the main contribution of the paper involves human subjects, then as much detail as possible should be included in the main paper.
According to the NeurIPS Code of Ethics, workers involved in data collection, curation, or 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 Sub- jects
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:
No crowdsourcing research with human subjects involved.
Guidelines:
The answer NA means that the paper does not involve crowdsourcing nor research with human subjects.
Depending on 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 the paper.
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 LLM usage
Question:
Does the paper describe the usage of LLMs if 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, or formatting purposes and does not impact the core methodology, scientific rigorousness, originality of the research, declaration is not required.
19

---

Answer: [NA]
Justification:
LLMs were not used as part of this work.
Guidelines:
The answer NA means that the core method development in this research does not involve LLMs as any important, original, or non-standard components.
Please refer to our LLM policy (https::/neerrips. /Dtp:///p.pcccccoeeeeee///////// what should or should not be described.
LLMM) for
20

---

A Proofs
A.1 Convergence proof of the CAD algorithm
Proof of Theorem 7 Let (:(k), (k) denote the CAD iterates with initial point I R" and m convex
Pi
constraints C₁ c R". We show that there exists two convex sets C1, C2 c Rmx= and an initial
point X Rmxx which, when used as inputs for the two-set Dykstra algorithm 2.1 produce iterates ((k), ((k), Q Rmxxa such that for all have
(A) (*)( ((k) T
(A.1)
P
where T is certain linear bijection. We do this by using sparse product space formulation to describe the component-averaging process as a linear projection. Using this idea, we can describe both the simultaneous projections and component averages computed in (2.3) projections convex sets. The convergence result for CAD iterates (k) then follows from the convergence of
(A) to Pe
Let C x C2 R be the Cartesian product of sets C, formulated For any c = (C1, (C1,..., Cm) C₁ we interpret row E C as variable values satisfying constraint C2. Hence, C₁ encodes the problem constraints, and is the set of all possible variable values satisfying the constraints. Let X; E I" denote the i-th row of matrix X E _=x"" Then, the projection of X onto C₁ is simply the projection of its rows onto the indivisual constaint sets:
of matrices.
Pc, (X)
[Po, (X1), Pc.. (Xm) X (X.),
where we use the cross symbol to denote the stacking row vectors into a matrix.
Define the set C2 R as the collection of matrices X such that:
Xij = 0 if variable & No, i.e., j is not involved in constraint C₁,
Xij X,y for all i, such that j N₂ o Nk, i.e., shared variable j has consistent values across all constraints they appear in.
In other words, C2 encodes the sparsity and variable-sharing structure the constraints C:: only relevant variables appear in each row, and any variable shared across constraints must take the same value in those rows.
We now show that the linear projection Pc, (X) is the component average over the columns of n),
X max" For this, we define basis matrices ( e Rmxx for each variable k where each V(A) encodes which constraints involve variable k. Specifically, the entries of (() are given by:
v:(k) (k)
1, if k and & N₁,
0, otherwise.
That is, the i-th of (4) is nonzero only if constraint C₁ depends variable ks and the nonzero entry appears in column k. The subspace C2 = Rmxn has an orthonormal basis given by the scaled matrices (k))) --11 where 4 denotes the number of constraints involving variable k. The linear projection then given as
V(1
Px, =
5
V(6), x
5-1 11
VO
Σ
X.j.
l
iL,,
21

---

There is a natural correspondence between elements X = 2/==1 74 VG) C2 and vectors =
ba 7aJ" e R", defined by the linear bijection T': C2 IR"
T
Σ 73 v() =
7nl".
This map identifies the representation of a matrix in C2 with its coefficient vector R"
For an initial input x E IR", let X = (x) E C2. The projection of X onto the constraint sets
followed by projection back onto C2. and mapped to R" by T, yields:
Σ V
X Pc. Px.()) )
To Pc, o Pc, (X) = Pc,
X Pc.
(x. (xi)
VO)
1
The inner product (vU), X x- Pc, (Xi)) simplifies to:
V X (X.( = (P),( Σ (Pc,
iEL,
iEL,
Hence,
To Pc. o Pc, (X) T
Σ
v
Σ (Pa Te
- 1
Σ (Pe. (x))
Σ iE
(A.2)
Note that the right-hand side of A.2 show two equivalences AA.Ip. For k 1 both equivalences hold by definition. For the induction step, we assume they hold iteration equivalence we have:
coincides with the component-avend -vverae in (2.3). We can now
and show they also hold iteration k. For the first
(k) = To Pc, o Pc, (x( + +P(R-1)
= To Pc, (x (x-1)) + p(k-1)) + Q(h-1))
= =(((())
where the second equality follows from the fact that Q (k-11) is orthogonal to C2. To show the second equivalence, recall that for all i E {l,..., m}:
(k) x (k-1) 1) + Ps (1-1) Po; (x-11) + Pt 1)
Pt
(k) and since 4.
(k-1) x(k-1) for all j E L, we have
For all $ Ls, we have Pij
(k)
x (k 1) ((k-1) Pc, (x(+-1) P)) = p
Pi
The convergence 2.1 implies the convergence of the CAD iterates (*) to the projection
T o Pc,nc (X) = T arg min \ECC,nC,
|| -X1
=T
arg min
Σ (Y,j - -x,i)²
yeC, j=1 lEL,
arg min Σ Σ x.y2 = (x)
yEC
Finally, we derive the limit of the CAD algorithm for input (zz/ (ny//[/) //y)sejnn By scaling and C. use algorithm (2.3) compute the orthogonal projection Pa(x) instead. Given any two vectors x, y E IR" with * 0 and set C = IR"" let z/y denote element-wiss wise division and let C/v- = {a/u : E C). Setting / = (h,,.... .....,,,)) 4n)) and using inputs x/VT, C/vī in (..3 converges to Po (s/Vi) arg min (vT, x2,2 (arg || - v|(3) / VI Po(x)/Vl. yEC/vi
Therefore, we obtain the projection Po(x) using the CAD algorithm 2.3 by scaling x and C by 1/vī, and scaling the output by vi.
22

---

A.2 Properties of the surrogate Jacobian
Proof. We prove each part of the proposition separately.
1. Rank guarantee. As r- - Po(x) = for x E C, this yields rank PH.) n. For x & C, we have Po(x) * implying that the direction vector d,, z- =-Pc(() [x-Pc(=)] is well-defined. Then PH, = II I -d.d7 is the orthogonal projection onto the hyperplane normal to dz. Since d_d is rank-1 matrix, PH, has rank n
It can be seen that there are no lower bounds on the rank of the exact gradient 2g Pc. We can demonstrate this using simple counter-example. Let C IR" be the unit hyper-cube defined by the inequality x & Il. For xt = 2 Il we have (c( = which a vertex of C. The tangent cone and supporting hyperplane at are given by
Tx {x 0)
H₂ = (z Σ rj 0}. These sets intersect only at the origin, therefore (a, Pc)
Exactness.
We prove each direction of the equivalence.
Case 1: x E iC) Then Pe(x) = x, so d5 = 0 and the surrogate projection becomes
PH. I. The projection map is locally the identity, the Jacobian o₂ Pc I as well.
Hence, PH_ = =d
Case 2: Pa(xx) E (C,) for all but one half-space C. In this case, Pa(a)) lies on the boundary of exactly one of the half-spaces defining C, and in the interior of the rest. Then, locally around x, the constraint set behaves like a single active hyperplane. Therefore, the projection Pe is locally the orthogonal projection onto that hyperplane. Since the Jacobian 0.PC in this case is the projection matrix onto the corresponding hyperplane, and Hg is defined the hyperplane orthogonal
Po(x), again have PH₃ ô.fc.
Converse: If PH. Pc Pr.o.e it must be that H Tg. This occurs when no constraints are active (i.e.. x E int(C')), since then H, = Tz. If æ & C' then then H₂ is hyperplane, which could only be contained in T3 if the boundary of C Pc(x) locally of dimension n - 1. Clearly, this only occurs if only one constraint set is active, i.c. when (x) (C; for all but one
Alignment with descent direction. Let v T.. Note that d. PePc(v) E Tx Hz, where Tn the tangent cone at (x),, and H2 is the hyperplane orthogonal to Pc(x). Since Ph. is the orthogonal projection onto Hz and 0, Pc(c( E H., applying the projection preserves the
vector:
PH. (H.Po)(( = PoPo).
Therefore, we have:
(O,Pe(v),Pn, - ( =
- (0,PP(())
4. Local equivalence of gradient steps: Since Pa is piece-wise affine and differentiable at x, there exists some such that Pc is affine we have
B-neighbourhood of r. this affine region
Pe(x + y) = Po(x) PracH, (v).
for any E IR" with ||y|| < 8. Moreover, since projections are non-expansive we have
|e < ||PH, (v)ll < || < 3.
Therefore, both surrogate and exact gradient steps are the affine region Po. We then compute
Po(x + (H)()) = Po(x) + PranH, PH, (y)
Po(x) + PranH, (v) = Po(x) + PranH, (r.nH
Polx + PgPoly)).
23

---

B Generating test problems
In this section, describe the distributions of test problems which were used to train and test each of the three application problems described in Section
B.1 Generating problem objectives
Components of the linear objective terms in (4.1 and (4.2) were sampled i.i.d uniformly as o4 U)- The sparse quadratic objective term (4.2) was generated by first determining its sparsity pattern, choosing each element Qij as non-zero with probability d/n. Each non-zero element was sampled uniformly as Q U 10, 10). Finally, for the transmit power task, the channel gain matrix was generated by sampling a geometric graph to determine the sparsity pattern H then setting non-zero entries to Hij (d, 1) - where d,j is the distance between nodes and in the geometric graph. This is standard distribution to use for transmit power problems [22)
B.2 Generating linear constraints
We used sets of sparse, randomly generated linear constraints of the form Ar which we now describe. Let d denote the average degree of the adjacency matrix A E Rmxn To determine the sparsity pattern of A, we use the method described in the pervious section with degree then setting one additional random non-zero element Aij in each row to ensure no constraints are redundant The value of d determines the average degree of the corresponding bi-partite constraint graph. Then, each column A, is sampled as random unit length vector satisfying the chosen sparsity pattern. Finally, = + As where u, (0.1,, 1) and 8; ~ (0,
and
Our random linear constraints produce non-empty bounded polytopes such that ball of radius 0.1 guaranteed to fit inside C. Consequently, although the test problem distribution is appears fairly diverse, all instances are well-scaled and feasible. We leave the handling of poor scaling and infeasibility to future work.
C
Stopping conditions
Stopping conditions are an essential part of any optimisation algorithm as they ensure that the com- puted solution is within specified tolerance the true solution required iterations. We discuss stopping conditions for the CAD algorithm and outline some issues with comparing solvers which different stopping conditions.
while minimising the number of
C.1
CAD conditions
Dykstra's algorithm be stopped by using certain monotone sequence. How- ever, we decided to use the feasibility condition max Az b} & c for a given precision > 0 which appears to work well empirically. We provide brief discussion this condition and leave more detailed theoretical analysis for future work.
As shown in (
For an initial point x and CAD iterates 20(0) [[2); is strictly increasing. We are not aware of this result the literature but, assuming it does hold, we may use the feasibility stopping condition max Ax b} & since the condition is then satisfied with
(f(() we observed numerically that the sequence (||æ
J; Pc( as
- Pc(x)ll) < || vll2.
Vy CM))(()).
Moreover, Hoffman's lemma
36 guarantees the existence of an error bound the form
di I C) min (|x(0) - ull2 & emax(Ar0 6}
PEC
24

---

for some constant e 0 which depends only on C. This allows us to place lower and upper bounds on the dual solution Pc(x) as
|| - x 200 & || Po(x)]l2
& | 12 + d(x < - | caxx Ax bj.
C.2 Comparing solvers
When comparing different optimisation algorithms at given numerical precision, is often difficult equate numerical tolerance parameters due to differences stopping conditions. We only use a primal feasibility condition for the CAD algorithm while most solvers use both primal and dual feasibility conditions. Gurobi appears only use absolute values of while PDLP uses both absolute and relative c's. The comparisons we offered are thus not completely precise but we would expect the general trends hold.
D Implementation details
In this section, we discuss some practical details regarding our numerical experiments. All experi- ments were conducted single machine with an NVIDIA GeForce RTX 4090 GPU and AMD Ryzen 9 7950X3D CPU.
D.1 Computing graph connected components
For given bipartite graph with n nodes of one type and m nodes of another, we describe a simple label propagation method for computing the set P of connected components. Initializing node labels as ((2) = (1,...,m) m) and iteratively updated them as
(k+1)
min (k) JEN.
Each label 1/t is the smallest node index in the k-hop neighbourhood of node i, therefore for any two connected nodes i n!} there exists some k' > such that l2,5 (x.j for all k > k'. This sequence labels converges to a vector Imm where constraints in the same connected com- ponent share the same label. This algorithm can be efficiently implemented using GPU scattering operations, described below.
D.2
CAD implementation
Our GPU implementation of the CAD algorithm relies on the scatter function which used to compute inhomogencous sums. Efficient implementations of seatter are provided by most popular ML frameworks including Pytorch (50) and TensorFlow [45 Given an index array No and
value array V e R", scattering computes
(scatter( 1, (,V));
Σ Vj. 5:1s-i
Scattering allows us to compute both the half-space projections (.4 and component inhomogencous
averages without directly resorting to sparse matrix methods, which found to be slower.
The full linear CAD algorithm is shown in Algorithm ! The algorithm takes a sparse representation of A as input where Auww Acal denote the row and column indices of the non-zero values of A with values Ay, and includes tolerance c. The matrix is assumed contain at least one non-zero element in each row, as otherwise that row is redundant. The while loop in lines 9-15 includes GPU-accelerated CAD iterates in (2.)p while lines 3, 5, and 16 perform change variables to ensure that the algorithm converges to the projection Pa(ar) opposed to the non-orthogonal projection, based on Theorem Algorithm (T) is executed concurrently for each independent set of constraints p P with independent stopping conditions.
25

---

Algorithm The linear CAD algorithm
1: Input: Arow, Acol, Av, b,
2: scatter Acol. 1)
3: Av < Av V[w
4: Anorm v Seatter(AAme At)
5: I < v/Vi
6: Av Av/Anorm
7: - b/Ann A/Aormm
8: p
9: while max Ax b) >
do IV 2 AAccol
11: - zV +p
min(b- scatter Arow; Av. 2),
13:
-v Anw
I scatter Acal, - p)
15: end while
16: x x- x-I.V V.
17: Output:
D.3
Computational efficiency experiments
The results in Figure 2 were obtained for randomly generated initial points z where components are sampled as 2; U (-5, 6), and random linear constraints Ar & b. Random constraints where sampled as described
B, but without the offset u. This ensures that the the origin-centred ball with radius 0.1 is feasible. The parameter 8 is hence a measure of the average distance of to C. Results were gathered for different values of and the tolerance parameter Five data points where sampled for each triple (6, e, n)) tested. Besides Gurobi, also tested the Clarabel [30] and OSPQ solvers but found them to be slower for this class of problems.
58
D.4
GNN architecture
GNN message passing layers on constrained graphs Ge were computed using three separate mes- sage passes hetween the two different node types in Gc- Using the notation X - to denote bipartite graph convolution which passes features of nodes X nodes Y, we computed the follow- ing message-passes order:
1. Variables Constraints
Variables Variables
3. Constraints Variables
Each of these message-passes was performed using Graph Attention Network with dynamic atten- tion
For linear programming, we trained ProjNet models with message passing layers for GNN and no SVC layers. We didn't use any SVC layers since LPs always have solutions the boundary of their feasible set, hence the boundary bias of projections doesn't it negatively impact performance. As result, the SVC baseline wasn't relevant for LPs. For all other application problems used 8 message passing layers and 3 SVC layers. As GNN and SVC components are fairly shallow, the CAD algorithm was the most expensive element of our ProjNet models Architectures of ablation baselines were adjusted to ensure that all models had roughly equal number of learned parameters.
D.5 Ablation experiments
We provide more details regarding the three ablation models considered. These models were ob- tained by replacing the CAD projection and sparse vector clipping (SVC) component in ProjNet. As a substitute for SVC, we used standard form of the vector clipping scheme which instead directly uses the full constraint scaling value ag and as a substitute for the CAD algorithm we instead took
26

---

2 to be point in the interior C' which was computed by solving the following linear program
2 arg max(min - -
#EC
Runtime results for these ablation models were omitted from the main paper as they generally not competitive to ProjNet in terms of their accuracy-speed trade-off. We provide more complete runtime results in Figure
101
Friject (our)) Guesi
DCP
trat--aasst
Timeee
Time
FOLD
tratie (cue))
POLP dPo
Prigeet (our))
12-4
IVC CAB
10
10
280
400
600
800
1000
Varinbles
Linear programming
Quadratic (ER)
Time lieeltel
Ttreen 14
10*
trabbcostrr
20-4
Frittel (qy1)
can
10,CAD
Pripte (ra))
Varibble
Variables
Quadratic (BA)
Transmit power
Figure 5: Plots showing method runtime function of number of number of variables for all baseline methods and application problems considered. Error bars denote upper/lower quartiles. Numerical tolerance was set to 1
D.6 Difference of convex programming
The DCP method for solving non-convex quadratics is based on heuristic method outlined in (42) For problem of the form (4.2). we write Q = D where D a diagonal matrix with entries 1. Choosing D this way ensures that C' is diagonally dominant, and
Du Σ, |Qiyl where hence CC is degative-definite Let æn I", by linearising the diagonal matrix D, we can obtain heuristic solution to (4.2) by solving the following sequence of convex problems:
x (n+1) arg arg(mx( (c+ Dz! If Cz Ax & b).
(D.I)
These sub-problems were solved using the PIQP solver [56]
E Surrogate gradients experiments
In Figure 6 we provide empirical results comparing the exact and surrogate gradients from Section
3.3 for training ProjNet models. The exact gradient computed using using the CAD algorithm
with subspace constraint sets C,, which in this case reduces to the component-averaging method
presented in 23
27

---

Surrogate
Exact
6000
Average abjactive
objective
300
4000
200
2000
100
Surrogate
1000
Exact
100
200
300
400
500
100
200
300
400
500
Iterations
Iterations
Lincar programming
Non-convex quadratic programming with ER graphs
Figure 6: Training curves comparing exact and surrogate gradients to train ProjNet models for solv- ing linear and non-convex quadratic programs. Error regions denote upper/llwer quartiles. Curves are smoothed for visibility.
28


Reviews:
[{"review_id": "NDxWqXp7cP", "reviewer": "Reviewer_dghc", "summary": "The authors propose a new framework for learning GNNs that satisfy complicated convex constraints. The key idea of the framework is to develop and integrate two kinds of layers: 1) a GPU-friendly iterative algorithm, the Component-Averaged Dykstra (CAD) algorithm and 2) sparse vector clipping layers, with a base GNN. The authors demonstrate convergence results for CAD and propose a numerical scheme to compute a surrogate gradient of the CAD derivative. The method is validated on a range of convex and nonconvex constrained optimization problems.", "questions": "How does your method compare to existing ml / GNN-based optimization frameworks?", "limitations": "yes", "rating": 4, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 3, "clarity": 3, "significance": 3, "originality": 3, "strengths_and_weaknesses": "**Strengths**\n\nThe paper is well-written and the relevant works are sufficiently described. The algorithm is novel and comprehensively analyzed- including runtime and convergence properties. Numerically efficient variants of the framework are introduced and evaluated empirically. Experimentally, the framework exhibits promising results on several benchmark problems.\n\n**Weaknesses**\n\nThe comparative experiments are lacking. Although some classical methods are compared to, no comparison is made to existing ML-based optimization methodologies, e.g. those outlined in the related work section (differentiable optimization layers or other iterative methods). The comparison is primarily based on ablative experiments for the various components of ProjNet. I think it would be interesting and important to demonstrate ProjNet’s performance in comparison to these kinds of algorithms.", "was_revised": true, "final_justification": "I maintain my support of this paper. They provide a novel analysis of a neural CAD algorithm, which outperforms existing methods. Scalability is benefitted by carefully designed sparse subroutines. A follow up response from the authors has provided more experiments comparing to dense differentiable optimization layers.", "rebuttal": "We thank the reviewer for their thoughtful and constructive comments. We are pleased that you found the paper well-written and appreciated the novelty, analysis, and empirical evaluation of our method. Regarding comparisons to existing ML-based optimization frameworks, we appreciate the suggestion and clarify our rationale and design choices below. We have also revised the manuscript to better emphasize these distinctions.\n\n$\\textbf{Weaknesses:} $\n\nWe appreciate the reviewer’s suggestion to include comparisons to ML-based optimization methods. In addition to the three ML baselines already included in the paper, we actively explored a number of other potential methods.\n\n$\\textbf{Differentiable optimization layers:}$ We experimented with the qpth solver [3], a differentiable QP layer used in several prior works. However, it does not support sparse constraints as CAD does, and relies on dense matrix operations. As a result, it was not possible to scale to the problem sizes we consider. Due to these scaling issues, it was significantly slower, and for larger instances, could not fit in memory. Given these limitations, we found it unsuitable as a practical baseline, particularly since ProjNet succeeds in both speed and scalability. \n\n$\\textbf{Penalty-based ML methods:}$ motivated by penalty-based ML methods, such as [21], we also tested a penalty-based variant by isolating and training only the GNN component of ProjNet, which was trained using a soft constraint loss of the form $c_h ||w - P_C(w)||^2$.  This allows us to assess how a purely learned, penalty-driven method performs in our specific setting. Linear programming results are shown in the table below where infeasibility measures $\\max(|Ax - b|)$, and problem sizes range from $10^4$ to $10^6$, consistent with Figure 2 in the manuscript. While this method is faster, it does not guarantee feasibility. Moreover, infeasibility errors grow rapidly with problem size, and solutions violate constraints by orders of magnitude more than ProjNet. For example, at $\\epsilon = 10^{-4}$, this method is roughly 6× faster than ProjNet but up to 10,000× less accurate in feasibility.\n\n$$\n\\begin{array}{c|c c c c c c c c c}\n    \\text{Problem size} & 1 & 2 & 3 & 4 & 5 & 6 & 7 & 8 & 9 \\newline \\hline\n    \\text{Infeasibility} & 0.63 & 0.65 & 0.72 & 0.7 & 0.77 & 0.81 & 0.86 & 0.97 & 1.01 \\newline\n    \\text{Time} (s) & 0.009 & 0.010 & 0.018 & 0.030 & 0.055 & 0.102 & 0.184 & 0.335 & 0.601\\newline \\hline\n    \\text{Optimality} &&&&& 1.09 \n\\end{array}\n$$\n\n$\\textbf{Other iterative methods:}$ We also considered simple iterative algorithms for handling constraints, such as Sinkhorn's algorithm [19, 50] and related methods [56]. However, these are typically designed for specific constraint types. Sinkhorn’s algorithm applies primarily to doubly stochastic matrix constraints and is not applicable to our setting, which involves sparse, input-dependent constraints.\n\n$\\textbf{Questions:}$\n\nExisting GNN- or ML-based optimization frameworks typically fall into two categories:\n\n1. Methods that use classical solvers internally, either during training (e.g., through differentiable optimization layers) or inference. Examples include convex optimization layers [1, 3, 12], which explicitly rely on classical solvers such as interior-point methods. Since we directly benchmark against a state-of-the-art classical solver (i.e., Gurobi), our results apply to any ML method that incorporates such solvers in its pipeline.\n\n2. Inexact or penalty-based ML methods which attempt to enforce constraints through soft penalties. These methods do not guarantee feasibility and often yield solutions that violate constraints, particularly at scale. As discussed in the weakness section above, we tested such an approach using the GNN component of ProjNet trained with a penalty loss. This method was faster but significantly less accurate, with constraint violations several orders of magnitude larger than our method (see table above).\nIn contrast, our approach ensures hard constraint satisfaction while maintaining competitive runtime, and does so without reliance on classical solvers during inference. To our knowledge, no existing ML-based method matches this combination of feasibility, scalability, and efficiency in our setting."}, {"review_id": "okUFnxMT5N", "reviewer": "Reviewer_SUYk", "summary": "The paper introduces ProjNet, a novel GNN framework for enforcing input-dependent convex constraints via a CAD algorithm and sparse vector clipping. A theoretical result on the convergence of the CAD algorithm is provided, and numerical comparisons with classical baselines are presented.", "questions": "No.", "limitations": "The authors present limitations using numerical experiments.", "rating": 4, "confidence": 2, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 3, "clarity": 3, "significance": 3, "originality": 2, "strengths_and_weaknesses": "__Strengths:__\n\n1. This paper is in general well-written and is easy to follow. The main ideas are smooth and reasonable.\n\n2. A theoretical result on the convergence of the CAD algorithm is provided.\n\n3. Surrogate gradients are proposed as a computationally efficient alternative to exact gradients, with some theoretical properties.\n\n4. The numerical results look nice compared to classical baselines.\n\n__Weaknesses:__\n\n1. I do not think that projecting the neural network output onto a feasible region is a completely novel idea. There have been some existing works, for example: Wang, R., Zhang, Y., Guo, Z., Chen, T., Yang, X., & Yan, J. (2023, July). LinSATNet: the positive linear satisfiability neural networks. In International Conference on Machine Learning (pp. 36605-36625). PMLR.\n\n2. There are only numerical comparisons with classical baselines. Given the claimed scope of this paper, I would expect to see how this approach can contribute to the \"learning-to-optimize\" community. In particular, I think the authors should include comparisons with existing GNN methods for optimization. This might be a major flaw.", "was_revised": true, "final_justification": "The authors addressed my concerns about related works and ML baselines. I increase the rating to 4.", "rebuttal": "We thank the reviewer for their thoughtful comments and for acknowledging the clarity and motivation of the paper, as well as the theoretical and empirical contributions. We also appreciate the points raised regarding novelty and comparisons to existing ML-based optimization methods. We address these concerns below and have revised the manuscript to better highlight our contributions relative to prior work.\n\n$\\textbf{Weaknesses:}$\n\n$\\textbf{Novelty:}$ We agree that the idea of projecting neural network outputs onto a feasible set is not novel, and we do not claim otherwise. As the reviewer correctly notes, this concept has been explored in prior work, including [56], which we cite on page 2 of the manuscript.\n\nOur contribution lies in the method used to compute these projections. Specifically, we propose a novel approach that leverages the CAD algorithm to perform projections in a scalable manner using GPU-acceleration. This contrasts with most existing methods, which rely on classical solvers that are less amenable to GPU implementations, especially for sparse problems. The effectiveness of our approach is demonstrated by the No CAD baseline in Table 1 in the manuscript, which shows that replacing CAD with a classical solver significantly degrades performance.\n\n$\\textbf{ML baselines:}$ We thank the reviewer for raising this point. In addition to classical baselines, we included three ML-based baselines in the form of ablation studies described in Appendix D.5. These include variations of our method with replacements for the proposed CAD and SVC layers, helping isolate the impact of each component. \n\nRegarding GNN-based methods for constrained optimization, we carefully considered relevant work, including [3, 21, 54]. However, such approaches typically fall into two categories:\n\n1. Methods which rely on classical solvers internally, such as differentiable convex optimization layers [1, 3, 12]. These often use interior-point solvers during training or inference. Since we benchmark directly against state-of-the-art solvers (e.g., Gurobi), we can confidently say that our method outperforms any ML approach that embeds such solvers. Moreover, unlike CAD, these methods often lack support for both sparse constraints and GPU acceleration, rendering them unsuitable for our problem setting.\n\n2. Inexact ML methods that use penalty terms or Lagrangian-style optimization [21, 54] are typically faster but lack feasibility guarantees, and their constraint satisfaction often degrades with problem scale. To demonstrate this, we implemented a penalty-based method obtained by isolating and training only the GNN component of ProjNet, which was trained using a soft constraint loss of the form $c_h \\|w - P_C(w)\\|^2$, with $c_h=1$. This allows us to assess how a purely learned, penalty-driven method performs in our specific setting. Linear programming results are shown in the table below, where infeasibility measures $\\max(|Ax - b|)$, and problem sizes range from $10^4$ to $10^6$, consistent with Figure 2 in the paper.. While this method is faster, it does not guarantee feasibility. Compared to ProjNet, the penalty-based GNN results in constraint violations that are several orders of magnitude larger, despite only modest gains in speed. For example, at $\\epsilon = 10^{-4}$, this method is roughly 6× faster than ProjNet but up to 10,000× less accurate in feasibility. Note that optimality $>1$ as constraints are not satisfied.\n\n    $$\n    \\begin{array}{c|c c c c c c c c c}\n    \\text{Problem size} & 1 & 2 & 3 & 4 & 5 & 6 & 7 & 8 & 9 \\newline \\hline\n    \\text{Infeasibility} & 0.63 & 0.65 & 0.72 & 0.7 & 0.77 & 0.81 & 0.86 & 0.97 & 1.01 \\newline\n    \\text{Time} (s) & 0.009 & 0.010 & 0.018 & 0.030 & 0.055 & 0.102 & 0.184 & 0.335 & 0.601\\newline \\hline\n    \\text{Optimality} &&&&& 1.09 \n    \\end{array}\n    $$\n\n    To directly address the reviewer’s concern about contributing to the learning-to-optimize community: our work presents a learned projection architecture that generalizes across input-dependent constraint sets and significantly outperforms baselines in both accuracy and efficiency. We believe this contributes a novel and practical direction for integrating deep learning with constrained optimization in a principled way. Furthermore, since ProjNet can be trained using any (potentially supervised) loss function, its applicability extends beyond the learning-to-optimize domain to broader constrained machine learning tasks—for instance, the task from [50]."}, {"review_id": "K9NRCgDPT7", "reviewer": "Reviewer_ev9W", "summary": "This paper introduces ProjNet, a GNN architecture designed to enforce input-dependent constraints. To implement this architecture, Sparse Vector Clipping and Component-Averaged Dykstra are proposed for refined solution qualities. A surrogate gradient for CAD enables end-to-end training without unrolling. The authors provide a theoretical convergence guarantee for CAD and implement it with GPU acceleration. ProjNet is evaluated on several synthetic datasets involving both convex and non-convex problems", "questions": "In addition to the points raised in the weaknesses section, another concern is the omission of a potentially relevant baseline. The overall structure of ProjNet is to predict an initial solution, followed by post-processing to improve feasibility and optimality, which bears resemblance to [1] that follows a similar approach for integer programs.\nWhile [1] focuses on discrete optimization, its method can be adapted to linear problems by relaxing the integrality constraints. A brief empirical or conceptual comparison would strengthen the paper. Moreover, this raises an important question about the individual contributions of the CAD and SVC components: could the predicted solutions from ProjNet serve as warm starts for established solvers like PDLP, potentially achieving comparable performance with less architectural complexity?\n\n[1] Qingyu Han, et. al., A gnn-guided predict-and-search framework for mixed-integer linear programming. arXiv preprint arXiv:2302.05636, 2023.", "limitations": "yes", "rating": 4, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 3, "clarity": 3, "significance": 2, "originality": 3, "strengths_and_weaknesses": "**Strengths**\n1. The paper provides a convergence guarantee for the CAD algorithm and introduces a well-motivated surrogate gradient for training.\n2. The GPU-accelerated implementation and sparse methods enable ProjNet to handle large-scale, input-dependent constrained problems efficiently.\n3. ProjNet achieves surprisingly good results on non-convex optimization tasks, showing its practical utility.\n\n**Weaknesses**\n1. The main results are established on randomly generated instances, which is as convincing as using public datasets, especially considering.\n2. The performance can vary significantly by adjusting the hyperparameter $c_h$, and not much details are provided.\n3. While the paper claims a GPU-accelerated implementation of the CAD algorithm, it does not explicitly report whether the experimental results are based on the GPU version. If so, comparing against CPU-only baselines like PDLP may be unfair, as PDLP also supports GPU acceleration. Clarifying this would improve the fairness and transparency of the comparisons.", "was_revised": true, "final_justification": "I appreciate the author’s thoughtful and thorough response. All of my concerns have been adequately addressed, and I will maintain my original score.", "rebuttal": "We thank the reviewer for their thoughtful feedback. We are glad you found the surrogate gradient, convergence guarantee, and GPU implementation of CAD to be valuable contributions, and that the performance on non-convex tasks was promising. We address the questions and concerns about baselines, dataset design, hyperparameters, and GPU usage in detail below, and we have clarified these aspects in the revised manuscript.\n\n$\\textbf{Weaknesses:}$\n\n$\\textbf{Public datasets:}$ While we agree that using public datasets would be ideal, to the best of our knowledge, existing publicly available optimization datasets are typically small and likely not suitable for training general-purpose ML models. To support the development of scalable neural solvers, we generated a broad and diverse set of feasible problem instances tailored to our setting. This allowed us to control key properties (e.g., sparsity, constraint structure) while ensuring sufficient scale for training. Details of the generation process are provided in Appendix B.\n\n$\\textbf{Penalty hyperparameter:}$ The hyperparameter $c_h$ controls the strength of the penalty term $c_h \\|w - P_C(w)\\|^2$ and encourages the GNN output $w$ to lie close to the constraint set $C$, thereby accelerating convergence of the CAD projection (see Figure 6).\n\nWe tested values $c_h \\in \\\\{0, 0.01, 0.1, 1\\\\}$ and observed a consistent trend: larger values of $c_h$ led to faster CAD convergence, but could interfere with task performance if set too high. Specifically, $c_h = 0.1$ reduced runtime compared to $c_h = 0.01$, but resulted in slightly worse solution quality. For our main experiments, we found $c_h = 0.01$ to offer a good trade-off between CAD speed and task accuracy. In practice, we recommend this as a robust default. The choice of $c_h$ will depend on the objective function's scale and can be adjusted depending on the desired trade-off between efficiency and performance.\n\n$\\textbf{GPU-acceleration:}$ We thank the reviewer for raising this important point. We confirm that the experimental results reported in the paper use our GPU-accelerated implementation of the CAD algorithm.\n\nRegarding fairness, we note that most classical solvers we compare against do not have accessible or mature GPU implementations. For example, while PDLP supports GPU acceleration, we were unable to find a public implementation with a Python interface we can use. The OR-Tools package only offered a CPU implementation of PDLP. Nvidia is developing a GPU implementation, but it was still in beta and not publicly available at the time of our experiments.\n\nMoreover, traditional methods like the simplex or interior-point algorithms are challenging to parallelize efficiently due to complex control flow and sparse matrix operations, which limits practical GPU acceleration -- though this is an active area of research. In contrast, our method was specifically designed to benefit from GPU acceleration, which is central to its performance gains.\n\nTo ensure transparency, we will clarify in the revised version that CAD was run on GPU, and we acknowledge this difference in platform when interpreting comparisons with CPU-only baselines.\n\n$\\textbf{Questions:}$ \n\nWe thank the reviewer for highlighting [1] as a related approach. While [1] focuses on discrete optimization (specifically, MILPs), its structure of GNN-guided prediction followed by solution refinement does share high-level similarities with our approach. However, there are important distinctions. [1] does not offer true feasibility guarantees and still relies on classical solvers for completing the MILP solution. In contrast, our method ensures feasibility through the CAD projection and is much more general. We view our work as offering a stronger guarantee of constraint satisfaction, which is essential in many real-world applications. For linear programming, SOTA methods include simplex, interior-point, and (more recently) PDLP methods, all of which were tested. We are not aware of any competitive ML algorithms.\n\nRegarding warm-starting: These methods can work only if established solver exist for the problem considered. Meanwhile, ProjNet can be directly trained using any objective. Also, most established solvers do not support GPU acceleration, so even with a warm-start, they would likely be slow for large-scale problems -- as observed in the classical solver results we reported.\n\nWe will add a short discussion of [1] in the revised related work section to acknowledge this connection."}, {"review_id": "QbhwVrUR5i", "reviewer": "Reviewer_RXdy", "summary": "This work proposed a projection and search framework, using CAD and SVC. The architecture first ouputs an unconstrained output, then projects it to feasible region, and performs a search in the feaible region. The authors provided convergence analysis, designed a surrogate gradient for the CAD and exhibited promising experiment results.", "questions": "- Is it a supervised learning loss to predict the objective value?\n- How is the SVC part trained? I mean how is the loss designed?\n- For LPs, the numbers in the table are ratio of reaching optimal value? But for continuous optimization, how can you reach optimality? There must be some numerical error. \n- What is the runtime bottleneck of the ProjNet? The CAD or MPNN or what?\n- Each CAD step calculation seems cheap, yes, but if it takes a lot of steps then GPU cannot accelerate it.", "limitations": "The authors discussed the limitation in the paper.  \nThe frameworks applies to linear constraints, which is a limitation.", "rating": 4, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 3, "clarity": 2, "significance": 3, "originality": 3, "strengths_and_weaknesses": "Strength:\n- The paper is well motivated. Even readers with little background knowledge can understand the research direction and target.\n- The paper is well written.\n- The project then search idea is interesting and clear, and seems easily implementable in practice.\n- The practical consideration for CAD convergence and SVC step size is nice and very convincing. \n- The experiment results are promising.\n\nWeakness:\n- Some implementation details are unclear, see questions below\n- There is a concurrent work https://arxiv.org/abs/2502.02446, which also uses projection and search method for feasible solutions. Maybe consider citing it and add comparisons or discussions to their method.", "was_revised": true, "final_justification": "I think this paper is novel, well written and exhibits good results. The main concern is the bottleneck of CAD component.  \nAfter consideration on the pros and cons, as well as other reviews, I personally believe borderline accept is a fair score.", "rebuttal": "We thank the reviewer for their thoughtful and encouraging feedback. We are glad that you found the paper well-written, clearly motivated, and the projection-and-search framework both interesting and practical. We also appreciate the specific questions and suggestions regarding implementation details, comparisons to concurrent work, and limitations. We address each point in detail below and have incorporated clarifications into the revised manuscript.\n\n$\\textbf{Weaknesses}$\n\n$\\textbf{Implementation details:}$ We address your questions regarding implementation details in the Question section below.\n\n$\\textbf{Related work:}$ We thank the reviewer for bringing this concurrent work to our attention. We have reviewed the cited paper [arXiv:2502.02446] and agree that it shares some similarities with our No CAD baseline (Appendix D.5). Their approach uses an interior-point (IP) solver to obtain a feasible initial point $x_0$, which is then refined using a GNN while maintaining feasibility—an idea aligned with projection-and-search strategies.\nWhile a direct empirical comparison is not currently available, we note that ProjNet is likely to be faster. As shown in our experiments, CAD achieves significantly lower runtime than IP solvers, and our method avoids the need for a full interior-point solve at inference time. Furthermore, unlike their approach, ProjNet is trained end-to-end, and includes a dedicated surrogate gradient to improve the backward pass.\n\nWe will cite and discuss this work in the revised manuscript.\n\n$\\textbf{Questions}$\n\n$\\textbf{Learning loss:}$ Our method is trained in an unsupervised manner by directly maximizing the objective function (plus a penalty term $c_h ||w - P_C(w)||^2$), without requiring ground-truth labels. For instance, in the linear programming task, we train ProjNet to maximize the known objective $c \\cdot x$, based solely on problem data.\n\n$\\textbf{SVC training:}$ We do not design a separate or explicit loss for the SVC component. The entire ProjNet architecture, including the SVC layer, is end-to-end differentiable and trained jointly using a single overall loss function (e.g., the objective function combined with a penalty term for constraint proximity). As a result, the SVC component learns naturally during training, in the same way layers in a standard neural network (e.g., MLP) are optimized via backpropagation.\n\n$\\textbf{LP numerical error:}$ All optimality ratios reported are measured relative to the solution obtained by PDLP, which we treat as a reference point. In our experiments, PDLP was run with both absolute and relative optimality tolerances set to various values of $\\epsilon$, as detailed in the paper. While this does introduce a small numerical error, the resulting solutions are sufficiently accurate to serve as a practical proxy for the ground truth. Moreover, the optimality gaps observed for other methods are typically much larger than this numerical margin, so the comparisons remain valid and informative.\n\n$\\textbf{Bottleneck:}$ The runtime bottleneck depends on several factors, including the size of the GNN, the number of SVC layers, and the value of the penalty weight $c_h$. In our setup, we used a relatively small GNN and at most 3 SVC layers as described as in Appendix D.4. Under these conditions, the CAD projection step was the primary bottleneck. For example, in the linear programming experiments with $c_h = 1$, the GNN was approximately 3$\\times$ faster than CAD, making CAD the dominant cost component.\n\n$\\textbf{GPU acceleration:}$ The reviewer is correct that the GPU does not reduce the number of CAD iterations. The GPU implementation only accelerates each iteration. To reduce the total number of steps, we use the penalty term $c_h ||w - P_C(w)||^2$ during training to encourage the initial point $w$ to lie closer to the feasible set $C$, thereby speeding up convergence. For instance, at $c_h=0, n=10^{4}, \\epsilon=10^{-3}$ CAD needed $\\sim$15k iterations, at $c_k=1$ that number dropped to about 300.\n\n$\\textbf{Limitations}$\n\nWhile our experiments focus on linear constraints, our method extends to a broader class of convex constraints. Specifically, Theorem 1 is stated for arbitrary convex constraint sets, and the vector clipping technique we use can be (and has been, e.g., in [54]) applied to more general settings. The primary challenge lies in representing such constraints in a way that is amenable to GNNs or ML models more generally. We will clarify in the revised manuscript that our framework is not limited to linear constraints and, in principle, applies to arbitrary convex constraint sets."}]

Author rebuttal:

Meta-review: This paper introduces ProjNet, a GNN architecture designed to enforce input-dependent constraints. To implement this architecture, Sparse Vector Clipping and Component-Averaged Dykstra are proposed for refined solution qualities. The paper provides a theoretical convergence guarantee for CAD and implements it with GPU acceleration. The main concern is the bottleneck of CAD component. In summary, all the reviewers agree that borderline accept would be reasonable and fair.
Decision: Accept (poster)