Decision-Focused Surrogate Modeling for Mixed-Integer Linear Optimization
Authors: Shivi Dixit
Rishabh Gupta
Qi Zhang
Venue: Accepted by TMLR ()

Abstract: Mixed-integer optimization is at the core of many online decision-making systems that demand frequent updates of decisions in real time. However, due to their combinatorial nature, mixed-integer linear programs (MILPs) can be difficult to solve, rendering them often unsuitable for time-critical online applications. To address this challenge, we develop a data-driven approach for constructing surrogate optimization models in the form of linear programs (LPs) that can be solved much more efficiently than the corresponding MILPs. We train these surrogate LPs in a decision-focused manner such that for different model inputs, they achieve the same or close to the same optimal solutions as the original MILPs. One key advantage of the proposed method is that it allows the incorporation of all of the original MILP’s linear constraints, which significantly increases the likelihood of obtaining feasible predicted solutions. Results from two computational case studies indicate that this decision-focused surrogate modeling approach is highly data-efficient and provides very accurate predictions of the optimal solutions. In these examples, the resulting surrogate LPs outperform state-of-the-art neural-network-based optimization proxies.

Published in Transactions on Machine Learning Research (03/2025)
Decision-Focused Surrogate Modeling for Mixed-Integer Linear Optimization
Shivi Dixit
dizit06, Gumm. edu
Deportment of Chemical Engineering and Materials Science, University of Minnesola Twin Cities, MN, USA.
Rishabh Gupta
rishabh, rishabh.gupta @eraonmobiil com
EczonMobil Technology and Engineering Comparey, Spring, TX, USA.
Qi Zhang
gizh Gumn edu
Department of Chemical Engineering and Materials Science, University of Minnesota Twin Cities, MN, USA.
Reviewed on OpenReview: https: // openreview. forun? id= A6t00kkkEZ
Abstract
Mixed-integer optimization is at the core of many online decision- makinn systems that demand frequent updates of decisions real time. However, due to their combinatorial nature, mixed-integer linear programs (MILPs) can be difficult solve, rendering them often unsuitable for time-critical online applications. To address this challenge, we develop data-driven approach for constructing surrogate optimization models in the form of lincar programs (LPs) that can be solved much more efficiently than the corresponding MILPs. We train these surrogate I.Ps in decision-focused manner such that for different model inputs, they achieve the same or close to the same optimal solutions as the original MILPs. One key advantage of the proposed method is that allows the incorporation all the original MILP's linear constraints, which significantly increases the likelihood of obtaining feasible predicted solutions. Results from two computational case studies indicate that this decision-focused surrogate modeling approach is highly data-efficient and provides very accurate predictions of the optimal solutions. In these examples, the resulting surrogate LPs outperform state-of-the-art neural-networke optimization proxies.
1 Introduction
Effectively operating complex systems in online environments, where input parameters constantly changing, requires efficient decision making in real time. Many online decision-making frameworks involve the solving of mathematical optimization problems; however, often the computational complexity of the given optimization problem presents a major challenge such that the long solution time renders it ineffective real-time applications. common approach to tackling this challenge to perform the online optimization with a surrogate model, which is an approximation of the original model that can be solved more efficiently (Bhosckar & Terapetritou, 2018).
Surrogate modeling for optimization typically involves three steps: (i) identify the complicating constraints in the original optimization model, (ii) generate surrogate models for the complicating constraint functions, and (iii) replace the complicating constraints with the surrogales to obtain a surrogate optimization model that can be solved more efficiently. This approach especially effective in applications where the computational complexity of the problem only due to small part of the model such that one can keep most of the original constraints. Here, a major underlying assumption is that a surrogate model that provides a good approximation of the complicating constraints will also, once incorporated into the optimization problem,

---

Published in Transactions on Machine Learning Research (03/2025)
lead to solutions that are close to the true optimal solutions. However, this assumption often does not hold. The original optimization problem and the surrogate optimization problem may achieve very different optimal solutions despite having a highly accurate (but not perfect) embedded surrogate model.
Recently, Gupta Zhang (2024) have developed a new surrogate modeling framework that explicitly aims to construct surrogate models that minimize the decision prediction error defined as the difference between the optimal solutions of the original and the surrogate optimization problems: it is hence referred to decision-focused surrogate modeling (DFSM). DFSM has been applied to construct convex decision-focused surrogate optindzation models (DFSOMs) for nonconvex nonlinear programs that difficult to solve to global optimality Gupta & Zhang (2024); however, these problems only involve continuous variables. In this work, we extend the DFSM framework to MILPs, which arise in many real-time decision-making problems such as online production scheduling and hybrid model predictive control. Due to their combinatorial nature, MILPs can take long times to solve. Here, we propose to construct DFSOMs in the form of LPs that are computationally significantly more efficient and achieve optimal solutions close the oncs of the original MIL.Ps.
The remainder of this paper is organized as follows. We first provide review of related works in Section 2, focusing on existing data-driven and parametric methods for solving online mixed-integer optimization problems. The proposed DFSM approach for MILPs is presented in Section 3, while its utility is highlighted in two computational case studies in Section 4. Finally, in Section 5, we close with some concluding remarks.
2
Related work
In this section, we provide brief review of related works. Note that we do not review the vast literature on surrogate modeling as the methods described therein are of less relevance to us given our specific focus on MILPs and decision-focusod modcling approaches.
2.1 ML-based optimization proxies
The line of research that is most closely related to this work is the constructing optimization proxies using machine learning (ML). Here, the goal is to learn direct mapping of the input parameters of an optimization problem to its optimal solution. In that sense, it has the same objective as DFSM, namely to obtain model that predicts optimal solutions close to the true ones. It uses the same offline generated data consisting of different model inputs and the corresponding true optimal solutions. The main difference is that while DFSM trains a model the form of surrogate optimization problem, optimization proxies are ML models, most commonly neural networks (NNs) (Sun al., 2018; Pan et al., 2019; Karg & Lucia, 2020; Kumar et al., 2021).
Optimization proxies take advantage of the expressive power of deep NNs; however, in deep learning, it is difficult to enforce constraints on the predictions, which often leads to predicted solutions that are infeasible in the original optimization problem. In their work on constructing NN-based optimization proxies for the AC optimal power flow (OPF) problem, Zamzam & Baker (2020) ensure feasibility by generating a training set of strictly feasible solutions through a modified AC OPF formulation and by using the natural bounds of the sigmoid activation function to enforce generation and voltage limits. When also addressing the AC OPF problem, Fioretto et al. (2020) consider constraints in their deep learning framework by augmenting the loss function with penalty terms derived from the Lagrangian dual of the original optimization problem and applying subgradient method to update the Lagrange multipliers. This approach has also been applicd to the becurity-constrained OPF problem (Velloso & Van Hentenryck, 2021) and job shop scheduling (Kotary et al., 2022). Most of these and similar approaches increase the likelihood of obtaining feasible predictions but cannot guarantce it; hence, in most cases, a fcasibility restoration step added obtain a feasible solution. One common approach correcting an infeasible prediction is project it onto a suitable set such that the projection is a feasible solution (Chen et al., 2018; Paulson & Mesbah, 2020).
2

---

