FMIP: Joint Continuous-Integer Flow For Mixed-Integer Linear Programming

## Abstract

Abstract Mixed-Integer Linear Programming (MILP) is a foundational tool for complex decision-making problems.
However, the NP-hard nature of MILP presents a significant computational challenge, motivating the development of machine learning-based heuristic solutions to accelerate downstream solvers.
While recent generative models have shown promise in learning powerful heuristics, they suffer from a critical limitation.
That is, they model the distribution of only the integer variables and fail to capture the intricate coupling between integer and continuous variables, creating an information bottleneck and ultimately leading to suboptimal solutions.
To this end, we propose Joint Continuous-Integer F low for M ixed- I nteger Linear P rogramming ( FMIP ), which is the first generative framework that models the joint distribution of both integer and continuous variables for MILP solutions.
Built upon the joint modeling paradigm, a holistic guidance mechanism is designed to steer the generative trajectory, actively refining solutions toward optimality and feasibility during the inference process.
Extensive experiments on eight standard MILP benchmarks demonstrate the superior performance of FMIP against existing baselines, reducing the primal gap by 41.34% on average.
Moreover, we show that FMIP is fully compatible with arbitrary backbone networks and various downstream solvers, making it well-suited for a broad range of real-world MILP applications.
Our code is available * * * https://github.com/Lhongpei/FMIP .

## 1 Introduction

Mixed-Integer Linear Programming (MILP) represents a cornerstone of mathematical optimization, providing a powerful framework for modeling complex decision-making problems that involve both discrete choices and continuous quantities Zhou et al. (2025); Kratica et al. (2014); He et al. (2015). Its ability to capture the intricate interplay between integer and continuous variables makes it indispensable across diverse domains, including combinatorial optimization (Della Croce and Paschos, 2014), energy systems (Miehling et al., 2023), and supply chain design (Ivanov et al., 2022).

Despite its versatility, solving a MILP instance remains a fundamental challenge due to its NP-hard property (Bertsimas and Tsitsiklis, 1997). Consequently, exact solvers often rely on heuristics to find high-quality solutions in a practical amount of time. This has motivated a surge of interest in using deep learning to predict powerful heuristics that can warm-start and guide the downstream solvers (Nair et al., 2020; Han et al., 2023; Liu et al., 2025; Chen et al., 2025).
Among them, generative modeling has emerged as a particularly promising approach. Modern techniques like diffusion models and flow matching reframe the difficult one-step task of predicting a complete solution into an iterative refinement process (Sun and Yang, 2023; Feng et al., 2024; Zeng et al., 2024). By learning to transform a simple noise distribution into a distribution of high-quality solutions, these models can effectively navigate the complex combinatorial search space, breaking down the challenge of satisfying intricate constraints into a series of more manageable steps.

However, existing generative methods for MILP suffer from a critical limitation: they model the distribution of only the integer variables. As illustrated in Figure [1](https://arxiv.org/html/2507.23390v2#S1.F1), this architectural choice ignores the strong, coupled relationship between integer and continuous variables, which is fundamental to the structure of MILP problems.
Since both variable types are often required to evaluate the objective function and check for constraint satisfaction, this incomplete modeling creates an information bottleneck, leading to suboptimal solution quality and hindering the effectiveness of any guidance mechanism.

Figure: Figure 1: The key advantages of our FMIP over existing works: 1) the joint distribution modeling on both integer and continuous variables and 2) the consequent holistic guidance during inference.
Refer to caption: pic/fmip_compare_exist.png

To this end, we propose Joint Continuous-Integer Flow for Mixed-Integer Linear Programming (FMIP), which is the first generative framework that models the joint distribution of both integer and continuous variables for MILP solutions. FMIP leverages a conditional flow matching process to progressively generate a complete solution, fully capturing the interdependence between all decision variables. This joint modeling approach unlocks a holistic guidance mechanism that can steer the generation process using complete, instance-wise feedback from both the objective function and constraint violations.
FMIP is fully compatible with arbitrary backbone networks and downstream solvers, showing 41.34% relative improvement on average across eight standard MILP benchmarks over SOTA baselines.
Our contributions are summarized as follows:

- •
We introduce FMIP, the first generative framework to jointly model the complete distribution of both integer and continuous variables for MILP solutions, effectively capturing their critical coupling relationship.
- •
We design a holistic guidance mechanism that leverages the jointly generated variables to steer the sampling process toward solutions with better objective values and constraint satisfaction.
- •
As a powerful generative learning paradigm, FMIP is agnostic to the choice of backbone networks and downstream solvers, making it well-suited for a broad range of real-world MILP applications.
- •
FMIP sets a new state-of-the-art for learning-based MILP heuristics, achieving superior results on eight benchmarks while showing compatibility with diverse backbones and solvers.

## 2 Related Works

Machine learning is widely used to accelerate Mixed-Integer Linear Programming (MILP) solvers. One major line of research focuses on enhancing the internal components of the Branch-and-Bound algorithm, such as learning policies for variable branching (Gasse et al., 2019a; Khalil et al., 2016; Gasse et al., 2019b; Gupta et al., 2020; Zarpellon et al., 2021; Gupta et al., 2022; Scavuzzo et al., 2022; Lin et al., 2024; Zhang et al., 2024) or cutting plane selection (Tang et al., 2020; Huang et al., 2022).
Another closely related direction, which our work belongs to, aims to predict high-quality heuristic solutions to warm-start the optimization process, thereby reducing the search space and speeding up convergence (Nair et al., 2020; Han et al., 2023; Huang et al., 2024b; Liu et al., 2025; Chen et al., 2025; Zeng et al., 2024).

Most existing solution-prediction methods are discriminative in nature, employing models like Graph Neural Networks (GNNs) to predict the variable assignments in a single forward pass (Huang et al., 2024b; Liu et al., 2025; Hu et al., 2024; Huang et al., 2024a; Ye et al., 2024). However, predicting a complete, feasible, and high-quality solution in one step is exceptionally challenging due to the complex combinatorial structure of MILPs. This difficulty has motivated the use of generative models, such as diffusion models and flow matching, for MILP heuristics (Feng et al., 2024; Zeng et al., 2024). These models reframe the solution prediction as an iterative refinement process, decomposing the hard single-step task into a more manageable multi-step generation, which is better suited for navigating complex constrained spaces (Zeng et al., 2024).
Despite these advances, prior approaches generally focus on partial distribution modeling on integer variables, leading to suboptimal guidance and inferior solutions.
In contrast, our proposed Joint Continuous-Integer Flow for MILP is the first generative framework that explicitly models the joint distribution of both continuous and integer variables, directly addressing this critical gap in the literature.

## 3 Preliminaries

### 3.1 MILP Definition and Graph Representation

Suppose we have $n$ decision variables denoted as $\boldsymbol{x}=(\boldsymbol{x}_{\mathcal{I}},\boldsymbol{x}_{\mathcal{C}})\in\mathbb{Z}^{q}\times\mathbb{R}^{n-q}$, where $\boldsymbol{x}_{\mathcal{I}}$ and $\boldsymbol{x}_{\mathcal{C}}$ represent the integer and continuous variables, respectively.
An MILP instance can be defined as:

$$ $\displaystyle\min_{\boldsymbol{x}}$ $\displaystyle\boldsymbol{w}^{\top}\boldsymbol{x}$ (1) s.t. $\displaystyle\boldsymbol{A}\boldsymbol{x}\leq\boldsymbol{b}$ $\displaystyle\boldsymbol{l}\leq\boldsymbol{x}\leq\boldsymbol{u}$ $\displaystyle\boldsymbol{x}_{\mathcal{I}}\in\{0,1,\dots,K\}^{q}$ $\displaystyle\boldsymbol{x}_{\mathcal{C}}\in\mathbb{R}^{n-q}$ $$

Here, $\boldsymbol{w}$ is the objective function coefficient vector. $\boldsymbol{l}$ and $\boldsymbol{u}$ are the lower bounds and upper bounds for decision variables.
$\boldsymbol{A}\in\mathbb{R}^{m\times n}$ and $\boldsymbol{b}\in\mathbb{R}^{m}$ denote the coefficient matrix and the right-hand-side vector of linear constraints.
Moreover, $K$ is a scalar integer defining the maximum feasible value of integer variables (i.e., $\boldsymbol{x}_{i}\in\{0,1,\dots,K\}$, $\forall i\in\mathcal{I}$).
An MILP instance with $m$ linear constraints and $n$ variables can thus be represented by the tuple $\boldsymbol{M}=(\boldsymbol{A},\boldsymbol{b,l,u,w})$. A solution $\boldsymbol{x}$ is considered feasible if all constraints are satisfied.
Note that, in this paper, we focus on MILP problems with bounded integer variables, as most real-world applications exhibit this characteristic. Specially, when $K=1$ such that $\boldsymbol{x}_{\mathcal{I}}\in\{0,1\}^{q}$, the problem reduces to Mixed Integer Binary Programming (MIBP), which is the most commonly studied subclass of MILP (Huang et al., 2024b; Liu et al., 2025).

