# Solution of a bilevel optimistic scheduling problem on parallel machines

Quentin Schau, Olivier Ploton, Vincent T'kindt, Han Hoogeveen, Federico

Della Croce, Jippe Hoogeveen

► To cite this version:

Quentin Schau, Olivier Ploton, Vincent T'kindt, Han Hoogeveen, Federico Della Croce, et al. Solution of a bilevel optimistic scheduling problem on parallel machines. 2025. hal-05536998

HAL Id: hal-05536998

https://hal.science/hal-05536998v1

Preprint submitted on 4 Mar 2026

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
# Solution of a bilevel optimistic scheduling problem on parallel machines

Quentin Schau¹,²,*, Olivier Ploton¹, Vincent T'kindt¹, Han Hoogeveen⁴, Federico Della Croce²,³, Jippe Hoogeveen⁴

¹University of Tours, LIFAT, France.

²Politecnico di Torino, DIGEP, Torino, Italy.

³CNR, IEIIT, Torino, Italy

⁴Utrecht University, Department of Information and Computing Sciences, Netherlands.

## Abstract

We consider the uniform parallel machines scheduling problem in the context of optimistic bilevel optimization, where two speed options are considered. In this scenario, the leader aims to minimize the weighted number of tardy jobs, while the follower seeks to minimize the total completion time on a set of uniform machines. This problem has practical applications in Industry 4.0. We show that this problem is $NP$-hard in the strong sense by providing a reduction from the Numerical 3-Dimensional Matching problem and we provide a moderately exponential-time dynamic programming algorithm. The problem is solved by means of a concise MIP formulation and a branch-and-bound algorithm that embeds a column generation approach for the lower bound computation. Computational experiments are presented for instances with up to 80 jobs and 4 machines while larger problems are out of reach for the proposed approaches.

**Keywords:** Bilevel scheduling, Exact algorithms, Complexity results.

## 1. Introduction

Bilevel scheduling problems model the situation in which two agents interact to schedule a set of jobs, but under a specific hierarchical setting. The first agent, called the leader, take part of the scheduling decisions to minimize its objective function $f^L$. Next, the second agent, called the follower, takes the remaining decisions to minimize its objective function $f^F$. The particularity of such problems is that these two objective functions can only be evaluated when all scheduling decisions have been taken. Bilevel scheduling is a challenging and emerging research topic that relates to bilevel optimization. The literature on such problems is rich and has focused on both continuous linear and non-linear problems (Bard, 1998; Dempe, 2002, 2020). It appears that, when the leader has taken his decisions, the follower's problem may have multiple optimal solutions: this leads to different classes of bilevel problems. In optimistic (*OPT*) bilevel problems, the follower returns the solution that leads to the smallest

<footnote>*Corresponding author.</footnote>

<footnote>E-mail addresses: quentin.schau@univ-tours.fr; quentin.schau@polito.it (Q. Schau)</footnote>
value of the leader's objective function among all optimal solutions. In pessimistic (PES) bilevel problems, the follower returns the optimal solution that is the worst for the leader's objective function. Finally, in adversarial (ADV) bilevel problems, the leader takes decisions to ensure that the optimal solution of the follower's problem is the worst possible.

To the best of our knowledge, the literature on bilevel scheduling is relatively limited. To solve a bilevel flow shop problem, Karlof and Wang (1996) proposed an exact two level branch-and-bound algorithm, where the leader assigns the operators to the machines in order to minimize the total flowtime, while the follower decides on the jobs' schedule in order to minimize the makespan. The considered approach was limited to instances with up to 10 jobs. Abass (2005) studied an extension of this problem by considering fuzzy processing times. Kasemset and Kachitvichyanukul (2010) proposed a bilevel scheduling multi-objective job shop problem that operates under the Theory of Constraints (TOC) policy. It was extended in Kasemset and Kachitvichyanukul (2012) who proposed a heuristic approach based on Particle Swarm Optimization. In this problem, the leader aims to minimize idle time on bottleneck machines before the follower, that schedules jobs in order to minimize the weighted aggregating function, i.e., the makespan, the maximum tardiness and the maximum earliness. Brown et al. (2005) considered a bilevel project scheduling problem where the leader takes decisions in order to maximize the completion time of the follower. We can see this problem as an adversarial setting. They investigated the computational complexity of variants of the problem. Lukač et al. (2008) proposed a heuristic based on Tabu search to solve a bilevel scheduling problem, where the objective of the leader is to assign the products to the machines to minimize the total sequence dependent setup time, while the objective of the follower is to minimize the production, storage and setup costs of the machine. Kovács and Kis (2011) used an exact constraint programming approach to model and solve discrete bilevel problems, and they applied it on optimistic bilevel single machine problems in which the leader selects the set of jobs the follower next schedules. Their approach failed finding the optimal solution already for several instances with 25 jobs. Following the three-field classification scheme (Graham et al. (1979) extended by T'kindt et al. (2024)), this problem is denoted by $1|OPT - n, r_j, \tilde{d}_j| \sum_j w_j^L R_j, \sum_j w_j^F C_j^F$, with OPT - n meaning that the leader selects n jobs in an optimistic setting. Moreover, the Boolean variable $R_j = 1$ means that job $J_j$ is rejected by the leader, who minimizes the total cost of the rejected jobs $\sum_j w_j^L R_j$, after which the follower schedules the jobs in order to minimize the weighted sum of completion times. This problem was shown weakly NP-hard (Kis and Kovács, 2012). Furthermore, Kis and Kovács (2012) studied the two versions optimistic and pessimistic of the bilevel scheduling problem: the $P|OPT - A_k| \sum_j C_j^L, \sum_j w_j^F C_j^F$ and $P|PES - A_k| \sum_j C_j^L, \sum_j w_j^F C_j^F$ problems, where $A_k$ is a set of jobs assigned to any machine k by the leader, while the follower sequences jobs on each machine. These problems are shown to be strongly NP-hard. Konur and Golias (2013) studied scheduling of inbound trucks at the inbound doors of a cross-dock facility under truck arrival time uncertainty. A Genetic Algorithm was proposed to solve both optimistic and pessimistic bilevel problems. The leader's problem is the minimization of the expected total service time of the cross-dock operator while the follower's problem is the minimization or maximization, regarding the setting, of the expected total waiting time by determining the expected arrival times. Bianco et al. (2015) proposed a bilevel scheduling problem occurring in grid computing. The leader assigns jobs to the follower aiming to minimize the total cost for rejecting or delaying tasks while the follower wants to maximize the computational resource usage efficiency. They also proposed an exact algorithm based on a bilevel programming formulation, as well as, a reform
mulation as a single-level Mixed Integer Programming (MIP) by using Karush-Kuhn-Tucker (KKT) optimality conditions. This MIP is often referred to as a Mathematical Program with Equilibrium Constraints (MPEC) (Dempe, 2003). Finally, the paper proposes a heuristic algorithm to solve large size instances. T'kindt et al. (2024) studied several bilevel scheduling problems, with one machine and the adversarial setting. The paper considers jobs selection and data modification for the leader's decision and several usual objective functions from the scheduling literature. They provided some complexity proofs as well as exact polynomial-time algorithms. At the current state of the art, these problems appear very challenging and exact solution approaches can just handle instances very small in size.

This study focuses on a scenario with one leader and one follower composed of uniform parallel machines operating at two speeds. We consider a periodic decision-making framework with a time horizon $T$ to account for dynamic effects. The leader minimizes the weighted number of tardy jobs, whereas the follower's objective is to minimize the sum of completion times. From a bilevel perspective, we consider an optimistic case where the follower returns a schedule that achieves the minimum number of weighted tardy jobs among all optimal schedules for the total completion time.

Such a problem may occur in the context of the Industry of the Future in which the manufacturing system is composed of autonomous cyber-physical systems (CPS) interacting with each other in a hierarchical way. Thus, the problem is to dispatch production orders to these CPS. In such an organization, the leader is the system and the followers are the CPS. The integration of maintenance in the scheduling process may lead to setup the speed to potentially defective machines to its lowest value while healthy machines work at their fastest speed.

Our contributions focus on two main aspects. First, we present several complexity results and, in particular, we establish that the considered bilevel scheduling problem is strongly $\mathcal{NP}$-hard. To this extent, we introduce a MIP formulation with a polynomial number of variables and constraints that rules out the possibility that the bilevel problem belongs to the second level of the polynomial hierarchy (Woeginger, 2021) and provide a polynomial-time reduction from the Numerical 3-Dimensional Matching problem which is known to be $\mathcal{NP}$-complete in the strong sense. Next, we present a moderately exponential-time exact dynamic programming procedure and a dedicated branch-and-bound algorithm embedding a column generation scheme for the computation of the lower bound. Experimental results reveal that the bilevel problem is extremely challenging to solve and only instances with up to 80 jobs and 4 machines can be attacked. The remainder of this paper is structured as follows. Section 2, presents structural properties and a MIP formulation as well as complexity results and a dynamic programming formulation. Section 3 focuses on the branch-and-bound algorithm. Section 4, reports on computational experiments conducted using randomly generated instances. Finally, Section 5 draws conclusions and outline potential directions for future research.

## 2. Properties, modeling and complexity issues

We assume that the leader has taken his maintenance decisions, i.e., changed the floor shop configuration: we have $m_1$ machines with high speed denoted by set $V_1$, and $m_0$ machines with low speed denoted by set $V_0$. Correspondingly, let have $m = m_1 + m_0$ machines. A set $\mathcal{J}$ of $N$ jobs is given among which $n \le N$ have to be scheduled, with $n$ a given value. Each job $J_j$ is characterized by its processing time $p_j$, weight $w_j$, and due date $d_j$. Thus, the leader selects a subset $I \subset \mathcal{J}$ of $n$ jobs and assigns them to the follower who, next, schedules them
on the machines. The leader aims to minimize the number of weighted tardy jobs $\sum_j w_j U_j^L$
whereas the follower aims to minimize the sum of completion times $\sum_j C_j^F$. Let $\mathcal{S}_I$ be the
set of machine schedules that contain $n$ jobs from $I$, and let $C_j^F(s)$ be the completion time of
job $J_j$ in machine schedule $s$.

Then, we define $U_j^L(s) = \begin{cases} 1 & \text{if } C_j^F(s) > d_j \\ 0 & \text{otherwise} \end{cases}$. In the follower's problem, the optimistic scenario implies that there is a lexicographical objective function, since we want the schedule that minimizes the leader's objective function among all optimal schedules for the follower's objective function. Consequently, the follower minimizes $Lex(\sum_j C_j^F, \sum_j w_j U_j^L)$. The bilevel problem can be formulated as:

$$
\begin{equation}
\begin{array}{@{}l@{}}
  \displaystyle \min_{\substack{\mathcal{I} \subset \mathcal{J} \\ |\mathcal{I}|=n}} \sum_{j \in \mathcal{I}} w_j U_j^L(s^*) \\
  \displaystyle \text{st.} \\
  \phantom{\text{st.}} s^* \in \operatorname*{argmin}_{s \in \mathcal{S}_{\mathcal{I}}} \left( \operatorname{Lex} \left( \sum_{j \in s} C_j^F(s), \sum_{j \in s} w_j U_j^L(s) \right) \right)
\end{array}
\tag{1}
\end{equation}
$$

Correspondingly, we denote this bilevel problem by $Q|V_i \in \{V_0, V_1\}, OPT-n|\sum_j C_j^F, \sum_j w_j U_j^L$.
The related follower's problem is denoted by $Q|V_i \in \{V_0, V_1\}|Lex(\sum_j C_j, \sum_j w_j U_j)$. For the
rest of the paper, we assume that all jobs are sorted by Shortest Processing Time (SPT) rule,
i.e., $p_1 \le p_2 \le \dots \le p_N$. In case of equal-size jobs, i.e., jobs with same processing times, we
sort them according to the Earliest Due Date (EDD) rule.

## 2.1. Block structure

In this section, we consider structural properties for the $Q||\sum_j C_j$ problem. We know from Conway et al. (1967) that the problem can be solved in $O(n \log n)$ time using the Shortest Processing Time job on the First Available Machine (STP-FAM) rule. The follower's problem is a special case where there are only two different speeds. This implies a special property for optimal schedules that is obtained by adapting the proof of Conway et al. (1967).

Given $n$ jobs and $m$ machines $M_i$ with a speed $V_i$, the processing time of a job $J_j$ on the machine $M_i$ is $p_j/V_i$. Our objective is to minimize $\sum_j C_j$. We represent a schedule by a non-rectangular and two-dimensional list that is a function of the machine number $i$ and of the chronological order $k$. We denote by $J_{[i,-\ell]}$ the $\ell$ th job on the machine $M_i$ from the end, with $p_{[i,-\ell]}$ its processing time. We also denote by $J_{[i,k]}$ the $k$ th job on the machine $M_i$ from the beginning, with $p_{[i,k]}$ its processing time.

A schedule $\sigma$ is given by : $\sigma = \begin{pmatrix} J_{[1,1]} & \cdots & J_{[1,\eta_1]} \\ \vdots & \ddots & \vdots \\ J_{[m,1]} & \cdots & J_{[m,\eta_m]} \end{pmatrix}$ where $\eta_i$ denotes the number of jobs assigned on machine $M_i$, with $\sum_{i=1}^m \eta_i = n$.

Since all $\eta_i$ depend on the schedule $\sigma$, we can have multiple values of $\eta_i$. To simplify the
notation, let us introduce some dummy jobs with a processing time and a weight equal to
zero. We will ensure that each machine is assigned exactly $n$ jobs, including the dummy ones.
It is clear that in an optimal schedule, these dummy jobs are scheduled on the left and can be removed without modifying the objective function nor the optimality of the schedule. To facilitate notation, we will denote all dummy jobs with the index 0, thus, $p_0 = 0$. Thus we add $(mn - n)$ dummy jobs. This results in the modified schedule $\sigma$ being defined as :