Published in Transactions on Machine Learning Research (03/2025)
2.2 Decision-focused learning
We now describe a framework that is closely related to the DFSM approach but is not motivated by the need for fast online optimization. Instead, it addresses the following problem: In traditional data-driven oplimization, we often follow a two-stage wredict-lhen-optimize approach, i.e. we first predict unknown input parameters from data with external features and then solve the optimization problem with those predicted inputs. Here, the learning step focuses on minimizing the parameter estimation error; however, this does not necessarily lead to the best decisions (evaluated with the true parameter values) in the optimization step. In contrast, decision-focused learning (Wilder et al., 2019), also known as smart wediccttteen-ppting (Elmachtoub & Grigas, 2022) or predict-and- and-optimize (Mandi et al., 2020), integrates the two steps to explicitly account for the quality of the optimization solution in the learning of the model parameters.
Elmachtoub & Grigas (2022) consider decision- objective function and a closed convex feasible region. They develop a convex surrogate loss function that is Fisher consistent with respect to the original loss function and allows the use of stochastic gradient descent. to efficiently solve the problem. The same algorithm has also been applied to combinatorial problems (Mandi et al., 2020). Decision-focused learning has further motivated rescarch on diffcrentiable optimization in deep learning (Amos & Kolter, 2017), which allows NNs to comprise differentiable layers that represent full optimization problems. Using this approach, quadratic programs (Donti et al., 2017), linear programs (Wilder et al., 2019), and mixed-integer linear programs (Ferber et al., 2020) have been considered in the training of decision-focusee NNs. This approach is applied by Ferber et al. (2023) to address a problem of similar nature our DFSM problem, namely to learn linear surrogate objective functions for discrete optimization problems with lincar constraints but nonlinear objective functions.
learning for optimization problems with a linear
2.3 ML-enhanced mixed-integer optimization
There is rapidly growing body of literature using ML to speed up the solution of mixed-integer optimization problems (Bengio et al., 2021). Our work well as the ML-based optimization proxies reviewed in Section 2.1 fall into this broad category. Other strategies include using ML to learn more effective primal heuristics (Bengio et al., 2020; Shen et al., 2021), branching rules (Lodi & Zarpellon, 2017; Khalil et al., 2016), the set of active constraints at the optimal solution (Bertsimas & Stellato, 2022), and how to warm-start MILPs (Jinénez-Corderoo et al., 2022). More closely related to our work are the contributions on learning how to improve the generation of cutting planes (Deza & Khalil, 2023; Tang et al., 2020; Huang al., 2022).
However, while these methods add cuts as part of branch-and-cut algorithm, our approach learns parametric cuts that can be added immediately to the LP relaxation for given model input, discussed in more detail in Section 3.
2.4 Multiparametric programming
While DFSM and optimization proxies can generally only provide approximate optimal solutions, multipara- metric programming is an exact approach that maps model inputs directly to the true optimal solutions through explicit functions (Oberdieck et al., 2016). The solution of a multiparametrie LP follows from the basic sensitivity theorem (Fiacco & Ishizuka, 1990), which states that the active set remains unchanged in the neighborhood of the given vector of input parameters and the corresponding optimal solution. This results in the partitioning of the fcasible parameter space into polytopic regions called critical regions. The goal of multiparameuric programming is to obtain these critical regions and their corresponding optimal policies offline, then the online optimization reduces to selecting and applying the right policy for given parameter vector that lies in the corresponding critical region.
When solving multiparametric MILPs, the same concept can be applied. One naive way is to enumerate all the possible discrete solutions and then treat the multiparametric MILP as a set of many multiparametric LPs (Roodman, 1972), one for each discrete solution. Thus, the resulting solution to the multiparametric MILP is again achieved by partitioning the parameter space into critical regions, each of which is associated with an optimal policy as an affine function of the model parameters. More advanced methods apply algorithms
3

---

Published in Transactions on Machine Learning Research (03/2025)
based branch-and-bound (Ohtake & Nishida, 1985; Acevedo & Pistikopoulos, 1999) and decomposition (Dua & Pistikopoulos, 2000).
Multiparametric programming has found several successful applications in explicit model predictive control (Alessio Bemporad, 2009; Pistikopoulos et al., 2015) with extensions to other frameworks such as integrated design and operations (Burnak et al., 2019) and optimization under uncertainty (Charitopoulos However, as the problem size grows with the number of decision variables, constraints, and model parameters, the construction of the full multiparametric solution becomes challenging the number of critical regions grows exponentially with the problem size. Moreover, given a parameter vector, finding the critical region for that point becomes a point location problem, which can also be computationally difficult problem for larger systems with many critical regions.
al., 2018).
3 DFSM for MILPs
We consider MILPs, which aim to solve efficiently in online setting, given in the following general form:
minimize
duu)aa
A(u) x & b(u) R" x zn-m
(1)
where is a vector of continuous and discrete variables. The cost vector, constraint matrix, and right-hand-side vector are denoted by c(u), A(u), and b(u), respectively, which generally all depend on the input parameters u.
In general, the complicating constraints in an MILP are the integrality constraints on the discrete varia barisllss. Here, the key idea is to replace them with simple linear inequalities. This is motivated by a fundamental property of MILPs that can be stated as follows (Schrijver, 1998): Given a mixed-integer polyhedral set S = fIr TM Zn-m; : Ar I b), there exist A and such that the convex hull of S is conv(S) = {x R": Ar b, Az < b}, illustrated in Figure 1. Then, loosely speaking, an MILP with a cost vector c and S as its feasible region will have the same optimal value as an LP with the same cost vector and constrained by Ax S b and Ar S b. Morcover, if that MILP has a unique optimal solution, the corresponding LP will have the same optimal solution, ie. arg min, {c x x E S} = arg min, {c x x E conv(S)).
X2
Ax S
conv(s)
Ax S h
S
Åx b
X1
Figure 1: Illustrative example in which for the given cost vector c, the MILP with integer variables i and 22 and constraints Ar I has the same optimal solution the LPs with conv(S) and S' as their feasible regions.
One does not have to recover the full convex hull to obtain an LP that achieves the same optimal solution as the MILP. As illustrated in Figure 1, for the given cost vector c, the LP with S' = {æ e R" : Ar < Ax < b) as its feasible region, where Ar < b are a smaller set of constraints, will equally lead the same optimal solution. This idea forms the basis for traditional exact satting-plane approaches, where linear incqualities are successively added to the linear relaxation of the MILP until the optimal solution is found (Marchand
4

---

Published in Transactions on Machine Learning Research (03/2025)
al., 2002). The downside of cutting-plane algorithm (or one implemented in branch-and-cut scheme) is that for every new instance of the MII.P, new cutting planes need to be generated in the same iterative fashion, which can be too time-intensive for real-time applications.
In our DFSM approach, we propose to learn LPs with added linear inequalities parameterized by the input parameters u as surrogate optimization models for the corresponding MILPs. As such, this approach can be interpreted learning parametric cutting planes that are generated a priori and change automatically with each new instance given by new values for the input parameters
3.1 General formulation
Given an MILP of the form (1), the surrogate LP that we aim to construct can be written as follows:
minimize 2ER"
c(u)
(2a)
A(u)ær S b(u)
(2b)
A(u, 0)z S b(u,
(2c)
where (2b) are the original inequality constraints of problem (1), and (2c) are the new set of linear constraints that we want to learn. The elements of both the constraint matrix Ā and right-hand-side vector b are functions of the imput parameters These functional relationships are parameterized by the parameters O and can in principle be of arbitrary complexity. In this work, we use low-order polynomials, which have proven to be sufficiently accurate in our case studies and also provide some computational advantages (see Section 3.2). In this case, are the coefficients of the polynomial terms. In Appendix A.1, we present an illustrative example that highlights how such parametric inequalities can help the surrogate LP recover the optimal solutions the original MILP.
The DFSM problem for learning DFSOM of the form (2) be formulated as follows:
minimize 8E #
Σ EE
| 2:ll1
(3a)
3, arg surgminea x A(u,)x S b(ūe), A(u,, 0)x s 0(ū,,0) ViEI, zER"
(3b)
where I = {(4, xt ) constitutes the training dataset with ot being the optimal solution to the original MILP (1) for input ū, As shown in (3a), we use the f--nrrm to measure the decision error, and the objective is choose from some appropriate set such that the sum of decision errors across all training data points is minimized. Problem (3) is a bilevel optimization problem with lower-level problems. As shown in (3b), each lower-level problem represents the surrogate LP for particular input Ui and T; constrained to be an optimal solution to that LP.
3.2 Solution strategy
To solve problem (3), we first reformulate it into single-level optimization problem by replacing the lower-level problems with their KKT optimality conditions, which results in the following formulation:
minimize 8E0,8, 8.
Σ 1ET
| - -
(4a)
c+ A(T,) X + Á(ū, 0) T
#4: =0 VIEI A(ū.);; < VIEI
(4b)
(4c)
A(ū,, -0:;0));;;,)l00 < Viez
(4d)
((A))(((7) = VIET
(4c)
((H:)(((T,,))) b(0, 0)) =0 ViEI
(4f)
A 0, Mi 2 0, EER r
(4g)
5

---

Published in Transactions on Machine Learning Research (03/2025)
where (4b) are the stationarity conditions, (4c)-(4d) are the primal feasibility conditions, (4f)- are the complementary slackness conditions, and (4g) are the dual feasibility conditions. The dual variables are denoted by À and μ. In (4f)- D-)) denotes a diagonal matrix formed from a given vector, R. represents the set of original constraints in the MILP, and V represents the set of additional constraints in the DFSOM.
Problem (4) is a nonconvex NL.P that violates constraint qualification (due to, for example, the complementary slackness conditions), which means that standard NLP solvers camnot guarantee convergence the optimal solution when directly applied to this problem. Here, we adopt the solution strategy proposed by Gupta & Zhang (2022) for data-driven inverse optimization problems, which turn out to have the same structure the DFSM problem. We apply an exact penalty reformulation that is achieved by penalizing the violation of constraints (4b)-(4f) in the objective function using the &m-orm. We can then iteratively update the penalty parameter until no more constraint violation is observed at which point the algorithm is guaranteed to have converged to a local solution of problem (4). The penalty reformulation takes the form of a. multiconvex optimization problem in terms of the three sets of variables 8, é, and (>, j)), i.c. the problem is convex with respect to each of these three sets of variables. This structure can be exploited using a block coordinate descent (BCD) algorithm that allows efficient decomposition of the problem, which has been crucial in solving large instances of the DFSM problem. For more details on the solution algorithm, see Appendix A.2.
4 Computational case studies
We conduct two computational case studies to assess the efficacy of the proposed DFSM approach. Both examples are representative of common real-time optimization problems involving both continuous and discrete decisions. All optimization problems were implemented in Julia v1.7.2 using the modeling environment JuMP v0.22.3 (Dunning al., 2017)- We applied Gurobi v10.3 to solve all LPs and MILPs, and all NLPs were solved using IPOPT v0.9.1. All DFSM instances were solved utilizing 24 cores and 60 GB memory on the Mesabi cluster of the Minnesota. Supercomputing Institute (MSI) equipped with Intel Haswell E5- processors. The codes developed for these case studies are available at https: //gtt///t.ttttttt/eeeee wttmm/ddolab/DSSMffr--tt
4.1 Hybrid vehicle control
We first consider a hybrid vehicle control problem adapted from Takapoui et al. (2015). Given a hybrid vehicle with a battery, an electric generator, and an engine, the objective to minimize the fuel consumption cost over given time horizon while also maintaining high terminal battery level. This problem can be formulated as the following MILP:
T-1
minimize E, phatt, poux,
(ot peng + B2x) + (E Er)
(5a)
Eg Finit
(5b)
S E. < Emaa Vt 0, ..... T
(5c)
E Es - TPP batt
Vt = 0, ..... T
(5d)
S peme cmg < 2 pmax Is
Vt = 0, ..... T -
(5e)
P. batt + pmng 2 Dr
Vt = 0, .....!?? T -
(5f)
2 E (0,1,.,, S} Vt T......T--
(5g)
Here, the decision variables given for each time period t, starting with the energy status of the battery E,, then the power to from the battery P. batt the power from the engine prm and the engine switch status The input parameters are the power demand D, over the T time periods and the initial battery state Emnt. The cost. function to be minimized, (5a), contains the fuel cost in each time period given by atPt Beas and a penalty on the deviation of the terminal battery charge from its maximum given by (yEmax - Er). The initial battery state is given by (5b). The bounds on the battery state are stated in (5c). The battery charge balance is given by (5d), where 7 denotes the length of each time interval. We assume that the engine can operate at different levels, modeled using the integer variable 2, The different levels represent different.
6

