Enforcing convex constraints in Graph Neural Networks

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.

$**Weaknesses**$

$**Implementation details:**$ We address your questions regarding implementation details in the Question section below.

$**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. While 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.

We will cite and discuss this work in the revised manuscript.

$**Questions**$

$**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.

$**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.

$**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.

$**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.

$**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.

$**Limitations**$

While 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.

Thank you so much for the clarification. I am curious how the parameter of CAD iterations would affect the precision. But considering the time limit, I would look forward to your final version. Besides, the projection method to ensure feasibility is not totally novel idea. I still believe the CAD part limits the applicability of this paper. As I said, I believe this paper definitely has its value, yet not conclusively should be accepted. I believe borderline accept would be a reasonable decision.

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.

$**Weaknesses:**$

$**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.

$**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).

We 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.

$**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.

Regarding 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.

Moreover, 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.

To 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.

$**Questions:**$

We 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.

Regarding 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.

We will add a short discussion of [1] in the revised related work section to acknowledge this connection.

Thank you for making clarifications to the raised concerns, but I still have the following questions:

1. I agree that existing datasets are still insufficient for training truly general-purpose models. However, the coefficient ranges generated by the current experimental setup appear to be quite limited. Suppose the ultimate goal is to develop a general-purpose neural solver. In that case, it should be capable of handling a wide variety of problems, including those with highly diverse numerical ranges, as seen in public benchmarks like QPLIB and Netlib. While I understand that this presents significant challenges for current methods, a short discussion on the proposed model’s ability or limitations in handling these cases would be interesting.
2. Regarding the fairness of comparison, I don't believe it's a sufficient justification to omit a comparison with the GPU-based PDLP. Integrating their Julia implementation with the current Python setup would require only minor modifications to the code. link: https://github.com/jinwen-yang/cuPDLP.jl

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.

$**Weaknesses:**$

$**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.

Our 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.

$**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.

Regarding GNN-based methods for constrained optimization, we carefully considered relevant work, including [3, 21, 54]. However, such approaches typically fall into two categories:

1. 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.

2. 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.

   $$\begin{array}{c|c c c c c c c c c} \text{Problem size} & 1 & 2 & 3 & 4 & 5 & 6 & 7 & 8 & 9 \newline \hline \text{Infeasibility} & 0.63 & 0.65 & 0.72 & 0.7 & 0.77 & 0.81 & 0.86 & 0.97 & 1.01 \newline \text{Time} (s) & 0.009 & 0.010 & 0.018 & 0.030 & 0.055 & 0.102 & 0.184 & 0.335 & 0.601\newline \hline \text{Optimality} &&&&& 1.09 \end{array}$$

   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].

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.

$**Weaknesses:**$

We 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.

$**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.

$**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.

$**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.

$**Questions:**$

Existing GNN- or ML-based optimization frameworks typically fall into two categories:

1. 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.

2. 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). In 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.