$$ \sigma = \left( \begin{array}{cccc} 0 & \overset{n-\eta_1 \text{ times}}{\cdots} & 0 & J_{[1,n-\eta_1+1]} & \cdots & J_{[1,n]} \\ & \ddots & \ddots & \ddots & \ddots & \vdots \\ \cdots & \cdots & \cdots & \cdots & \cdots & \vdots \\ 0 & \overset{n-\eta_m \text{ times}}{\cdots} & 0 & J_{[m,n-\eta_m+1]} & \cdots & J_{[m,n]} \end{array} \right) = \left( \begin{array}{ccc} J_{[1,1]} & \cdots & J_{[1,n]} \\ & \ddots & \ddots \\ J_{[m,1]} & \cdots & J_{[m,n]} \end{array} \right) $$

We focus on a ranking from right to left, i.e., in anti-chronological order $\ell = n + 1 - k$ (with $\ell = 1$ for the last, 2 for the second-to-last job, and so forth). Since the completion time of the job on machine $M_i$ at position $k$ can be expressed as $C_{[i,k]} = \sum_{h=1}^{k} p_{[i,h]}/V_i$, we see that $\sum_{k=1}^{n} C_{[i,k]} = \sum_{k=1}^{n} \sum_{h=1}^{k} p_{[i,h]}/V_i = \sum_{h=1}^{n} (n+1-h)p_{[i,h]}/V_i$. Hence, we find that $\sum_{j=1}^{n} C_j = \sum_{\substack{1 \le i \le m \\ 1 \le \ell \le n}} \ell/V_i \times p_{[i,-\ell]}$

Correspondingly, we must minimize the weighted sum of processing times $p_{[i,-\ell]}$ with positional weights $\omega_{[i,\ell]} = \ell/V_i$, for position $\ell$ ($\ell = 1, \dots, n$) from the end on machine $M_i$. The $\omega_{[i,\ell]}$'s do not depend on the jobs or on the schedule. It is well known that such a sum of pairwise products of two sequences of reals is achieved if one sequence is arranged in increasing order and the other in decreasing order (Hardy et al., 1934). To get the optimal solution, we have to sort the $\omega_{[i,\ell]}$'s in non-decreasing order and the jobs in the reverse order, i.e. Longest Processing Time (LPT) rule. Specifically, the $k^{th}$ smallest coefficient $\omega_{[k]}$ is associated with the $k^{th}$ largest processing time $p_{[k]}$. Consequently, we have the following theorem.

**Theorem 2.1.** The $Q|V_i \in \{V_0, V_1\}| \sum_j C_j$ problem is solvable in $O(n \log n)$ time.

□

**Remark.** Firstly, the proof of the general case, i.e., $Q \parallel \sum_j C_j$, provides the algorithm by constructing the schedule from the right to the left, i.e., from the end of the machine schedule to the beginning. Alternatively, we can equivalently consider an SPT-FAM rule, thus scheduling jobs from the left to the right.

Secondly, sets with identical positional weights correspond to batches of jobs. Furthermore, when there are only two machine speeds involved, an optimal schedule is characterized by a repetition of blocks. Figure 1 illustrates the blocks structure, where jobs within the dashed group signify their membership to a block and can subsequently undergo permutation without changing the value of $\sum_j C_j$. We denote such a block as $B_b$, which is a set composed of pairs $(i, k)$ where $i$ is the machine index and $k$ is the position on the machine schedule.

This is an important result because it enables us to characterize the set of optimal schedules for the first criterion ($\sum_j C_j$). Thus, the problem is to determine which permutation of jobs within a block should be applied to obtain an optimal schedule for the second criterion ($\sum_j w_j U_j$).
<figure><img src="image_2.png" /><figcaption>Figure 1: Scheduling pattern with block structure for the solution of $Q|V_i \in \{V_1, V_0\}| \sum_j C_j$</figcaption></figure>

## 2.2. Mixed Integer Programming

In this subsection we show that the bilevel problem can be modeled by means of Mixed Integer Programming. Before presenting the model, we introduce a few notations.

We have $b_{\max}$ blocks $\mathcal{B}_b$, $1 \le b \le b_{\max}$, given by the block structure described in Section 2.1. Let $E = \{(i, k)\}$ be a set of locations (just like $\mathcal{B}_b$) where $i$ is the machine index and $k$ the position on the machine schedule in chronological order; we have $E = \mathcal{B}_1 \cup \dots \cup \mathcal{B}_{b_{\max}}$. We need to assign a job to each position in the blocks $\mathcal{B}_2, \dots, \mathcal{B}_{b_{\max}}$, but the first block $\mathcal{B}_1$ may not be completely filled. However, we can compute the number of jobs $Q$ that should go into $\mathcal{B}_1$ as $Q = n - \sum_{b=2}^{b_{\max}} |\mathcal{B}_b|$.

We use binary variables $x_{i,j,k}$, such that $x_{i,j,k} = 1$ if job $J_j$ is scheduled at location $(i, k)$; and zero otherwise. We also use binary variables $U_{i,j,k}$, with $U_{i,j,k} = 1$ if job $J_j$ is scheduled at location $(i, k)$ and is penalized. The equations only state that tardiness implies penalty. So, a job $J_j$ may be a priori unduly penalized, but this situation is, at the end, avoided by the minimization of penalties. So, we have $U_j \le \sum_{(i,k) \in E} U_{i,j,k}$, and the minimization of $\sum_j \sum_{(i,k) \in E} w_j U_{i,j,k}$ leads to minimizing $\sum_j w_j U_j$. We use variables $C_{i,k}$ to express the completion time of the job in position $k$ on machine $M_i$.

The following MIP model holds:
$$ \min \sum_{j=1}^{N} \sum_{(i,k) \in E} w_j U_{i,j,k} \qquad (2) $$

$$ \text{s.t.} \qquad \sum_{(i,k) \in E} x_{i,j,k} \le 1 \qquad \forall j \in \{1, \dots, N\} \qquad (3) $$

$$ \sum_{j=1}^{N} x_{i,j,k} \le 1 \qquad \forall (i,k) \in E \qquad (4) $$

$$ \sum_{j=1}^{N} \sum_{(i,k) \in E \setminus B_1} x_{i,j,k} = n - Q \qquad (5) $$

$$ \sum_{j=1}^{N} \sum_{(i,k) \in B_1} x_{i,j,k} = Q \qquad (6) $$

$$ \sum_{j=1}^{N} x_{i,j,k} p_j \le \sum_{j=1}^{N} x_{i',j,k'} p_j \qquad \forall 1 \le b \le b_{\max} - 1, \forall (i,k) \in B_b, \forall (i',k') \in B_{b+1} \qquad (7) $$

$$ C_{i,k} = C_{i,k-1} + \sum_{j=1}^{N} x_{i,j,k} p_j / V_i \qquad \forall (i,k) \in E \qquad (8) $$

$$ C_{i,0} = 0 \qquad \forall i \in \{1, \dots, m\} \qquad (9) $$

$$ U_{i,j,k} \le x_{i,j,k} \qquad \forall j \in \{1, \dots, N\}, \forall (i,k) \in E \qquad (10) $$

$$ C_{i,k} \le \sum_{j=1}^{N} x_{i,j,k} d_j + \left( \sum_{j=1}^{N} U_{i,j,k} \right) H(i,k) \quad \forall (i,k) \in E \qquad (11) $$

Here $H(i,k) = k/V_i \times \max_j p_j$, which is an upper bound on the completion time of the $k$th job on machine $M_i$.

The objective function (2) minimizes the weighted number of tardy jobs. Constraints (3) ensure that no jobs are duplicated. The disjunctive constraints (4) ensure that each location has at most one job. Constraints (5) and (6) guarantee that all blocks, except maybe the first one, must be completely filled: these constraints imply that we select exactly $n$ jobs. Constraints (7) enforce that the processing times of assigned jobs between consecutive blocks do not decrease. Constraints (8) allow to compute the completion times. Constraints (10) ensure that only a scheduled job can be penalized. Finally, constraints (11) guarantee that a job that is not penalized is necessarily early, and so, we need to check that its completion time is smaller than its due date.

This mathematical formulation has $O(N(n+m))$ variables and $O(N^2+Nm^2)$ constraints.

## 2.3. Complexity results

### 2.3.1. The follower's problem

Firstly, we demonstrate that in the case of equal-size jobs, the follower's problem, with an arbitrary number of machine speeds can be solved in polynomial time.
**Theorem 2.2.** The $Q|p_j = p|Lex(\sum_j C_j, \sum_j w_j U_j)$ problem is polynomially solvable.

*Proof.* We know that all optimal schedules for the $Q||\sum_j C_j$ problem are characterized by the block structure, and the SPT-FAM rule gives an optimal solution (Conway et al., 1967) in $O(n \log n)$ time. Thus, we can compute all completion times, as jobs' processing times are equal. Moreover, by extending the algorithm mentioned by Dessouky et al. (1990), we can solve this problem. The idea is to determine which job is on time if assigned to a completion time, and then choose among them the one with the maximum weight. To do this, the jobs are initially sorted according to the reverse EDD order, in $O(n \log n)$ time. Then, we start with the largest completion time and assign the early job with the largest weight. If no early job exists then this completion time is temporarily left unassigned, and we move to the second largest completion time, repeating the above operation. At the end of this first pass, when all completion times have been considered, then the remaining jobs can be randomly assigned to the remaining unassigned completion times: they are all tardy and with minimum weight. This can be implemented to run in $O(n \log n)$ time. □

Now, we focus on the general case with arbitrary processing times. In the case of 2 identical machines, we show that the follower's problem is $\mathcal{NP}$-hard in the weak sense.

**Theorem 2.3.** The $P2||Lex(\sum_j C_j, \sum_j U_j)$ problem is $\mathcal{NP}$-hard in the weak sense.

*Proof.* Let us consider the following problem which is known to be $\mathcal{NP}$-complete in the weak sense (Garey et al., 1988):

**Even-Odd Partition (EOD)**

**Data:** $2n$ positive integers $a_1 < a_2 < \dots < a_{2n}$ where $A = 1/2 \times \sum_{j=1}^{2n} a_j$.

**Question:** Does there exist a partition into two sets $A_1$ and $A_2$, each containing exactly one element from the pair $\{a_{2j-1}, a_{2j}\}$ for $j = 1, \dots, n$, such that $\sum_{j \in A_1} a_j = \sum_{j \in A_2} a_j = A$?

We show that this problem reduces to the decision version of problem $P2||Lex(\sum_j C_j, \sum_j U_j)$:

**DEC**

**Data:** A set of jobs $J_j$ with a processing time $p_j$, a weight $w_j$ and a due date $d_j$, to be processed on 2 identical machines. Let $C^*$ be the value of an optimal solution of $P2||\sum_j C_j$, and $U$ a given value.

**Question:** Does there exist a schedule $\sigma$ such that $\sum_j C_j(\sigma) = C^*$ and $\sum_j U_j(\sigma) \le U$?

We show that there exists a reduction $\propto$ from EOD to DEC. Let $I$ be a given instance of EOD and $\propto(I)$ is the corresponding instance of DEC. This instance $\propto(I)$ contains $2n$ jobs that must be scheduled on 2 identical machines. Job $J_j$ corresponds to $a_j$. The processing times are defined as follows: $\forall k = 0, \dots, n-1$ we put $p_{2k+1} = 2kA + a_{2k+1}$, $p_{2k+2} = 2kA + a_{2k+2}$. Furthermore, each job $J_j$ ($j = 1, \dots, 2n$) has due date $d_j = D = 1/2 \times \sum_{i=1}^{2n} p_i = \sum_{k=0}^{n-1} (p_{2k+1}/2 + p_{2k+2}/2) = \sum_{k=0}^{n-1} 2A + A = n(n-1)A + A$. Finally, we put $U = 0$.

We demonstrate that the EOD problem's answer is true if and only if there exists a schedule $\sigma$ that has $\sum_j C_j(\sigma) = C^*$ and in which all jobs are on time, i.e., $\sum_j U_j(\sigma) = 0$.
We show the forward implication (⇒). Let us consider the scenario where the answer to the EOD problem is true. This implies that there exists a partition in two sets, $A_1$ and $A_2$, each containing exactly one element from the pair $\{a_{2j-1}, a_{2j}\}$ for $j = 1, \dots, n$, such that $\sum_{j \in A_1} a_j = \sum_{j \in A_2} a_j = A$. We create our schedule by putting $J_j$ on $M_1$ if $j \in A_1$, and on $M_2$, otherwise. Consequently, the completion times of the last job on each machine are $C_{M_1} = \sum_{j \in A_1} p_j = \sum_{j \in A_1} 2[(j-1)/2]A + a_j = 2A \sum_{k=0}^{n-1} k + A = n(n-1)A + A$ and similarly $C_{M_2} = \sum_{j \in A_2} p_j = n(n-1)A + A$. Since the common due date equals $n(n-1)A + A$, this indicates that all jobs can, indeed, be scheduled on time and $\sum_j C_j$ is minimum as this schedule follows the SPT-FAM rule (see Figure 2).

<figure><img src="image_3.png" /><figcaption>Figure 2: Machine schedules for the problem $P2||Lex(\sum_j C_j, \sum_j U_j)$</figcaption></figure>

Now suppose that there exists a schedule $\sigma$ with $\sum_j C_j(\sigma) = C^*$ and $\sum_j U_j(\sigma) \le U = 0$. Because of the choice of the processing times, we know that the $k$th position (counted from the start) on machines $M_1$ and $M_2$ is occupied by the jobs $J_{2k-1}$ and $J_{2k}$. Moreover, since the total processing time of all jobs is equal to $2D$ and each job is on time, both machines must finish exactly at time $D$. It is easily verified that the desired sets $A_1$ and $A_2$ can be obtained as the sets containing the indices of the jobs executed by machines $M_1$ and $M_2$ in $\sigma$, respectively.

<figure><img src="image_4.png" /></figure>