---

Published in Transactions on Machine Learning Research (03/2025)
fractions of the maximum engine power available for use in given time period. As per constraints (5e), when = O, no power can be derived from the engine, and when 2n S, the maximum amount of power can be drawn. Finally, the combined power output from the battery and the engine must at least match the demand Dj as stated in (5f).
4.1.1 Data generation and surrogate design
To generate different problem instances, we choose " (225,,5.5), ar ~ (6, Br ~ N/0.5, 2), and =
Emit ~ (900, 97), and set == 5, pmax = 1, and A = 100. Moreover, the value of S is also varied between (1.2.3) 3) across these instances. This way we generate 10 different instances, each representing a hybrid vehicle with different attributes. For each instance, the training dataset is generated by varying the power demand profile D, which is the model input, and solving the original MILP for that demand profile to obtain the corresponding optimal solution (E", peng*, =- Repeating this process provides a set of data points, each consisting of an input-solution pair.
The additional linear incqualitics in the postulated surrogate LP take the form of Á(D, 0] penR =|_ S b(D, 0), where the number of constraints [V] is a hyperparameter that can be adjusted. Each element of À and is cubic polynomial in D with 0 being the coefficients. In this case, the time required to train the surrggate model for 100 training data points is around 10,000 seconds. This surrogate design remains the same across the 10 different instances; however, for each instance, different linear inequality parameters 0 have to be derived.
We also compare the performance of the DFSM approach with that of three NN-based optimization proxics, one that applies a simple feedforward NN architecture, second one where the loss function contains penalty terms derived from the Lagrangian dual of the optimization problem to impose feasibility on the predicted solutions (Fioretto et al., 2020), and a third one that is based on a graph NN (GNN) architecture. To represent the MILP as a GNN, we derive a bipartite graph whose two sets of nodes correspond to the problem's variables and constraints, respectively. An edge connects a variable node with a constraint node if that variable appears in the constraint. Also, fcatures such as variable bounds and right-hand-side values are assigned to the nodes to fully describe the MILP. These NN-based optimization proxies are similarly trained in decision-focused manner using the same datasets. See Appendix A.3 for more details on the design and training of the NN-based optimization proxies.
4.1.2 Surrogate model performance
All reported values are averaged over the 10 different randomly generated instances. For the DFSOM and the NN-based optimization proxies, if a predicted solution is not directly integer-feasible, an inexpensive projection-based fcasibility restoration problem is solved to recover a feasible solution. More details on the feasibility restoration step can be found in Appendix A.4.
The decision prediction errors and optimality gaps for each instance are determined using 100 unseen test data points, cach given by different power demand profile. The relative prediction error is defined the fi-nrrm of the difference between the predicted solution (after feasibility restoration) and true optimal solution divided by the 4i-nrrm of the true optimal solution. Similarly, the optimality gap is the difference between the costs of the predicted and true optimal solutions. The size of the problem increases with T'. In the following, we present results for the case with T = 30, which has 61 continuous variables, 30 integer variables, and 122 constraints in the resulting MTL.; additional results for 10 and = 20 are shown Appendix A.5.
Figure 2 compares the performance of the DFSOM and the NN-based optimization proxies. For these methods, it shows how the prediction errors and optimality gaps change with the number of training data points. Note that while the maximum number of training data points we use to construct each DFSOM is 100, we use up to 500 data points to (rain each of the three NN-based optimization proxies. Also, to aid the interpretation of the results, we show the prediction errors with respect to the discrete (Figure 2(a)) and continuous (Figure 2(b)) decisions separately. We observe that the Lagrangian-based NN (NN+LD) achieves slightly better solutions than the standard NN, likely because its predictions are closer to being integer-feasible The GNN performs significantly better. Yet. the DFSOM outperforms all three NN-based optimization proxies in terms
7

---

Published in Transactions on Machine Learning Research (03/2025)
of prediction accuracy, where the prediction errors plateau at much lower values as the size of the training dataset increases. DFSM also proves to be very data-efficient as it achieves high prediction accuracy with only 50 data points while the NNs seem to require much more data.
(9)
Discrete decisions
[b)
Continuous decisions
(0)
% decciion prediction anror
NOHLD
ddecision prediction
% optimalitt gap
29 111
30
90 70
90
190
200 300 489 ea0
180
200
100 383
of training data points
tralninn data points
training data points
Figure 2: DFSOM performance compared to NN-based optimization proxies in the hybrid vehicle control problem.
Notice that the relative prediction errors that the DFSOM achieves are about 20% and 8% in terms of the discrete and continuous decisions, respectively, which seem relatively high. The corresponding optimality gaps, however, are close to zero, as shown in Figure 2(c). While surprising at first, upon closer investigation, we can attribute this bohavior to the specific nature of the hybrid vehicle control problem. Given a fixed solution in terms of the discrete variables z for problem (5), due to the flexibility that the battery provides, there are many solutions in terms of the continuous variables that are optimal or close to optimal. Moreover, as the cost parameter a has much larger value than 3,, the contribution of the continuous variables to the total cost outweighs that of the discrete variables. As a result, from a cost standpoint, it is more important to predict the continuous variables accurately, hence the difference in the relative prediction errors.
We see that it is important to predict feaegbre solutions, which is where DFSM has an advantage as it keeps all linear inequalities from the original MILP in the DFSOM. Here, the difference between the solutions directly predicted by the DFSOM and the integer-feasible solutions obtained after feasibility restoration is average 6.6%. The standard NN tends to predict solutions that are far from being integer-feasible such that the obtained solutions after feasibility restoration are highly suboptimal, as indicated in Figure 2(c). This is also reflected in the on average 26.8% difference between the solutions predicted by the NN and those post fensibility restoration. This difference decreases to 23.1% in the case of the Lagrangian- baaod NN, which is still much higher than what is achieved by the DFSOM. In the case of the GNN, the difference further decreases to 12.8%, highlighting the advantage of using GNNs as NN-based optimization proxies for MILPs in problems of this format.
We further investigate the impact of the number of added constraints, [V], on the performance of the resulting DFSOM. Figure 3 shows the decision prediction errors and optimality gaps for different (V] Here, we can see that when 7/3 linear inequalities are added, the DFSOM cannot reduce the prediction error to a reasonable level even with a large number of training data points. As we increase [V] to T and then to 2T, the prediction accuracy increases significantly. However, there is only a marginal improvement from | = 27 to |V] = 37, which indicates that [V] = 2T likely a good choice that provides a good balance between prediction accuracy and computational complexity.
Finally, we consider all 1,000 test data points across the 10 different instances to compare the computation times of solving the original MILPs, the corresponding DFSOMs, and the original MILPs but only to the same objective function values as achieved by the DFSOMs. The third method serves as a benchmark heuristic solution method that allows a fair comparison with the DFSOM. Figure 4 shows the number of instances solved by cach of the three mothods as the solution time increases. One can see that the surrogate LP clearly outperforms the MILP, even when the MILP only solved to achieve a solution of the same quality as the DFSOM.
8

---

Published in Transactions on Machine Learning Research (03/2025)
(N)
Discrete decislons
(h)
Continuuus decisions
(4)
prediction BITER
pmntta & pyet D
Motoa d Net
Nenta de Nert a MoeT Nedt
decision prediction
3 MriT
% optimaltty gmu
8
# training date points
P training darte points
8 training data points
Figure 3: DFSOM performance for varying number of constraints learned ((2)) in the hybrid vehicle control problem.
Model predictive control
100
% of Instances solved
80
60
DFSOM
20
MILP
MILP solved DFSON
objective value
10
10
10
10 time (s)
102
10a
10
Figure 4: Computational performance of DFSOMs evaluated across 1,000 random instances in the hybrid vehicle control case compared directly solving the original MILPs and solving the MILPs the same objective values achieved by the corresponding DFSOMs
4.2 Production scheduling
In the second case study, we consider a single-stage production scheduling problem with parallel units, which can be formulated as the following MILP (Maravelias, 2021):
minimize S,T,
ΣΣΣ My Xije JET T3 T3--1 VjEJ,tET
(6a)
(6b)
ΣΣ XIj== =1 ViEI
(6c)
JET IET
Σxit Xijt 21 VjEJ,LET
(6d)
Sy M Σ Xijt ViEI,jEJ
(Ge)
Sy 2 Ty.-1 - 7/11 - X(x) Vi,
(6f)
Sy + STy, + 7/1 -ii)) ViEt.j
(6g)
Σ S4 S522 2p: ViEI
(6h)