Graph Representation for MILP. When applying deep learning to MILP, it is common to represent an MILP instance as a bipartite graph (Gasse et al., 2019a).
The bipartite graph consists of two sets of nodes: one representing the $n$ decision variables and the other representing the $m$ constraints.
Edges in this graph connect variables to the constraints in which they appear, effectively encoding the sparsity pattern of the constraint matrix.
Building on this foundation, as shown in the top-left part of Figure [2](https://arxiv.org/html/2507.23390v2#S3.F2), we employ a more granular tripartite graph structure to explicitly distinguish between variable types, which is crucial for our joint modeling approach.
Specifically, an MILP instance is denoted as a graph $\boldsymbol{G}=(\mathcal{V}_{\mathcal{I}},\mathcal{V}_{\mathcal{C}},\mathcal{V}_{\text{con}},\mathcal{E})$, where the nodes are partitioned into three distinct sets:

- •
Integer Variable Nodes ($\mathcal{V}_{\mathcal{I}}$): One node for each of the $q$ integer variables.
- •
Continuous Variable Nodes ($\mathcal{V}_{\mathcal{C}}$): One node for each of the $n-q$ continuous variables.
- •
Constraint Nodes ($\mathcal{V}_{\text{con}}$): One node for each of the $m$ linear constraints.

An edge is drawn between a variable node $j\in\mathcal{V}_{\mathcal{I}}\cup\mathcal{V}_{\mathcal{C}}$ and a constraint node $i\in\mathcal{V}_{con}$ if the variable $x_{j}$ has a non-zero coefficient in the $i$-th constraint.
Based on the rich, relational MLIP graph, various graph neural networks (e.g., GCN (Kipf, 2016) or GAT (Veličković et al., 2017)) can be employed to extract the rich topological representations for later discriminative or generative learning processes.
More details can be found in Appendix [A](https://arxiv.org/html/2507.23390v2#A1) and Appendix [B](https://arxiv.org/html/2507.23390v2#A2).

### 3.2 Flow Matching

Flow Matching (FM) is a powerful framework for training generative models (Lipman et al., 2023; Albergo and Vanden-Eijnden, 2023; Liu et al., 2023). The core idea is to learn a time-dependent vector field $v_{t}$ that transforms samples from a simple prior distribution $p_{0}$ (e.g., Gaussian noise) into samples from a target data distribution $p_{1}$.
This transformation is defined by a probability path $(p_{t})_{t\in[0,1]}$ interpolating between $p_{0}$ and $p_{1}$, which is induced by an Ordinary Differential Equation (ODE): $\frac{d\boldsymbol{c}_{t}}{dt}=\bm{v}_{t}(\boldsymbol{c}_{t})$.
FM allows the vector field $\bm{v}_{t}$ to be learned directly via a simple regression objective, given a defined path between noise and data samples.
FM is flexible and can model both continuous and discrete data, which is essential for our work.

- •
For continuous data, the prior $p_{0}$ is typically a standard Gaussian $\mathcal{N}(\mathbf{0},\mathbf{I})$. The interpolating path between data point $\boldsymbol{c}_{1}$ and noise $\boldsymbol{c}_{0}$ is defined as $\boldsymbol{c}_{t}|\boldsymbol{c}_{1}\sim\mathcal{N}(t\boldsymbol{c}_{1},(1-t)^{2}\mathbf{I})$. This yields a simple, closed-form target vector field for model training, i.e., $\bm{v}_{t}(\boldsymbol{c}_{t}|\boldsymbol{c}_{1})=\frac{\boldsymbol{c}_{1}-\boldsymbol{c}_{t}}{1-t}$.
- •
For discrete data, the prior $p_{0}$ is a uniform categorical distribution over $K$ categories. The conditional path for each component $i$ of a data point $\boldsymbol{d}_{1}$ is defined as $\boldsymbol{d}_{t}^{(i)}|\boldsymbol{d}_{1}^{(i)}\sim\operatorname{Cat}(t\delta(\boldsymbol{d}_{t}^{(i)},\boldsymbol{d}_{1}^{(i)})+(1-t)/K)$, where $\delta$ is the Kronecker delta and $\operatorname{Cat}(\cdot)$ denotes the categorical distribution. This path implies a target conditional rate matrix $R_{t|1}(\cdot,\cdot|\boldsymbol{d}_{1}^{(i)})$ for model training:
$R_{t|1}(\boldsymbol{d}_{t}^{(i)},j\;\,|\;\boldsymbol{d}_{1}^{(i)})=\frac{\delta(\boldsymbol{d}_{1}^{(i)},j)}{1-t}\left(1-\delta(\boldsymbol{d}_{1}^{(i)},\boldsymbol{d}_{t}^{(i)})\right).$
(2)

During inference, these targets are replaced by the model’s predictions to generate new samples by evolving the state over time.
At each step, the learned model provides an estimate of the vector field $\hat{\bm{v}}$ or rate matrix $\hat{R}$, which is used to update the current state:

$$ $\displaystyle\boldsymbol{c}_{t+\Delta t}=\boldsymbol{c}_{t}+\hat{\bm{v}}_{t|1}(\boldsymbol{c}_{t}|\boldsymbol{c}_{1})\Delta t,$ $\displaystyle\text{for continuous data},$ (3) $\displaystyle\boldsymbol{d}_{t+\Delta t}^{(i)}\sim\operatorname{Cat}\left(\delta(\boldsymbol{d}_{t+\Delta t}^{(i)},\boldsymbol{d}_{t}^{(i)})+\hat{R}_{t|1}(\boldsymbol{d}_{t}^{(i)},\boldsymbol{d}_{t+\Delta t}^{(i)})\Delta t\right),$ $\displaystyle\text{for discrete data}.$ $$

Figure: Figure 2: The overall framework of our proposed FMIP.
Refer to caption: pic/new_banner.png

## 4 FMIP: Joint Continuous-Integer Flow for MILP

In this section, we introduce FMIP, a generative framework that learns the joint distribution of integer and continuous variables to find high-quality solutions for MILP instances.
As shown in Figure [2](https://arxiv.org/html/2507.23390v2#S3.F2), FMIP is built on a time-dependent graph representation, a joint continuous-integer flow model, a holistic guidance mechanism for inference, and seamless integration with downstream solvers.

### 4.1 Time-Dependent Graph Representation

FMIP operates on a time-dependent solution graph, which captures both the static structure of an MILP instance and its dynamic state during the generative process. As defined in Section [3.1](https://arxiv.org/html/2507.23390v2#S3.SS1), the static structure is encoded in a tripartite graph $\boldsymbol{G}=(\mathcal{V}_{\mathcal{I}},\mathcal{V}_{\mathcal{C}},\mathcal{V}_{\text{con}},\mathcal{E})$.

During the flow matching process, we define the solution graph at time-step $t\in[0,1]$ as the tuple $\boldsymbol{G}_{t}=(\boldsymbol{G},\boldsymbol{d}_{t},\boldsymbol{c}_{t})$, where $\boldsymbol{d}_{t}$ and $\boldsymbol{c}_{t}$ are the values of the integer and continuous variables, respectively. This dynamic state information is incorporated directly into the graph by augmenting the static features of each variable node.
We define $\boldsymbol{x}_{t}=(\boldsymbol{d}_{t},\boldsymbol{c}_{t})$ as the vector of solution values for all variables. $\boldsymbol{x}_{t}^{(i)}$ denotes the state of the $i$-th variable at time step $t$.

### 4.2 Joint Continuous-Integer Flow

To overcome the limitations of integer-only models, we design a flow matching process over the mixed solution space of both continuous and integer variables. The process constructs a probability path $p_{t|1}(\boldsymbol{G}_{t}|\boldsymbol{G}_{1})$ for $t\in[0,1]$, which transforms a simple noise distribution $p_{0}(\boldsymbol{G}_{0})$ into the target distribution of high-quality solutions $p_{1}(\boldsymbol{G}_{1})$.

#### Joint Training Objective

Our model, parameterized by $\theta$, learns this transformation by predicting the target solution $(\boldsymbol{d}_{1},\boldsymbol{c}_{1})$ from a noisy state $\boldsymbol{G}_{t}$.
It features two prediction heads over a shared learnable graph backbone network (e.g., GCN or GAT).
One head outputs the denoised continuous values $\hat{\boldsymbol{c}}_{1}(\boldsymbol{G}_{t})$ and another outputs the probability distribution of the integer solution $\hat{p}(\boldsymbol{d}_{1}|\boldsymbol{G}_{t})$. The model is trained by minimizing a joint loss function that combines a regression loss for the continuous variables and a cross-entropy loss for the integer variables:

$$ $\mathcal{L}(\theta)=\mathbb{E}_{t,\boldsymbol{G}_{1}}\left[\frac{\left\|\hat{\boldsymbol{c}}_{1}(\boldsymbol{G}_{t})-\boldsymbol{c}_{1}\right\|_{2}^{2}}{1-t}-\omega\log\hat{p}(\boldsymbol{d}_{1}|\boldsymbol{G}_{t})\right],$ (4) $$

where $t\sim\mathcal{U}(0,1)$, $\boldsymbol{G}_{1}\sim p_{\text{data}}$, and $\omega$ is a hyperparameter balancing the two terms.

#### Sampling Process

At inference, we generate a solution by simulating the forward ODE from $t=0$ to $t=1$. At each time step $t$, the model takes the current noisy graph $\boldsymbol{G}_{t}$ as input and produces its predictions, i.e., $\hat{\boldsymbol{c}}_{1}(\boldsymbol{G}_{t})$ and $\hat{p}(\boldsymbol{d}_{1}|\boldsymbol{G}_{t})$.
As illustrated in Eq. [3](https://arxiv.org/html/2507.23390v2#S3.E3), these outputs are used to estimate the conditional vector field $\bm{v}_{t}$ and rate matrix $R_{t|1}$, enabling the next step of the simulation. We use a cosine schedule for time discretization, which allocates more steps to the low-noise region near $t=1$ for better quality (Nichol and Dhariwal, 2021). This generative formulation is a key advantage over traditional one-step discriminative predictors, as it enables the multi-step iterative refinement with a holistic guidance mechanism during sampling, which we detail next.

### 4.3 Holistic Guidance Mechanism

A key benefit of FMIP is that the model predicts a complete solution candidate $(\hat{\boldsymbol{d}},\hat{\boldsymbol{c}})$ at any step $t$ during generation.
Hence, holistic guidance can be directly derived using instance-specific information from the MILP formulation itself (both the objective function and constraints), which steers the generation process towards solutions that are both feasible and optimal.
Such a capability is absent in previous works that only consider integer variables.

Guidance Target Function. We define a target function $f(\boldsymbol{x})$ that captures the two primary goals of an MILP solution: minimizing the objective value and satisfying the constraints. The function is a weighted sum of the objective and a penalty for constraint violations:

$$ $f(\boldsymbol{d},\boldsymbol{c})=f(\boldsymbol{x})=\boldsymbol{w}^{\top}\boldsymbol{x}+\gamma\sum_{i=1}^{m}\left[\max(0,\boldsymbol{A}_{i,*}\boldsymbol{x}-\boldsymbol{b}_{i})\right]^{2},$ (5) $$

where $\boldsymbol{x}=(\boldsymbol{d},\boldsymbol{c})$ is the set of all variables, $\gamma$ is a balancing hyperparameter, and $\boldsymbol{A}_{i,*}$ is the $i$-th row of the constraint matrix. The model is guided to generate solutions that minimize this function. Following Lin et al. (2025), we implement guidance by modifying the update steps for the continuous and integer variables separately.

Guidance on Continuous Variables.
For the continuous variables, we employ gradient-based guidance. At each step $t$, we perform gradient descent on the target function $f$ w.r.t. the predicted continuous values $\hat{\boldsymbol{c}}_{1}(\boldsymbol{G}_{t})$, steering them towards regions of lower objective value and higher feasibility:

$$ $\boldsymbol{c}_{t+\Delta t}\leftarrow\text{Project}_{\boldsymbol{l},\boldsymbol{u}}\left(\boldsymbol{c}_{t}-\rho_{t}\nabla_{\boldsymbol{c}_{t}}f\left(\hat{\boldsymbol{d}}_{1|t},\hat{\boldsymbol{c}}_{1}(\boldsymbol{G}_{t})\right)\right),$ (6) $$

where $\hat{\boldsymbol{d}}_{1|t}$ is sampled from the predicted integer distribution $\hat{p}(\boldsymbol{d}_{1}|\boldsymbol{G}_{t})$ and $\rho_{t}$ is the step size.
The projection function $\text{Project}_{\boldsymbol{l},\boldsymbol{u}}(\cdot)$ constrains the continuous variables within their specified bounds.

Guidance on Integer Variables.
For the integer variables, we adopt an effective sampling-and-reweighting scheme motivated by recent work (Lin et al., 2025). Instead of computing gradients, at time step $t$, we sample a batch of $B$ candidate integer solutions $\{\boldsymbol{d}_{1|t,r}\}_{r=1}^{B}$ from the model’s current predicted integer distribution $\hat{p}(\boldsymbol{d}_{1}|\boldsymbol{G}_{t})$.
We then reweight the transition probabilities in the rate matrix $\hat{R}$ based on the quality of each integer candidate, as evaluated by $f(\boldsymbol{d}_{1|t,r},\hat{\boldsymbol{c}}_{1}(\boldsymbol{G}_{t}))$, steering the categorical distribution towards more promising integer assignments:

$$ $\hat{R}(\boldsymbol{d}_{t}^{(i)},\cdot)\leftarrow\frac{\sum_{r=1}^{B}\operatorname{exp}\left(f\left(\boldsymbol{d}_{1|t,r},\hat{\boldsymbol{c}}_{1}(\boldsymbol{G}_{t})\right)/\psi\right)\cdot R_{t\mid 1}(\boldsymbol{d}_{t}^{(i)},\cdot\mid\boldsymbol{d}_{1\mid t,r}^{(i)})}{\sum_{r=1}^{B}\operatorname{exp}\left(f\left(\boldsymbol{d}_{1|t,r},\hat{\boldsymbol{c}}_{1}(\boldsymbol{G}_{t})\right)/\psi\right)}$ (7) $$

where $\psi$ is a temperature parameter. This updated rate matrix $\hat{R}$ is then used to sample the next integer state followed by the projection function for feasibility constraints, i.e., $\boldsymbol{d}_{t+\Delta t}\leftarrow\text{Project}_{\boldsymbol{l},\boldsymbol{u}}(\boldsymbol{d}_{t+\Delta t})$.

### 4.4 Integration with Downstream Solvers

The output of the guided sampling process of FMIP is a high-quality candidate solution $(\boldsymbol{d}_{1},\boldsymbol{c}_{1})$ and a probability distribution over the integer variables.
This provides a powerful warm-start for downstream solvers to efficiently search for better solutions.
FMIP is fully compatible with a wide range of downstream solvers and can provide a powerful warm start for them to efficiently search for optimal solutions.
For example, Predict-and-Search (Han et al., 2023)) employs the integer distribution to define a promising search region and uses the continuous values as an initial feasible point for the solver’s LP relaxations.

## 5 Experiments

### 5.1 Experiment Settings

#### Datasets & Benchmarks

We evaluate FMIP on a comprehensive suite of eight MILP benchmarks.
Five of them focus on classic combinatorial optimization problems: Combinatorial Auctions (CA) (Gasse et al., 2019a), Generalized Independent Set (GIS) (Colombi et al., 2017), Maximum Independent Set (MIS) (Gasse et al., 2019a), Fixed-Charge Multi-Commodity Network Flow (FCMNF) (Hewitt et al., 2010), and Set Covering (SC) (Gasse et al., 2019a).
We also include two real-world MILP datasets with both binary and continuous variables from the NeurIPS ML4CO 2021 competition (Gasse et al., 2022): Load Balancing (LB) and Item Placement (IP). Finally, we adopt the standard MIPLIB2017 benchmark (MIPLIB) which is one of the most common standards for MILP solver evaluation (Gleixner et al., 2021a).

#### Baselines

We compare FMIP against three representative learning-based methods and the commercial solver Gurobi (Gurobi Optimization, LLC, 2024). The learning-based baselines are:

- •
SL: Standard supervised learning serves as a strong discriminative baseline. The model is trained to directly predict the variable assignments in a one-step manner.
- •
DIFUSCO: A state-of-the-art diffusion model for generating solutions to integer linear programs (Sun and Yang, 2023; Feng et al., 2024).
- •
IP-Guided-Diff: A guided discrete diffusion framework designed for integer programs that incorporates problem-specific guidance (Zeng et al., 2024).

Since DIFUSCO and IP-Guided-Diff are designed for integer-only problems, we adapt them to the MILP setting by including continuous variables and their constraints in the graph structure for the backbone graph encoder.
In this way, we allow the baselines to perceive both the continuous and integer variables for fair comparison.

FMIP and three baselines are all learning methods that are agnostic to backbone graph neural networks and downstream solvers.
In later experiments, we employ four different backbone graph encoders, i.e., Tri-GCN (Gasse et al., 2019b), Bi-GCN (Gasse et al., 2019b), GAT (Brody et al., 2021), and ClusterGCN (Chiang et al., 2019), with Tri-GCN as the default choice.
As for downstream solvers, we adopt Neural Diving (ND) (Nair et al., 2020), Predict-and-Search (PS) (Han et al., 2023), PMVB (Chen et al., 2025) and Apollo-MILP (Apollo) (Liu et al., 2025).

#### Metrics

Following previous works (Han et al., 2023; Liu et al., 2025), we report the objective value (OBJ) found by each method within a fixed time limit.
To compare performance, we first determine the best-known solution (BKS), defined as the best objective value found across all methods, including a long run of Gurobi.
We then calculate the absolute primal gap: $\text{GAP}=|\text{OBJ}-\text{BKS}|$.
Since MIPLIB contains instances from diverse domains with very different optimal values, we report the relative primal gap for MIPLIB instead: $\text{Rel.GAP}=|\text{OBJ}-\text{BKS}|/(|\text{BKS}|+1)$. A lower absolute/relative primal gap indicates better performance.

#### Implementation Details

Due to the page limitation, we provide implementation details in Appendix [C](https://arxiv.org/html/2507.23390v2#A3), including the backbone model architecture, training & inference hyperparameters of FMIP and baselines, as well as the downstream solver configuration.

**Table 1: The overall performance of FMIP and other baselines combined with different downstream solvers. Tri-GCN is chosen as the graph encoder. Best results are given in bold, and the second-best values are underlined. The relative improvement is computed as $\text{Rel.Imprv.}=(\text{GAP}_{\text{best-bsl}}-\text{GAP}_{\text{FMIP}})/({\text{GAP}_{\text{best-bsl}}}+1e-6)$, where $\text{GAP}_{\text{FMIP}}$ is the gap of the FMIP method, and $\text{GAP}_{\text{best-bsl}}$ is the gap of the best-performing baseline method. For MIPLIB, since we already give the relative primal gap (Rel.GAP), the relative improvement is computed as $\text{Rel.Imprv.}=(\text{Rel.GAP}_{\text{best-bsl}}-\text{Rel.GAP}_{\text{FMIP}})$**
| Downstream<br>Solver | Training<br>Method | CA | GIS | MIS | FCMNF |  |  |  |  |
| --- | --- | --- | --- | --- | --- | --- | --- | --- | --- |
| OBJ $\uparrow$ | GAP $\downarrow$ | OBJ $\uparrow$ | GAP $\downarrow$ | OBJ $\uparrow$ | GAP $\downarrow$ | OBJ $\downarrow$ | GAP $\downarrow$ |  |  |
| ND(400s) | SL | 14691.07 | 1074.38 | 1250.10 | 573.60 | 440.30 | 9.70 | 1240317.40 | 160200.00 |
| IP-Guided-Diff | 15132.65 | 632.80 | 1542.70 | 281.00 | 440.30 | 9.70 | 1226793.60 | 146676.20 |  |
| DIFUSCO | 15042.35 | 723.10 | 1487.50 | 336.20 | 439.20 | 10.80 | 1250323.40 | 170206.00 |  |
| FMIP (Ours) | 15452.75 | 312.70 | 1554.20 | 269.50 | 449.60 | 0.40 | 1221710.50 | 141593.10 |  |
|  | Rel.Imprv. | - | 50.58% | - | 4.09% | - | 95.88% | - | 3.47% |
| PS(600s) | SL | 15765.45 | 0.00 | 1783.00 | 40.70 | 450.00 | 0.00 | 1080117.40 | 0.00 |
| IP-Guided-Diff | 15765.45 | 0.00 | 1815.30 | 8.40 | 449.60 | 0.40 | 1080117.40 | 0.00 |  |
| DIFUSCO | 15761.85 | 3.60 | 1732.60 | 91.10 | 448.50 | 1.50 | 1120110.30 | 39992.90 |  |
| FMIP (Ours) | 15765.45 | 0.00 | 1816.50 | 7.20 | 450.00 | 0.00 | 1080117.40 | 0.00 |  |
|  | Rel.Imprv. | - | 0.00% | - | 14.29% | - | 0.00% | - | 0.00% |
| PMVB(600s) | SL | 15765.45 | 0.00 | 1690.90 | 132.80 | 449.90 | 0.10 | 1080117.40 | 0.00 |
| IP-Guided-Diff | 15765.45 | 0.00 | 1695.90 | 127.80 | 448.30 | 1.70 | 1080117.40 | 0.00 |  |
| DIFUSCO | 15745.25 | 20.20 | 1673.80 | 149.90 | 443.20 | 6.80 | 1114069.00 | 33951.60 |  |
| FMIP (Ours) | 15765.45 | 0.00 | 1695.90 | 127.80 | 450.00 | 0.00 | 1080117.40 | 0.00 |  |
|  | Rel.Imprv. | - | 0.00% | - | 0.00% | - | 100.00% | - | 0.00% |
| Apollo(800s) | SL | 15765.45 | 0.00 | 1823.70 | 0.00 | 448.20 | 1.80 | 1080117.40 | 0.00 |
| IP-Guided-Diff | 15765.45 | 0.00 | 1823.70 | 0.00 | 449.10 | 0.90 | 1080117.40 | 0.00 |  |
| DIFUSCO | 15431.43 | 334.02 | 1763.24 | 60.46 | 447.50 | 2.50 | 1087517.20 | 7399.80 |  |
| FMIP (Ours) | 15765.45 | 0.00 | 1823.70 | 0.00 | 450.00 | 0.00 | 1080117.40 | 0.00 |  |
|  | Rel.Imprv. | - | 0.00% | - | 0.00% | - | 100.00% | - | 0.00% |
| Gurobi(3600s) | — | 15765.45 | 0.00 | 1823.70 | 0.00 | 450.00 | 0.00 | 1080117.40 | 0.00 |
| Downstream<br>Solver | Training<br>Method | SC | LB | IP | MIPLIB |  |  |  |  |
| OBJ $\downarrow$ | GAP $\downarrow$ | OBJ $\downarrow$ | GAP $\downarrow$ | OBJ $\downarrow$ | GAP $\downarrow$ | OBJ $\downarrow$ | Rel.GAP $\downarrow$ |  |  |
| ND(400s) | SL | 401.00 | 0.30 | 719.67 | 13.77 | 14.62 | 0.88 | 1619541.12 | 6.21% |
| IP-Guided-Diff | 400.95 | 0.25 | 712.34 | 6.44 | 14.63 | 0.89 | 1609725.33 | 5.85% |  |
| DIFUSCO | 402.10 | 1.40 | 717.73 | 11.83 | 15.41 | 1.67 | 1625188.49 | 6.48% |  |
| FMIP (Ours) | 400.80 | 0.10 | 706.30 | 0.40 | 13.99 | 0.25 | 1604917.80 | 5.42% |  |
|  | Rel.Imprv. | - | 60.00% | - | 93.79% | - | 71.59% | - | 0.43% |
| PS(600s) | SL | 400.80 | 0.10 | 749.60 | 43.70 | 15.34 | 1.60 | 1612105.77 | 5.67% |
| IP-Guided-Diff | 400.75 | 0.05 | 725.34 | 19.44 | 15.21 | 1.47 | 1601226.94 | 5.31% |  |
| DIFUSCO | 401.30 | 0.60 | 714.11 | 8.21 | 14.73 | 0.99 | 1611839.42 | 5.65% |  |
| FMIP (Ours) | 400.70 | 0.00 | 705.90 | 0.00 | 13.92 | 0.18 | 1595308.15 | 4.84% |  |
|  | Rel.Imprv. | - | 100.00% | - | 100.00% | - | 74.86% | - | 0.47% |
| PMVB(600s) | SL | 422.10 | 21.40 | 706.10 | 0.20 | 15.39 | 1.65 | 1620951.88 | 6.34% |
| IP-Guided-Diff | 405.50 | 4.80 | 706.30 | 0.40 | 15.17 | 1.43 | 1614667.01 | 5.88% |  |
| DIFUSCO | 412.60 | 11.90 | 714.21 | 8.31 | 15.03 | 1.29 | 1628433.01 | 6.62% |  |
| FMIP (Ours) | 400.80 | 0.10 | 706.10 | 0.20 | 14.41 | 0.67 | 1608125.64 | 5.55% |  |
|  | Rel.Imprv. | - | 97.92% | - | 0.00% | - | 48.06% | - | 0.57% |
| Apollo(800s) | SL | 400.80 | 0.10 | 706.00 | 0.10 | 14.16 | 0.42 | 1604557.34 | 4.98% |
| IP-Guided-Diff | 405.50 | 4.80 | 705.95 | 0.05 | 14.33 | 0.59 | 1593821.50 | 4.55% |  |
| DIFUSCO | 403.30 | 2.60 | 709.32 | 3.42 | 14.01 | 0.27 | 1591488.19 | 4.43% |  |
| FMIP (Ours) | 400.70 | 0.00 | 705.90 | 0.00 | 13.74 | 0.00 | 1586142.71 | 4.15% |  |
|  | Rel.Imprv. | - | 100.00% | - | 100.00% | - | 100.00% | - | 0.28% |
| Gurobi(3600s) | — | 400.80 | 0.10 | 706.00 | 0.10 | 17.73 | 3.99 | 1522661.56 | 0.00% |

### 5.2 Overall Performance

We choose Tri-GCN as the backbone graph encoder and evaluate FMIP against baselines based on four downstream solvers (ND, PS, PMVB, and Apollo), with respective time limits of 400, 600, 600, and 800 seconds.
The time limit for Gurobi is set to 3600 seconds.
The results are reported in Table [1](https://arxiv.org/html/2507.23390v2#S5.T1). We observe that FMIP generally achieves the best performance among learning-based methods, substantially reducing the primal gap across all benchmarks and downstream solvers.
On simpler benchmarks (i.e., CA, GIS) where multiple methods find the optimal solution, FMIP performs on par with the best baselines.
However, on more challenging benchmarks (i.e., LB, IP), the benefits of our joint modeling and holistic guidance are clear, culminating in the significant performance gains in terms of both objective value and primal gap.

**Table 2: The model compatibility study, where we apply FMIP and baselines to different backbone graph neural networks. We report the objective value (OBJ) metric. The best results are given in bold, and the second-best values are underlined.**
| Downstream<br>Solver | Training<br>Method | Item Placement (IP) | Load Balancing (LB) |  |  |  |  |  |  |
| --- | --- | --- | --- | --- | --- | --- | --- | --- | --- |
| Bi-GCN | Tri-GCN | GAT | ClusterGCN | Bi-GCN | Tri-GCN | GAT | ClusterGCN |  |  |
| ND(400s) | SL | 14.85 | 14.62 | 14.65 | 14.69 | 719.54 | 719.67 | 719.20 | 720.15 |
| IP-Guided-Diff | 14.82 | 14.63 | 14.63 | 14.63 | 713.32 | 712.34 | 712.20 | 713.10 |  |
| DIFUSCO | 15.44 | 15.41 | 15.39 | 15.50 | 717.75 | 717.73 | 717.68 | 717.75 |  |
| FMIP | 14.11 | 13.99 | 13.95 | 14.06 | 707.24 | 706.30 | 706.50 | 706.50 |  |
| PS(600s) | SL | 15.38 | 15.34 | 15.33 | 15.32 | 750.34 | 749.60 | 749.26 | 749.69 |
| IP-Guided-Diff | 15.24 | 15.21 | 15.21 | 15.21 | 726.04 | 725.34 | 725.46 | 725.72 |  |
| DIFUSCO | 14.73 | 14.73 | 14.76 | 14.81 | 714.01 | 714.11 | 714.77 | 714.90 |  |
| FMIP | 13.92 | 13.92 | 13.94 | 13.92 | 705.88 | 705.90 | 706.10 | 706.00 |  |
| PMVB(600s) | SL | 15.36 | 15.39 | 15.36 | 15.44 | 706.50 | 706.10 | 706.30 | 706.30 |
| IP-Guided-Diff | 15.25 | 15.17 | 15.15 | 15.17 | 706.62 | 706.30 | 706.20 | 706.30 |  |
| DIFUSCO | 15.07 | 15.03 | 15.03 | 15.09 | 714.89 | 714.21 | 714.01 | 714.55 |  |
| FMIP | 14.45 | 14.41 | 14.42 | 14.39 | 706.32 | 706.10 | 706.10 | 706.20 |  |
| Apollo(800s) | SL | 14.19 | 14.16 | 14.11 | 14.19 | 706.30 | 706.25 | 706.20 | 706.30 |
| IP-Guided-Diff | 14.36 | 14.03 | 14.33 | 14.40 | 706.00 | 705.95 | 706.10 | 706.10 |  |
| DIFUSCO | 14.01 | 14.01 | 13.99 | 14.05 | 709.92 | 709.32 | 709.40 | 709.40 |  |
| FMIP | 13.72 | 13.74 | 13.74 | 13.75 | 705.90 | 705.90 | 705.90 | 706.00 |  |

### 5.3 Compatibility Analysis

The downstream-solver compatibility of FMIP is already validated in Table [1](https://arxiv.org/html/2507.23390v2#S5.T1) in Section [5.2](https://arxiv.org/html/2507.23390v2#S5.SS2), where FMIP consistently achieves the best performance across four different solvers.
To further study the model compatibility of FMIP, we apply FMIP and baselines to four different backbone graph encoders on IP and LB benchmarks: Bi-GCN (Gasse et al., 2019b), Tri-GCN (Gasse et al., 2019b), GAT (Brody et al., 2021), and ClusterGCN (Chiang et al., 2019)
We report the objective value (OBJ) metric in Table [2](https://arxiv.org/html/2507.23390v2#S5.T2), where FMIP maintains its superior performance across all backbone graph encoders.
This backbone independence demonstrates that FMIP acts as a powerful and general framework that can enhance a wide range of GNN-based models for various downstream solvers, highlighting the fundamental value of its joint distribution modeling on both continuous and integer variables.

**Table 3: The inference time (s) per MILP instance of FMIP and baselines. We only compare the neural model inference here and exclude the time that downstream solvers take.**
| Downstream Solver | Training Method | CA | GIS | MIS | FCMNF | SC | LB | IP | MIPLIB |
| --- | --- | --- | --- | --- | --- | --- | --- | --- | --- |
| ND, PS, PMVB | SL | 0.047 | 0.089 | 0.138 | 0.065 | 0.087 | 0.088 | 0.023 | 0.115 |
|  | DIFUSCO | 0.187 | 0.652 | 0.340 | 0.238 | 0.412 | 1.154 | 0.093 | 1.561 |
|  | IP-Guided-Diff | 0.213 | 0.712 | 0.452 | 0.374 | 0.533 | 1.781 | 0.117 | 2.153 |
|  | FMIP | 0.281 | 0.761 | 0.476 | 0.310 | 0.623 | 1.294 | 0.176 | 2.043 |
| Apollo | SL | 0.082 | 0.273 | 0.184 | 0.310 | 0.209 | 0.213 | 0.093 | 0.327 |
|  | DIFUSCO | 0.512 | 2.125 | 1.341 | 0.545 | 1.112 | 3.129 | 0.153 | 3.231 |
|  | IP-Guided-Diff | 0.612 | 2.173 | 1.549 | 0.634 | 1.268 | 3.151 | 0.167 | 4.185 |
|  | FMIP | 0.718 | 2.060 | 1.805 | 0.546 | 1.782 | 3.229 | 0.121 | 4.812 |

### 5.4 Inference Efficiency

We report the averaged inference time of FMIP and baselines to generate a heuristic solution for one MILP problem on different benchmarks.
Tri-GCN is adopted as the default backbone graph encoder.
Note that most downstream solvers (ND, PS, and PMVB) only need to obtain the heuristic solution based on one inference process from the model.
However, Apollo is a special solver that requires the model to infer multiple times for each single MILP instance, with the problem size gradually reduced and the solution iteratively refined.
In this experiment, we iterate 3 times for Apollo.

As shown in Table [3](https://arxiv.org/html/2507.23390v2#S5.T3), we report the inference time per MILP instance separately for ND/PS/PMVB and Apollo.
Despite the superior performance of FMIP, its generative process is also highly efficient.
The inference time of FMIP is comparable to that of other generative baselines (i.e., DIFUSCO and IP-Guided-Diff).
Crucially, this time represents a negligible fraction (often less than 1%) of the total time spent by the downstream solver, which can be hundreds or thousands of seconds.
FMIP thus provides its substantial performance improvements with minimal computational overhead, striking an excellent balance between solution quality and efficiency.

**Table 4: Ablation study of FMIP and its variants. The best results are given in bold, and the second-best values are underlined.**
|  | ND | PS | PMVB | Apollo |  |  |  |  |
| --- | --- | --- | --- | --- | --- | --- | --- | --- |
| Variant | OBJ $\downarrow$ | GAP $\downarrow$ | OBJ $\downarrow$ | GAP $\downarrow$ | OBJ $\downarrow$ | GAP $\downarrow$ | OBJ $\downarrow$ | GAP $\downarrow$ |
| Full FMIP | 13.99 | 0.25 | 13.92 | 0.18 | 14.41 | 0.67 | 13.74 | 0.00 |
| w/o Feasibility Guidance | 14.75 | 1.01 | 14.67 | 0.93 | 15.47 | 1.73 | 14.97 | 1.23 |
| w/o Objective Guidance | 14.65 | 0.91 | 14.76 | 1.02 | 22.32 | 8.58 | 14.53 | 0.79 |
| w/o Guidance | 14.58 | 0.84 | 13.98 | 0.24 | 14.45 | 0.71 | 14.11 | 0.37 |
| w/o Continuous | 14.43 | 0.69 | 14.28 | 0.54 | 14.95 | 1.21 | 14.45 | 0.71 |

### 5.5 In-Depth Analysis

Ablation Study.
We derive four variants of our proposed FMIP: (i) removing the objective guidance in Eq. [5](https://arxiv.org/html/2507.23390v2#S4.E5) (w/o Feasibility Guidance), (ii) removing the objective guidance in Eq. [5](https://arxiv.org/html/2507.23390v2#S4.E5) (w/o Objective Guidance), (iii) removing the entire holistic guidance mechanism (w/o Guidance), and (iv) reverting to an integer-only model by disabling continuous variable generation (w/o Continuous).
As shown in Table [4](https://arxiv.org/html/2507.23390v2#S5.T4), the full version of FMIP significantly outperforms all ablated variants, confirming that both the joint generation of all variables and the holistic guidance mechanism are critical to the performance gain.
Notably, the “w/o Continuous” variant suffers from the worst performance, providing direct evidence that capturing the coupled relationship between integer and continuous variables is essential for high-quality solution prediction in MILP.

Impact of Sampling Steps.
We analyze the effect of the number of sampling steps on the heuristic solution quality.
Specifically, we feed the heuristic solutions from different sampling steps into the downstream solvers, and report the primal gap metric (i.e., GAP) in Figure [3](https://arxiv.org/html/2507.23390v2#S5.F3).
We can observe that, at the beginning, more sampling steps lead to better solutions as the model has more opportunities for refinement.
However, the quality of heuristic solution starts to worsen when a plethora of steps are performed.
We attribute this phenomenon to the trade-off between generative refinement and search diversity.
An excessive number ofsampling steps might cause the guidance to make the predicted variable distribution overly sharp, which reduces the diversity of the solutions being explored and can prematurely trap the search in a local optimum.
This insight reveals a nuanced interplay between generative refinement and search diversity, which is critical for optimizing performance.

Figure: Figure 3: The performance of FMIP with downstream solvers w.r.t. different sampling steps. We report the absolute primal gap (GAP) as the metric.
Refer to caption: pic/iter_compare.png

## 6 Conclusion

In this paper, we propose Joint Continuous-Integer Flow for Mixed Integer Linear Programming (i.e., FMIP), which is the first generative framework to model the joint distribution of both integer and continuous variables.
Based on this, we further design the holistic guidance mechanism that uses instance-specific objective and constraint information to steer the generative process toward higher-quality heuristic solutions.
As a powerful generative learning framework, FMIP is fully compatible with arbitrary backbone graph encoders and various downstream solvers, demonstrating its great potential for real-world MILP applications.
Extensive experiments show that FMIP sets new state-of-the-art performance for learning-based MILP heuristics, reducing the primal gap by 41.34% on average compared to the best baseline method. Two key directions for future work are: enhancing the graph representation to better capture task-specific features, and developing a customized solver tailored to our FMIP framework to optimize solution quality and computational efficiency.

## References

- Albergo and Vanden-Eijnden [2023]
Michael S Albergo and Eric Vanden-Eijnden.
Building normalizing flows with stochastic interpolants.
*International Conference on Learning Representations*, 2023.
- Bertsimas and Tsitsiklis [1997]
Dimitris Bertsimas and John N. Tsitsiklis.
*Introduction to Linear Optimization*.
Athena Scientific, 1997.
- Brody et al. [2021]
Shaked Brody, Uri Alon, and Eran Yahav.
How attentive are graph attention networks?
*ArXiv*, abs/2105.14491, 2021.
URL [https://api.semanticscholar.org/CorpusID:235254358](https://api.semanticscholar.org/CorpusID:235254358).
- Chen et al. [2025]
Yanguang Chen, Wenzhi Gao, Wanyu Zhang, Dongdong Ge, Huikang Liu, and Yinyu Ye.
Data-driven mixed integer optimization through probabilistic multi-variable branching, 2025.
URL [https://arxiv.org/abs/2305.12352](https://arxiv.org/abs/2305.12352).
- Chiang et al. [2019]
Wei-Lin Chiang, Xuanqing Liu, Si Si, Yang Li, Samy Bengio, and Cho-Jui Hsieh.
Cluster-gcn: An efficient algorithm for training deep and large graph convolutional networks.
In *Proceedings of the 25th ACM SIGKDD International Conference on Knowledge Discovery & Data Mining*, KDD ’19, page 257–266, New York, NY, USA, 2019. Association for Computing Machinery.
ISBN 9781450362016.
doi: 10.1145/3292500.3330925.
URL [https://doi.org/10.1145/3292500.3330925](https://doi.org/10.1145/3292500.3330925).
- Colombi et al. [2017]
Marco Colombi, Renata Mansini, and Martin Savelsbergh.
The generalized independent set problem: Polyhedral analysis and solution approaches.
*European Journal of Operational Research*, 260(1):41–55, 2017.
- Della Croce and Paschos [2014]
Fabio Della Croce and Vangelis Th. Paschos.
Mixed integer linear programming models for combinatorial optimization problems.
In *Concepts of Combinatorial Optimization*, pages 101–133. John Wiley & Sons, Inc., 2014.
doi: 10.1002/9781119005216.ch5.
- Feng et al. [2024]
Shengyu Feng, Zhiqing Sun, and Yiming Yang.
DIFUSCO-LNS: Diffusion-guided large neighbourhood search for integer linear programming, 2024.
URL [https://openreview.net/forum?id=9QV7Q9gKl9](https://openreview.net/forum?id=9QV7Q9gKl9).
- Fey and Lenssen [2019]
Matthias Fey and Jan Eric Lenssen.
Fast graph representation learning with pytorch geometric.
*arXiv preprint arXiv:1903.02428*, 2019.
- Gasse et al. [2019a]
Maxime Gasse, Didier Chételat, Nicola Ferroni, Laurent Charlin, and Andrea Lodi.
Exact combinatorial optimization with graph convolutional neural networks.
*Advances in neural information processing systems*, 32, 2019a.
- Gasse et al. [2019b]
Maxime Gasse, Didier Chételat, Nicola Ferroni, Laurent Charlin, and Andrea Lodi.
Exact combinatorial optimization with graph convolutional neural networks.
*Advances in neural information processing systems*, 32, 2019b.
- Gasse et al. [2022]
Maxime Gasse, Simon Bowly, Quentin Cappart, Jonas Charfreitag, Laurent Charlin, Didier Chételat, Antonia Chmiela, Justin Dumouchelle, Ambros Gleixner, Aleksandr M. Kazachkov, Elias Khalil, Pawel Lichocki, Andrea Lodi, Miles Lubin, Chris J. Maddison, Morris Christopher, Dimitri J. Papageorgiou, Augustin Parjadis, Sebastian Pokutta, Antoine Prouvost, Lara Scavuzzo, Giulia Zarpellon, Linxin Yang, Sha Lai, Akang Wang, Xiaodong Luo, Xiang Zhou, Haohan Huang, Shengcheng Shao, Yuanming Zhu, Dong Zhang, Tao Quan, Zixuan Cao, Yang Xu, Zhewei Huang, Shuchang Zhou, Chen Binbin, He Minggui, Hao Hao, Zhang Zhiyu, An Zhiwu, and Mao Kun.
The machine learning for combinatorial optimization competition (ml4co): Results and insights.
In Douwe Kiela, Marco Ciccone, and Barbara Caputo, editors, *Proceedings of the NeurIPS 2021 Competitions and Demonstrations Track*, volume 176 of *Proceedings of Machine Learning Research*, pages 220–231. PMLR, 06–14 Dec 2022.
URL [https://proceedings.mlr.press/v176/gasse22a.html](https://proceedings.mlr.press/v176/gasse22a.html).
- Ge et al. [2024]
Dongdong Ge, Qi Huangfu, Zizhuo Wang, Jian Wu, and Yinyu Ye.
Cardinal optimizer (copt) user guide, 2024.
URL [https://arxiv.org/abs/2208.14314](https://arxiv.org/abs/2208.14314).
- Gleixner et al. [2021a]
Ambros Gleixner, Gregor Hendel, Gerald Gamrath, Tobias Achterberg, Michael Bastubbe, Timo Berthold, Philipp Christophel, Kati Jarck, Thorsten Koch, Jeff Linderoth, et al.
Miplib 2017: data-driven compilation of the 6th mixed-integer programming library.
*Mathematical Programming Computation*, 13(3):443–490, 2021a.
- Gleixner et al. [2021b]
Ambros M. Gleixner, Gregor Hendel, Gerald Gamrath, Tobias Achterberg, Michael Bastubbe, Timo Berthold, Philipp Christophel, Kati Jarck, Thorsten Koch, Jeff T. Linderoth, Marco E. Lübbecke, Hans D. Mittelmann, Derya B. Özyurt, Ted K. Ralphs, Domenico Salvagnin, and Yuji Shinano.
Miplib 2017: data-driven compilation of the 6th mixed-integer programming library.
*Math. Program. Comput.*, 13(3):443–490, 2021b.
URL [https://doi.org/10.1007/s12532-020-00194-3](https://doi.org/10.1007/s12532-020-00194-3).
- Gupta et al. [2020]
Prateek Gupta, Maxime Gasse, Elias Khalil, Pawan Mudigonda, Andrea Lodi, and Yoshua Bengio.
Hybrid models for learning to branch.
*Advances in neural information processing systems*, 33:18087–18097, 2020.
- Gupta et al. [2022]
Prateek Gupta, Elias B Khalil, Didier Chetélat, Maxime Gasse, Yoshua Bengio, Andrea Lodi, and M Pawan Kumar.
Lookback for learning to branch.
*arXiv preprint arXiv:2206.14987*, 2022.
- Gurobi Optimization, LLC [2024]
Gurobi Optimization, LLC.
Gurobi Optimizer Reference Manual, 2024.
URL [https://www.gurobi.com](https://www.gurobi.com).
- Han et al. [2023]
Qingyu Han, Linxin Yang, Qian Chen, Xiang Zhou, Dong Zhang, Akang Wang, Ruoyu Sun, and Xiaodong Luo.
A gnn-guided predict-and-search framework for mixed-integer linear programming, 2023.
URL [https://arxiv.org/abs/2302.05636](https://arxiv.org/abs/2302.05636).
- He et al. [2015]
Xiaozheng He, Hong Zheng, and Srinivas Peeta.
Model and a solution algorithm for the dynamic resource allocation problem for large-scale transportation network evacuation.
*Transportation Research Procedia*, 7:441–458, 2015.
ISSN 2352-1465.
doi: https://doi.org/10.1016/j.trpro.2015.06.023.
URL [https://www.sciencedirect.com/science/article/pii/S2352146515000915](https://www.sciencedirect.com/science/article/pii/S2352146515000915).
21st International Symposium on Transportation and Traffic Theory Kobe, Japan, 5-7 August, 2015.
- Hewitt et al. [2010]
Mike Hewitt, George L Nemhauser, and Martin WP Savelsbergh.
Combining exact and heuristic approaches for the capacitated fixed-charge network flow problem.
*INFORMS journal on computing*, 22(2):314–325, 2010.
- Hu et al. [2024]
Xinyi Hu, Jasper Lee, Jimmy Lee, and Peter Stuckey.
Multi-stage predict+ optimize for (mixed integer) linear programs.
*Advances in Neural Information Processing Systems*, 37:64794–64827, 2024.
- Huang et al. [2024a]
Taoan Huang, Aaron M Ferber, Arman Zharmagambetov, Yuandong Tian, and Bistra Dilkina.
Contrastive predict-and-search for mixed integer linear programs.
In *Forty-first International Conference on Machine Learning*, 2024a.
- Huang et al. [2024b]
Taoan Huang, Aaron M Ferber, Arman Zharmagambetov, Yuandong Tian, and Bistra Dilkina.
Contrastive predict-and-search for mixed integer linear programs.
In Ruslan Salakhutdinov, Zico Kolter, Katherine Heller, Adrian Weller, Nuria Oliver, Jonathan Scarlett, and Felix Berkenkamp, editors, *Proceedings of the 41st International Conference on Machine Learning*, volume 235 of *Proceedings of Machine Learning Research*, pages 19757–19771. PMLR, 21–27 Jul 2024b.
URL [https://proceedings.mlr.press/v235/huang24f.html](https://proceedings.mlr.press/v235/huang24f.html).
- Huang et al. [2022]
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.
- Ivanov et al. [2022]
Boyan Ivanov, Desislava Nikolova, Elisaveta Kirilova, and Rayka Vladova.
A milp approach of optimal design of a sustainable combined dairy and biodiesel supply chain using dairy waste scum generated from dairy production.
*Computers & Chemical Engineering*, 166:107976, 2022.
- Khalil et al. [2016]
Elias B. Khalil, Pierre Le Bodic, Le Song, George Nemhauser, and Bistra Dilkina.
Learning to branch in mixed integer programming.
In Dale Schuurmans and Michael Wellman, editors, *Proceedings of the Thirtieth AAAI Conference on Artificial Intelligence*, pages 724–731, United States of America, 2016. Association for the Advancement of Artificial Intelligence (AAAI).
URL [http://www.aaai.org/Conferences/AAAI/aaai16.php](http://www.aaai.org/Conferences/AAAI/aaai16.php).
AAAI Conference on Artificial Intelligence 2016, AAAI 2016 ; Conference date: 12-02-2016 Through 17-02-2016.
- Kipf [2016]
TN Kipf.
Semi-supervised classification with graph convolutional networks.
*arXiv preprint arXiv:1609.02907*, 2016.
- Kratica et al. [2014]
Jozef Kratica, Djordje Dugošija, and Aleksandar Savić.
A new mixed integer linear programming model for the multi level uncapacitated facility location problem.
*Applied Mathematical Modelling*, 38(7):2118–2129, 2014.
ISSN 0307-904X.
doi: https://doi.org/10.1016/j.apm.2013.10.012.
URL [https://www.sciencedirect.com/science/article/pii/S0307904X13006240](https://www.sciencedirect.com/science/article/pii/S0307904X13006240).
- Lin et al. [2025]
Haowei Lin, Shanda Li, Haotian Ye, Yiming Yang, Stefano Ermon, Yitao Liang, and Jianzhu Ma.
Tfg-flow: Training-free guidance in multimodal generative flow, 2025.
URL [https://arxiv.org/abs/2501.14216](https://arxiv.org/abs/2501.14216).
- Lin et al. [2024]
Jiacheng Lin, Meng XU, Zhihua Xiong, and Huangang Wang.
CAMBranch: Contrastive learning with augmented MILPs for branching.
In *The Twelfth International Conference on Learning Representations*, 2024.
URL [https://openreview.net/forum?id=K6kt50zAiG](https://openreview.net/forum?id=K6kt50zAiG).
- Lipman et al. [2023]
Yaron Lipman, Ricky TQ Chen, Heli Ben-Hamu, Maximilian Nickel, and Matt Le.
Flow matching for generative modeling.
*International Conference on Learning Representations*, 2023.
- Liu et al. [2025]
Haoyang Liu, Jie Wang, Zijie Geng, Xijun Li, Yuxuan Zong, Fangzhou Zhu, Jianye Hao, and Feng Wu.
Apollo-milp: An alternating prediction-correction neural solving framework for mixed-integer linear programming, 2025.
URL [https://arxiv.org/abs/2503.01129](https://arxiv.org/abs/2503.01129).
- Liu et al. [2023]
Xingchao Liu, Chengyue Gong, and Qiang Liu.
Flow straight and fast: Learning to generate and transfer data with rectified flow.
*International Conference on Learning Representations*, 2023.
- Miehling et al. [2023]
Sebastian Miehling, Andreas Hanel, Jerry Lambert, Sebastian Fendt, and Hartmut Spliethoff.
Energy system optimisation using (mixed integer) linear programming, 2023.
URL [https://arxiv.org/abs/2308.01882](https://arxiv.org/abs/2308.01882).
- Nair et al. [2020]
Vinod Nair, Sergey Bartunov, Felix Gimeno, Ingrid Von Glehn, Pawel Lichocki, Ivan Lobov, Brendan O’Donoghue, Nicolas Sonnerat, Christian Tjandraatmadja, Pengming Wang, et al.
Solving mixed integer programs using neural networks.
*arXiv preprint arXiv:2012.13349*, 2020.
- Nichol and Dhariwal [2021]
Alexander Quinn Nichol and Prafulla Dhariwal.
Improved denoising diffusion probabilistic models.
In *International Conference on Machine Learning*, pages 8162–8171, 2021.
- Paszke et al. [2019]
Adam Paszke, Sam Gross, Francisco Massa, Adam Lerer, James Bradbury, Gregory Chanan, Trevor Killeen, Zeming Lin, Natalia Gimelshein, Luca Antiga, et al.
Pytorch: An imperative style, high-performance deep learning library.
*Advances in neural information processing systems*, 32, 2019.
- Scavuzzo et al. [2022]
Lara Scavuzzo, Feng Chen, Didier Chetelat, Maxime Gasse, Andrea Lodi, Neil Yorke-Smith, and Karen Aardal.
Learning to branch with tree mdps.
In S. Koyejo, S. Mohamed, A. Agarwal, D. Belgrave, K. Cho, and A. Oh, editors, *Advances in Neural Information Processing Systems*, volume 35, pages 18514–18526. Curran Associates, Inc., 2022.
URL [https://proceedings.neurips.cc/paper_files/paper/2022/file/756d74cd58592849c904421e3b2ec7a4-Paper-Conference.pdf](https://proceedings.neurips.cc/paper_files/paper/2022/file/756d74cd58592849c904421e3b2ec7a4-Paper-Conference.pdf).
- Sun and Yang [2023]
Zhiqing Sun and Yiming Yang.
DIFUSCO: Graph-based diffusion solvers for combinatorial optimization.
In *Thirty-seventh Conference on Neural Information Processing Systems*, 2023.
URL [https://openreview.net/forum?id=JV8Ff0lgVV](https://openreview.net/forum?id=JV8Ff0lgVV).
- Tang et al. [2020]
Yunhao Tang, Shipra Agrawal, and Yuri Faenza.
Reinforcement learning for integer programming: Learning to cut.
In *International conference on machine learning*, pages 9367–9376. PMLR, 2020.
- Vaswani et al. [2017]
Ashish Vaswani, Noam Shazeer, Niki Parmar, Jakob Uszkoreit, Llion Jones, Aidan N Gomez, Łukasz Kaiser, and Illia Polosukhin.
Attention is all you need.
*Advances in neural information processing systems*, 2017.
- Veličković et al. [2017]
Petar Veličković, Guillem Cucurull, Arantxa Casanova, Adriana Romero, Pietro Lio, and Yoshua Bengio.
Graph attention networks.
*arXiv preprint arXiv:1710.10903*, 2017.
- Ye et al. [2024]
Huigen Ye, Hua Xu, and Hongyan Wang.
Light-milpopt: Solving large-scale mixed integer linear programs with lightweight optimizer and small-scale training dataset.
In *The Twelfth International Conference on Learning Representations*, 2024.
- Zarpellon et al. [2021]
Giulia Zarpellon, Jason Jo, Andrea Lodi, and Yoshua Bengio.
Parameterizing branch-and-bound search trees to learn branching policies.
In *Thirty-Fifth AAAI Conference on Artificial Intelligence, AAAI 2021, Thirty-Third Conference on Innovative Applications of Artificial Intelligence, IAAI 2021, The Eleventh Symposium on Educational Advances in Artificial Intelligence, EAAI 2021, Virtual Event, February 2-9, 2021*, pages 3931–3939. AAAI Press, 2021.
doi: 10.1609/AAAI.V35I5.16512.
URL [https://doi.org/10.1609/aaai.v35i5.16512](https://doi.org/10.1609/aaai.v35i5.16512).
- Zeng et al. [2024]
Hao Zeng, Jiaqi Wang, Avirup Das, Junying He, Kunpeng Han, Haoyuan Hu, and Mingfei Sun.
Effective generation of feasible solutions for integer programming via guided diffusion, 2024.
URL [https://arxiv.org/abs/2406.12349](https://arxiv.org/abs/2406.12349).
- Zhang et al. [2024]
Changwen Zhang, Wenli Ouyang, Hao Yuan, Liming Gong, Yong Sun, Ziao Guo, Zhichen Dong, and Junchi Yan.
Towards imitation learning to branch for MIP: A hybrid reinforcement learning based sample augmentation approach.
In *The Twelfth International Conference on Learning Representations*, 2024.
URL [https://openreview.net/forum?id=NdcQQ82mfy](https://openreview.net/forum?id=NdcQQ82mfy).
- Zhou et al. [2025]
Chenyu Zhou, Tianyi Xu, Jianghao Lin, and Dongdong Ge.
Steporlm: A self-evolving framework with generative process supervision for operations research language models, 2025.
URL [https://arxiv.org/abs/2509.22558](https://arxiv.org/abs/2509.22558).

## Appendix A Graph Representation Details

Building on previous work, we select a similar graph structure. To more naturally accommodates our joint generation purpose, we explicitly distinguishing three types of nodes: Integer Variables (Ivar), Continuous Variables (Cvar), and Constraints (con).Edges, denoted as $e_{ij}$, are drawn between the variable nodes (both Ivar and Cvar) and the constraint nodes to capture their algebraic relationships within the MILP formulation.The edge weights $e_{ij}$
correspond to the coefficients of variables in their associated constraints. The features of nodes
and edges of a MILP graph $\boldsymbol{G}=(\mathcal{V}_{\text{Ivar}},\mathcal{V}_{\text{Cvar}},\mathcal{V}_{\text{con}},\mathcal{E})$ are:

- •
Ivar/Cvar Nodes ($\mathcal{V}_{\text{Ivar}}/\mathcal{V}_{\text{Cvar}}$): The feature vector for variable $v_{i}\in\mathcal{V}_{\text{Ivar}}\cup\mathcal{V}_{\text{Cvar}}$ is defined as:
$\boldsymbol{h}_{v_{i}}=\left[\underbrace{\boldsymbol{w}_{i}}_{\text{objective coefficient}},\,\underbrace{\boldsymbol{l}_{i}}_{\text{lower bound}},\,\underbrace{\boldsymbol{u}_{i}}_{\text{upper bound}},\,\underbrace{\boldsymbol{hlb}_{i}}_{\text{if }l_{i}>-\infty},\,\underbrace{\boldsymbol{hub}_{i}}_{\text{if }u_{i}<+\infty}\right]\in\mathbb{R}^{5},$
where $\boldsymbol{w},\boldsymbol{l},\boldsymbol{u}$ are defined in Eq. [1](https://arxiv.org/html/2507.23390v2#S3.E1).
- •
Constraint Nodes ($\mathcal{V}_{\text{con}}$): The feature for constraint $con_{j}$ is scalar $b_{j}\in\mathbb{R}$, representing its right-hand side constant.
- •
Edges ($\mathcal{E}$): We construct a sparse bipartite graph based on the nonzero entries of the constraint matrix $\mathbf{A}$. Specifically, an edge is added between variable node $j$ and constraint node $i$ if and only if $\mathbf{A}[i,j]>0$. The edge feature is set as the corresponding coefficient, i.e., the edge from $j$ to $i$ carries weight $\mathbf{A}[i,j]$.

When processing with a graph neural network, the state in generative model $\boldsymbol{x}_{t}=(\boldsymbol{d}_{t},\boldsymbol{c}_{t})$ is incorporated into the features of the corresponding variable nodes in $\boldsymbol{G}$. Specifically, for each variable node $v_{i}\in\mathcal{V}_{\text{Ivar}}\cup\mathcal{V}_{\text{Cvar}}$, its static feature vector $\boldsymbol{h}_{v_{i}}$ from the MILP graph is augmented with its current value from the solution vector $\boldsymbol{x}_{t}$, forming the time-dependent feature vector $\boldsymbol{h}_{v_{i},t}=[\boldsymbol{h}_{v_{i}},\boldsymbol{x}_{t}^{(i)}]\in\mathbb{R}^{6}$. This implies that the static variable features $\boldsymbol{h}_{v_{i}}$ are 5-dimensional. Constraint node features are typically static.
Figure [4](https://arxiv.org/html/2507.23390v2#A1.F4) provides a simple example illustrating the structure of a solution graph derived from an MILP instance.

Figure: (a) Problem
Refer to caption: pic/example.png

## Appendix B Graph Neural Networks

The input to our neural network is a tuple
$(\boldsymbol{G}_{t},t),$
where $\boldsymbol{G}_{t}$ represents the solution graph as described in Section [A](https://arxiv.org/html/2507.23390v2#A1). This graph encodes the relationships between variables and constraints, effectively capturing the structure of the MILP problem. The scalars $t$ and $t$ are temporal parameters associated with the integer and continuous variables, respectively, with both $t\in[0,1]$. These time parameters are introduced to enable a dynamic, progressive refinement of the solution space during the training process, aligning with the flow matching paradigm where time plays a critical role in guiding the evolution of the generated solution. We’ll details our neural networks structure in the following subsection.

### B.1 TripartiteGCNConv Layer

We adopt the following mathematical notations for the following parts: $\text{Linear}_{\cdot}(\cdot)$ denotes a learnable affine transformation (linear layer), $\text{LN}_{\cdot}(\cdot)$ represents layer normalization, and $\text{GELU}(\cdot)$ indicates Gaussian Error Linear Unit activation. The element-wise sigmoid function is written as $\sigma(\cdot)$, while $\odot$ signifies element-wise multiplication. Vector concatenation is denoted by $[\cdot\ ;\ \cdot]$, and $\mathbb{I}_{\text{res}}\in\{0,1\}$ serves as an indicator variable for residual connection usage.

The output feature $h_{v}^{\prime}$ of a target node $v$ is computed as:

$$ $h_{v}^{\prime}=\text{LN}_{\text{out}}\left(\text{Linear}_{o2}\left(\text{GELU}\left(\text{Linear}_{o1}\left([\text{LN}_{\text{PC}}(g_{v}\odot a_{1}+(1-g_{v})\odot a_{2})\ ;\ h_{v}]\right)\right)\right)\right)+\mathbb{I}_{\text{res}}\cdot h_{v},$ $$

where the aggregated messages from source node types $k=1,2$ are defined as:

$$ $a_{k}=\sum_{u\in\mathcal{N}_{k}(v)}\text{Linear}_{\text{final}}\left(\text{GELU}\left(\text{LN}_{\text{PN}}\left(\text{Linear}_{t}(h_{v})+\text{Linear}_{s_{k}}(h_{u}^{(k)})+\text{Linear}_{e}(e_{uv}^{(k)})\right)\right)\right),$ $$

and the gating vector is computed by:
$g_{v}=\sigma\left(\text{Linear}_{fg}([a_{1};a_{2}])\right)$.

### B.2 Model Structure

Overview  We propose a graph-based neural architecture tailored for mixed-variable optimization problems. The model captures the interactions among integer variables, continuous variables, and constraints using iterative message passing with temporal conditioning.

Initial Embeddings  We define three types of nodes:

- •
Integer variables ($\mathcal{V}_{\text{Ivar}}$): $n_{\text{Ivar}}$ nodes, with features $\mathbf{X}_{\text{Ivar}}\in\mathbb{R}^{n_{\text{Ivar}}\times d_{\text{Ivar}}}$
- •
Continuous variables ($\mathcal{V}_{\text{Cvar}}$): $n_{\text{Cvar}}$ nodes, with features $\mathbf{X}_{\text{Cvar}}\in\mathbb{R}^{n_{\text{Cvar}}\times d_{\text{Cvar}}}$
- •
Constraints ($\mathcal{C}_{\text{con}}$): $n_{\text{con}}$ nodes, with features $\mathbf{X}_{\text{con}}\in\mathbb{R}^{n_{\text{con}}\times d_{\text{con}}}$

Initial node embeddings are computed via type-specific MLPs:

$$ $\mathbf{h}_{\text{Ivar}}^{(0)}=\text{MLP}_{\text{Ivar}}(\mathbf{X}_{\text{Ivar}}),\quad\mathbf{h}_{\text{Cvar}}^{(0)}=\text{MLP}_{\text{Cvar}}(\mathbf{X}_{\text{Cvar}}),\quad\mathbf{h}_{\text{con}}^{(0)}=\text{MLP}_{\text{con}}(\mathbf{X}_{\text{con}}),$ $$

where each embedding lies in $\mathbb{R}^{\cdot\times h}$ with $h$ denoting the shared hidden dimension.

The time embedding $\mathbf{e}_{t}\in\mathbb{R}^{h}$ is initialized using positional encoding[Vaswani et al., 2017]. This embedding is broadcasted and added to the variable node features at every iteration.

Message Passing Layers  For each layer $\ell=1,\dots,L$, we update the node representations via a message passing scheme over a tripartite graph. The relevant edge sets are defined as:

- •
$\mathcal{E}_{\text{Ivar2con}}\subseteq\mathcal{V}_{\text{Ivar}}\times\mathcal{C}_{\text{con}}$: integer-to-constraint edges
- •
$\mathcal{E}_{\text{Cvar2con}}\subseteq\mathcal{V}_{\text{Cvar}}\times\mathcal{C}_{\text{con}}$: continuous-to-constraint edges
- •
Reverse edges $\mathcal{E}_{\text{con2Ivar}},\mathcal{E}_{\text{con2Cvar}}$: constraint-to-variable edges for backward message flow

The node updates at layer $\ell$ are:

$$ $\displaystyle\mathbf{h}_{t}^{(\ell)}$ $\displaystyle=\text{MLP}_{t}^{(\ell)}\big(\mathbf{e}_{t}),$ $\displaystyle\mathbf{h}_{\text{con}}^{(\ell)}$ $\displaystyle=\mathbf{h}_{\text{con}}^{(\ell-1)}+\mathbf{h}_{t}^{(\ell)}+\text{MLP}_{\text{con}}^{(\ell)}\big(\text{TriConv}\big(\mathbf{h}_{\text{Ivar}}^{(\ell-1)},\mathbf{h}_{\text{Cvar}}^{(\ell-1)},\mathbf{h}_{\text{con}}^{(\ell-1)},\mathcal{E}_{\text{Ivar2con}},\mathcal{E}_{\text{Cvar2con}}\big)\big),$ $\displaystyle\mathbf{h}_{\ast}^{(\ell)}$ $\displaystyle=\mathbf{h}_{\ast}^{(\ell-1)}+\mathbf{h}_{t}^{(\ell)}+\text{MLP}_{\ast}^{(\ell)}\big(\text{BiConv}_{\ast}(\mathbf{h}_{\text{con}}^{(\ell)},\mathbf{h}_{\ast}^{(\ell-1)},\mathcal{E}_{\text{con2}\ast})\big),\quad\ast\in\{\text{Ivar},\text{Cvar}\}.$ $$

TriConv is a tripartite message aggregation module that aggregates signals from both integer and continuous variable nodes into constraint nodes. BiConv has the same architecture as TriConv, but applies bipartite message passing from constraints back to one variable type.

Output Layer  After $L$ layers of message passing, we apply type-specific MLP heads to predict outputs from the final node embeddings:

$$ $\mathbf{o}_{\text{Ivar}}=\text{MLP}_{\text{out}}^{\text{Ivar}}(\mathbf{h}_{\text{Ivar}}^{(L)}),\quad\mathbf{o}_{\text{Cvar}}=\text{MLP}_{\text{out}}^{\text{Cvar}}(\mathbf{h}_{\text{Cvar}}^{(L)}),$ $$

where $\mathbf{o}_{\text{Ivar}}\in\mathbb{R}^{n_{\text{Ivar}}\times d_{\text{out}}^{\text{Ivar}}}$ and $\mathbf{o}_{\text{Cvar}}\in\mathbb{R}^{n_{\text{Cvar}}\times d_{\text{out}}^{\text{Cvar}}}$ are the prediction outputs for integer and continuous variables, respectively.

## Appendix C Implementation Details

We’ll present the implementation details in this section, which contains infrastructure and hyperparameter in our proposed methods.

### C.1 Training Details

Training labels for FMIP and other baselines were generated by solving the MILP instances using COPT [Ge et al., 2024], an outperforming MILP solver [Gleixner et al., 2021b], with a time limit of 3600 seconds per instance and 12 threads. The training loss for the GNN baseline is binary cross entropy for the prediction of binary variables, and the training loss for FMIP is described in Eq. [4](https://arxiv.org/html/2507.23390v2#S4.E4).
For MIPLIB, the training set is collected from all other benchmarks due to its diverse problem types. For the other datasets, the training set is split from the entire data with a 9:1 ratio.

### C.2 Infrastructure

Our flow model is implemented in PyTorch[Paszke et al., 2019] and PyTorch Geometric[Fey and Lenssen, 2019] and trained on a single NVIDIA H100 GPU.For CPU, we we use 12 cores of
an Intel Xeon Platinum 8469C at 2.60 GHz CPU with 512 GB RAM. We select Gurobi[Gurobi Optimization, LLC, 2024]
, which is the most famous and well-implement MILP solver, in Decoding Strategy and baseline methods due to its high efficiency. For all experiments, we set 12 threads for backend solver, Gurobi, which is a standard setting [Gleixner et al., 2021b].

**Table 5: Experimental Parameters Comparison**
| Method | ND | PS | PMVB | Apollo |
| --- | --- | --- | --- | --- |
| (Parameters) | $[K,\alpha]$ | $[k_{0},k_{1},\Delta]$ | $[\delta,\tau]$ | $[k_{0},k_{1},\Delta,K]$ |
| Cauctions | $[50,0.1]$ | $[0.3,0.06,0.3]$ | $[0.7,0.9]$ | $[0.3,0.06,0.3,2]$ |
| GISP | $[50,0.1]$ | $[0.2,0.02,0.2]$ | $[0.7,0.9]$ | $[0.2,0.02,0.2,2]$ |
| Independent Set | $[50,0.1]$ | $[0.3,0.2,0.3]$ | $[0.7,0.9]$ | $[0.3,0.2,0.3,2]$ |
| FCMNF | $[50,0.1]$ | $[0.3,0.03,0.2]$ | $[0.7,0.9]$ | $[0.3,0.03,0.2,2]$ |
| Item Placement | $[50,0.1]$ | $[0.3,0.08,0.4]$ | $[0.7,0.9]$ | $[0.3,0.08,0.4,2]$ |
| Load Balancing | $[50,0.1]$ | $[0.2,0.2,0.2]$ | $[0.7,0.9]$ | $[0.2,0.2,0.2,2]$ |
| Set Covering | $[50,0.1]$ | $[0.3,0.04,0.2]$ | $[0.7,0.9]$ | $[0.3,0.04,0.2,2]$ |

### C.3 Downstream Solvers

Table [5](https://arxiv.org/html/2507.23390v2#A3.T5) summarizes the hyperparameter settings used by each downstream solver across different problem instances. The meaning of each parameter (e.g., $K$, $\alpha$, $k_{0}$, $\Delta$, etc.) is consistent with the notation defined in Appendix [D](https://arxiv.org/html/2507.23390v2#A4), where detailed descriptions of downstream solvers and their hyperparameters are provided. For a fair comparison, we use the same parameters for the downstream solvers when evaluating FMIP and all baselines.

### C.4 FMIP

**Table 6: Training and Inference Configuration for FMIP**
| Parameter Description | Value |
| --- | --- |
| Training epochs | 300 |
| Learning rate | 2e-4 |
| Weight decay | 1e-4 |
| Learning rate scheduler | cosine-decay |
| GNN layers | 12 |
| Hidden dimension | 64 |
| Inference schedule | cosine |
| Inference Steps | 12 |
| integer guidance temperature ($\psi$) | 0.1, 1, 10 |
| continuous guidance stepsize ($\rho$) | 0.1 |

This section outlines the key hyperparameters employed during model training and inference, covering: (1) general training configurations (e.g., learning rate, batch size), (2) neural network architecture specifications, and (3) inference-specific settings such as sampling strategies and guidance control. The fixed default parameters are provided in Table [6](https://arxiv.org/html/2507.23390v2#A3.T6), while two critical parameters, Batch Size and Guidance Temperature, are dynamically adjusted based on dataset characteristics, detailed in following parts.

Choice of Batch Size  In practice, Batch size is dynamically adjusted according to the problem scale and available GPU memory. We adopt the largest feasible batch size that fits into memory to maximize hardware utilization and training efficiency.

Guidance Temperature  We temperature coeffient ($\psi$) is introduced to control the strength of guidance during inference, where lower temperature represents higher strength. This parameter is selected based on validation performance. Specifically, we set the temperature to 0.1 for the problem Load Balancing and Item Placement, 10 for the Set Covering problem, and 1 for all other problem types.

## Appendix D Downstream Solver Algorithm

In this section, we introduction detailed algorithms employed in our experiments, including Neural Diving, Predict-and-Search, PMVB and Apollo-MILP.

### D.1 Neural Diving

For a new instance, the neural network first generates multiple partial variable assignments based on predicted probability distributions. Let $\mathcal{V}$ denote the set of variables. The algorithm fixes a subset $\mathcal{S}\subset\mathcal{V}$ where:

$$ $|\mathcal{S}|=\alpha\cdot|\mathcal{V}|\quad(\alpha\in(0,1))$ $$

and generates $K$ parallel sub-MIPs through $\{\mathcal{S}_{i}\}_{i=1}^{K}$ assignments. Through a selective prediction mechanism, it decides which variables to fix and samples their specific values, forming various assignment combinations that only involve subsets of variables. Subsequently, each partial assignment is transformed into a smaller sub-MIP problem by fixing these assigned variables. These sub-MIP problems are then solved in parallel using traditional solvers to generate multiple complete feasible solutions. Finally, the solution with the optimal objective function value is selected as the final result. This process combines the predictive capabilities of neural networks with the efficiency of traditional solvers, achieving an end-to-end workflow from partial assignments to optimal solutions.

### D.2 Predict-and-Search

In the predict phase, a graph neural network estimates binary variable assignments $\hat{x}_{j}\in\{0,1\}$. During the search phase, the algorithm constructs a trust region based on the predictions, restricting the original problem to a feasible solution space within the neighborhood of the predicted solutions. The trust region is defined as follows:

$$ $\displaystyle\mathcal{T}_{0}=\{x_{j}\mid\hat{x}_{j}\leq k_{0}\}\quad\text{(variables near 0)}$ $\displaystyle\mathcal{T}_{1}=\{x_{j}\mid\hat{x}_{j}\geq 1-k_{1}\}\quad\text{(variables near 1)}$ $\displaystyle\sum_{j\in\mathcal{T}_{0}}\mathbb{I}(x_{j}^{*}=1)+\sum_{j\in\mathcal{T}_{1}}\mathbb{I}(x_{j}^{*}=0)\leq\Delta\cdot(|\mathcal{T}_{0}|+|\mathcal{T}_{1}|)$ $$

where $k_{0},k_{1}\in(0,1)$ control selection thresholds and $\delta$ defines permissible deviation. Specifically, by adding constraints (e.g., limiting the number of variables with scores close to 0 or 1 and controlling deviations from predicted values), the problem is transformed into a smaller subproblem, which is then solved using traditional solvers to efficiently search for high-quality feasible solutions. Compared to methods that directly fix variables, the trust region strategy retains flexibility to avoid infeasibility while improving solution quality through localized search.

### D.3 PMVB

The PMVB (Probabilistic Multi-Variable Cardinality Branching) method leverages predictions from a graph neural network to construct statistically grounded branching constraints. Specifically, the binary variables are partitioned into two disjoint sets based on their predicted probabilities:

$$ $\mathcal{U}=\{j:\hat{y}_{j}\geq\tau\},\quad\mathcal{L}=\{j:\hat{y}_{j}\leq 1-\tau\}$ $$

where $\tau\in(0.5,1]$ is a confidence threshold indicating how certain the model is about a variable being 1 or 0, respectively. Using the principles of risk pooling and concentration inequalities, the algorithm constructs two soft constraints that serve as branching hyperplanes:

$$ $\displaystyle\mathcal{C}_{\mathcal{U}}:\sum_{j\in\mathcal{U}}y_{j}\geq\left\lceil(1-\delta)\sum_{j\in\mathcal{U}}\mathbb{E}[\hat{y}_{j}]-\gamma\right\rceil,\quad\mathcal{C}_{\mathcal{L}}:\sum_{j\in\mathcal{L}}y_{j}\leq\left\lfloor\delta\sum_{j\in\mathcal{L}}\mathbb{E}[\hat{y}_{j}]+\gamma\right\rfloor$ $$

where $\gamma=\sqrt{\frac{|\mathcal{S}|\ln(2/\delta)}{2}}$ for $\mathcal{S}\in\{\mathcal{U},\mathcal{L}\}$, and $\delta$ is a tunable confidence parameter.

These hyperplanes probabilistically constrain the variable assignments, forming a high-confidence trust region that filters out unlikely configurations while preserving feasible and high-quality candidates. By intersecting the feasible region with both $\mathcal{C}_{\mathcal{U}}$ and $\mathcal{C}_{\mathcal{L}}$, the algorithm defines a subproblem $\mathcal{P}_{\mathcal{C}_{\mathcal{U}},\mathcal{C}_{\mathcal{L}}}$ that is both reduced in size and refined in quality.

This subproblem is then passed to a traditional MILP solver for efficient resolution. Compared to methods that directly fix variables, PMVB offers greater robustness by retaining feasible flexibility while introducing statistically grounded directional guidance. As a result, it effectively balances confidence-driven pruning with solution diversity, improving both computational efficiency and solution quality in challenging problem instances.

### D.4 Apollo-MILP

The algorithm begins by employing a graph neural network to predict values for currently unfixed variables, yielding a predicted solution. Subsequently, the Predict-and-Search method is applied to solve subproblems derived from this prediction, generating a reference solution. By comparing the predicted and reference solutions, variables with identical values in both solutions are fixed, thereby constructing a reduced problem. At each iteration $t\in\{1,...,K\}$, the algorithm applies parameters $[k_{0}^{(t)},k_{1}^{(t)},\Delta^{(t)}]$ in Predict-and-Search algorithm,
where $K$ controls total iterations. This process iteratively repeats, progressively identifying high-quality partial solutions and expanding the subset of fixed variables to reduce problem dimensionality.

Throughout the iterations, Apollo-MILP alternates between prediction and correction steps, continuously refining the predicted solution and reducing the complexity of the MILP problem. This approach ensures optimality while significantly enhancing both solution quality and computational efficiency.

## Appendix E Pseudo code for FMIP Inference

Figure: Algorithm 1 Inference Phase of FMIP

## Peer Reviews

confidence: 2
contribution: 3
initial_rating: 6
presentation: 3
questions: None
soundness: 3
strengths: The claimed empirical performance is good, achieving decent relative improvement in primal gap over state-of-the-art baselines on some benchmarks. The framework's compatibility with multiple GNN architectures and solvers, as shown in ablation studies and compatibility analyses, highlights its flexibility and generalizability. Additionally, the holistic guidance, supported by gradient-based updates for continuous variables and sampling-reweighting for integers, effectively leverages complete solution candidates at each step, leading to higher-quality heuristics.
summary: The paper introduces FMIP (Joint Continuous-Integer Flow for Mixed-Integer Linear Programming), a generative framework that addresses limitations in existing MILP heuristics by modeling the joint distribution of both integer and continuous variables using conditional flow matching. This approach captures the interdependencies between variable types, enabling a holistic guidance mechanism during inference that steers solutions toward optimality and feasibility via objective function and constraint feedback. FMIP, integrated with various graph neural network backbones and downstream solvers, is evaluated on some benchmarks.
weaknesses: FMIP's focus on bounded integer variables limits its applicability to general MILP problems with unbounded or non-binary integers. The reliance on flow matching introduces higher inference times compared to discriminative baselines like supervised learning. While ablations confirm the value of joint modeling and guidance, the paper could benefit from deeper analysis on failure cases or scalability to very large-scale instances beyond the tested benchmarks.

---

confidence: 3
contribution: 2
initial_rating: 6
presentation: 4
questions: 1. Does the inference time in Table 3 include the time cost of the guidance process?
2. In Table 4, how does the "w/o Continuous" variant perform the guidance? Without continuous variables, how are constraint violations and objective value calculated?
soundness: 4
strengths: 1. The paper is well-organized with a clear logical flow.
2. Extending the application of generative methods from ILP to MILP is a natural progression.
3. Empirical experiments are comprehensive, and the results show improvements over the selected baselines.
summary: This manuscript is an incremental extension of [1], with the following core enhancements: i) an extention from integer linear programs (ILP) to mixed-integer linear programs (MILP), and ii) replacing diffusion models with flow matching models. The authors conduct extensive experiments, demonstrating improvements over existing methods.


[1] Zeng, Hao, et al. "Effective generation of feasible solutions for integer programming via guided diffusion." ACM SIGKDD 2024.
weaknesses: I have reviewed this paper at NeurIPS 2025, and I am satisfied with the revisions made by the authors. I do not have any major concerns, but I have a minor suggestion:

The paper states that "existing generative methods for MILP suffer from a critical limitation: they model the distribution of only the integer variables." However, the transition from ILP to MILP seems relatively straightforward, and this limitation may not constitute a major challenge. I recommend that the authors emphasize the more difficult aspects and make the challenges clearer to the reader.

---

confidence: 3
contribution: 2
initial_rating: 4
presentation: 3
questions: 1. Where is the 41.34% improvement mentioned in the abstract and the conclusion found in the experiments? Maybe it's lost 
2. In Table 4, how is the system "w/o Continuous" implemented? How does it handle continuous variables? This is important since, as mentioned by authors, this performance gap gives a direct evidence that modelling joint distributions as advocated in FMIP is important.
3. Why do you call your graph a tripartite graph instead of a bipartite graph with continuous/discrete binary features for variable nodes? It overcomplicates the model and makes the reader thinks edges between discrete and continuous nodes are allowed.
   Besides, line 726 mentions a bipartite graph
4.The discussion lines 468-474 is interesting. Could an additional reason to explain why the number of sampling steps may have a detrimental effect is that the coefficient $\gamma$ given to the violation (and possibly the measure of the violation itself based on the max operator) may lead to a situation where the function $f$ favors constraints over objective?
soundness: 3
strengths: This paper identifies and  highlights one of the blind posts from previous generative models for MILPs: the lack of a proper treatment  of continuous variables. I believe this is the main contribution of the paper
summary: This paper introduces FMIP a flow model for MILP solving, cast at sampling from a generative model, that defines a method to sample from joint distributions over discrete and continuous variable assignments. This can help capturing interactions between these two sets of variables, which is often missing in previous works on generative models for MILPs.
The model is based on a 'tripartite-graph'  a variation of the standard bipartite graph where integer and continuous variables are mapped to two disjoint set of nodes.

However, while the distribution is joint, it is assumed to be factorized and amounts to one independent prediction per variable (discrete or continuous). Only the conditioning is from a complete, and noisy, assignment.

Additionally, authors introduce guidance mechanisms for the two types of
variables helping steer the model's trajectories towards feasible assignments.
For the continuous variables the guidance is based on projected gradient descent to enforce box constraints (the other constraints are taken into account via the violation penalty)
For integer variables, authors use a sample-and-reweight method where the denoiser is used to generate several candidates which are evaluated. The rates are updated to promote (ie, increase the probability) of candidates which combine good objective value and low constraint violation. The final samples are also projected to enforce box constraints.



Eventually, FMIP returns an assignment for the continuous relaxation of the
input MILP, and must be processed further by a solver. In that sense, FMIP
remains a tool to warm-start a solver.


Typos:
- The two first references to Gasse et al. are the same
- The two first references to Taoan Huang et al. are the same
- l. 712 missing space before Edges
- l. 728 why $> 0$ and not $\neq 0$ ?
- l. 753 $t$ and $t$ ??
weaknesses: - The model is straightforward and amounts to the concatenation of a continuous flow and a discrete flows that can both be conditioned on predicted assignments.
- The tri-partite representation is surprising and questionable (see question below)

---

confidence: 4
contribution: 2
initial_rating: 4
presentation: 3
questions: 1. The guidance combines the objective with squared constraint violations via a hyperparameter $\gamma$. Please clarify how it balances objective vs. violation terms, its sensitivity, and the tuning/selection procedure used.
2. The reweighting step appears to require evaluating continuous variables, which is extra computational efforts. If many continuous variables are functionally dependent or cheaply recoverable, why is explicit prediction necessary? Could a GNN/GAT that aggregates global information suffice without predicting continuous components? Is it necessary because of the generative model? Please clarify the additional computational cost or justify the benefit of jointly predicting discrete and continuous variables (e.g., tighter guidance, reduced repair time, better warm starts).
soundness: 3
strengths: 1. The joint modeling of integer variables and continuous variables is a notable contribution, where previous works usually only focus on discrete value predictions, which ignores some global relationships.
2. The author proposed the holistic guidance mechanism integrated into FMIP’s sampling. During inference, FMIP uses the MILP’s objective and constraint violations to guide the generative trajectory, which is expected to produce better predictions.
summary: FMIP uses a multimodal flow matching generative model to learn the distribution of high-quality MILP solutions in a combined continuous–discrete space. The author proposes a **“holistic” guidance mechanism** that steers the generation process using the MILP’s objective function and constraint satisfaction feedback, refining candidate solutions toward feasibility and optimality during sampling. This approach is expected to addresses the limitation of previous graph neural network (GNN) predictors that predicted only integer variables and left continuous ones to a solver, thereby missing the intricate coupling between them .
weaknesses: 1. Given the context where integer–continuous coupling is important, the motivation for a generative approach is under-developed. Why is a generative model preferable to discriminative predictors under identical solver budgets? Why is the Flow Matching used in the proposed framework?
2. In Tables 1–2 the downstream solver time limits should be standardized (e.g., a flat 1000 s for all methods) to ensure fairness; the mixed 400/600 s budgets are short for some instances and confound comparisons. A unified budget would make the FMIP gains more persuasive.

---

confidence: 4
contribution: 3
initial_rating: 6
presentation: 4
questions: - Why don't you use any exact downstream solvers as Gurobi, SCIP etc.?
- Looking at the experiments, why should I use your (very complex) model involving data collection and training instead of just running Gurobi?
- How good is your model for hard instances which are harder than the instances in the training set?
- How much does your model improve the Gurobi process compared to the standard heuristics implemented in Gurobi?
soundness: 3
strengths: The paper is very clearly written and all concepts are presented in a detailed but concise way. The methodology is very sophisticated and the idea of incorporating the continuous variables into the generative model is interesting and hence the whole paper provides an important contribution to the knowledge of the field. The experiments are extensive, although there are some subtle issues which should be improved.
summary: The authors study mixed-integer optimization problems (MILPS) which are computationally hard to solve to optimality in general. Solvers for MILPs benefit from good heuristic start solutions, which can be generated by using state-of-the-art learning methods. The authors propose a generative model which is able to generate good feasible solutions by predicting (compared to previous models) both, integer and continuous variables. The performance of the model is tested on a wide variety of instances and compared to several benchmarks.
weaknesses: 1. The goal of the paper is to provide (exact) solvers with good feasible solutions. Hence I don't understand why there is no state-of-the-art exact MILP solver used as downstream solver (e.g. Gurobi, SCIP or HiGHS)? It would be important to see on very hard instances, how much the predicted solutions improve the performance of the exact solver (e.g. in terms of optimality gap after 1 hour). Solvers as Gurobi have sophisticated methods implemented to search for good feasible solutions. If your extremely complex model (which has to be trained etc.) cannot provide better solutions then it is useless for exact solvers.

2. In the experiments Gurobi is able to provide the same objective values as your model for nearly all instances (except LB and IP, although on LB it nearly achieves the same value). This indicates that the instances which are solved here are quite easy for exact solvers as Gurobi and at the same time Gurobi provides you with an optimality gap (which you should report). So looking at the results as a user I would prefer using Gurobi, which maybe has to run for 3600 seconds but on the other hand it is extremely simple to use: there is no training involved and it provides an optimality gap. As mentioned above, it would be interesting to see how good your predicted solutions are for harder instances.

3. It is not clear from the main paper, what data the model is trained on. Since you have to provide good solutions for the training instances these can only be instances which are easy to solve by the MILP solvers. Your model should be evaluated on instances which are of larger dimension/harder than the training set.

Minor comments:
- The paragraph "Metrics" is not clear enough. What exactly do you mean hear by OBJ and what is BKS? Is BKS the best known solution after running the downstream solver? Why do you compare it to OBJ which I understand is the objective value of the solution predicted by the generative models?
- The statement in abstract and conclusion that the primal gap is reduced by 41.34% is misleading. Gurobi is able to provide primal solutions nearly as good as your model in one hour.


Decision: Accept (Poster)