We also show that when the number of machines is a part of the instance, then the follower's problem is $\mathcal{NP}$-hard in the strong sense.

**Theorem 2.4.** The $P||Lex(\sum_j C_j, \sum_j U_j)$ problem is $\mathcal{NP}$-hard in the strong sense.

*Proof.* Let us consider the following problem which is known to be $\mathcal{NP}$-complete in the strong sense (Garey and Johnson, 1979):

### Numerical 3-Dimensional Matching (NUM-3DM)

**Data:** Three disjoint sets of integers $X = \{x_1, \dots, x_n\}$, $Y = \{y_1, \dots, y_n\}$, and $Z = \{z_1, \dots, z_n\}$ along with a bound $b$.

**Question:** Can $X \cup Y \cup Z$ be partitioned into $n$ disjoint sets $A_1, A_2, \dots, A_n$ such that each $A_h$ ($h = 1, \dots, n$) contains exactly one element from each of $X$, $Y$, and $Z$, i.e., $A_h = \{x_{[h]}, y_{[h]}, z_{[h]}\}$ and for $1 \le h \le n$, $x_{[h]} + y_{[h]} + z_{[h]} = b$?

We show that this problem reduces to the decision version of problem $P||Lex(\sum_j C_j, \sum_j U_j)$:

**DEC**
**Data:** A set of jobs $J_j$ with a processing time $p_j$ and a due date $d_j$, to be processed on $m$ identical machines. Let $C^*$ be the value of an optimal solution of $P \parallel \sum_j C_j$ and $U$ a given value.

**Question:** Does there exist a schedule $\sigma$ such that $\sum_j C_j(\sigma) = C^*$ and $\sum_j U_j(\sigma) \le U$?

We show that there exists a reduction $\propto$ from NUM-3DM to DEC. Let $I$ be a given instance of NUM-3DM and let $\propto (I)$ be the corresponding instance of DEC. We have: $U = 0$ and given three multisets of integers $X$, $Y$, and $Z$, each containing $n$ elements, we create $3n$ jobs with the following processing times $\forall i = 1, \dots, n$:

$p_i = x_i$, $p_{n+i} = y_i + \sum_{x \in X} x$, $p_{2n+i} = z_i + \sum_{y \in Y} y + \sum_{x \in X} x$ and a common due date $D = b + \sum_{y \in Y} y + 2 \sum_{x \in X} x$. Furthermore, we introduce $m = n$ machines.

Since we have $\forall i, j, k = 1, \dots, n$, $p_i < p_{n+j} < p_{2n+k}$, $\sum_j C_j$ is minimized by scheduling exactly one triplet $(x_{[h]}, y_{[h]}, z_{[h]}) = (x_i, y_j, z_k)$, with $i, j, k = 1, \dots, n$, on each machine $M_h$ ($h = 1, \dots, n$) (see Figure 3).

<figure><img src="image_5.png" /><figcaption>Figure 3: Machine schedules for the $P \parallel Lex (\sum_j C_j, \sum_j U_j)$</figcaption></figure>

We demonstrate that the NUM-3DM problem's answer is true if and only if all jobs can be scheduled in time.

Let us demonstrate the forward implication ($\Rightarrow$). Assuming the solution to the NUM-3DM problem is true implies the existence of $n$ disjoint triples $A_h = \{x_i, y_j, z_k\}$ such that $x_i + y_j + z_k = b$. Consequently, we schedule the jobs $\{p_i, p_{n+j}, p_{2n+k}\}$ on the machine $M_h$. The completion time of the last job on machine $M_h$ is:

$C_{M_h} = x_i + y_j + z_k + \sum_{y \in Y} y + 2 \sum_{x \in X} x = b + \sum_{y \in Y} y + 2 \sum_{x \in X} x = D$ ensuring all jobs are scheduled in time. Now, we establish the reverse implication ($\Leftarrow$). Suppose that there exists a schedule $\sigma$ such that $\sum_j C_j(\sigma) = C^*$ and $\sum_j U_j(\sigma) \le U = 0$. Because of the choice of the processing times, we see that each machine $M_h$ ($h = 1, \dots, n$) executes exactly one job from the sets $\{J_1, \dots, J_n\}$, $\{J_{n+1}, \dots, J_{2n}\}$, and $\{J_{2n+1}, \dots, J_{3n}\}$. Furthermore, since the total processing time of all jobs is equal to $nD$ and all jobs are on time, all machine must finish exactly at time $D$. A feasible solution to the NUM-3DM problem is immediately obtained by constructing triplet $h$ ($h = 1, \dots, n$) according to the indices of the jobs executed by machine $M_h$. $\square$

We can also derive the following result.
**Corollary 2.4.1.** The $P||Lex(\sum_j C_j, \sum_j U_j)$ problem is not approximable in polynomial or pseudo-polynomial time unless $\mathcal{P} = \mathcal{NP}$.

*Proof.* This result directly follows from the proof of Theorem 2.4 where it is shown that the particular problem when $\sum_j U_j^* = 0$ is strongly $\mathcal{NP}$-hard. This rules out the existence of an approximation algorithm for ratio $\rho = \frac{\sum_j U_j}{\sum_j U_j^*}$ $\square$

From Theorem 2.4 and Corollary 2.4.1, the following result immediately holds (as the considered problem is a generalization of the $P||Lex(\sum_j C_j, \sum_j U_j)$ problem).

**Corollary 2.4.2.** The $Q|V_i \in \{V_0, V_1\}|Lex(\sum_j C_j, \sum_j w_j U_j)$ problem is $\mathcal{NP}$-hard in the strong sense and is not approximable in polynomial or pseudo-polynomial time unless $\mathcal{P} = \mathcal{NP}$.

### 2.3.2. The bilevel problem

From the reductions of the previous section, we can derive that the bilevel problem is at least $\mathcal{NP}$-hard in the strong sense. With respect to the polynomial hierarchy in complexity theory, the bilevel problem cannot be complete with respect to class $\Sigma_2^p$ of level 3 due to the existence of the MIP formulation in Section 2.2 (Woeginger, 2021): consequently, the bilevel is $\mathcal{NP}$-hard in the strong sense.

Here, we show that when $m$ is fixed, the problem is weakly $\mathcal{NP}$-hard by providing a moderately exponential-time dynamic programming algorithm based on the block structure of an optimal schedule as described in Section 2.1. The main idea of the dynamic programming approach is to fill the blocks one by one with the correct set of jobs. We use state variables that keep track of each possible set of completion times, the set of machines that are available, the current block, and the currently considered job. Thus, we try to assign the job to each available machine in the current block; once it is assigned, we can update the corresponding machine's completion time and proceed to the next job. When the block is full, we update the state variables for the next block and continue until all blocks are fulfilled.

We need to select and schedule $n$ jobs from a set $\mathcal{J}$ containing $N$ jobs. Initially, we consider the case where there are no equal-size jobs, i.e., $p_j \neq p_\ell$ for all $j, \ell = 1, \dots, N$, $j \neq \ell$. We index blocks $B_b$, which represent sets of pairs (machine index, position on machine schedule), in chronological order with $1 \le b \le b_{\max}$, where block $B_{b_{\max}}$ is the final block in the schedule, which contains the final jobs on the machines with speed $V_0$. Additionally, let $\vec{C} = (C_{M_i})_{1 \le i \le m}$ denote the vector of machine completion times. To facilitate reading, we define $\omega(\vec{C}, i, j) = \begin{cases} w_j & \text{if } C_{M_i} + p_j/V_i > d_j \\ 0 & \text{otherwise} \end{cases}$,

and $M(B_b) = \{i \mid \exists (i, k) \in B_b\}$ as the set of machine indices in block $B_b$. We also define $\vec{C}' = \vec{C} \oplus (i, p_j)$ as $C'_{M_i} = C_{M_i} + p_j/V_i$ and $C'_{M_\ell} = C_{M_\ell}$ for $\ell \neq i$.

In our dynamic programming algorithm we work from behind. We assume that we have a schedule that follows the block structure in which each machine $M_i$ ($i = 1, \dots, m$) has a given completion time $C_{M_i}$; the number of such possible completion times amounts to $O(\sum p_j)$ per machine. In each step we select one job to fill a position on one of the machines in the current block (starting with block $b_{\max}$), until this block is full, after which we continue with the next
block. To facilitate notation, we denote the set of available machines in the current step of the dynamic program by $\mathcal{A} \subset \mathcal{V}_0 \cup \mathcal{V}_1$. Here a machine is available, if it has a position in the currently considered block $\mathcal{B}_b$ to which we have not yet assigned a selected job.

We define $OPT[\vec{C}, b, \mathcal{A}, g]$ as the optimal value of the total weighted number of tardy jobs in the bilevel problem, given that all machine completion times are equal to $\vec{C}$, the next job will be scheduled on the current block $\mathcal{B}_b$, with the available machines $\mathcal{A}$ and where $J_g$ is the smallest indexed job that the follower can schedule. Remark here that $OPT[\vec{C}, b, \mathcal{A}, g]$ only accounts for the cost of the jobs that have been selected. Since, all blocks, except maybe the first one, must be completely filled, we have that the block $\mathcal{B}_b$ with the set of available machines $\mathcal{A}$ is completely filled when $\text{full}(b, \mathcal{A})$ is true with

$$ \text{full}(b, \mathcal{A}) \iff \begin{cases} |\mathcal{A}| = n - \sum_{2 \le b' \le b_{\max}} |B_{b'}| & \text{if } b=1 \\ |\mathcal{A}| = 0 \ (\text{i.e. } \mathcal{A} = \emptyset) & \text{if } b > 1. \end{cases} $$

We initialize the recursion with $OPT[\vec{C}, b, \mathcal{A}, g] = +\infty$ whenever $g > N$, and $OPT[\vec{C}, b, \mathcal{A}, g] = 0$ if $b = b_{\max}$ and $\text{full}(b, \mathcal{A})$ is true; this enforces that we have a schedule containing $n$ (yet unselected) jobs. Our recurrence relation is then

$$ OPT[\vec{C}, b, \mathcal{A}, g] = \begin{cases} \min_{\substack{i \in \mathcal{A} \\ g \le j \le N}} (w(\vec{C}, i, j) + OPT[\vec{C} \oplus (i, p_j), b, \mathcal{A} \setminus \{i\}, j+1]) & \text{if full}(b, \mathcal{A}) \text{ is false} \\ OPT[\vec{C}, b+1, M(\mathcal{B}_{b+1}), g] & \text{otherwise} \end{cases} \quad (12) $$

The first line corresponds to the case that the current block is not completely full. Hence, we select the best job $J_g$ with $g \le j \le N$ to be scheduled on the best available machine $M_i \in \mathcal{A}$. In the second line, the current block is completely filled, and we need to pass to the next block.

The value of the optimal solution is given by $OPT[\vec{0}, 1, M(\mathcal{B}_1), 1]$.

Next, we focus on the general case with arbitrary processing times. We apply the formula described above and when we face $k$ equal-size jobs, we apply an algorithm that selects $\ell$ jobs, with $\ell = 0, \dots, k$. The main idea of the algorithm is that these $\ell$ jobs can completely fill the current block $\mathcal{B}_b$ and some more, and partially fill a last block. Consequently, we need to choose which machines on the last block will be selected; we enumerate all these possibilities. For each choice, we can determine the machine completion times since we have $\vec{C}$ and $\ell$ equal-size jobs with a processing time $p_j$. Then, we can find the best set $\ell$ out of these $k$ equal-size jobs by applying the algorithm that we used to prove Theorem 2.2 in Section 2.3.1. Note that we can update the dynamic programming state when we use this algorithm for any value of $\ell$ and any choice of selected machines on the last block. This general case multiplies the time complexity by a factor of $2^m$ since we have to make $\mathcal{O}(2^m)$ choices and by a factor of $n \log n$ for the time complexity of the jobs' assignment.

Therefore, this dynamic programming algorithm requires $\mathcal{O}(2^m (\sum_{j=1}^N p_j)^m N n^2 \log n)$ time and space. Correspondingly, the following theorem holds.

**Theorem 2.5.** The $Q|V_i \in \{V_0, V_1\}, OPT - n|Lex(\sum_j C_j^F, \sum_j w_j U_j^L)$ problem is $\mathcal{NP}$-hard in the strong sense and when the number of machines $m$ is fixed it turns to be $\mathcal{NP}$-hard in the
ordinary sense.

*Proof.* We have shown these results for the follower's problem, which is a special case of the bilevel problem, where the leader selects $n = N$ jobs. In the scenario where the number of machines is fixed, the moderate-exponential worst-case time complexity dynamic programming runs in pseudo-polynomial time. □

### 2.3.3. A related polynomially solvable bilevel scheduling problem

As a side result, we consider the adversarial bilevel optimization problem, where the follower still minimizes $\sum_j C_j^F$. In this scenario, the leader takes decisions so that the optimal solution of the follower's problem is the worst possible. For the rest of the section, we assume that the jobs are sorted according to the LPT rule.

**Theorem 2.6.** The $Q|ADV-n| \sum_j C_j^F$ problem can be solved in $O(N \log(N))$ time as follows:

(i) The leader selects the $n \le N$ jobs with the largest processing times,

(ii) The follower applies the SPT-FAM rule, on these $n$ jobs.

*Proof.* The objective value of any schedule $\sigma_s$ with $n$ jobs determined by the follower can be rewritten as $\sum_{k=1}^n \omega_{[k]} p_{[k]}$, where $\omega_{[k]}$ and $p_{[k]}$ denote the $k$th smallest positional weight and the $k$th largest processing time of the selected jobs, for $k = 1, \dots, n$.