---

Published in Transactions on Machine Learning Research (03/2025)
Σs EΣ
Tay Tuy Se, VieI 1 SET
(6i)
S 2 0, T3 20, Xiji E (0.1) VieI,jEJ,tET
(6j)
where I denotes the set of units, I is the set of batches to be processed, T is the set of time slots, and ff is the horizon length. For each batch i, we are given a release time Pi, a due time E,, and processing time Tij if it is processed on unit j. The binary variable Xija equals 1 if batch i is processed which is associated with cost denoted by 7aj; Siy denotes the start time of batch i on unit, j, and T3 is the end time of time slot t for unit j. The objective is to minimize the total production cost given in (6a). Constraints (6b) ensure the non-negative length of each time slot, equations (6c) ensure that every batch is processed, (6d) allow at most one batch to be processed on a given unit at a given time, (6e) set the start time of a batch on a unit. to 0 if the batch is not being processed at any time on that unit and otherwise provide a bound on the start time, (6f) and (6g) combined determine the end time of time slot on each unit, (6h) only allow batch to start processing after its release, and (6i) ensure that the processing of each batch is completed before its due time.
unit j in time slot t,
4.2.1 Data generation and surrogate design
In all the 10 randomly generated problem instances, we use a horizon length of n 40. We choose Tij (),,,)), P ~ 9),, e ~ (10,399) Ti ~ (300, 100), and we set M max;ez(ex). Here, the model input parameters are p, e, and T, which we vary to generate the training datasets. We also set T = {1,... 8),, which is the smallest set that is still sufficiently large for all considered instances. With |Z] 13, | = 4, and [7] = 8, the MILP has 84 continuous variables, 416 binary variables, and 987 constraints.
In the hybrid vehicle control case study, we learned additional linear inequalities that involve all decision variables. If we were to take the same approach here, each new constraint would involve []]]]||]| + [ZI + [J][T] decision variables; this would lead to a very large number of parameters to be learned, which would significantly increase the computational complexity of the DFSM problem. Hence, instead, we try to exploit the structure of the scheduling problem in designing the additional constraints. Specifically, we notice that most of the constraints in (6) are written for all E I, and it also makes intuitive sense that the relationships are strongest among the variables associated with the same batch. Using this rationale, we propose to add for each batch é linear incqualitics that only involve variables associated with that batch. The constraint coefficients and right-hand-sides are given as quadratic functions of the input parameters. In this case, the time taken to train the surrogate model for 100 training data points is around 3,000 seconds.
4.2.2 Surrogate model performance
Like in the first case study, we assess the performance of the DFSOM and the three NN-based optimization proxies. Note that the MILP solutions that serve as the basis for comparison are obtained at 1% optimality gap since many problem instances cannot be solved to full optimality within several hours. The results in Figure 5 show that the DFSOM clearly outperforms the NN, NN+LD, and GNN. Here, unlike in the first case study, the GNN's performance does not come close to that of the DFSOM, which could be attributed to the significant increase in problem size with this production scheduling problem having 500 decision variables and close to 1000 constraints. Moreover, the training time for GNN here is an order of magnitude higher compared to the DFSM approach as there are many more parameters to be learned in the GNN. This further showcases the advantage of DFSM where it is straightforward to leverage domain knowledge to reduce the size of the learning problem by choosing which variables to include in the additional constraints.
The DFSM approach is again very data-cfficient, achieving decision prediction error with respect to binary variables of less than 10% with only 50 training data points (Figure 5(a)). In this case, the prediction error with respect to the continuous variables is much higher and docs not improve with incrcasing number of training data points (Figure 5(b)), yet. the DFSOM still achieves very small optimality gaps (Figure 5(c)). The reason is that the objective function in problem (6) only involves the binary variables, which is why the optimality gap will be small as long as the solution in terms of the binary variables is prodicted accurately. In addition, for a fixed assignment of batches to units and time slots (discrete decisions), there are many start times (continuous variables) that are feasible and hence optimal since they do not affect the cost. Due
10

---

Published in Transactions on Machine Learning Research (03/2025)
to this multiplicity of optimal solutions, it is not surprising that the surrogate LP does not achieve the same solutions in terms of the continuous variables. But as long as the predicted solutions in terms of the binary variables are close to the true ones, the DFSOM will perform well. The difference between the surrogate solutions and the solutions obtained post fcasibility restoration here is 8.3% when DFSOM is used, while 50.4% and 56.6% are achieved by the standard and Lagrangian-based NNs, respectively. The GNN shows some improvement with the difference reaching 38.1%.
(0)
Binary decisions
(b)
Continuous decisions
(c)
N6 decision prediction error
% dectoion prediction atror
WHHLD
% optimallt gap
8
5.10
training
points
training data points
training points
Figure 5: DFSOM performance compared to NN-based optimization proxies in the production scheduling problem.
Figure 6 shows the impact of the number of added constraints on the DFSOM's performance. Here, [ denotes the number of constraints added per batch. Interestingly, we observe that it is sufficient to add just one constraint for each batch. In fact, the DFSOMs for | and | 10 have much larger prediction errors with respect to the binary variables and optimality gaps for smaller training datasets, which indicates overfitting.
decisions
Continuous dectsions
(8)
(4)
+ Mrr
1/441
O MITT
9 Mis Morts
0144
tive o prrte
ddecciion prediction 2
predtction
11p110
optimmiity gpp.
training
points
Figure 6: DFSOM performance for varying number of constraints learned (VV) in the production scheduling problem.
We again solve the 1,000 problem instances taken from the 10 different instances to compare the computation times of solving the original MILPs, the corresponding DFSOMs, and the original MILPs but only to the same objective function values as achieved by the DFSOMs. In this case, the original MILPs were solved to 1% optimality gap as it would have taken excessively long to solve them to full optimality. Even with 1% optimality, as shown Figure 7, some instances of the MILP take multiple hours solve while the DFSOM solves within a few seconds all instances.
In addition, we test the performance of the DESOM on some hard problem instances for which the corre- sponding MILPs could not be solved to 1% optimality within the given time limit of 3,000 seconds. We use 30 such hard instances to form the test dataset, while the DFSOM is trained only using instances that could be solved to optimality. The results are shown in the performance plot in Figure 8, where we observe that the DFSOM obtains integer-feasee feasible solutions significantly faster when compared to solving the MILP to the same objective function values. We should note that our hope here was that the DFSOM might occasionally achieve even better solutions than the MILP (at the time limit) since it was trained using optimal solutions. While this does not seem to be true in this case, it is an aspect we plan to explore more in our future work.
11

---

Published in Transactions on Machine Learning Research (03/2025)
Production scheduling
100
% of Instances solved
80
60
DFSOM
20
DFSOM ozjective value
10 100
10
10
10
10a
time (s)
Figure 7: Computational performance of DFSOMs evaluated across 1,000 random instances in the production scheduling case compared to directly solving the original MILPs and solving the MILPs to the same objective function values achieved by the corresponding DFSOMs.
100
% of Instances solved
80
40
DFSOM
20
MILP solved DFSOM
objective value
10
10
102
10
time (s)
Figure 8: Computational performance of DFSOMs evaluated across 30 hard instances in the production scheduling case compared to solving the MILPs to the same objective function values achieved by the corresponding DFSOMs.
5 Conclusions
Motivated by the need solve difficult MILPs in real-time settings, we developed a data-driven approach for constructing efficient surrogate LPs that can replace the MILPs. When generating these surrogate LPs, explicitly try to minimize the decision prediction error defined as the difference between the optimal solutions of the surrogate and the original optimization problems. The resulting docision-focused surrogate modeling problem is large-scale bilevel program that we solve using a penalty-based block coordinate descent algorithm. In our computational case sludies, we demonstrated the elficacy of the proposed approach both in terms of prediction accuracy and data efficiency. It also allows the use of problem-specific knowledge design surrogate model structures that enhance the training and/or the performance of the surrogate optimization models. The results further show that the surrogate LPs trained using the proposed approach outperfomm, often significantly, state-of-the-art NN-based optimization proxies.
12

---

Published in Transactions on Machine Learning Research (03/2025)
Acknowledgments
The authors gratefully acknowledge financial support from the National Science Foundation under Grant #2044077 and from the State of Minnesota through an appropriation from the Renewable Development Account to the University of Minnesota West Central Research and Outreach Center. They would also like to acknowledge the Minnesota. Supercomputing Institute (MST) the University of Minnesota for providing resources that contributed the rescarch results reported in this paper.
References
JoaquN Acevodo and Efstratios Pistikopoulos. An algorithm for multiparametric mixed-integer lincar programming problems. Operalions research letlers, 24(3):139-144 148, 1999.
Alessandro Alessio and Alberto Bemporad. A survey on explicit model predictive control. Lecture Noles in Control and Information Sciences, 384:345-3663 2009.
Brandon Amos and J. Zico Kolter. Optnet: Differentiable optimization as a layer in neural networks. 2017.
Yoshua Bengio, Emma Frejinger, Andrea Lodi, Rahul Patel, and Sriram Sankaranarayanan. A learning-based algorithm to quickly compute good primal solutions for stochastic integer programs. In Integration of Constraint Programming, Artificial Intelligence, and Operations Research: 17th International Conference, CPAIOR 2020, Vienna, Austria, September 21 24, 2020, Proceedings 17, pp. 99 111. Springer, 2020.
Yoshua Bengio, Andrea Lodi, and Antoine Prouvost. Machine learning for combinatorial optimization: a methodological tour d' d'horizon. European Journal of Operational Research, 290(2):405 421, 2021.
Dimitris Bertsimas and Bartolomeo Stellato. Online mixed-intege optimization in milliseconds. INFORMS
Journal on Computing, 34(4):2229-2245 2022.
Athary Bhosckar and Marianthi Terapetritou. Advances in surrogate based modeling, feasibility analysis, and optimization: A review. Computers and Chemical Engineering, 108:250 267, 2018.
Baris Burnak, Nikolaos A. Diangelakis, Justin Katz, and Efstratios N. Pistikopoulos. Integrated process design, scheduling, and control using multiparametric programming. Compulers Chemical Engineering, 125:164-184, 6 2019.
Vassilis M Charitopoulos, Lazaros G Papageorgiou, and Vivek Dua. parameeareeetric mixed integer linear programming under global uncertainty. Computers & Chemical Engineering, 116:279-222 00:22992995 2018.
Steven Chen, Kelsey Saulnier, Nikolay Atanasov, Daniel D. Lee, Vijay Kumar, George J. Pappas, and Manfred Morari. Approximating Explicit Model Predictive Control Using Constrained Neural Networks. In Proceedings of the American Control Conference, pp. 1520 1527. AACC, 2018.
Arnaud Deza and Elias B Khalil. Machine learning for cutting planes in integer programming: survey. arXiv preprint Xrxvx22302,0 09166, 2023.
Priya L Donti, Brandon Amos, and J Zico Kolter. Task-based End-to-End Model Learning in Stochastic Optimization. In Advances in Neural Information Processing Systems 30, 2017.
Vivek Dua and Efstratios N Pistikopoulos. An algorithm for the solution of multiparametric mixed integer linear programming problems. Annals of operalions research, 99:123 139, 2000.
lain Dumning, Joey Huchette, and Miles Lubin. Jump: A modeling language for mathematical optimization.
SIAM Review, 59(2):295-33) 320, 2017.
Adam N. Elmachtoub and Paul Grigas. Smart 'Predict, then Optimize'. Management Science, 68(1):9-26, 2022.
Aaron Ferber, Bryan Wilder, Bistra Dilkina, and Milind Tambe. MIPaaL: Mixed integer program as a layer. Tn AAAT 2020 34th AAAT Conference on Artificial Intelligence, pp. 1504-1511, 2020.
13

---

Published in Transactions on Machine Learning Research (03/2025)
Aaron M Ferber, Thoan Huang, Daochen Zha, Martin Schubert, Benoit Steiner, Bistra Dilkina, and Yuandong Tian. Surco: Learning linear surrogates for combinatorial nonlinear optimization problems. In International Conference on Machine Learning, pp. 10034 10052.. PMLR, 2023.
Anthony V. Fiacco and Yo Ishizuka. Sensitivity and stability analysis for nonlinear programming. Annals of Operations Research, 27:215-2583 235, 12 1990.
Ferdinando Fioretto, Terrence W. K. Mak, and Pascal van Hentenryck. Predicting AC optimal power flows: Combining doep learning and Lagrangian dual methods. In Proceedings of the AAAI Conference on Arlificial Intelligence, pp. 630-637, 2020.
Maxime Gasse, Didier Chételat, Nicola Ferroni, Laurent Charlin, and Andrea Lodi. Exact combinatorial optimization with graph convolutional neural networks. Advances in Newral Information Processing Systems, 32, 2019.
Rishabh Gupta and Qi Zhang. Efficient learning of decision-making models: A penalty block coordinate descent algorithm for data-driven inverse optimization. Compulers & Chemical Engineering, pp. 108123, 2022.
Rishabh Gupta and Qi Zhang. Data-driven decision-focused surrogate modeling. AIChE Journal, pp. e18338, 2024.
Zeren Huang, Kerong Wang, Furui Liu, Hui-Ling Zhen, Weinan Zhang, Mingxuan Yuan, Jianye Hao, Yong Yu, and Jun Wang. Learning to select cuts for efficient mixed-integer programming. Pattern Recognition, 123:108353, 2022.
Asunción Jiménez- Jiménez-Cordero, Juan Miguel Morales, and Salvador Pineda. Warm-starting constraint generation for mixed-integer optimization: machine learning approach. Knowledge-Band Based Systems, 253:109570, 2022.
Benjamin Karg and Sergio Lucia. Efficient Representation and Approximation of Model Predictive Control Laws via Deep Learning. IEEE Transactions on Cybernetics, 50(9):3866-33 3878, 2020.
Elias Khalil, Pierre Le Bodic, Le Song, George Nemhauser, and Bistra Dilkina. Learning to branch in mixed integer programming. In Proceedings of the AAAT Conference on Artificial Intelligence, volume 30, 2016.
James Kotary, Ferdinando Fiorelto, and Pascal Van Hentenryck. Fast Approximations for Job Shop Scheduling: A Lagrangian Dual Deep Learning Method. In Proceedings of the AAAI Conference on Artificial Intelligence, 2022.
Pratyush Kumar, James B. Rawlings, and Stephon J. Wright. Industrial, largo-scale model prodictive control with structured neural networks. Computers and Chemical Engineering, 150: 107291, 2021.
Andrea Lodi and Giulia Zarpellon. On learning and branching: a survey. Top, 25(2):207-22) 2017.
Jayanta Mandi, Emir Demirović, Peter J. Stuckey, and Tias Guns. Smart predict-and-optimize for hard combinatorial optimization problems. In AAAI 2020 34th AAAI Conference on Artificial Intelligence, pp- 1603 1610,, 2020.
Christos T Maravelias. Chemicol Production Scheduling: Skired-Inteeer Programming Models and Methods. Cambridge University Press, 2021.
Hugues Marchand, Alexander Martin, Robert Weismantel, and Laurence Wolsey, Cutting planes in integer and mixed integer programming. Discrete Applied Mathematics, 123(1-3):397-446, 2002.
Richard Oberdieck, Nikolaos A. Diangelakis, Toana Nascu, Maria M. Papathanasion, Muxin Sun, Styliani Avraamidou, and Efstratios N. Pistikopoulos. On multi-parametric programming and its applications in process systems engineering. Chemical Engineering Research and Design, 116:61-82, 12 2016.
Yoshiaki Ohtake and Naonori Nishida. A branch-and-a bound algorithm for 0-1 parametric mixed integer programming. Operations Research Letters, 4(1):41-44 1985.
14