Let us assume that the leader selects an optimal set $\mathcal{I}_1$ of $n$ jobs such that they are not the ones with the largest processing times. Suppose that the leader selects job $J_u$ instead of job $J_v$, whereas $p_u < p_v$. Let $\mathcal{I}_2$ denote the set of $n$ jobs that we obtain by replacing job $J_u$ with $J_v$ in $\mathcal{I}_1$. Let $p_{[k]}^1$ and $p_{[k]}^2$ denote the $k$th largest processing time of the jobs in the set $\mathcal{I}_1$ and $\mathcal{I}_2$, respectively. Then we see that $p_{[k]}^1 \le p_{[k]}^2$ for all $k = 1, \dots, n$, where at least one of the inequalities is strict. Since the positional weights $w_{[k]}$ are positive and independent of the selected jobs, we see for the resulting schedules $\sigma_s^1$ and $\sigma_s^2$ created by the follower that

$$ \sum_{j=1}^{n} C_j(\sigma_s^1) - \sum_{j=1}^{n} C_j(\sigma_s^2) = \sum_{k=1}^{n} w_{[k]}(p_{[k]}^1 - p_{[k]}^2) < 0, $$

which contradicts the optimality of the set $\mathcal{I}_1$. □

### 2.4. Properties of the solutions

In this section, we describe some dominance conditions between two partial solutions. Let $s$ and $s'$ be two partial solutions. We have $s$, that is composed of a partial schedule $\sigma_s$, a set of remaining jobs $\Omega_s$, and a number of jobs $k_s^1$ (respectively $k_s^0$) that must be scheduled on high-speed (respectively low-speed) machines to complete the solution, i.e. $|\sigma_s| + k_s^1 + k_s^0 = n$. Similarly, we have $s'$, that is composed of a partial schedule $\sigma_{s'}$, a set of remaining jobs $\Omega_{s'}$, and a number of jobs $k_{s'}^1$ (respectively $k_{s'}^0$) that must be scheduled on high-speed (respectively low-speed) machines to complete the solution, i.e. $|\sigma_{s'}| + k_{s'}^1 + k_{s'}^0 = n$. Schedules $\sigma_s$ and $\sigma_{s'}$ start at time zero without any unnecessary idle time. We denote by $C_{M_i}(\sigma_s)$ the completion time of machine $M_i$ in any partial solution $s$.

**Theorem 2.7.** If the following conditions hold:

(i) $C_{M_i}(\sigma_s) \le C_{M_i}(\sigma_{s'}) \quad \forall i \in \{1, \dots, m\}$,
$$ (ii) \sum_{j \in \sigma_s} w_j U_j(\sigma_s) \le \sum_{j \in \sigma_{s'}} w_j U_j(\sigma_{s'}), $$

$$ (iii) k_s^1 \le k_{s'}^1, $$

$$ (iv) k_s^0 \le k_{s'}^0, $$

$$ (v) \Omega_{s'} \subseteq \Omega_s, $$

then the best solution $s^*$ reachable from the partial solution $s$ is better than the best solution $s^*$ reachable from the partial solution $s'$. We say that $s$ dominates $s'$.

**Proof.** Let $(\sigma_s||\alpha_s^1 \cup \alpha_s^0)$ be the best schedule that begins with $\sigma_s$ where $\alpha_s^1$ is a partial schedule on high-speed machines and $\alpha_s^0$ is a partial schedule on low-speed machines. Both are composed only of jobs coming from $\Omega_s$ and $|\alpha_s^1| = k_s^1$ and $|\alpha_s^0| = k_s^0$. We say that $\alpha_s^1 \cup \alpha_s^0$ is the best complementary schedule from $\Omega_s$ containing $k_s^1 + k_s^0$ jobs.

We denote by $\tilde{\alpha}_s^1 \cup \tilde{\alpha}_s^0$ the best complementary schedule to $\sigma_s$ with $k_s^1 + k_s^0$ jobs from $\Omega_{s'}$. As we have $\Omega_{s'} \subseteq \Omega_s$, we obtain:

$$ \sum_{j \in (\sigma_s || \alpha_s^1 \cup \alpha_s^0)} w_j U_j(\sigma_s || \alpha_s^1 \cup \alpha_s^0) \le \sum_{j \in (\sigma_s || \tilde{\alpha}_s^1 \cup \tilde{\alpha}_s^0)} w_j U_j(\sigma_s || \tilde{\alpha}_s^1 \cup \tilde{\alpha}_s^0) $$

Let $\alpha_{s'}^1 \cup \alpha_{s'}^0$ be the best complementary schedule to $\sigma_{s'}$ with $k_{s'}^1 + k_{s'}^0$ jobs from $\Omega_{s'}$, and as we have $k_s^1 \le k_{s'}^1$ and $k_s^0 \le k_{s'}^0$, scheduling $\alpha_{s'}^1 \cup \alpha_{s'}^0$ after $\sigma_s$ guarantees that jobs in $\alpha_{s'}^1$ are executed on high-speed machines, while those in $\alpha_{s'}^0$ are executed on low-speed machines. We have:

$$ \sum_{j \in (\sigma_s || \alpha_{s'}^1 \cup \alpha_{s'}^0)} w_j U_j(\sigma_s || \alpha_{s'}^1 \cup \alpha_{s'}^0) \ge \sum_{j \in (\sigma_s || \tilde{\alpha}_s^1 \cup \tilde{\alpha}_s^0)} w_j U_j(\sigma_s || \tilde{\alpha}_s^1 \cup \tilde{\alpha}_s^0) \ge \sum_{j \in (\sigma_s || \alpha_s^1 \cup \alpha_s^0)} w_j U_j(\sigma_s || \alpha_s^1 \cup \alpha_s^0) $$

Since the completion times of the machines in $s$ are smaller than those in $s'$, i.e. $\forall i \in \{1, \dots, m\}, C_i(s) \le C_i(s')$, we have:

$$ \sum_{j \in (\sigma_{s'} || \alpha_{s'}^1 \cup \alpha_{s'}^0)} w_j U_j(\sigma_{s'} || \alpha_{s'}^1 \cup \alpha_{s'}^0) \ge \sum_{j \in (\sigma_s || \alpha_{s'}^1 \cup \alpha_{s'}^0)} w_j U_j(\sigma_s || \alpha_{s'}^1 \cup \alpha_{s'}^0) \ge \sum_{j \in (\sigma_s || \alpha_s^1 \cup \alpha_s^0)} w_j U_j(\sigma_s || \alpha_s^1 \cup \alpha_s^0) $$

So, $s$ dominates $s'$.

<figure><img src="image_6.png" /></figure>

**Remark.** For the leader's problem which is a selection problem, we could think about Knapsack-like dominance conditions, e.g. given two jobs $J_j$ and $J_\ell$, $J_j$ dominates $J_\ell$ if and only if $p_j \le p_\ell$ and $d_j \ge d_\ell$ and $w_j \le w_\ell$. Such a condition would possibly hold in a single level optimization problem. However, in our bilevel setting, the following counterexample holds, where the leader has to select two jobs, with one machine of speed $V = 1$ and three jobs defined as follow:

Here, the Knapsack-like dominance condition would indicate that $J_1$ dominates $J_3$ and, correspondingly, the schedule $(J_2, J_3)$ should be worse than $(J_1, J_2)$. However, $s = (J_2, J_3)$ is better than $s = (J_1, J_2)$ with a weighted number of tardy jobs equal to 0 against 3. Indeed, the
<table><thead><tr><th>Job</th><th>p<sub>j</sub></th><th>d<sub>j</sub></th><th>w<sub>j</sub></th></tr></thead><tbody><tr><td>1</td><td>1</td><td>7</td><td>4</td></tr><tr><td>2</td><td>2</td><td>2</td><td>3</td></tr><tr><td>3</td><td>3</td><td>5</td><td>5</td></tr></tbody></table>

follower's block structure imposes a specific assignment for the jobs that may invalidate the above dominance condition.

# 3. A branch-and-bound algorithm

## 3.1. Branching scheme

We propose a branch-and-bound algorithm that relies on a lower bound based on column generation. We use a branching strategy based on the block structure property of the follower's problem, i.e. we branch on the blocks in the order $B_1, B_2, \dots, B_{b_{max}}$. At a given node, we have taken decisions on jobs $J_1, \dots, J_{j-1}$, so we need to take a decision for job $J_j$. We have the following decisions:

* We select job $J_j$. Then we iterate over all available locations in the first available block from the parent node: let us denote by $(i, k)$ such location. Then, we assign job $J_j$ on the machine $M_i$ at the position $k$. Thus, we create as many children nodes as available locations in this block.

* We do not select job $J_j$, so we remove it in the corresponding child node.

Consequently, if $N_a$ is the number of available locations in the first available block, we create exactly $(N_a + 1)$ children nodes.

The nice aspect of this branching scheme is that we know the completion time of each machine. If we have $k$ equal-size jobs of which we want to select and schedule $\ell$ ($\ell = 1, \dots, k$) jobs, then we can use the algorithm that we used to prove Theorem 2.2 to find these $\ell$ jobs. We create a node for each $\ell$ value and for each machine's selection on the last block, resulting in $O(2^m n)$ nodes.

## 3.2. Memorization

Shang et al. (2021) introduce the Branch-and-memorize paradigm, describing how to apply memorization within search tree algorithms for sequencing problems. The authors present several memorization schemes that can be used to prune the search tree and avoid solving identical or dominated subproblems.

The first scheme is solution memorization, which involves storing in memory the best solution reachable from each node. This approach can be used when the problem is decomposable, meaning that an optimal solution of a subproblem does not depend on the decisions taken to generate the subproblem. Before branching on a node, we check if there exists an optimal solution of the subproblem that has already been explored in memory. If so, we do not need to branch on this node.

The second scheme is called passive node memorization and involves storing information about the explored nodes in memory. This scheme requires a function `check(s, s')` that returns `true`
if partial solution $s$ dominates partial solution $s'$. The idea is to check before branching on a node whether there exists a node in memory that dominates it.

The third and final scheme is predictive node memorization. Like passive node memorization, it shares the same idea: storing information about explored nodes in memory. Before branching on a node, we check whether there exists a node in memory that dominates it. If no such node exists, we try to create a partial solution that has not been explored in the tree and dominates the current branching node. If such a better partial solution is found, we store it in memory, and we do not branch on the current node since it is dominated.

Our problem is non-decomposable because the decision made during the partial schedule has a significant impact on the optimal solution for the remaining jobs. Moreover, each node $\mathcal{N}$ of the search tree represents a partial solution $s$, that is composed of a partial schedule $\sigma_s$, a set of remaining jobs $\Omega_s$, and a number of jobs $k_s^1$ (respectively $k_s^0$) that must be scheduled on high-speed (respectively low-speed) machines to complete the solution.

In this context, we employ passive node memorization based on the function $check(s, s')$. To that end, we associate a data structure with each node. This structure consists of: a bit set $\{0,1\}^N$ encoding the set of remaining jobs; an array with two parts representing the total completion time of high-speed (respectively low-speed) machines, both sorted in non-increasing order; the partial solution value $\sum_{j \in \sigma_s} w_j U_j(\sigma_s)$; and integers for $k_s^1$ and $k_s^0$. Thus, given two nodes $\mathcal{N}_1$ and $\mathcal{N}_2$, their data structures enable implementation of the dominance condition from Theorem 2.7 as a function $check(s_1, s_2)$.

We store the data structure of each explored node in a database $Db$. To accelerate the search for a dominating node, we group data structures in the database by the size $|\Omega_s|$. Thus, $Db[k]$ is the set of data structures associated with explored nodes for which $|\Omega_s| = k$. Consequently, given a node $\mathcal{N}_1$, we search for a dominating node $\mathcal{N}_2$ in $Db[k]$ (for $k = 1, \dots, |\Omega_{s_1}|$).

However, the database can become full due to the exponential number of explored nodes. Therefore, we need to define a policy to clear it: at any time the database is full, we remove all nodes that have never dominated any other node.

As authored by Shang et al. (2021), the memorization mechanism has an impact on the choice of the most efficient search strategy. Experiments on that matter are reported in Section 4.

### 3.3. Initial upper bound computation

To compute a good upper bound at the root node of the search tree, we first compute a solution with the MIP by setting a time limit to 20 seconds (experimental results showed that 20 seconds is enough to often get near-optimal solutions).

We also use a heuristic that improves this solution. Let $s = (\sigma_s, \Omega_s)$ be a solution, with $\sigma_s$ the schedule, $\mathcal{I}_s$ the set of jobs that have been scheduled, and $\Omega_s$ the set of jobs that have not been selected. The main idea, is to loop over each block $B_b$, for $b = 1, \dots, b_{max}$. For a given block, we make all jobs in it available. Afterwards, we solve an assignment problem with all the jobs that can be scheduled in $B_b$, i.e., jobs $J_j \in \mathcal{I}_s \cup \Omega_s$ with the processing time $p_j$ such as $\max_{(i,k) \in B_{b-1}} (p_{[i,k]}) \le p_j \le \min_{(i,k) \in B_{b+1}} (p_{[i,k]})$ where $p_{[i,k]}$ is the processing time of the job in position $k$ on machine $M_i$. Finally, we optimize all equal-size jobs using the same algorithm as described in the proof of Theorem 2.2 in Section 2.3.1.
### 3.4. Lower bound computation by Column Generation

We propose an extension of the column generation approach presented in Van den Akker et al. (1999) in order to compute a lower bound at each node of the search tree.

The master problem is formulated as a set-covering problem with an exponential number of binary variables and a linear number of constraints. We use machine schedules, which are defined as a sequence of jobs that can be assigned to any machine. In the master problem, we cannot enforce the block structure, because we do not know the selected jobs. We do know, however, that the machine schedules must be part of a potential schedule constructed by the follower that follows the block structure. Hence, we can impose several constraints on the machine schedules that we consider. First of all, we know that the jobs must be in SPT order. Second, since the block structure is independent of the selection of the jobs, we know that we must fill the blocks $B_b$, for $b = 1, \dots, b_{max}$, which provides information on the number of jobs per machine. Since the blocks $B_b$ are full for $b = 2, \dots, b_{max}$, we know the resulting number of jobs for each fast and slow machine. The only freedom with respect to the number of jobs per machine schedule that we have concerns the $Q$ jobs in the first block, $B_1$. If this block is not filled completely, then we must distribute these $Q$ jobs over the relevant machines. This leads to machine schedules with a minimum or maximum number of jobs, which corresponds to a difference of one job. Finally, we can limit the maximum number of jobs with equal processing time that can occur in a machine schedule. We will work this out later.

Each machine schedule $s$ is characterized by the following set of parameters:

$$a_{j,s} = \begin{cases} 1 & \text{if the job } J_j \text{ is in the schedule } s \\ 0 & \text{otherwise} \end{cases}$$

$$b_{s\gamma+} = \begin{cases} 1 & \text{if the machine schedule } s \text{ is to be scheduled on a machine with speed } V_{\gamma} \\ & \text{and contains the maximum number of jobs} \\ 0 & \text{otherwise} \end{cases}$$

$$b_{s\gamma-} = \begin{cases} 1 & \text{if the machine schedule } s \text{ is to be scheduled on a machine with speed } V_{\gamma} \\ & \text{and contains the minimum number of jobs} \\ 0 & \text{otherwise} \end{cases}$$

$$t_{j,s} = \begin{cases} 1 & \text{if the job } J_j \text{ is late in the schedule } s \\ 0 & \text{otherwise} \end{cases}$$

These constants can be computed for any given schedule $s$.

In the column generation approach, we aim to compute a lower bound for $\sum_j w_j U_j$. Notice that we cannot compute it as in the bilevel scenario, because we cannot impose a constraint that ensures that, when $n$ jobs are selected, they must be scheduled to minimize $Lex(\sum_j C_j, \sum_j w_j U_j)$: this is a consequence of the block structure that cannot be enforced
before we know the $n$ selected jobs. As a result, we focus on minimizing the single-level $\sum_j w_j U_j$ problem in a scenario where we select $n$ jobs. So, the optimal solution of our restricted master problem is a lower bound to the optimal solution of the master problem that is itself a lower bound to the optimal solution of the bilevel problem. For a schedule $s$, we associate the cost $c_s$ defined as follows: $c_s = \sum_{j=1}^{N} w_j t_{j,s}$.

Let $S$ be the set containing all feasible machine schedules. We introduce variables $x_s \in \{0,1\}$, with $s \in S$, and $x_s = \begin{cases} 1 & \text{if the machine schedule } s \text{ is selected} \\ 0 & \text{otherwise} \end{cases}$

Now, our master problem consists of selecting $m$ machine schedules that contain exactly $n$ distinct jobs corresponding to the leader's decision, and must minimize the total cost. These $m$ machine schedules are partitioned in $m_0$ and $m_1$ machine schedules for the machines with speed $V_0$ and $V_1$, respectively. Furthermore, these $m_0$ machine schedules are partitioned into $m_{0+}$ and $m_{0-}$ machine schedules with a maximum and minimum number of jobs; a similar partition holds for the $m_1$ machine schedules. If $B_1$ corresponds to positions on machines with speed $V_0$ (or $V_1$) only, then we can compute the corresponding values of $m_{0+}, m_{0-}, m_{1+}, m_{1-}$, and we find the following master problem formulation:

$$ \min \sum_{s \in S} c_s x_s \qquad (13) $$

**s.t.**

$$ \sum_{s \in S} a_{j,s} x_s \le 1 \quad \forall j = 1, \dots, N \qquad (14) $$

$$ \sum_{s \in S} b_{s0+}x_s = m_{0+} \qquad (15) $$

$$ \sum_{s \in S} b_{s0-}x_s = m_{0-} \qquad (16) $$

$$ \sum_{s \in S} b_{s1+}x_s = m_{1+} \qquad (17) $$

$$ \sum_{s \in S} b_{s1-}x_s = m_{1-} \qquad (18) $$

$$ x_s \in \{0,1\} \quad \forall s \in S $$

In this formulation Constraints (14) ensure that each job is processed at most once. Constraints (15) and (16) ensure that we select the correct number of machine schedules for fast machines with the maximum and minimum number of scheduled jobs. Constraints (17) and (18) ensure the same for the slow machines. Constraints (15), (16), (17) and (18) ensure that we select exactly $n$ jobs, reflecting the leader's decisions.

If $B_1$ is processed on both high and low speed machines, then we don't know the exact values $m_{0+}, m_{0-}, m_{1+}, m_{1-}$. Nevertheless, we know that there are $Q$ jobs to be scheduled in $B_1$. Since there are $m_0$ fast machines, we know that we need a total of $m_0$ machine schedules for fast machines. Similarly, we need $m_1$ machine schedules for slow machines. Furthermore, we need a total of exactly $Q$ machines schedules with a maximum number of jobs.

Given these considerations, we can rewrite the last four constraints as follows:

* Constraints (15) and (16) become: $\sum_{s \in S} (b_{s0+} + b_{s0-})x_s = m_0$. This constraint ensures
that the total number of machine schedules for fast machines with minimal and maximal number of scheduled jobs selected is equal to the number of high speed machines.

* Constraints (17) and (18) become: $\sum_{s \in S} (b_{s1+} + b_{s1-})x_s = m_1$. This constraint ensures that the total number of machine schedules for slow machines with minimal and maximal number of scheduled jobs selected is equal to the number of low speed machines.

* We add a constraint $\sum_{s \in S} (b_{s0+} + b_{s1+})x_s = Q$ that ensures that we select exactly the number of machines with the maximum number of scheduled jobs to be equal to the number of jobs in the first block.

We obtain the linear programming relaxation by replacing the definition set of variables $x_s \in \{0, 1\}$ by $x_s \in [0, 1]$.

*Generating Columns.* In the analysis below we refer to the model with constraints (15), (16), (17) and (18). The analysis can easily be adjusted to deal with the second model. We compute the reduced cost $r_s$ of any machine schedule $s$ as follows:

$$r_s = c_{s\gamma} - \sum_{j=1}^{N} \lambda_j a_{j,s\gamma} - \tilde{\lambda}_{\gamma} \quad (19)$$

$$\text{where } \tilde{\lambda}_{\gamma} = \begin{cases} \lambda_{0+} & \text{if } s \text{ is a fast machine with the maximum number of jobs} \\ \lambda_{0-} & \text{if } s \text{ is a fast machine with the minimum number of jobs} \\ \lambda_{1+} & \text{if } s \text{ is a slow machine with the maximum number of jobs} \\ \lambda_{1-} & \text{if } s \text{ is a slow machine with the minimum number of jobs} \end{cases}$$

and respectively $\lambda_{0+}, \lambda_{0-}, \lambda_{1+}$ and $\lambda_{1-}$ represent the given values of the dual variables corresponding to Constraints (15), (16), (17) and (18). Notice that $\lambda_1, \dots, \lambda_N$ constitute the given values of the dual variables corresponding to constraints (14).

It is well-known from the linear programming theory that a solution to a minimization problem is optimal if the reduced cost of each non basic variable is nonnegative. Therefore, we need to determine if there exists a machine schedule $s$ with such a negative cost. We solve the *pricing problem*, that involves finding a machine schedule in $S$ with the minimum reduced cost. Removing the constant, we need to minimize:

$$r_s = c_s - \sum_{j=1}^{N} \lambda_j a_{j,s} = \sum_{j=1}^{N} (w_j t_{j,s} - \lambda_j a_{j,s})$$

We achieve the calculation of the schedule with a minimal reduced cost by dynamic programming. It is applied for each speed $V_0$ and $V_1$ and for the maximum and minimum number of scheduled jobs.

Define the job groups $G_1, \dots, G_c$ such that jobs within each group have equal processing time. We assume that $p(G_1) < p(G_2) < \dots < p(G_c)$, where $p(G_g)$ denotes the processing time of any job in $G_g$. Thus, the set of jobs $\mathcal{J}$ is partitioned as $\mathcal{J} = G_1 \cup \dots \cup G_c$. Assume that we solve the dynamic programming for speed $V_\gamma$. Define $F_k^\gamma(t, g)$ as the minimum reduced cost of a machine schedule with speed $V_\gamma$ such that the first job from sets $G_g, \dots, G_c$ starts
at time $t$ and where exactly $k$ jobs have been selected from the job sets $G_g, \dots, G_c$. The initialization is $F_k^\gamma(t, g) = 0$ for $k = 0$, and is $+\infty$ otherwise; this implies that we must start with an empty schedule. We find the best machine schedule for a minimum and maximum number of scheduled jobs in one run of the dynamic programming, adding the multiplier $\tilde{\lambda}_\gamma$ for the corresponding constraint at the end. Given the relevant values of $F_k^\gamma(t, g)$, we compute the values of $F_k^\gamma(t, g-1)$ by taking the jobs from $G_{g-1}$ into account. If we select $\ell$ jobs from set $G_{g-1}$, then the first of these jobs will start at time $t$ and the $\ell$-th job will end at time $t + \ell p/V_\gamma$, where $p$ is the processing time of all jobs in set $G_{g-1}$. We can use a subroutine to compute the cost of completing these $\ell$ jobs, which is a polynomial-time algorithm and is described at the end of the section. We denote $S(\ell, g-1, t)$ as the cost for completing $\ell$ jobs from group $G_{g-1}$ that start at time $t$.

The advantage of recording the number of jobs already scheduled is two-fold. First, we only select machine schedules with the right number of jobs. Secondly, we can reduce the number of jobs that we select from a group $G_g$ by setting $\ell_{max}$ as the maximum number of jobs allowed from that group. Suppose that we want to determine $F_k^\gamma(t, g)$, including the jobs from $G_g$; we need to schedule $k$ jobs on the machine with speed $V_\gamma$. Thus, using the block structure, we can identify which block $B_h$ contains the $k^{th}$ job from the end on the machine. If the leader selects all jobs from $G_g$ with at least one job in block $B_h$, then this will fill subsequent blocks, and we denote the last filled block as $B_u$. Thus, $\ell_{max}$ is the maximum number of jobs that can be assigned to one machine with speed $V_\gamma$ from the relevant blocks $B_h$ up to and including $B_u$. These computations are performed as preprocessing – prior to solving any pricing problem – since the block structure remains unchanged throughout the column generation process.

The recurrence relation then becomes equal to:

$$F_k^\gamma(t, g-1) = \min_{0 \le l \le l_{max}} \{F_{k-l}^\gamma(t + \ell p/V_\gamma, g) + S(\ell, g-1, t)\}$$

In this formula $\ell$ corresponds to the number of jobs from $G_{g-1}$ that we add; here $\ell = 0$ implies that no jobs from $G_{g-1}$ are added.

The optimal value $F_k^\gamma$ for the machine with the speed $V_\gamma$ is obtained by computing $F_{\eta_\gamma+}^\gamma(0, 1)$, where $\eta_\gamma+$ is the maximum number of jobs that can be scheduled on machine with speed $V_\gamma$.

Now, let us return to the subroutine $S(\ell, g, t)$. First, if $G_g$ contains only one job denoted $J_g$, then the cost of completing this single job is given by $S(\ell, g, t) = \begin{cases} -\lambda_g & \text{if } t + p_g/V_\gamma \le d_g \\ w_g - \lambda_g & \text{if } t + p_g/V_\gamma > d_g \end{cases}$

Secondly, we consider a set of equal-size jobs that start at time $t$ with $n_g$ jobs $J_1, \dots, J_{n_g}$ in $G_g$, where $p_1 = p_2 = \dots = p_{n_g} = p$. We aim to minimize the reduced cost by selecting and scheduling some of these jobs. We consider all the possible scenarios where 1 job is selected, 2 jobs are selected, ..., $\ell_{max}$ jobs are selected. We solve the equivalent maximization reduced cost problem, where each coefficient is multiplied by -1. Thus, each job has a contribution defined by $(\lambda_j - w_j U_j)$, with $U_j$ equals 1 if the job is late 0 otherwise.

Let $v_j = \lambda_j$ be the contribution of job $J_j$ when it is executed on time, and let $q_j = \lambda_j - w_j$ be its contribution of $J_j$ when it executed tardy. There are three possibilities for job $J_j$: either it is selected as an on-time job, or as a tardy job, or not at all. A naive approach to
solving this problem involves solving $\ell_{max}$ assignment problems, one for each scenario where we select $\ell = 1, \dots, \ell_{max}$ jobs. The total time complexity of this naive algorithm that computes all possible schedules of negative reduced cost is $O(n_g^4)$. We propose below a dedicated algorithm with a total time complexity of $O(n_g^3)$.

Since we do not know the number of jobs that need to be selected beforehand, this may force us to select a job even if its contribution is negative when tardy. To do so, we subtract some value $M$ from each job's contribution. This step aims to reduce to two the number of choices for each job. If the job is on-time or tardy but its contribution is positive, i.e. $v_j - M \ge 0$ or $q_j - M \ge 0$, then we select it. Otherwise, the job's contribution is negative, and we do not select it at all. We can adjust the number of jobs that are selected by modifying the value of $M$.

Assume that $M$ is given, and define, for each job $J_j$, the value $w'_j$ as the difference between being on time with adjusted value ($v_j - M$) and the second option where either the job is executed tardy with the adjusted value ($q_j - M$) or it is not selected. Thus, the value $w'_j$ becomes $w'_j = (v_j - M) - \max(0, q_j - M) = \min(v_j - M, v_j - q_j)$. The problem now essentially becomes a $1 \parallel \sum w'_j U_j$ problem, which can be solved using Dynamic Programming which assumes that the jobs are in EDD order. Let $f_j(k)$ be the total cost of the on-time jobs from $\{J_1, \dots, J_j\}$, where we have exactly $k$ jobs on time. The recurrence relation is

$$f_{j+1}(k) = \max\{f_j(k), f_j(k-1) + w'_{j+1}\} \quad \text{if } d_{j+1} \ge t + kp,$$

and $f_j(k)$, otherwise. The running time of this algorithm is $O(n_g^2)$, since $k$ and $j$ run from 0 to $n_g$. We can reduce the running time of the above dynamic programming algorithm to $O(n_g)$ by renumbering the jobs in order of non-increasing $w'_j$ value and using a clever data structure; we describe this algorithm in Appendix A.

We focus on choosing the value of $M$. First, we select a small perturbation $\epsilon > 0$; since scheduling or removing a job is allowed by the differences $v_j - M$ and $q_j - M$, we have at least four possible values of $M$: $M = v_j - \epsilon$, $M = v_j + \epsilon$, $M = q_j - \epsilon$, and $M = q_j + \epsilon$. Furthermore, the value of $M$ impacts $w'_j = \min(v_j - M, v_j - q_j)$; for all $M < q_j$, $w'_j$ is constant and equal to $v_j - q_j$, while it becomes linear for $M \ge q_j$. Since $w'_j$ determines the order in which jobs are carried out, we need to add values of $M$ where the linear part of $w'_j$ intersects the constant part of $w'_k$ for all distinct jobs $j, k = 1, \dots, n_g$, specifically at $M = v_j - v_k + q_k + \epsilon$ and $M = v_j - v_k + q_k - \epsilon$. If we must select all jobs, the value of $M$ is set to $\min_{1 \le j \le n_g} \{q_j\} - \epsilon$.

If there are equal contributions for on-time jobs, we sort them in non-increasing order of weight $w_j$, since we prefer to select the job with the greatest weight first. Conversely, if there are equal contributions for tardy jobs, we sort them in non-decreasing order of weight $w_j$, since we prefer to select the job with the smallest weight first. After sorting, we apply a slight perturbation to the contributions to ensure that all $M$ values become distinct, leading to the computation of $O(n_g^2)$ distinct values of $M$.

Finally, as the value of $M$ decreases, the number of selected jobs increases. Therefore, we sort the $M$ values in non-increasing order. If, for a given value of $M$, the number of selected jobs exceeds $\ell_{max}$, then we can terminate the subroutine.

So, the total complexity of the algorithm that computes all possible schedules of negative
reduced cost is $O(n_g^3)$ time.

Notably, in many cases, $\ell_{max}$ is small, typically less than or equal to three. In these situations, we employ a brute-force approach that exhaustively tests all possible combinations of job selection. For the cases where $\ell_{max} \ge 4$ jobs, we use the algorithm described above. For implementation details on how to optimize the column generation approach we refer to Appendix B. Also, we attempted to use heuristics to accelerate the generation of improved columns; unfortunately, experimental results indicate that this approach does not enhance the method. Conversely, generating multiple columns during the backtracking phase in dynamic programming has shown to improve the method, as evidenced by our experiments where generating three columns during this phase yields best performances.

To resume, the pricing problem is solved by a dynamic programming that requires

$O(\sum_{j=1}^{N} p_j \times N \times n^3)$ time and space in the worst case.

## 4. Experimental results

In this section we computationally compared CPLEX solver version 22.1.1. applied to the MIP formulation (denoted as MIP) and the branch-and-bound algorithm (denoted as BaB). All algorithms were coded in language C++, and the tests were run on a PC 2GHz AMD Epyc 7702. The source code can be found on the GitHub repository at https://github.com/QuentinSchau/Exact_Algorithms_Bilevel_Scheduling_Industry_Futur. We randomly generated instances where we have $m_1$ high-speed machines with their speed $V_1 = 2$ and $m_0$ low-speed machines with their speed $V_0 = 1$. The processing times were drawn at random from the uniform distribution $U[1, 100]$ and weights from the uniform distribution $U[1, 10]$. Due dates were generated using the uniform distribution $U[P \times (1 - tf - rdd/2), P \times (1 - tf + rdd/2)]$ with $tf, rdd \in \{0.2, 0.4, 0.6, 0.8, 1.0\}$ and $P = \frac{\sum_{j=1}^{N} p_j}{m_1 \times V_1 + m_0 \times V_0}$. Here, $tf$ and $rdd$ define a class of instances and we have 25 classes of instances. For each class, we randomly generated 10 instances of size $N = 40, 50, 60, 70, 80$. For each value of $N$ we set $n = N/4, N/2, 3N/4$, where values were rounded down in case they were fractional. Moreover, we set the total number of machines to $m = 2, 4$, where we have the same number of high-speed and low-speed machines.

Primary experimental results indicated us that using memorization yields the best results for the BaB algorithm when a depth-first search strategy is used. In the remainder, we consider this version of BaB.

We divide these instances into "easy" instances for which, all exact methods solve them to optimality in at most 300 seconds and "hard" instances for which, at least one method does not solve them to optimality within this time budget.

Table 1 presents the results on easy instances. In this table, we report the number of instances solved to optimality (#OPT) out of the 250 per class by both algorithms. For each algorithm, we provide the minimum, average, and maximum computation times, in seconds, for each value of $N$ and $n$.

Table 2 presents the results for "hard" instances. There, #OPT shows the number of instances solved to optimality by each method, as well as the minimum, average, and maximum computation times.
Tables 3 and 4 present the minimum, average, and maximum number of nodes explored by both methods for easy and "hard" instances, respectively.

Tables 1 and 2 demonstrate that BaB enables solving more instances to optimality with a significant improvement in computation time compared to MIP. MIP can solve to optimality all instances with $n = N/2$ and $N = 40$ on two and four machines, whereas BaB can solve all instances for $n = N/2$ until $N = 60$ for two machines and until $N = 50$ for four machines. Additionally, for $n = 3N/4$ and for two machines, MIP does not solve to optimality all instances for $N = 40$, whereas BaB can solve these instances for $N = 50$. Furthermore, although MIP can solve more instances for $n = 3N/4$ and $N = 80$ for four machines on hard instances, on average on easy instances, BaB is faster. Nevertheless, in terms of computation time, BaB is generally faster than MIP.

Moreover, Tables 3 and 4 demonstrate that for instances solved optimally with $n = N/4$ and $N = 40, 50, 60$ on two machines and for $N = 40, 50, 60, 70$ on four machines, the number of nodes explored by BaB is 0. Since we run MIP for several seconds to compute an upper bound, it follows that MIP has identified and validated the optimal solution. For hard instances, BaB explores significantly fewer nodes (at least a factor of 2 reduction on average until $N = 80$ and $n = N/2$) compared to MIP for two machines. However, for instances with four machines, MIP appears to explore fewer nodes. This disparity can be attributed to the branching schemes used by each algorithm, specifically $O((m + 1)^N)$ for BaB. Although this scheme remains competitive with MIP for two machines, its exponential growth with respect to the number of machines has a negative impact on BaB's performance due to the increased number of lower bound computations, which are relatively time-consuming.

Since the leader must select $n$ out of $N$ jobs, leading to $\binom{N}{n}$ possibilities, the combinatorics of selection is greatest when $n = N/2$. Furthermore, we observe that this value does not necessarily lead to the most challenging scenario for the bilevel problem, whereas it does for the leader's problem. This difference can be attributed to the block structure of the follower's problem, as having more jobs allows for more permutations within this block compared to the selection options of the leader.
<table><caption>Table 1: Computation times on easy instances</caption><thead><tr><th rowspan="2">m</th><th rowspan="2">N</th><th rowspan="2">n</th><th rowspan="2">#OPT</th><th colspan="3">MIP</th><th colspan="3">BaB</th></tr><tr><th>t<sub>min</sub></th><th>t<sub>avg</sub></th><th>t<sub>max</sub></th><th>t<sub>min</sub></th><th>t<sub>avg</sub></th><th>t<sub>max</sub></th></tr></thead><tbody><tr><td rowspan="3">40</td><td></td><td>10</td><td>250/250</td><td>0.05</td><td>0.08</td><td>0.33</td><td>0.05</td><td>0.09</td><td>0.32</td></tr><tr><td></td><td>20</td><td>250/250</td><td>0.16</td><td>4.41</td><td>238.83</td><td>0.16</td><td><strong>1.63</strong></td><td><strong>23.17</strong></td></tr><tr><td></td><td>30</td><td>222/250</td><td>0.38</td><td>17.08</td><td>281.96</td><td>0.38</td><td><strong>8.55</strong></td><td><strong>45.97</strong></td></tr><tr><td rowspan="3">50</td><td></td><td>12</td><td>250/250</td><td>0.08</td><td>0.14</td><td>3.18</td><td>0.08</td><td>0.15</td><td>3.13</td></tr><tr><td></td><td>25</td><td>240/250</td><td>0.30</td><td>10.02</td><td>289.95</td><td>0.29</td><td><strong>3.91</strong></td><td><strong>27.90</strong></td></tr><tr><td></td><td>37</td><td>157/250</td><td>0.76</td><td>48.83</td><td>271.72</td><td>0.76</td><td><strong>28.60</strong></td><td><strong>247.08</strong></td></tr><tr><td rowspan="3">60</td><td></td><td>15</td><td>250/250</td><td>0.15</td><td>0.30</td><td>3.23</td><td>0.15</td><td>0.30</td><td>3.23</td></tr><tr><td></td><td>30</td><td>219/250</td><td>0.49</td><td>15.57</td><td>208.50</td><td>0.51</td><td><strong>7.90</strong></td><td><strong>56.44</strong></td></tr><tr><td></td><td>45</td><td>75/250</td><td>1.26</td><td>37.04</td><td>279.53</td><td>1.27</td><td><strong>27.22</strong></td><td><strong>234.12</strong></td></tr><tr><td rowspan="3">70</td><td></td><td>17</td><td>250/250</td><td>0.20</td><td>0.72</td><td>43.22</td><td>0.20</td><td><strong>0.63</strong></td><td><strong>20.46</strong></td></tr><tr><td></td><td>35</td><td>166/250</td><td>0.76</td><td>26.36</td><td>284.26</td><td>0.78</td><td><strong>14.30</strong></td><td><strong>250.45</strong></td></tr><tr><td></td><td>52</td><td>55/250</td><td>1.82</td><td>25.96</td><td><strong>225.60</strong></td><td>1.84</td><td><strong>21.44</strong></td><td>278.99</td></tr><tr><td rowspan="3">80</td><td></td><td>20</td><td>249/250</td><td>0.25</td><td>4.05</td><td>280.34</td><td>0.24</td><td><strong>1.81</strong></td><td><strong>24.32</strong></td></tr><tr><td></td><td>40</td><td>136/250</td><td>1.05</td><td>15.06</td><td>235.28</td><td>1.08</td><td><strong>9.08</strong></td><td><strong>172.46</strong></td></tr><tr><td></td><td>60</td><td>44/250</td><td>2.46</td><td>12.08</td><td>256.70</td><td>2.47</td><td><strong>8.23</strong></td><td><strong>96.26</strong></td></tr><tr><td rowspan="3">40</td><td></td><td>10</td><td>250/250</td><td>0.06</td><td>0.10</td><td>0.20</td><td>0.06</td><td>0.10</td><td>0.22</td></tr><tr><td></td><td>20</td><td>250/250</td><td>0.17</td><td>2.63</td><td>69.37</td><td>0.17</td><td><strong>2.36</strong></td><td><strong>45.78</strong></td></tr><tr><td></td><td>30</td><td>132/250</td><td>0.41</td><td>40.36</td><td>292.22</td><td>0.41</td><td><strong>33.97</strong></td><td><strong>285.52</strong></td></tr><tr><td rowspan="3">50</td><td></td><td>12</td><td>250/250</td><td>0.08</td><td>0.14</td><td>0.44</td><td>0.08</td><td>0.14</td><td>0.43</td></tr><tr><td></td><td>25</td><td>232/250</td><td>0.35</td><td>29.90</td><td>294.48</td><td>0.35</td><td><strong>15.98</strong></td><td><strong>293.98</strong></td></tr><tr><td></td><td>37</td><td>78/250</td><td>0.74</td><td>23.06</td><td>286.44</td><td>0.73</td><td><strong>13.44</strong></td><td><strong>171.93</strong></td></tr><tr><td rowspan="3">60</td><td></td><td>15</td><td>250/250</td><td>0.19</td><td>0.39</td><td>4.18</td><td>0.20</td><td>0.39</td><td>4.17</td></tr><tr><td></td><td>30</td><td>169/250</td><td>0.46</td><td>20.27</td><td>264.98</td><td>0.47</td><td><strong>10.48</strong></td><td><strong>157.46</strong></td></tr><tr><td></td><td>45</td><td>66/250</td><td>1.22</td><td>10.87</td><td>155.13</td><td>1.28</td><td><strong>5.05</strong></td><td><strong>30.88</strong></td></tr><tr><td rowspan="3">70</td><td></td><td>17</td><td>250/250</td><td>0.15</td><td>0.84</td><td>17.11</td><td>0.15</td><td>0.85</td><td>17.07</td></tr><tr><td></td><td>35</td><td>144/250</td><td>0.67</td><td>9.37</td><td>292.39</td><td>0.67</td><td><strong>8.02</strong></td><td><strong>285.59</strong></td></tr><tr><td></td><td>52</td><td>55/250</td><td>1.99</td><td>7.38</td><td><strong>137.67</strong></td><td>1.99</td><td><strong>6.31</strong></td><td>140.04</td></tr><tr><td rowspan="3">80</td><td></td><td>20</td><td>250/250</td><td>0.24</td><td>1.70</td><td>45.36</td><td>0.24</td><td><strong>1.52</strong></td><td><strong>23.20</strong></td></tr><tr><td></td><td>40</td><td>132/250</td><td>1.01</td><td><strong>5.43</strong></td><td><strong>106.35</strong></td><td>1.00</td><td>7.28</td><td>274.19</td></tr><tr><td></td><td>60</td><td>53/250</td><td>2.71</td><td>9.39</td><td>153.82</td><td>2.71</td><td><strong>5.28</strong></td><td><strong>20.43</strong></td></tr></tbody></table>

<table><caption>Table 2: Computation times on hard instances</caption><thead><tr><th rowspan="2">m</th><th rowspan="2">N</th><th rowspan="2">n</th><th rowspan="2">#OPT</th><th colspan="3">MIP</th><th colspan="3">BaB</th></tr><tr><th>t<sub>min</sub></th><th>t<sub>avg</sub></th><th>t<sub>max</sub></th><th>t<sub>min</sub></th><th>t<sub>avg</sub></th><th>t<sub>max</sub></th></tr></thead><tbody><tr><td rowspan="3">40</td><td></td><td>30</td><td>0/28</td><td>300.00</td><td>300.00</td><td>300.01</td><td><strong>28/28</strong></td><td>21.01</td><td>28.69</td><td>39.76</td></tr><tr><td></td><td>25</td><td>0/10</td><td>300.00</td><td>300.00</td><td>300.01</td><td><strong>10/10</strong></td><td>20.76</td><td>22.85</td><td>25.01</td></tr><tr><td></td><td>37</td><td>0/93</td><td>300.00</td><td>300.01</td><td>300.01</td><td><strong>93/93</strong></td><td>29.29</td><td>60.14</td><td>185.26</td></tr><tr><td rowspan="3">60</td><td></td><td>30</td><td>0/31</td><td>300.00</td><td>300.01</td><td>300.01</td><td><strong>31/31</strong></td><td>22.08</td><td>33.58</td><td>72.06</td></tr><tr><td></td><td>45</td><td>1/175</td><td>90.91</td><td>298.81</td><td>300.01</td><td><strong>155/175</strong></td><td>28.31</td><td>153.69</td><td>300.58</td></tr><tr><td></td><td>35</td><td>0/84</td><td>300.01</td><td>300.01</td><td>300.01</td><td><strong>83/84</strong></td><td>35.30</td><td>77.21</td><td>300.04</td></tr><tr><td rowspan="3">70</td><td></td><td>52</td><td>11/195</td><td>58.21</td><td>290.34</td><td>300.02</td><td><strong>59/195</strong></td><td>73.46</td><td>268.91</td><td>301.30</td></tr><tr><td></td><td>20</td><td>0/1</td><td>300.01</td><td>300.01</td><td>300.01</td><td><strong>1/1</strong></td><td>21.31</td><td>21.31</td><td>21.31</td></tr><tr><td></td><td>40</td><td>1/114</td><td>261.60</td><td>299.67</td><td>300.02</td><td><strong>107/114</strong></td><td>33.22</td><td>117.99</td><td>302.06</td></tr><tr><td rowspan="3">80</td><td></td><td>60</td><td>9/206</td><td>26.99</td><td>294.52</td><td>300.02</td><td><strong>16/206</strong></td><td>58.59</td><td>293.56</td><td>301.63</td></tr><tr><td></td><td>40</td><td>19/118</td><td>31.58</td><td>274.93</td><td>300.01</td><td><strong>48/118</strong></td><td>22.40</td><td>225.61</td><td>300.02</td></tr><tr><td></td><td>25</td><td>0/18</td><td>300.00</td><td>300.01</td><td>300.01</td><td><strong>18/18</strong></td><td>21.40</td><td>49.84</td><td>144.24</td></tr><tr><td rowspan="3">60</td><td></td><td>37</td><td>8/172</td><td>38.49</td><td>294.23</td><td>300.01</td><td><strong>18/172</strong></td><td>31.21</td><td>284.66</td><td>300.18</td></tr><tr><td></td><td>30</td><td>2/81</td><td>122.81</td><td>296.64</td><td>300.01</td><td><strong>61/81</strong></td><td>21.15</td><td>149.51</td><td>300.04</td></tr><tr><td></td><td>45</td><td>6/184</td><td>34.65</td><td>293.26</td><td>300.02</td><td><strong>8/184</strong></td><td>31.36</td><td>291.16</td><td>300.29</td></tr><tr><td rowspan="3">70</td><td></td><td>35</td><td>3/106</td><td>52.62</td><td>294.08</td><td>300.02</td><td><strong>39/106</strong></td><td>24.85</td><td>240.49</td><td>300.48</td></tr><tr><td></td><td>52</td><td>4/195</td><td>62.19</td><td>296.77</td><td>300.02</td><td><strong>7/195</strong></td><td>118.10</td><td>295.80</td><td>300.86</td></tr><tr><td></td><td>40</td><td>6/118</td><td>18.87</td><td>288.99</td><td>300.02</td><td><strong>15/118</strong></td><td>36.42</td><td>276.95</td><td>300.49</td></tr><tr><td rowspan="3">80</td><td></td><td>60</td><td><strong>11/197</strong></td><td>30.08</td><td>290.52</td><td>300.02</td><td>0/197</td><td>300.00</td><td>300.03</td><td>301.21</td></tr></tbody></table>
<table><caption>Table 3: Number of explored nodes on easy instances</caption><thead><tr><th rowspan="2">m</th><th rowspan="2">N</th><th rowspan="2">n</th><th colspan="2">MIP</th><th colspan="4">BaB</th></tr><tr><th>N<sub>min</sub></th><th>N<sub>avg</sub></th><th>N<sub>max</sub></th><th>N<sub>min</sub></th><th>N<sub>avg</sub></th><th>N<sub>max</sub></th></tr></thead><tbody><tr><td></td><td>10</td><td>0</td><td>11.12</td><td>1223</td><td>0</td><td>0</td><td>0</td><td>0</td></tr><tr><td rowspan="3">40</td><td>20</td><td>0</td><td>1750.36</td><td>97543</td><td>0</td><td>14.42</td><td>1300</td><td></td></tr><tr><td>30</td><td>0</td><td>4743.21</td><td>125273</td><td>0</td><td>280.86</td><td>5593</td><td></td></tr><tr><td>12</td><td>0</td><td>2.19</td><td>158</td><td>0</td><td>0</td><td>0</td><td></td></tr><tr><td rowspan="3">50</td><td>25</td><td>0</td><td>2815.82</td><td>59516</td><td>0</td><td>44.50</td><td>1510</td><td></td></tr><tr><td>37</td><td>0</td><td>9595.18</td><td>94478</td><td>0</td><td>1356.18</td><td>17518</td><td></td></tr><tr><td>15</td><td>0</td><td>73.30</td><td>2637</td><td>0</td><td>0</td><td>0</td><td></td></tr><tr><td rowspan="6">2</td><td rowspan="3">60</td><td>30</td><td>0</td><td>3955.66</td><td>62184</td><td>0</td><td>218.45</td><td>6404</td><td></td></tr><tr><td>45</td><td>0</td><td>3319.24</td><td>20040</td><td>0</td><td>909.21</td><td>20398</td><td></td></tr><tr><td>17</td><td>0</td><td>255.98</td><td>23073</td><td>0</td><td>0.09</td><td>22</td><td></td></tr><tr><td rowspan="3">70</td><td>35</td><td>0</td><td>3949.63</td><td>47013</td><td>0</td><td>556.61</td><td>31449</td><td></td></tr><tr><td>52</td><td>0</td><td>2019.62</td><td>15805</td><td>0</td><td>792.40</td><td>12706</td><td></td></tr><tr><td>20</td><td>0</td><td>1257.62</td><td>87631</td><td>0</td><td>9.67</td><td>442</td><td></td></tr><tr><td rowspan="2">80</td><td>40</td><td>0</td><td>1696.01</td><td>27584</td><td>0</td><td>337.74</td><td>15842</td><td></td></tr><tr><td>60</td><td>0</td><td>301.95</td><td>6455</td><td>0</td><td>178.77</td><td>4303</td><td></td></tr><tr><td></td><td>10</td><td>0</td><td>11.53</td><td>651</td><td>0</td><td>0</td><td>0</td><td></td></tr><tr><td rowspan="3">40</td><td>20</td><td>0</td><td>1678.14</td><td>29660</td><td>0</td><td>229.43</td><td>47437</td><td></td></tr><tr><td>30</td><td>0</td><td>12026.83</td><td>142538</td><td>0</td><td>14963.45</td><td>170468</td><td></td></tr><tr><td>12</td><td>0</td><td>25.62</td><td>1939</td><td>0</td><td>0</td><td>0</td><td></td></tr><tr><td rowspan="3">50</td><td>25</td><td>0</td><td>10522.53</td><td>115001</td><td>0</td><td>7275.95</td><td>178482</td><td></td></tr><tr><td>37</td><td>0</td><td>3196.86</td><td>45725</td><td>0</td><td>1060.41</td><td>25233</td><td></td></tr><tr><td>15</td><td>0</td><td>120.92</td><td>3581</td><td>0</td><td>0</td><td>0</td><td></td></tr><tr><td rowspan="6">4</td><td rowspan="3">60</td><td>30</td><td>0</td><td>4412.47</td><td>59521</td><td>0</td><td>3662.15</td><td>104169</td><td></td></tr><tr><td>45</td><td>0</td><td>910.39</td><td>11873</td><td>0</td><td>30.05</td><td>770</td><td></td></tr><tr><td>17</td><td>0</td><td>486.33</td><td>15517</td><td>0</td><td>0</td><td>0</td><td></td></tr><tr><td rowspan="3">70</td><td>35</td><td>0</td><td>1193.62</td><td>53012</td><td>0</td><td>2880.17</td><td>201848</td><td></td></tr><tr><td>52</td><td>0</td><td>270.25</td><td>5441</td><td>0</td><td>2923.82</td><td>160809</td><td></td></tr><tr><td>20</td><td>0</td><td>765.20</td><td>29961</td><td>0</td><td>27.23</td><td>2811</td><td></td></tr><tr><td rowspan="2">80</td><td>40</td><td>0</td><td>550.84</td><td>12639</td><td>0</td><td>3044.73</td><td>221898</td><td></td></tr><tr><td>60</td><td>0</td><td>260.09</td><td>10272</td><td>0</td><td>0.04</td><td>1</td><td></td></tr></tbody></table>

<table><caption>Table 4: Number of explored nodes on hard instances</caption><thead><tr><th rowspan="2">m</th><th rowspan="2">N</th><th rowspan="2">n</th><th colspan="2">MIP</th><th colspan="4">BaB</th></tr><tr><th>N<sub>min</sub></th><th>N<sub>avg</sub></th><th>N<sub>max</sub></th><th>N<sub>min</sub></th><th>N<sub>avg</sub></th><th>N<sub>max</sub></th></tr></thead><tbody><tr><td rowspan="6">2</td><td>40</td><td>30</td><td>13370</td><td>32334.43</td><td>107382</td><td>79</td><td>1566.07</td><td>5127</td></tr><tr><td>50</td><td>25</td><td>31952</td><td>63335.10</td><td>120191</td><td>84</td><td>498.80</td><td>1294</td></tr><tr><td></td><td>37</td><td>3510</td><td>28523.70</td><td>89151</td><td>192</td><td>3338.08</td><td>18599</td></tr><tr><td>60</td><td>30</td><td>14447</td><td>39547.10</td><td>151066</td><td>126</td><td>1260.06</td><td>5371</td></tr><tr><td>70</td><td>45</td><td>1062</td><td>17580.82</td><td>65564</td><td>150</td><td>5918.02</td><td>36704</td></tr><tr><td>80</td><td>35</td><td>5318</td><td>31884.17</td><td>115348</td><td>254</td><td>2313.39</td><td>14788</td></tr><tr><td rowspan="6">4</td><td>52</td><td>905</td><td>11075.91</td><td>31534</td><td>374</td><td>6791.83</td><td>33511</td></tr><tr><td>20</td><td>87389</td><td>87389</td><td>87389</td><td>44</td><td>44</td><td>44</td></tr><tr><td>40</td><td>4078</td><td>22124.56</td><td>46693</td><td>204</td><td>3810.09</td><td>23907</td></tr><tr><td>60</td><td>651</td><td>8956.86</td><td>23618</td><td>17</td><td>7103.34</td><td>30108</td></tr><tr><td>40</td><td>30</td><td>7650</td><td>51219.26</td><td>119710</td><td>2171</td><td>134821.12</td><td>260359</td></tr><tr><td>50</td><td>25</td><td>33910</td><td>61079.11</td><td>104747</td><td>505</td><td>23021.56</td><td>120436</td></tr><tr><td rowspan="6">60</td><td>37</td><td>5323</td><td>27257.72</td><td>56562</td><td>1683</td><td>119503.60</td><td>287996</td></tr><tr><td>30</td><td>20591</td><td>46902.17</td><td>108168</td><td>270</td><td>78885.80</td><td>263596</td></tr><tr><td>45</td><td>1428</td><td>11963.22</td><td>34792</td><td>716</td><td>105314.18</td><td>252623</td></tr><tr><td>35</td><td>10431</td><td>25746.80</td><td>62082</td><td>462</td><td>97975.81</td><td>232447</td></tr><tr><td>52</td><td>668</td><td>8647.37</td><td>27620</td><td>2688</td><td>98116.84</td><td>281695</td></tr><tr><td>80</td><td>40</td><td>1951</td><td>18563.06</td><td>64581</td><td>1523</td><td>94636.77</td><td>309197</td></tr><tr><td></td><td>60</td><td>532</td><td>6536.15</td><td>21029</td><td>7254</td><td>106843.60</td><td>294649</td></tr></tbody></table>
## 5. Conclusion

In this paper, we tackled a challenging bilevel parallel machine scheduling problem in an optimistic setting, i.e. when the follower returns among all its optimal solutions, the one that leads to the leader's optimal value. We have shown the $NP$-hardness in a strong sense of the bilevel problem. We also developed several exact approaches, namely, a moderately exponential-time dynamic programming algorithm, a MIP formulation and a branch-and-bound algorithm using memorization mechanism in which we use column generation to compute a lower bound.

While the moderately exponential-time dynamic programming algorithm is only relevant for theoretical purposes, our experimental results demonstrate that the branch-and-bound algorithm outperforms the MIP on all instances, except for $N = 80$ and $n = 3N/4$ on four machines, but its performance deteriorates when dealing with more machines. Both approaches can efficiently solve all instances with up to 40 jobs and 2 or 4 machines.

These results show that the considered bilevel scheduling problem is very challenging and remains difficult to solve in practice, even for a dedicated exact algorithm. At first glance, this seems surprising as optimizing only the follower's objective function can be done in polynomial time. Moreover, the hardest scenario for the leader's problem is not the hardest for the bilevel problem. As a future research direction, we intend to develop heuristic approaches to tackle large instances. These heuristics could exploit the BaB algorithm as a beam-search heuristic or leverage machine learning techniques to approximate the follower's decisions.

## Appendix A. Optimization of the pricing algorithm

In this section, we describe an algorithm with improved computational complexity compared to the dynamic programming approach for computing all possible schedules with negative reduced cost in the case of $n_g$ equal-sized jobs. Recall that this problem boils down to $1|p_j = p| \sum w'_j U_j$, where $w'_j = \min\{v_j - M, q_j - M\}$. To ease notation, we assume without loss of generality that $p = 1$ and that we start at time 0 on the machine instead of starting at time $t$. To compensate for this, we put the due dates equal to

$$d_j \leftarrow \lfloor \frac{d_j - t'}{p} \rfloor;$$

rounding down the due dates does not change the problem, since the jobs must finish at fixed times $t$, which are limited to values within the set $\{1, 2, \dots, n_g\}$. For the remainder of this section, whenever we refer to $t$, it will represent any possible value from this predetermined set. We renumber the jobs in order of non-increasing $w'_j$ value; initially the set of on-time jobs is empty. We add the jobs one by one in order of index to the current set of on time jobs. If the job can be on time together with the already selected on-time jobs, then we put it on time. Otherwise, we mark the job as tardy and execute it at the end of the schedule if $q_j - M > 0$, and discard it if $q_j \le M$. To test whether the next job $J_j$ can be set as on time, we need some an appropriate data structure. Any job with a due date smaller than or equal to zero will be marked as tardy immediately and not be considered as a candidate for the on time set.

We first focus on an efficient way to decide whether or not a job $J_j$ can be scheduled on time together with the already selected on-time jobs. For a given set of on-time jobs, we keep
track of the number of jobs $n_t$ in the on time set that have a due date that is smaller than or equal to $t$, for all $t = 1, \dots, n_g$. It is easily verified that the jobs in the on time set can only be scheduled on time together if $n_t \le t$, for all $t$. If we add job $J_j$ to the on-time set, then this implies that $n_t$ will increase by 1 for all $t$ such as $t \ge d_j$ and remain equal for all $t$ with $t < d_j$. If all new $n_t$ values remain smaller than or equal to $t$, then job $J_j$ is added to the on-time set. Updating this in a straightforward manner takes $O(n_g)$ time, and leads to an $O(n_g^2)$ algorithm. To speed up the process of updating the $n_t$ values, we compute the values $k_t = t - n_t$. The feasibility of the set of on-time jobs is guaranteed as long as $k_t \ge 0, \forall t$.

Determine for each $t = 1, 2, \dots, n_g$ the value $m_t$ as the minimum $k_{t^*}$ value for all $t^* \ge t$, and $t^* \in \{1, 2, \dots, n_g\}$ that is, $m_t = \min_{t^* \ge t} k_{t^*}$. If we plot the $m_t$ values as a function of $t$, then we would see a number of horizontal line segments that move up. We are only interested in the jump points $(m_t, t)$. Therefore, we store the $t$ values corresponding to a jump point and add a pointer to the predecessor and successor $t$ value for each jump point.

When we consider a new job $J_j$ with due date $d_j$, we first determine the largest value $t_1$ such that $t_1 \le d_j < t_2$, where $t_1$ and $t_2$ correspond to jump points. Let $t_0$ be the value of the jump point that precedes $t_1$. Adding job $J_j$ to the set of on time jobs implies that $m_t$ is decreased by 1, $\forall t \ge t_1$, which is feasible if they all remain nonnegative. To avoid useless updates, we rather store the difference between two consecutive jump points, i.e. $(m_{t_{k+1}} - m_{t_k})$. Consequently, adding job $J_j$ only requires to decrease by 1 the value $(m_{t_1} - m_{t_0})$. If this difference becomes 0, then this means that $t_1$ does not correspond to a jump point anymore. This is easily taken care of by adjusting the pointers: the successor of $t_0$ becomes $t_2$ and the predecessor of $t_2$ becomes $t_0$. We further mark job $J_j$ as on-time and add it to the on time set. If $t_1$ is the origin point, i.e., $t_0$ does not exist then it means that the job cannot be on time because $m_{t_1} = 0$. At the end of the algorithm, we then process the jobs in the on-time set by increasing order of their due date.

The remaining point is how to efficiently compute the value $t_1$ for each job $J_j$. This can be done using a disjoint set union find. This data structure supports taking the union of two segments and finding the segment to which a certain value belongs. Its running time is $O(a(n))$, where $a(n)$ is the inverse of the Ackermann function (Tarjan and van Leeuwen, 1984). This is less than 5 for reasonable $n$ values.

We now focus on how to efficiently update the $w'_j$ when the value $M$ changes. We initially create two lists of the jobset: one has the jobs sorted by non-increasing order of the $v_j$'s while, in the other, the jobs are sorted by non-increasing order of the $(v_j - q_j)$'s. As $w'_j = \min(v_j - M, v_j - q_j)$ we compute them by a simple pass through the two lists. Notice that, due to the sorting of the two lists, the $w'_j$'s are computed in non-increasing order of their computed value. Consequently, the running time of computing and sorting the $w'_j$ decreases from $O(n_g \log(n_g))$ to $O(n_g)$.

So, the time complexity of this algorithm to compute a machine schedule, for a given value $M$, is $O(n_g)$.
## Appendix B. Optimization of the Column Generation within the Branch and bound

Since Column Generation is used to provide a lower bound, and solving the master problem does not necessarily lead to obtaining a feasible solution of the bilevel problem, we can stop the column generation process earlier. Thereto, we compute at each iteration $k$ of the column generation a lower bound $LB(k)$ on the optimal value $z_{MP}^*$ of the master problem. Let $r_{\gamma+}$ and $r_{\gamma-}$ denote the outcome values of the pricing problem for a machine with speed $V_{\gamma}$ that contains a maximum or minimum number of jobs. Then it is well known that $LB(k) = z_{RMP}(k) + m_{0+}r_{0+} + m_{0-}r_{0-} + m_{1+}r_{1+} + m_{1-}r_{1-} \le z_{MP}^*$, where $z_{RMP}(k)$ is the optimal value of the restricted master problem. Thus, we have a global lower bound $LB = \max_k LB(k)$. We stop the column generation at iteration $k$ when $z_{RMP}(k) - LB < 1$, in this case the lower bound of a node is given by the weighted number of tardy jobs of the restricted master problem.

Then, when exploring a node in the branch-and-bound algorithm, the scheduling decisions for the jobs have already been made. Consequently, at a given node $N_{child}$, the dynamic programming procedure begins with the last undecided set of jobs, $G_{g'}$, and the current total completion time, $t'$, for machine $M_i$, reflecting the previously made scheduling decisions. We then use these parameters, $(t', g')$, to compute the minimal reduced cost for machine $M_i$ at node $N_{child}$. Furthermore, if $B_1$ is processed on both high and low speed machines, then for the given node $N_{child}$, we might have already scheduled jobs on $B_1$. Let $\bar{m}_{0+}$ denote the number of high speed machines with a maximum of jobs in the node $N_{child}$. Thus we add the constraint $\sum_{s \in S} b_{s0+}x_s \ge \bar{m}_{0+}$, which helps the master problem and gives better experimental results in practice. Remark, at the root node, if the number of jobs to schedule in the first block is greater than the number of low-speed machines, then we can have at least $\max(Q - m_1, 0) = \bar{m}_{0+}$ high-speed machines with a maximum number of jobs.

Moreover, we reuse columns from parent nodes; however, as we explore the search tree using branch-and-bound, a subset of these columns becomes infeasible. Rather than removing each infeasible column immediately, we count the number of times this occurs and define a threshold $\epsilon$. If the proportion of all infeasible columns exceeds $\epsilon$, then we remove them.

Finally, in a search tree, when facing the subproblem of filling the last block, instead of exploring the tree further, we solve it as an assignment problem. Indeed, by the time we need to fill the last block, we already know the completion times of the machines and the jobs that can be scheduled in this last block. Therefore, we can formulate and solve the assignment problem, where our goal is to assign the remaining jobs to the last block in order to minimize the weighted number of tardy jobs.

## Declarations

**Conflict of interest** The authors declare that they have no conflict of interest.

**Ethical approval** This article does not contain any studies with human participants or animals performed by any of the authors.
References

Abass, S.A., 2005. Bilevel programming approach applied to the flow shop scheduling problem under fuzziness. Computational Management Science 2, 279–293. doi:10.1007/s10287-005-0035-z.

Bard, J.F., 1998. Practical Bilevel Optimization. volume 30 of *Nonconvex Optimization and Its Applications*. Springer US. doi:10.1007/978-1-4757-2836-1.

Bianco, L., Caramia, M., Giordani, S., Mari, R., 2015. Grid scheduling by bilevel programming: A heuristic approach. European Journal of Industrial Engineering 9, 101–125. doi:10.1504/EJIE.2015.067450.

Brown, G.G., Carlyle, W.M., Royset, J.O., Wood, R.K., 2005. On the Complexity of Delaying an Adversary's Project, in: The Next Wave in Computing, Optimization, and Decision Technologies, Springer US. pp. 3–17. doi:10.1007/0-387-23529-9_1.

Conway, R.W., Maxwell, W.L., Miller, L.W., 1967. Theory of Scheduling. Addison-Wesley Publishing Company.

Dempe, S., 2002. Foundations of Bilevel Programming. volume 61 of *Nonconvex Optimization and Its Applications*. Kluwer Academic Publishers. doi:10.1007/b101970.

Dempe, S., 2003. Annotated Bibliography on Bilevel Programming and Mathematical Programs with Equilibrium Constraints. Optimization 52, 333–359. doi:10.1080/0233193031000149894.

Dempe, S., 2020. Bilevel Optimization: Theory, Algorithms, Applications and a Bibliography, in: Bilevel Optimization: Advances and Next Challenges. Springer International Publishing, pp. 581–672. doi:10.1007/978-3-030-52119-6_20.

Dessouky, M.I., Lageweg, B.J., Lenstra, J.K., Van de Velde, S.L., 1990. Scheduling identical jobs on uniform parallel machines. Statistica Neerlandica 44, 115–123. doi:10.1111/j.1467-9574.1990.tb01276.x.

Garey, M.R., Johnson, D.S., 1979. Computers and Intractability: A guide to the theory of NP-Completeness. W.H. Freeman and Co.

Garey, M.R., Tarjan, R.E., Wilfgong, G.T., 1988. One-processor scheduling with symmetric earliness and tardiness penalties. Mathematics of Operations Research 13, 330–348.

Graham, R.L., Lawler, E.L., Lenstra, J.K., Kan, A.H.G.R., 1979. Optimization and Approximation in Deterministic Sequencing and Scheduling: A Survey, in: Annals of Discrete Mathematics. Elsevier. volume 5 of *Discrete Optimization II*, pp. 287–326. doi:10.1016/S0167-5060(08)70356-X.

Hardy, G.H., Littlewood, J.E., Pólya, G., 1934. Inequalities (Cambridge Mathematical Library). Cambridge University Press.

Karlof, J.K., Wang, W., 1996. Bilevel programming applied to the flow shop scheduling problem. Computers & Operations Research 23, 443–451. doi:10.1016/0305-0548(95)00034-8.
Kasemset, C., Kachitvichyanukul, V., 2010. Bi-level multi-objective mathematical model for job-shop scheduling: The application of Theory of Constraints. International Journal of Production Research 48, 6137-6154. doi:10.1080/00207540903176705.

Kasemset, C., Kachitvichyanukul, V., 2012. A PSO-based procedure for a bi-level multi-objective TOC-based job-shop scheduling problem. International Journal of Operational Research 14, 50. doi:10.1504/IJOR.2012.046343.

Kis, T., Kovács, A., 2012. On bilevel machine scheduling problems. OR Spectrum 34, 43-68. doi:10.1007/s00291-010-0219-y.

Konur, D., Golias, M.M., 2013. Analysis of different approaches to cross-dock truck scheduling with truck arrival time uncertainty. Computers & Industrial Engineering 65, 663-672. doi:10.1016/j.cie.2013.05.009.

Kovács, A., Kis, T., 2011. Constraint programming approach to a bilevel scheduling problem. Constraints 16, 317-340. doi:10.1007/s10601-010-9102-3.

Lukač, Z., Šorić, K., Rosenzweig, V.V., 2008. Production planning problem with sequence dependent setups as a bilevel programming problem. European Journal of Operational Research 187, 1504-1512. doi:10.1016/j.ejor.2006.09.029.

Shang, L., T'Kindt, V., Della Croce, F., 2021. Branch & Memorize exact algorithms for sequencing problems: Efficient embedding of memorization into search trees. Computers & Operations Research 128, 105171. doi:10.1016/j.cor.2020.105171.

Tarjan, R.E., van Leeuwen, J., 1984. Worst-case Analysis of Set Union Algorithms. J. ACM 31, 245-281. doi:10.1145/62.2160.

T'kindt, V., Della Croce, F., Agnetis, A., 2024. Single machine adversarial bilevel scheduling problems. European Journal of Operational Research 315, 63-72. doi:10.1016/j.ejor.2023.11.018.

Van den Akker, J.M., Hoogeveen, J.A., Van de Velde, S.L., 1999. Parallel Machine Scheduling by Column Generation. Operations Research 47, 862-872. doi:10.1287/opre.47.6.862.

Woeginger, G.J., 2021. The trouble with the second quantifier. 4OR 19, 157-181. doi:10.1007/s10288-021-00477-y.