---

Published in Transactions on Machine Learning Research (03/2025)
Xiang Pan, Tianyu Zhao, and Minghua Chen. DeepOPF: Deep neural network for DC optimal power flow. In Proceedings of the 2019 IEEE International Conference on Communicalions, Control, and Compuling Technologies for Smart Grids (SmartGridCommm) IEEE, 2019.
Jool A. Paulson and Ali Mcsbah. Approximate Closed-Loop Robust Model Predictive Control with Guaranteod Stability and Constraint Satisfaction. IEEE Control Systems Letters, 4(3):719 724, 2020.
Efstratios N. Pistikopoulos, Nikolaos A. Diangelakis, Richard Oberdieck, Maria M. Papathanasion, Ioana Nascu, and Muxin Sun. PAROC- An integrated framework and software platform for the optimisation and advanced model-based control of process systems. Chemical Engineering Science, 136:115 138, 2015.
Gary M Roodman. Postoptimality analysis in zero-one programming by implicit enumcration. Noval Research Logislies Quarterly, 19(3):435 447, 1972.
Alexander Schrijver. Theory of linear and integer programming. John Wiley & Sons, 1998.
Yunzhuang Shen, Yuan Sun, Andrew Eberhard, and Xiaodong Li. Learning primal heuristics for mixed integer programs. In 2021 international joint conference on neural networks (tjenn), pp. 1 8. LEEE, 2021.
Haoran Sun, Xiangyi Chen, Qingjiang Shi, Mingyi Hong, Xiao Fu, and Nicholas D. Sidiropoulos. Learning to
Optimize: Training Deep Neural Networks for Interference Management. IEEE Transactions Signal Processing, 66(20):5438-56) 54438,5533311 2018.
Reza Takapoui, Nicholas Moehle, Stephen Boyd, and Alberto Bemporad. A simple effective heuristic for embedded mixed-integer quadratic programming. 9 2015.
Yunhao Tang, Shipra Agrawal, and Yuri Faenza. Reinforcement learning for integer programming: Learning cut. In International conference on machine learning, pp. 9367 9376. PMLR, 2020.
Pascal Van Hentenryck. Optimization learning. arXiv preprint arXit:2501.04443 2025.
Alexandre Velloso and Pascal Van Hentenryck. Combining Deep Learning and Optimization for Preventive Sensriy- DC Optimal Power Flow. IEEE Transactions on Power Systems, 36(4):3618 3628, 2021.
Bryan Wilder, Bistra Dilkina, and Milind Tambe. Melding the data-decisions pipeline: Decision-focused learning for combinatorial optimization. In 33rd AAAI Conference on Artificial Intelligence, AAAI 2019, 31st Innovative Applications of Artificial Intelligence Conference, IAAI 2019 and the 9th AAAI Sympostum Educational Advances in Artificial Intelligence, EAAI 2019, pp. 1658-16666 2019.
Ahmed S. Zamzan and Kyri Baker. Learning optimal solutions for extremely fast optimal power flow. Proceedings of the 2020 IEEE International Conference on Communications, Control, and Computing Technologies for Smart Grids (SmartGridComm), 2020.
A Appendix
A.1 Illustrative example
To illustrate the working of the surrogate LP model, we present a 2-dimensional multi-constrain. knapsack problem with one of its constraints changing function of the input parameter
maximize #y 22, 2268
4.8s, + 6%2
4r + 3z2 70
(7)
100kz + 85.g < 800u + 680
Φ Σ 17, 32 S 17.
15

---

Published in Transactions on Machine Learning Research (03/2025)
We obtain surrogate LP for problem (7) by adding linear inequalities to the original problem whose parameters are trained by solving the DFSM problem (3). The resulting DFSOM is the following L.P: maximize #1 ER, #268
1.801 + 6x2
40 3x² S 70
100wr; 85 - 800u + 680
(17.9- 56.9u + 45.1u² a=1 + (-4.5 25.3u 14.3u Joz S 100
(8)
(7.4 + 6.4u - 7.9n²) 21 v 100
11 (-0.2 + 14.5u - 9.42 ) + (12 9.7u + 3.7m S 100
21 - 17, 22 - 17.
In Figure we show for different values of the feasible regions of the LP relaxation of problem (7), which we denote S, in blue. The feasible region of problem (8), denoted by 8', is shown in green. The integer-feasible points for problem (7) constitute the set S. From the figure, one can see that the solutions obtained by solving the surrogate LP (8) for different values of u are the same as those obtained by solving the original MIL.P (7).
20
20
20
15
15
10
10
C
5
5
10
15
20
10 31
15
20
5
10
15
20
T1
T
(a) 1.45
(b) u
(c) 0.61
Figure Comparison of feasible regions and optimal solutions of (7) and (8). For the three instances of u (1.45; 0.61), the corresponding optimal solutions are x² = (2,, 17); (12, 7); (12,5)).
Note that the added inequalities A(u)ær S b(u) may not be valid for all integer-feasible points of the original problem as DFSM only aims minimize the prediction error with respect to the optimal solutions. Indeed, for the three cases of u-values shown in Figure 9, conv(S) g S in the first two instances as shown in Figures 9(a) and 9(b).
A.2 The BCD algorithm
To describe the BCD approach, we first re-write the formulation (4) in the following form:
minimize A.
Σ 1 3:lli
subject to
6 A/ dM =0 ViEI
A/a: - b, 50 ViEI
(9)
(A/) - b)) =0 ViEI
d,2 0, Z2 E R. ViEI
4/= = A(ū,))
b(ū,)
[X [t]
b = b(u, 0)
d
ViEI.
Á(u, 0)
Problem (9) a nonconvex NLP which becomes intractable for large instances. Furthermore, some constraints, such as the complementarity constraints, the problem cause it violate constraint qualification conditions leading to convergence difficulties for standard NLP solvers.
16

---

Published in Transactions on Machine Learning Research (03/2025)
We address the lack of regularity in (9) by considering its penalty reformulation in (10), which is standard strategy employed in the NL.P literature for degenerate NLPs such as (9):
Eiez |cc A/" dill
minimize B
Σ || -ill1 + q" EE
Σier max (0,4 - -
(10)
EiEz Ild (AfT: - 61)1
subject to 220 ViEI,
where we penalize the violation of the constraints using a nonsmooth (--orrm-base penalty function. This particular penalty function is known to be exact in the sense that, for values of the penally parameters q larger than certain threshold, every optimal solution of (9) also minimizes (10). Therefore, instead of (9), one can solve (10) which is amenable to solution via standard NLP solvers if the description of the set θ satisfies necessary constraint qualification conditions.
However, the threshold value for q above which the exactness of the penalty reformulation holds generally hard to determine a priori. Therefore, we employ an itcrative approach where we start with small values for q and successively increase them by factor in every iteration until we arrive at solution that is feasible for (9) (Algorithm 1).
The reformulated problem (10) is still nonconvex NLP that is difficult to solve when the mumber of data points ([ZI) is large.. We alleviate this issue by observing that (10) is multiconvex optimization problem if the set θ is convex. This feature of (10) allows us to solve it via an efficient block-cordinate--desenn (BCD) approach. While implementing BCD, we tackle the nonsmoothness of the objective functions of BCD subproblems by augmenting them with the proximal operator.
Algorithm 1 The penalty-based BCD algorithm.
initialize (#, 0, d') with a feasible solution from IPOPT
initialize g + qº, fix p+= p while convergence criteria. for blocks are not satisfied do solve (10) with BCD or IPOPT,
Σier ||c + A," di|l
q= -
EiEz max (0, A/T: - b3)
Eier [Id (A/T: k0Il
end while
A.3 Design of NN-based optimization proxies
In case studies, contrast the performance of the proposed DFSOM approach with three different NN-based optimization proxies as surrogates for the MILPs. In the following, we briefly outline the design and training of these NNs. Note that there are no hard constraints embedded in these NN models; hence, the feasibility of the predicted optimal solutions with respect to the constraints of the original model is not guaranteed. Therefore, we use the same feasibility restoration strategy that is applied to the outputs of the DFSOM obtaín integer-feasible solutions from the NN-based optimization proxies; discussion of this strategy provided in the next section.
A.3.1 Feedforward NN architecture
To train the NN models, we split the data containing input parameters and optimal solutions into training and test datasets. We use the mcan squared error between the true and prodicted values as the loss function and apply the ADAM optimizer to estimate the NN parameters. Similar to DFSOM, for given input, the NN model is also trained to predict optimal solution for the original problem.
17

---

Published in Transactions on Machine Learning Research (03/2025)
For the feedforward NN, the numbers of hidden layers and neurons in each layer are determined through trial-and-error approach aiming to minimize the training loss over a fixed number of epochs. The number of hidden layers is based on the empirical observation that increasing the number of hidden layers resulted in little improvements in the model's accuracy. The design of the specific NN model used the hybrid vehicle control case study is as follows: the input layer contains T nodes, the same as the dimension of the input vector u. The two hidden layers have 37 and 67 nodes, and finally, the output layer consists of 37 nodes aligning with the dimensionality of the decision space. We also account for the effect of different learning rates (namely 0.1, 0.01, and 0.001) the training loss. Based analysis, determine 0.01 to be the most suitable learning rate for our purpose.
A.3.2 Augmented Lagrangian-basse based NN architecture
Here, the NN has the same feedforward architeclure as described above, but the loss function is augmented with penalty terms associated with the Lagrangian dual of the optimization problem. The Lagrange multipliers are updated using subgradient approach in each iteration of the NN training process as discussed in Van Hentenryck (2025). Assume the given MILP is of the following form:
minimize
c(u) Tx
(u( x S b1(4) VIEI
(11)
ge (u) T x = ha(u)
EEE
E Th x zn-m
Then, the augmented loss function can be written as follows:
C (V,
ly
Ø1 + Σ Aelgo(u) he(u)) + Σ NET
Ms m0,(0) b;,)),,
(12)
cE
where 9 are the NN model predictions given by 1 = Mo(u), and the model parameters 0 are estimated in each epoch by keeping the Lagrange multipliers fixed from the last iteralion.
After these model parameters are tuned from 0f - gt+1, the multipliers are adjusted using the constraint violations as follows:
x+ No p n
ge (ux) Ma (ux) he (ux)l.
(13)
1+ p + p
max((0,a;(ux))
be(ux))
A.3.3 GNN architecture
Here, we leverage the strong representation power of GNNs to capture the structure of MILPs. We use a weighted bipartite graph consisting of variables nodes, constraints nodes, and edges, where an edge connects a variable node with a constraint node if the variable appears in that constraint. For example, problem (7) of the illustrative example can be represented as the graph shown in Figure 10.
(4.8, -00 17,1) 81
(70.
100u
000 17,1)
(800u 680,
85
Figure 10: Bipartite graph representation of the MILP (7) with variable nodes on the left and constraint. nodes on the right.
The GNN incorporates three sets of features:
18

---

Published in Transactions on Machine Learning Research (03/2025)
1. The features of each variable node contain information about the correspponding coefficient in the objective function, the lower and upper bounds of that variable, and whether the variable is discrete (1) or contimuous (0).
The features of each constraint node contain information about the corresponding right-hand-siee value and whether the constraint is inequality (5) or equality (=).
3. The edge weights are the coefficients of the variables in the constraints that they appear in.
Using this graph representation of the MILP as input, we train a graph convolutional neural network as described in Gasse et al. (2019) to predict the optimal solutions.
A.4 Feasibility restoration
The predicted decisions from the DFSOM may not always satisfy the integrality constraints. For such instances, we employ an inexpensive feasibility restoration strategy to obtain integer-feasible solution from the predicted solution. Let ž; be the optimal solution of the DFSOM for ú, Suppose the discrete decisions associated with é,, denoted by 2; are not integer. To restore integer feasibility, we solve the following problem:
minimize DE C x
11 - 17:11
(14)
A(ū,)x S b(ū.).
We then solve the original MILP with the discrete decision variables fixed to the values of the optimal solution to problem (14) to obtain the corresponding optimal values for the contimous variables. Figure 11(a) provides illustration of an infeasible prediction from the DFSOM. By using the feasibility restoration strategy described above, we obtain the closest integer-feasible point as shown in Figure 11(b).
X2
X2
Ax S
Ax S
conv
conv S
Ax = h
Åx S h
S'
Åx h
X1
(a) Infeasible prediction
(b) Feasible projection
Figure 11: The predicted decision variables that are not integor-fcasible are used a projection problem to return discrete decision values.
A.5 Additional results
For the hybrid vehicle control problem, Figures 12 and 13 compare the performance of the DFSOM and the NN-based optimization proxy for 10 and T = 20, respectively. The DFSOMs were trained with | = 2T. The results indicate DFSOM's superior performance even for these smaller problems, emphasizing the benefit. of incorporating original constraints DFSM. Notably, with only 50 training data points, DFSOM reliably yields solutions with optimality gap close to zero in both cases.
19

---

Published in Transactions on Machine Learning Research (03/2025)
Discrete decisions
Continuous decisions
(X)
55 preditonttdittton OITOr
DESCM
DFSOM
OFGON
dedrsion prediction
training data points
training data points
#oatnann data points
Figure 12: DESOM performance compared to NN-based optimization proxy for T = 10.
decisions
Continuous decisions
(6)
(6)
predictin arror
AM DFSON
100
% deciston predrcton
9 B
% optimaity
12
1>
training points
training points
fraining points
Figure 13: DFSOM performance compared NN-based optimization proxy for T = 20.
20


Reviews:
[{"review_id": "XHKmBHVdRV", "reviewer": "Reviewer_RhjW", "summary": "", "questions": "", "limitations": "", "rating": null, "confidence": null, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "Strengths: \n1. The idea of interpreting the learning of LP surrogates as learning parametric cutting planes is interesting and novel.\n2. The paper is generally well-written and clear.\n\nWeaknesses: \n1. While the paper compares DFSOM to NN-based optimization proxies, it does not evaluate their performance against SOTA, i.e., other ML-based optimization methods, limiting the broader impact of the findings and their positioning within the existing body of literature.\n2. Conducting experiments on realistic datasets would be beneficial.\n3. While the penalty-based BCD algorithm provides an efficient solution for training the surrogate models, its reliance on iterative optimization could become computationally expensive for very large MILPs. Additionally, the feasibility restoration process, used to correct infeasible predictions, adds computational overhead that may not be trivial in time-sensitive applications.", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "Yes", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "paper_topic_and_main_contributions": "", "excitement": null, "reproducibility": null, "ethical_concerns": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "summary_of_contributions": "This paper proposes a data-driven approach for online mixed-integer linear programs (MILPs), which learns surrogate linear programs (LPs) by adding cutting planes, i.e., linear inequalities, that can be solved more efficiently and yield similar optimal solutions compared to the original MILPs. Instead of learning new cutting planes whenever new parameters of the MILP are provided, the approach learns the parameters of the cutting planes as functions of the MILP parameters. Thus, when given a set of new MILP parameters, the learned functions automatically provide a new set of parameters for the cutting planes. The authors use the l1-norm to measure the difference between the LP solution and the MILP solution, first formulating the learning task as a bilevel optimization problem and then reformulating it into a single-level optimization problem by replacing the lower-level problems with their KKT optimality conditions. Computational case studies on hybrid vehicle control and production scheduling demonstrate that DFSOM outperforms NN-based optimization proxies in terms of prediction accuracy, data efficiency, and computational performance.", "claims_and_evidence_yn": "Yes", "claims_explanation": "", "audience_yn": "Yes", "audience_explanation": "", "requested_changes": "Please refer to \"Weaknesses\".", "broader_impact_concerns": "N/A", "additional_comments": "", "extra_scores": {}, "extra_text": {}}, {"review_id": "iOggpfVmX0", "reviewer": "Reviewer_ujd8", "summary": "", "questions": "", "limitations": "", "rating": null, "confidence": null, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "**Strengths**\n1. The proposed method is data-efficient, requiring relatively fewer instances compared to neural network proxies. However, further exploration of additional capabilities achievable with the method could enhance its impact (see weaknesses).\n\n2. The proposed method shows potential for application in online settings, which is a common requirement in real-world scenarios. However, this potential should be further clarified (see weaknesses).\n\n3. The method effectively solves large instances of DFSM by leveraging the problem's multiconvexity and employing block coordinate descent (BCD), ensuring computational tractability for complex problems.\n\n**Weaknesses**\n\n1. While the authors suggest that the method is applicable to online settings, it is unclear what this entails. From the experiments, some problems are solved in approximately 100 seconds, which, depending on the domain, might align more with offline settings. To manage expectations about where the method can be effectively used, it would be beneficial to provide a clearer definition of \"online setting\" (e.g., in the introduction) and its relevance to the application context. It is worth acknowledging that solving within 100 seconds could still be considered \"online\" in certain domains. Highlighting such domains explicitly would strengthen the claim. However, a key concern is whether the cost of solving an exact MILP or another expensive optimization approach in these domains makes the proposed method sufficiently competitive. Exploring a comparison of solving times and costs across these settings would provide a more comprehensive evaluation of the method's practical utility.\n\n\n2. A discussion of scenarios where the method might fail or perform suboptimally would strengthen the paper. \nIt would be interesting to see how the proposed method could benefit from using NNs to predict constraint cofficients instead of polynomial functions. For instance, if NNs were used to predict parameter coefficients, could the method learn better constraints and achieve improved generalizability? \n\n3. The benchmarking results leave some room for improvement.\n    - Graph Neural Networks (GNNs), which are well-suited for MILPs, should serve as a baseline for prediction tasks. While using MLPs is a reasonable choice, it is less common in the community for such tasks. Reassessing benchmarks using GNNs will provide a stronger evaluation and validation of the claims. However, it is important to note that GNNs are less interpretable and data-efficient than the proposed method, so the proposed method still holds relevance to the community.  \n\n    - While the related works section discusses ML-based optimization proxies, the paper does not directly compare its method with these approaches. Such comparisons would provide valuable insights into the method’s relative performance and unique advantages.\n\n\n4. A standard practice is to test a method’s robustness and generalizability by evaluating performance on harder instances. This would help identify the boundaries of the method's applicability. For example, displaying geometric mean solving times across problem categories would offer a clearer picture of its performance across varying complexities.", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "No", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "paper_topic_and_main_contributions": "", "excitement": null, "reproducibility": null, "ethical_concerns": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "summary_of_contributions": "The authors propose a novel decision-focused modeling strategy to approximate the solution of MILPs, offering a method that is interpretable, data-efficient, and computationally inexpensive during inference. By leveraging surrogate linear programs (LPs) to approximate the structure of MILPs, the approach achieves near-optimal solutions while significantly reducing computation time.", "claims_and_evidence_yn": "No", "claims_explanation": "", "audience_yn": "Yes", "audience_explanation": "", "requested_changes": "Some benchmarks as mentioned in the weaknesses section. \n\n**Clarifications:**\n- In Section 4.1.1, Does the training dataset include several D for a single draw of parameters? Will this violate the i.i.d assumption?\n- DFSOM & NLP are not expanded", "broader_impact_concerns": "N/A", "additional_comments": "", "extra_scores": {}, "extra_text": {}}, {"review_id": "VdYNfYcjt5", "reviewer": "Reviewer_NpQz", "summary": "", "questions": "", "limitations": "", "rating": null, "confidence": null, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "### Strengths\n - The paper is well-written and clear\n - The idea of \"learning cuts\" for turning MILPs into LPs is intuitive, yet novel in this context. \n - The experiments are good. I appreciate the use of problems grounded in real world scenarios.\n\n### Weaknesses\n - I think the learning problem/solution strategy this paper hones in on (i.e., $\\ell_1$ loss function, use KKT conditions to turn a bilevel problem into a single level, constrained problem) is a bit narrow. I think this paper could benefit from further discussion of alternate approaches (e.g., using the $\\ell_2$ loss: $\\\\|x_i^* - \\hat{x}_i\\\\|_2^2$).\n \n### Minor points/questions\n 1. Throughout the paper, the abbreviations DFSM and DFSOM are used, but only DFSM is defined. Is this a typo, do both abbreviations reference the same thing, or is there a subtle difference?\n 2. As mentioned on pg. 5, you use low-order polynomials to parametrize the added constraints, and I agree that they work well in your case studies. Did you try any other parameterization schemes? Would you suggest other parameterization schemes to practitioners interested in implementing your framework?\n 3. (minor) Equation 2 is split over two pages. If you can fiddle with the formatting to get it onto one page, that might be nice.\n 4. The use of hats to denote LP solutions ($\\hat{x}$) as well as added constraints  ($\\hat{A},\\hat{b}$) is confusing. I suggest using $\\tilde{A},\\tilde{b}$ instead.\n 5. Why do you use $\\ell_1$ error in equation (3), as opposed to, e.g., $\\ell_2$?\n 6.  For both DFSM and the NN approach, a feasibility restoration step is applied (as discussed in Section 4.1.2). Could the authors comment on how close to integral solutions from both methods tend to be?", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "Yes", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "paper_topic_and_main_contributions": "", "excitement": null, "reproducibility": null, "ethical_concerns": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "summary_of_contributions": "Consider a mixed-integer linear program (MILP) in the form:\n$$\\begin{split} \\operatorname*{minimize}_x c(u)^{\\top}x \\\\\\\\ A(u)x \\leq b(u) \\\\\\\\ x \\in \\mathbb{R}^m\\times \\mathbb{Z}^{n-m} \\end{split} \\tag{1}$$\nwhere $u$ is an input parameter. Suppose that this MILP needs to be solved repeatedly, for different values of $u$. As is well-known, solving an MILP can be computationally expensive, whereas solving a Linear Program (LP) is fairly cheap. Thus, this paper considers the problem of _learning_ LP surrogates, parametrized by $u$, for a parametrized MILP as above. \n\nThe surrogate LPs are constructed by dropping the integrality constraints in (1) and adding additional, learned linear inequality constraints: $\\hat{A}(u,\\theta)x \\leq \\hat{b}(u,\\theta)$. Here $\\theta$ represents parameters encoding the functional form of the learned constraints.  The parameters $\\theta$ are learned in a data-driven way, from a training set of parameter-solution pairs, $(u_i,x_i^*)$, where $x_i^*$ represents the solution to (1) given $u = \\bar{u}_i$.  An innovative optimization scheme is proposed to do so.\n\nThe paper is rounded out by two interesting case studies, both of which involve challenging MILPs. The results of both case studies indicate that the proposed approach works well---providing fast solutions with small optimality gap.", "claims_and_evidence_yn": "Yes", "claims_explanation": "", "audience_yn": "Yes", "audience_explanation": "", "requested_changes": "See above.", "broader_impact_concerns": "NA", "additional_comments": "", "extra_scores": {}, "extra_text": {}}]

Author rebuttal:

Meta-review: The paper presents a method to learn and compute surrogate LPs for MILPs that take less time to solve, yet still has optimum close to that of the original MILP. Specific attention is given to the goal of achieving feasibility of the surrogate optimum in the original MILP. This is a problem of broad interest in discrete optimization, while the problem itself is a machine learning problem requiring development of new/better ML techniques. After the author-reviewer discussion and paper revisions, reviewers are satisfied with the presentation and experimental evaluation of the proposed method. I recommend acceptance of this paper to TMLR, given that both the "Claims and Evidence" and "Audience" criteria are satisfied.
Decision: Accept (poster)