# Bin Packing with Thresholds: Mathematical Models and Theoretical Results

Ernest Foussard, John Martinovic, Marie-Laure Espinouse, Grégory Mounié, Margaux Nattaf

► To cite this version:

Ernest Foussard, John Martinovic, Marie-Laure Espinouse, Grégory Mounié, Margaux Nattaf. Bin Packing with Thresholds: Mathematical Models and Theoretical Results. 2025. hal-05230011

HAL Id: hal-05230011

https://hal.univ-grenoble-alpes.fr/hal-05230011v1

Preprint submitted on 29 Aug 2025

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 Attribution 4.0 International License
# Bin Packing with Thresholds: Mathematical Models and Theoretical Results

Ernest Foussard¹, John Martinovic², Marie-Laure Espinouse¹, Grégory Mounié³, and Margaux Nattaf¹

¹Univ. Grenoble Alpes, CNRS, Grenoble INP, G-SCOP, 38000 Grenoble, France

²Technische Universität Dresden, Institute of Numerical Mathematics, 01062 Dresden,
Germany

³Univ. Grenoble Alpes, Inria, CNRS, Grenoble INP, LIG, 38000 Grenoble, France

**Abstract**

We present a new variant of one-dimensional bin packing, called the bin packing with thresholds problem. In this scenario, any item is equipped with a size and an additional threshold, and packing an item into a bin is possible if and only if the remaining capacity is larger than both these item-specific parameters. Thus, the feasibility of a subset of items typically depends on the packing order. This problem originates from a joint production/maintenance scheduling application in which processing jobs deteriorate the condition of the machine, meaning that the threshold constraint represents requirements on the machine's condition to process a job. In this article, this novel variant of bin packing is investigated from a theoretical and numerical point of view. First, some important structural properties of an optimal solution are derived. Based on these observations, three integer linear programming (ILP) formulations are presented: a basic assignment model, an exponential-size pattern-based approach, and a pseudo-polynomial arcflow formulation. For these approaches, relations between the associated LP bounds are established mathematically. Finally, all three formulations are compared based on extensive numerical experiments involving differently characterized benchmark sets.

**Keywords** — Packing, Bin Packing Problem, Threshold Constraint, Integer Linear Programming, Column Generation

# 1 Introduction

The *Bin Packing Problem (BPP)* is a well-known strongly *NP*-hard combinatorial optimization problem, in which items of various weights must be packed in bins of finite capacity, such that the number of occupied bins is minimized. In this article, we introduce an original “threshold” constraint, stating that items assigned to each bin must be ordered so that whenever an item is packed, the remaining capacity before packing must be larger than the item’s threshold (and, of course, the item’s size as well). This new problem will be referred to as the *Bin Packing with Thresholds Problem (BPTP)*. Note that choosing the threshold parameter (of every item) identical to the item size leads to a classic BPP, meaning that BPTP inherits BPP’s complexity and negative approximation results, including the strong *NP*-hardness [15].

In terms of applications, the problem under consideration originates from the field of maintenance optimization. Due to recent innovations in sensors, it has become possible to gather precise data on the evolution of the equipment's condition and to derive degradation models from it. Notably,
a health index has been developed and used in recent research studies to derivate maintenance policies [22, 25]. For the sake of exposition, let us consider a scheduling problem on a single machine with a deterministic health index. A health index of zero means that the machine is inoperable for production (e.g., breakdown or risk of failure). Each job processed consumes a fixed amount of health, and may require a minimum health index to be processed, e.g., for risk prevention or safety reasons. Identical maintenance activities may be performed to restore the health index of the machine. Under the assumption that maintenance restores the machine to the good-as-new condition (perfect maintenance), the makespan minimization problem consists in packing jobs between consecutive maintenance activities such that the number of maintenance activities is minimized. This joint production/maintenance scheduling problem is equivalent to the BPTP, and the threshold constraint corresponds to a criterion on the state of health of the equipment, enabling it to perform a task under satisfactory conditions.

A similar constraint was studied by Chen et al. [5] and Chung et al. [6], who considered a flowtime minimization job scheduling problem with jobs causing contamination on equipment, with a contamination limit and the possibility of performing maintenance to reduce it. To the best of our knowledge, there are no other occurrences of such constraint in the literature outside of maintenance. Nevertheless, we believe that this type of constraint is likely to appear in other contexts as well to model criteria on risk, level of service or temporal limits for example. For the sake of exposition, let us briefly mention two such examples:

* In online bin packing problems with uncertainties, the item size is typically only known in a distributional sense, e.g., by its mean value and its variance. In this scenario, the thresholds may act as an implicit capacity buffer to ensure the feasibility of a pattern (at least with high probability) after the item is added and the true item sizes are revealed.

* Let us consider a set of jobs which needs to be assigned to the shifts of a worker or, alternatively, to several workers. It is generally known that employee's resilience decreases proportionally to the duration of work performed during the course of a day. Since, for instance, creative problem-solving requires a different freshness of mind than mundane bureaucratic tasks, the threshold in this setting could be interpreted as the minimum level of mental energy to perform a task efficiently and correctly.

In this article, we study the structural properties of this constraint, derive mathematical models for the BPTP, and compare these approaches based on theoretical and numerical investigations.

Many variants of the classic BPP have been studied in the literature. Coffman et al. [7] provide a comprehensive survey and a classification method for these problems. In the following, we review BPP variants that involve some sort of ordering of the items. Among others, partial orders in bin packing were considered to model specific assembly line balancing problems [31] and resource-constrained scheduling problems [14]. In these variants, an order is imposed on the bins, and the items must be packed in such a way that some precedence constraints are respected. Some online variants also involve the packing order of the items: in the *Ordered Open-End Bin Packing Problem (OOEBPP)*, items are released one by one and the last item of a bin is allowed to exceed the bin capacity. Thus, the feasibility of a bin depends on the order of the items [1]. On the other hand, the threshold constraint implies that depending on its packing order, a subset of items may or may not be feasible. To the best of our knowledge, there are no offline bin packing variants that involve an ordering of the items within the bins. In that regard, the threshold constraint we present in this article is a new take on the ordering of items in cutting and packing.

Many exact methods and mathematical programming models have been developed over the years for
BPP and its equivalent variant known as the *Cutting Stock Problem (CSP)*. Most of these methods have been covered and compared in an in-depth survey in 2016, see [12]. However, in recent years, a few more sophisticated exact approaches, typically involving problem-specific branching schemes and other tailored methods, have been proposed in [8, 10]. These solution methods show very convincing performance on a wide variety of classic and challenging benchmark sets and can be considered state-of-the-art. Although the field of exact algorithms is extensive, many rely on an initial mathematical formulation of the problem – typically an *Integer Linear Programming (ILP)* model – which is then solved using a procedure whose complexity can vary considerably. Hence, ILP-based approaches are still among the most important solution methods studied in the literature. To this end, in the following review we focus on existing ILP formulations for the classic BPP (and CSP):

* The assignment-based ILP formulation is mostly attributed to Kantorovich [19] and is often cited as a formal definition of the BPP and CSP. Despite its polynomial amount of variables and constraints, it is known to perform poorly in practice due to its weak continuous relaxation and symmetric solution space. For this reason, pattern-based formulations and network-flow formulations, involving respectively an exponential or a pseudo-polynomial amount of variables, are preferred for solving bin packing problems.

* The set covering (or equivalently set partitioning) model is a pattern-based formulation for the BPP, and was introduced by Gilmore and Gomory [16]. It is based on the (hypothetical) enumeration of all possible bin fillings, called *patterns*, and a feasible solution is described as a subset of patterns covering every item. The model involves an exponential number of variables and is, therefore, intractable for direct solving. However, its continuous relaxation can be dealt with efficiently using column generation techniques [16] and provides very tight bounds on the optimal value. Indeed, let *I* be an instance of an ILP, $z^*$ be the optimum of *I*, and $z_c^*$ be the optimum of its continuous relaxation. The instance *I* is said to verify the *Integer Property (IP)* if $z^* = z_c^*$, the *Integer Round-Up Property (IRUP)* if $z^* = \lceil z_c^* \rceil$, and the *Modified Integer Round-Up Property (MIRUP)* if $z^* \le \lceil z_c^* \rceil + 1$. It is conjectured that all BPP instances verify the MIRUP property. This conjecture still holds as of today [12, 27]. In practice, when working on randomly generated benchmarks, most instances verify IRUP and some even verify IP. For these reasons, computing the continuous relaxations of such models is of theoretical and practical interest. While column generation does not necessarily yield an integer solution, it is a core component of branch-and-price algorithms [2, 30] and diving heuristics [26], that are respectively among the most effective exact and heuristic approaches for BPP and CSP.

* The one-cut formulation, was proposed independently by Rao [24] and Dyckhoff [13] for the CSP. It relies on the idea that a bin can be separated into a left part and a right part, the left part being an item and the right part a smaller bin known as residual. Contrary to the set covering formulation, it only involves a pseudo-polynomial number of variables and thus can be dealt with directly by ILP solvers. Since the initial publication, the one-cut model has received very little interest in the literature, probably because being less descriptive than other approaches. However, recent works suggest that it could be competitive with the other much more extensively studied ILP models for BPP and CSP [11, 21].

* These other pseudo-polynomial formulations fall under the category of arcflow formulations. Arcflow formulations have been widely used for modeling and solving various *NP*-hard cutting and packing problems, as well as routing [9, 23] and scheduling [17] problems. They prove to
be very efficient in practice and do not require the implementation of complex algorithms, as they can be tackled directly by an ILP solver. In this approach, the problem is reformulated as a network flow minimization problem with decision variables representing the flow through each arc of the network. In fact, the bin is discretized (with the possible states represented by nodes), and the positioning of the items is encoded by arcs linking suitable states. While the theoretical foundations of such approaches for cutting and packing problems are due to (author?) [32], the first formulation of this type was actually proposed by Valério de Carvalho [29]. Later, Brandão and Pedroso provided a different formulation involving a more compact graph construction that proved to be very efficient in practice [3]. More recently, Delorme and Iori developed a new reflect arcflow formulation for the CSP involving a smaller amount of variables and constraints by considering only half of the bin capacity in the network [11]. For completeness, we also mention the DP-flow formulation proposed by Cambazard and O'Sullivan, which is based on the bin packing dynamic programming decision graph. While it does not compete with the other arcflow formulations in practice, it still has an application in constraint programming for the propagation of the knapsack constraint, see [4].

Set covering, one-cut and arcflow formulations have been proven to be equivalent with respect to their LP bounds [11, 21, 29]. As a last note, we highlight that column generation based procedures based on arcflow models have also been investigated in the literature. Valério de Carvalho [28] developed a first branch-and-price algorithm based on such model. Moreover, very recently, de Lima et al. proposed a framework for exactly solving network flow models with a strong linear relaxation. It involves a column generation procedure applying Dantzig-Wolfe decomposition combined with a strong branching scheme. These methods prove to be very efficient when used with Valério de Carvalho's model [10].

In this article, we specifically focus on Valério de Carvalho's arcflow formulation and demonstrate, how it can be adapted to include the threshold constraint, as the more recent formulations involve more restrictive hypotheses that are not compatible with our constraint. In particular, the threshold constraint cannot be directly enforced in the reflect model, as the folding operation makes it impossible to distinguish the first half of the bin from its second half, and may lead to infeasible patterns.

The remainder of the article is organized as follows. The problem, its notations and main properties are formalized in Section 2. In Section 3, we present the three ILP-based solution methods for BPTP: a polynomial assignment model, a column generation procedure based on an exponential set-partitioning formulation, and a pseudo-polynomial arcflow formulation. All three methods are theoretically assessed and numerically compared on benchmark instances in Section 4. Finally, Section 5 concludes this work and states various research perspectives.

# 2 Preliminaries and Structural Properties

In this section, we introduce the problem under consideration from a mathematical point of view. To this end, we first present the problem statement together with the most important notations and a detailed example. After that, we investigate some first theoretical properties of this new bin packing variant. These findings will later be helpful in deriving suitable ILP formulations.

## 2.1 Problem Statement and Notation

Formally, we are given a set of items $i \in I := \{1, \dots, n\}$, each of which specified by a size (or weight) $w_i$ and a threshold $t_i$, and a sufficiently large number of bins of capacity $C$. The goal of BPTP
is to pack the given items into as few bins as possible subject to a classic capacity condition and the following new constraint: When adding an item $i \in I$ of size $w_i$ to a (partially filled) bin, the remaining capacity of that bin has to be made up of at least $t_i$ units (so called threshold condition). The latter implicitly imposes a constraint on the order of packing items into the bins, separating the problem under considerations from traditional bin packing. Without loss of generality, we assume all input data to be positive integers. Moreover, a finite set of bins $K$ with $|K| \le n$ representing an upper bound for the number of bins required in an optimal solution, can be assumed for modeling purposes.

**Definition 1.** By grouping these input data into vectors $\mathbf{w} = (w_i)_{i \in I}$ and $\mathbf{t} = (t_i)_{i \in I}$, we can simply refer to a particular BPTP by an instance $E = (n, C, \mathbf{w}, \mathbf{t})$.

**Definition 2.** Let $E$ be an instance of the BPTP. Then, any feasible items-to-bin assignment is called a pattern.

Remember that the feasibility of a subset of items may depend on the specific packing order. Among other things, this fact is illustrated in the following example.

**Example 1.** Let us consider an instance $E$ with $C = 11$ and $n = 5$ items having the following specifications:

<table><tr><td>Item i</td><td>1</td><td>2</td><td>3</td><td>4</td><td>5</td></tr><tr><td>Size w<sub>i</sub></td><td>3</td><td>4</td><td>8</td><td>3</td><td>3</td></tr><tr><td>Threshold t<sub>i</sub></td><td>9</td><td>6</td><td>8</td><td>9</td><td>4</td></tr></table>

Two items-to-bins assignments are presented in Figures 1a and 1b. The weights of the items are represented as the height of the dark gray rectangles, the thresholds as the height of both dark and light gray rectangles. In both figures, the items are assigned to three bins, so that the capacity constraint is satisfied. However, in terms of the threshold constraint for bin $k = 2$, only the packing shown in Figure 1a is feasible. Indeed, after packing item $i = 4$, the remaining capacity is $C - w_4 = 8 \ge t_2$, and therefore item $i = 2$ can be packed. In contrast, the constellation represented in Figure 1b, obtained by exchanging the packing order of items $i = 4$ and $i = 2$ in bin $k = 2$, is no longer feasible, because the remaining capacity (after packing $i = 2$) is given by $C - w_2 = 7 < t_4$, and thus item $i = 4$ no longer fits in the bin.

## 2.2 Structural Properties and First Theoretical Results

Let us start by establishing two fundamental structural properties of the BPTP. These results will also help to derive some of the models presented later in the article. More precisely, we first prove that without loss of generality, we can assume the threshold $t_i$ of an item $i$ to be at least as large as its weight $w_i$. Afterwards, we establish that the items can be sorted in non-increasing order of $t_i - w_i$ within each bin (or pattern).

**Proposition 1.** Let $E = (n, C, \mathbf{w}, \mathbf{t})$ be an instance of the BPTP. Then, a pattern is feasible for $E$ if and only if it is feasible for $E' = (n, C, \mathbf{w}, \mathbf{t}')$, where $E'$ is obtained from $E$ by setting $t'_i := \max\{w_i, t_i\}$ for all $i \in I$.

*Proof.* Let us consider a feasible pattern belonging to $E$. At the moment, an item is added to that pattern, the remaining capacity has to be $\ge w_i$ (capacity condition) and $\ge t_i$ (threshold condition). Hence, it has to be $\ge \max\{w_i, t_i\}$, and the pattern is also feasible for $E'$. In the opposite direction, let us consider a pattern belonging to $E'$. Obviously, a pattern cannot become infeasible by reducing $t'_i$ (to $t_i$). Thus, the claim is proved. $\square$
<figure><img src="image_2.png" /><figcaption>(a) Configuration with only feasible patterns</figcaption></figure>

<figure><img src="image_3.png" /><figcaption>(b) Configuration with unfeasible pattern in bin $k=2$</figcaption></figure>

**Figure 1:** Two items-to-bins assignments for Example 1

This property is also reasonable from a practical perspective, since the threshold condition would be irrelevant in the case $t_i < w_i$. Hence, without loss of generality $t_i \ge w_i$ for all $i \in I$ can be assumed for any given instance $E$ from now on.

**Proposition 2.** Let $E = (n, C, \mathbf{w}, \mathbf{t})$ be an instance of the BP TP. Then, the items in a feasible pattern can always be sorted with respect to non-increasing values of $t_i - w_i$.

*Proof.* Let us assume that two successive items, say $i=1$ and $i=2$, appearing in a pattern satisfy $t_1 - w_1 < t_2 - w_2$. Moreover, let $C_r \le C$ denote the remaining capacity of the bin before adding item $i=1$. Hence, when packing item $i=1$, we certainly have $C_r - t_1 \ge 0$. In addition, when adding item $i=2$, we have $C_r - w_1 - t_2 \ge 0$. We now show that it is also possible to swap the order of these two items in the packing process. To this end, we note that $C_r - t_2 \ge C_r - w_1 - t_2 \ge 0$ holds, so item $i=2$ can be packed first. After that, we have $C_r - w_2 - t_1 > C_r - w_1 - t_2 \ge 0$ due to the assumption from the beginning of the proof. Hence, item $i=1$ can also be placed second, and the claim is proved. $\square$

**Remark 1.** Note that other classic sorting criteria might fail. For instance, ordering the items with respect to non-increasing values of $t_i$ (or $w_i$) typically excludes feasible patterns. In fact, let us consider an instance $E = (n, C, \mathbf{w}, \mathbf{t})$ with $n=2$ items having $w_1=1, t_1=3, w_2=3, t_2=4$ and bin size $C=5$. Placing item $i=1$ before item $i=2$ is allowed, but the opposite is not possible.

As a consequence of the property shown in Proposition 2, the items of a given instance $E$ can always be assumed to be ordered with respect to non-increasing values of $t_i - w_i$. From a mathematical point of view, this observation also helps to find a tractable characterization of feasible patterns:
**Lemma 1.** Let $E = (n, C, \mathbf{w}, \mathbf{t})$ be an instance of the BPPT in which the items are sorted with respect to non-increasing values of $t_i - w_i$. Then, a vector $\mathbf{p} \in \{0, 1\}^n$ is a feasible pattern if and only if

$$ \sum_{j=1}^{i-1} w_j p_j + t_i p_i \leq C $$

holds for all $i \in I$.

*Proof.* Since the items of any pattern can always be sorted with respect to non-increasing values of $t_i - w_i$, the items of a pattern can be processed in their canonical order (indicated by the item index). Hence, when adding the next (higher-indexed) item, the remaining capacity has to be able to satisfy the threshold constraint. The latter is sufficient, because $w_i \le t_i$ can be assumed without loss of generality. □

Hence, for any given instance $E$, the pattern set $\mathcal{P} := \mathcal{P}(E)$ can now be described by $n$ linear constraints in the following way:

$$ \mathcal{P}(E) := \left\{ \mathbf{p} \in \{0, 1\}^n : \sum_{j=1}^{i-1} w_j p_j + t_i p_i \leq C \text{ for all } i \in I \right\}. $$

This representation will be used in Section 3 to derive some of the integer linear programming formulations for the BPPT.

Prior to the examination of these exact solution approaches, we would like to conclude the current section with an important observation related to the approximation guarantee of common heuristics. To this end, let us introduce the worst-case performance ratio

$$ \sigma(ALG) := \sup_{E} \frac{ALG(E)}{OPT(E)}, $$

where $OPT(E)$ denotes the optimal value of instance $E$, and $ALG(E)$ is the heuristic value obtained by applying an algorithm (that is, a heuristic) to $E$.

**Proposition 3.** We have $\sigma(ALG) = \infty$, whenever ALG refers to first-fit, best-fit, or worst-fit (without any presorting of the items).

*Proof.* Let $(F_k)_{k \in \mathbb{Z}_+}$ be the $k$-th term of the Fibonacci sequence, defined as follows: $F_0 = 0$, $F_1 = 1$, and $F_{n+2} = F_{n+1} + F_n$ for all $n \in \mathbb{Z}_+$. Let $E_n$ denote an instance with bin capacity $C = \sum_{i=1}^n w_i$ and $n$ items with $w_i = F_{n-i+1}$ and $t_i = \sum_{j=1}^i w_j$ for $i \in I$. By packing the items in reversed order of indices, an optimal solution with a single bin is obtained, as $\sum_{i=1}^n w_i = C$ and $C - \sum_{j=k+1}^n w_j = t_k$ for all $k \in \{1, \dots, n\}$. Hence, we have $OPT(E_n) = 1$.

We then prove by induction that any of the heuristics mentioned above yields $ALG(E_n) = n$ bins, when the items are released one by one in increasing order of indices. Let $P_k$ be the predicate: “ALG uses $k$ bins to pack the $k$ first items”. $P_1$ is easily verified, as one bin is required to accommodate a single item. Let us now prove that $P_k \implies P_{k+1}$ for all $k \in \{1, \dots, n\}$. Assuming $P_k$, we show that item $i = k + 1$ does not fit in any of the first $k$ bins, and therefore a $(k + 1)$-th bin has to be opened for it. To this end, let us consider a bin index $j \le k$. We now compute the remaining capacity in bin $j$:

$$ C - w_j = \sum_{i=1}^{n} w_i - w_j = \sum_{i=1}^{k+1} w_i + \sum_{i=k+2}^{n} w_i - w_j = t_{k+1} + \sum_{i=1}^{n-k-1} F_i - F_{n-j+1} $$
Since the Fibonacci sequence is strictly increasing and since

$$ \sum_{i=1}^{n} F_i = \sum_{i=1}^{n} (F_{i+2} - F_{i+1}) = F_{n+2} - F_2 = F_{n+2} - 1 $$

holds for all $n \in \mathbb{Z}_+$, we further have $C - w_j = t_{k+1} + F_{n-k+1} - 1 - F_{n-j+1} \le t_{k+1} - 1 < t_{k+1}$. Hence, the item $i = k + 1$ does not fit in any of the first $k$ bins. Therefore, a new one is opened leading to $P_{k+1}$. Thus $P_k \implies P_{k+1}$ and we are done, since the performance ratio for ALG is at least $O(n)$, leading to an unbounded worst-case for $n \to \infty$. □

**Remark 2.** The same proof can be applied to show that any online heuristic has an unbounded worst-case for the BPTP. This is an important difference to the classic BPP, where also these traditional heuristics (like first fit, best fit, or worst fit) lead to feasible solutions of reasonable quality, see for instance [18].

# 3 Integer Linear Programming Formulations for the BPTP

In this section, we present three different types of ILP formulations to tackle the problem under consideration: an assignment model, a pattern-based approach, and a flow formulation. Besides describing the models, some reduction methods and theoretical results will also be proposed. In particular, the LP bounds obtained by the three formulations will be investigated.

## 3.1 An Assignment-based Formulation

Let us assume the items to be sorted by non-increasing values of $t_i - w_i$. Our first model involves two sets of decision variables.

* $z_k \in \{0, 1\}$ denotes whether bin $k \in K$ is used ($z_k = 1$) or not ($z_k = 0$).

* $x_{ik} \in \{0, 1\}$ denotes whether item $i \in I$ is assigned to bin $k \in K$ ($x_{ik} = 1$) or not ($x_{ik} = 0$).

Based on these definitions, an assignment model can be built as follows:

### Assignment Model for the BPTP (Model 1)

$$ z^1 = \sum_{k \in K} z_k \to \min \\ \text{s.t.} \qquad \sum_{k \in K} x_{ik} = 1, \qquad i \in I, \qquad (1) $$

$$ \sum_{j=1}^{i-1} w_j x_{jk} + t_i x_{ik} \le C z_k, \qquad i \in I, k \in K, \qquad (2) $$

$$ z_k \in \{0, 1\}, \qquad k \in K, \qquad (3) $$

$$ x_{ik} \in \{0, 1\}, \qquad (i, k) \in I \times K. \qquad (4) $$

Conditions (1) state that any item is packed precisely once. In addition, Conditions (2) ensure that there are enough capacity units left for the threshold when packing an item into a used bin ($z_k = 1$). Moreover, no item can be packed into an unused bin ($z_k = 0$). In this version, linking the two variable types is done in the capacity dimension.
**Remark 3.** An alternative way of linking the two types of constraints is given by removing the z-variables from Constraints (2), but to attach them to a new set of cardinality constraints. These constraints are given by

$$ \sum_{i \in I} x_{ik} \leq \mu \cdot z_k $$

for all $k \in K$. In that restriction, $\mu$ is an upper bound for the maximum number of items possible in a bin. In this manuscript, we do not study this alternative variant, because the model size becomes larger, and the modified capacity constraint is mostly weaker than Constraint (2) since the z-variables on the right-hand side do not appear anymore.

In classic bin packing, assignment-based formulations typically suffer from two bad features, see [12]: a relatively weak LP bound, and a highly symmetric solution space. Also for the BPPTP, the server indices of a given solution $(x, z)$ can be permuted without changing the objective value. To investigate the LP bound, we use an auxiliary relaxation of the above model.

**Lemma 2.** Let (LP1) denote the LP relaxation of Model 1. Then, the optimal value of (LP1) is identical to the optimal value of (LP1'), defined by:

**Aggregate Relaxation of Model 1**

$$ \begin{array}{ll} & Z \to \min \\ \text{s.t.} & \displaystyle\sum_{j=1}^{i-1} w_j + t_i \leq C \cdot Z, \quad i \in I, \qquad (5) \\ & Z \in [0, |K|]. \end{array} $$

*Proof.* At first, we note that the optimal value of (LP1) cannot increase by aggregating some constraints. Hence, (LP1') has an optimal value less than or equal to that of (LP1). On the other hand, let $Z$ be an optimal solution of (LP1'). Then, by setting

$$ z_k := \frac{Z}{|K|}, \quad x_{ik} := \frac{1}{|K|} $$

for all $i \in I$ and $k \in K$, we obtain a feasible solution of (LP1) with the same objective value. Indeed, for any $i \in I$, we have

$$ \sum_{k \in K} x_{ik} = |K| \cdot \frac{1}{|K|} = 1. $$

In addition, the validity of (2) follows from scaling (5) with a factor of $1/|K|$. Hence, the optimal value of (LP1) is less than or equal to that of (LP1'), meaning that both are identical. $\square$

With these insights, we can prove the following:

**Proposition 4.** The optimal value of (LP1) is given by

$$ z_{LP}^{1,*} = \frac{1}{C} \cdot \max_{i \in I} \left\{ \sum_{j=1}^{i-1} w_j + t_i \right\}. $$

*Proof.* This result can be obtained by considering the aggregate relaxation, having the same optimal value as (LP1). Conditions (5) state that

$$ Z \geq \frac{1}{C} \cdot \left( \sum_{j=1}^{i-1} w_j + t_i \right). $$
holds for all $i \in I$. Hence, $Z$ as the only variable appearing in the objective function can take the lowest possible value allowed by these conditions. This proves the claim. $\square$

Hence, we have obtained a closed-form term for the optimal value of (LP1).

**Remark 4.** Note that, due to $t_i \ge w_i$ for all $i \in I$, the LP bound of Model 1 is typically larger than the material bound $\sum_{i \in I} w_i/C$ known from classic bin packing.

Despite the observation made in the previous remark, this lower bound has an unbounded worst-case behavior for the BPTP.

**Lemma 3.** We have

$$ \sup_E \frac{z^{1,*}(E)}{z_{LP}^{1,*}(E)} = \infty. $$

*Proof.* Let us consider an instance $E$ with bin capacity $C$ and $n = C$ items having $w_i = w = 1$ and $t_i = t = C$ for all $i \in I$. Then, every item requires a separate bin leading to $z^{1,*}(E) = n = C$. On the other hand, we have

$$ z_{LP}^{1,*}(E) = \frac{1}{C} \cdot \max_{i \in I} \left\{ \sum_{j=1}^{i-1} w_j + t_i \right\} = \frac{\sum_{j=1}^{n-1} w_j + t_n}{C} = \frac{n-1+C}{C} = \frac{2C-1}{C} = 2 - \frac{1}{C}. $$

Hence, for $C \to \infty$ the lower bound becomes arbitrarily weak, and the claim is proved. $\square$

Note that, for classic bin packing, the material bound $\sum_{i \in I} w_i/C$ is at least a 2-approximation of the optimal value. Hence, for the BPTP the quality of the LP bound of Model 1 is even worse. As a consequence of that, this type of formulation cannot be expected to perform well in computational tests.

However, at least the second problem (that is, the symmetric solution space) can be addressed to some extent. To this end, we apply a smart numbering of the bins. By that, we mean that item $i=1$ is always placed in bin $k=1$. Then, item $i=2$ can either be positioned in the same bin or in another bin, which is then called $k=2$, and so on. This strategy implies that only index pairs $(i,k)$ with $k \le i$ are required for a correct model. By that, we eliminate roughly half of the assignment variables and also many constraints of type (2), making the model size much smaller. So, by defining

$$ \Delta := \{(i,k) \in I \times K : k \le i\}, $$

we end up with the following reduced version of the model presented before:

**Reduced Assignment Model for the BPTP (Reduced Model 1)**

$$ z = \sum_{k \in K} z_k \to \min \\ \text{s.t.} \qquad \sum_{k \le i} x_{ik} = 1, \qquad i \in I, \qquad (7) $$

$$ \sum_{j=k}^{i-1} w_j x_{jk} + t_i x_{ik} \le C z_k, \qquad (i, k) \in \Delta, \qquad (8) $$

$$ z_k \in \{0, 1\}, \qquad k \in K, \qquad (9) $$

$$ x_{ik} \in \{0, 1\}, \qquad (i, k) \in \Delta. \qquad (10) $$
Note that also some constraints were slightly modified, because a bin $k \in K$ can only contain items with index $i$ larger than or equal to $k$. As a consequence of that, the summation in (8) now has to start at $j = k$. Obviously, the LP value cannot become worse by renumbering the bins.

**Remark 5.** Another way of reducing the amount of symmetry is to apply a different type of assignment variables. By that, we mean that $x_{ik}$ with $i \ge k$ states whether item $i$ has been assigned to a bin, that was “initialized” (or started) by item $k$. In that interpretation, $x_{kk}$ would indicate that a new bin (that is, bin $k$) was opened, so these variables would have to appear in the objective function. In fact, this approach works without further types of variables.

## 3.2 A Pattern-based Approach and a Column Generation Scheme

Contrary to the assignment-based ILP, the pattern-based (or set partitioning) formulation involves an exponential number of variables, each representing a feasible packing pattern. To be more precise, let us define a variable $\xi_p \in \{0, 1\}$ indicating whether pattern $p \in P$ is used ($\xi_p = 1$) or not ($\xi_p = 0$). Then, an exponential-size formulation is given as follows:

**Pattern-based Model for the BPTP (Model 2)**

$$z^2 = \sum_{p \in \mathcal{P}} \xi_p \rightarrow \min \\ \text{s.t.} \qquad \sum_{p \in \mathcal{P}} \xi_p \cdot p_i = 1, \qquad i \in I, \qquad (11)$$

$$\xi_p \in \{0, 1\}, \qquad p \in \mathcal{P}. \qquad (12)$$

In this model, the feasible patterns are referred to as vectors $p \in \{0, 1\}^n$, and $p_i \in \{0, 1\}$, $i \in I$, denotes the $i$-th component of pattern $p$. Hence, Constraints (11) demand that every item has to appear precisely once in the patterns used. Since the patterns form the columns of the coefficient matrix in the constraints, the words pattern and column will be used interchangeably in the remainder of this article.

One way of dealing with the exponential-size model is Column Generation (CG). The main idea behind CG is that rather than solving the model directly, only a dynamically generated subset of relevant columns (and, thus, variables) is considered. To be more precise, the term column generation only refers to the process of solving the linear relaxation of the pattern-based model. Afterwards, the columns identified within this process can be used in the integer formulation, which is then solved to proven optimality with a subset of columns only. Hence, this method, referred to as the Restricted Master Heuristic (RMH), typically leads to sub-optimal heuristic solutions. However, these feasible points can be expected to be of reasonable quality, especially since the linear relaxation of Model 2 (whose solution is a preparatory step for RMH) is stronger than that of Model 1.

**Proposition 5.** The LP bound of Model 2 is at least as good as the LP bound of Model 1.

*Proof.* Let $\xi = (\xi_p)_{p \in \mathcal{P}}$ be a feasible solution to the LP relaxation of Model 2. For any pattern $p \in \mathcal{P}$, we have $\xi_p \ge 0$ and $\sum_{j=1}^{i-1} w_j p_j + t_i p_i \le C$, implying that

$$0 \le \sum_{p \in \mathcal{P}} \left( C - \sum_{j=1}^{i-1} w_j p_j - t_i p_i \right) \cdot \xi_p = \sum_{p \in \mathcal{P}} C \cdot \xi_p - \sum_{p \in \mathcal{P}} \sum_{j=1}^{i-1} w_j p_j \xi_p - \sum_{p \in \mathcal{P}} t_i p_i \xi_p$$
$$ = C \cdot \sum_{\mathbf{p} \in \mathcal{P}} \xi_{\mathbf{p}} - \sum_{j=1}^{i-1} w_j \underbrace{\sum_{\mathbf{p} \in \mathcal{P}} \xi_{\mathbf{p}} \cdot p_j}_{=1} - t_i \underbrace{\sum_{\mathbf{p} \in \mathcal{P}} \xi_{\mathbf{p}} \cdot p_i}_{=1} = C \cdot \sum_{\mathbf{p} \in \mathcal{P}} \xi_{\mathbf{p}} - \sum_{j=1}^{i-1} w_j - t_i $$

By rearranging the terms, we can conclude

$$ \frac{1}{C} \cdot \left( \sum_{j=1}^{i-1} w_j + t_i \right) \leq \sum_{\mathbf{p} \in \mathcal{P}} \xi_{\mathbf{p}} $$

for all $i \in I$. This proves the claim.

□

**Remark 6.** The previous proof also reveals that the two LP values are the same if and only if

$$ \sum_{\mathbf{p} \in \mathcal{P}} \left( C - \sum_{j=1}^{i-1} w_j p_j - t_i p_i \right) \cdot \xi_{\mathbf{p}} = 0 $$

holds for at least one index $i \in I$. In other terms, there has to be an index $i \in I$, so that

$$ \sum_{j=1}^{i-1} w_j p_j - t_i p_i = C $$

holds for any pattern $\mathbf{p}$ used in an optimal solution of the LP relaxation of Model 2. Note that this condition is very restricting, and hence the LP bound associated with Model 2 is typically better than that of Model 1.

Let us now provide some more details on the application of CG for the problem under consideration. Since the general structure of Model 2 is basically the same as in the set partitioning model of the BPP, the CG method presented in this subsection is a direct generalization of the one proposed by Gilmore and Gomory [16]. We note that (compared to the classic BPP) the set of columns is typically smaller, as some items may not be packed together due to thresholds. For the purpose of column generation, the integrality of the decision variables is relaxed, and $\mathcal{P}$ is replaced by a subset $\mathcal{P}'$ dynamically generated using a specific procedure detailed further in the article. The resulting model is known as the *restricted master problem (RMP)*. At each iteration, a new column that improves the optimal solution of the RMP is computed and added to the set $\mathcal{P}'$. Then, the RMP is solved again, until no column can improve the solution anymore (or a suitably defined stopping criterion is met). In this case, the solution obtained in the last iteration is the optimal solution of the linear relaxation of the master problem (that is, the linear relaxation of Model 2). The procedure is detailed in Algorithm 1.

**Algorithm 1 Column Generation Procedure**

**Step 1:** Initialize $\mathcal{P}'$ such that a feasible solution using a subset of columns of $\mathcal{P}'$ exists.

**Step 2:** Solve the RMP.

**Step 3:** Search for a new column $\mathbf{p} \in \mathcal{P} \setminus \mathcal{P}'$ of negative reduced cost.

**Step 4:** If such a column $\mathbf{p}$ exists, add $\mathbf{p}$ to $\mathcal{P}'$ and go to Step 2. Otherwise, the solution found in Step 2 is optimal.

To find a new pattern in Step 3 of Algorithm 1, the column $\mathbf{p}$ minimizing the reduced cost $r(\mathbf{p}) = 1 - \sum_{i \in I} \alpha_i p_i$ is computed, where $\alpha_i$ are the dual values associated with Constraints (11). This
problem is known as the pricing problem, and it is equivalent to a Knapsack with Thresholds Problem (KTP), that is, a maximization problem that can be seen as a restriction of the BPTP to a single bin. From an abstract point of view, it is defined by the following ILP.

## Knapsack with Thresholds Problem

$$z = \sum_{i \in I} \alpha_i p_i \rightarrow \max \\ \text{s.t.} \qquad \sum_{j=1}^{i-1} w_j p_j + t_i p_i \leq C, \qquad i \in I, \qquad (13)$$

$$p_i \in \{0, 1\} \qquad i \in I. \qquad (14)$$

In the objective function, the total utility of the chosen items is maximized. Constraints (13) display the characterization of a feasible pattern proposed in Section 2.

**Remark 7.** In fact, the KTP is equivalent to the well-known single machine throughput maximization scheduling problem, denoted $1 \parallel \sum w_i U_i$: Given a single machine and a set of jobs with runtime (or duration) $r_i$, deadlines $d_i$, and a penalty $\omega_i$ for late jobs, find an ordering of the jobs that minimizes the total penalty $\sum \omega_i U_i$. Indeed, through a change of variables, thresholds can be interpreted as health deadlines in the knapsack. Each item $i \in I$ is converted into a job such that $r_i := w_i$, $d_i := c - t_i + w_i$, and $\omega_i = v_i$. This conversion is illustrated in Figure 2.

<figure><img src="image_4.png" /><figcaption>Figure 2: Conversion from KTP to throughput maximization scheduling</figcaption></figure>

**Proposition 6.** The KTP can be solved in pseudo-polynomial time by dynamic programming.

*Proof.* The KTP is equivalent to $1 \parallel \sum w_i U_i$. This problem is well known to be weakly $\mathcal{NP}$-hard and solvable by dynamic programming in pseudo-polynomial time [20]. Under the assumption that the items are now sorted by non-increasing order of $t_i - w_i$, let $V(l, q)$ be the optimal value of the pricing problem when only the first $l$ items are considered and the remaining capacity is given by $q$. The dynamic program is expressed as follows:

$$V(l, q) = \begin{cases} 0 & \text{if } l = 0 \text{ or } q = 0, \\ V(l-1, q) & \text{if } l \geq 1 \text{ and } t_l > q, \\ \max\{V(l-1, q), v_l + V(l-1, q - w_l)\} & \text{otherwise.} \end{cases}$$

This recurrence equation allows to solve the problem in $O(nC)$ time. $\square$
**Remark 8.** Note that the sortation applied in the proof is chosen intentionally, as one can think of dynamic programming as a look-up table, and reconstructing the solution then begins in the very last row of this table. Hence, we require the item with largest value $t_i - w_i$ at that position to make sure that all feasible patterns can be obtained.

After completing the column generation procedure, constructing an integer optimal solution is very challenging, as it typically requires tailored branch-and-price approaches. Hence, in this article, we only consider the restricted master heuristic introduced before.

## 3.3 A Flow-based Formulation

Arcflow methods consist in expressing packing problems as flow minimization problems in digraphs that can be solved efficiently using integer linear programming. In this section, we adapt the arcflow formulation proposed by Valério de Carvalho [29] to include the threshold constraint.

Let $G = (V, A)$ be a directed acyclic graph with $C + 1$ vertices representing all possible remaining capacities in a bin, i.e., we have $V = \{C, \dots, 0\}$. Each arc $(u, v) \in V \times V$ with $u > v$ describes a portion of a bin that can possibly be occupied by an item of suitable size, that is, an arc $(u, v)$ exists between two vertices if there exists an item $i \in I$ with $w_i = u - v$. Additionally, dummy arcs, connecting each vertex to the sink 0, are introduced to model unused capacity units. In practice, the size and encoding of the digraph have a strong impact on the efficiency of the method. It is therefore necessary to prune unnecessary or redundant vertices and arcs, and to exploit the dominance rules at our disposal to reduce its size. Algorithm 2 exploits the dominance rule established in Proposition 2 to greedily build this reduced digraph while retaining the feasibility of solutions verifying the rule.

<caption>**Algorithm 2** Graph construction for the arcflow model</caption>

**Require:** Instance $E$ with item list sorted by non-increasing order of $t_i - w_i$
1: $V \leftarrow \{C\}; A \leftarrow \emptyset$ \quad ▷ initialize a graph with only vertex $C$
2: **for** $i \in I$ **do** \quad ▷ for every item
3: $V' \leftarrow V$ \quad ▷ make a copy of the vertex set
4: **for** $u \in V'$ **do** \quad ▷ for every possible starting point
5: **if** $u - t_i \ge 0$ **then** \quad ▷ if the threshold is respected
6: $V \leftarrow V \cup \{u - w_i\}$ \quad ▷ add the ending point to the original vertex set
7: $A \leftarrow A \cup \{(u, u - w_i)\}$ \quad ▷ create the item arc
8: **end if**
9: **end for**
10: **end for**
11: $V \leftarrow V \cup \{0\}$ \quad ▷ add 0 to the vertex set as it may not have been done before
12: **for** $u \in V \setminus \{0, C\}$ **do** \quad ▷ for every vertex in the set besides 0 and $C$
13: $A \leftarrow A \cup \{(u, 0)\}$ \quad ▷ create the loss arc
14: **end for**
**Ensure:** $G = (V, A)$

**Example 2 (Continuation of Example 1).** For the sake of exposition, the network $G$ corresponding to the instance from Example 1 is illustrated in Figure 3.

Note that two of these items share the same input data (that is, the same size and the same threshold), so they can be grouped to effectively deal with four item classes only. For clarity, these item classes are referred to by colored item arcs in the figure: red arcs are caused by an item of size $w = 3$ and threshold $t = 9$ (these are items $i = 1$ and $i = 4$), blue arcs are generated by an item of size $w = 4$ and threshold $t = 6$ (that is, item $i = 2$), green arcs are caused by an item of size $w = 3$ and threshold $t = 4$ (that is, item $i = 5$), and brown arcs are generated by an item with $t = w = 8$ (that is, item $i = 3$). Moreover, dummy arcs (dashed in the figure) are added to connect vertices
<figure><img src="image_5.png" /><figcaption>**Figure 3:** Illustration of the graph *G* belonging to the instance from Example 1</figcaption></figure>

7, 5, 3 and 1 to the sink. Notice that Algorithm 2 allowed to prune the vertices corresponding to the remaining capacities 10, 9, 6 and 2 as well as many arcs, effectively reducing the size of the network. In addition, we highlight that long arcs can typically follow shorter arcs, because the items are not processed according to their sizes.

As a direct consequence of this construction, we obtain the following result.

**Proposition 7.** Let $E$ be an instance of the BPTP. Then, any feasible pattern $\mathbf{p} \in \mathcal{P}(E)$ can be represented as a path from the source $C = 0$ to the sink 0.

Thus, a (C, 0)-flow, as a combination of directed paths, represents a feasible solution of BPTP. The value of the flow then equals the number of bins in the corresponding solution. The other direction of Proposition 7 is not true, as the network allows to build improper patterns. For the sake of exposition, note that two green arcs (representing item $i = 5$) can be part of one path in Figure 3.

To properly state an ILP model, items of same weights and thresholds are grouped together, as in classic bin packing. By that, we obtain the set of all *item classes* (or *item types*), denoted by $\bar{I} = \{1, \dots, \bar{n}\}$, as it was already sketched in the previous example. Note that, by construction, the items in one and the same class all share the same value $t - w$, hence sorting the classes is still possible. To remove any ambiguity, we will use the macron to differentiate the attributes of the item classes, by opposition to the attributes of the items. Hence, the size and the threshold of a given class $i \in \bar{I}$ are referred to by the symbols $\bar{w}_i$ and $\bar{t}_i$, respectively.

An important aspect not covered so far is given by the fact that arcs of the same length can refer to different item classes. By way of example, the red and the green arcs in Figure 3 both refer to items of size 3, not belonging to the same class. Hence, when thinking about classic demand constraints in arcflow models, using such arcs would contribute to satisfying the demand for several item types at once. To avoid this, we introduce the concept of potential demands $\bar{b}_i = |\{j \in I | w_j = \bar{w}_i, t_j \le \bar{t}_i\}|$ for any class $i \in \bar{I}$. By that, we artificially raise the demand of some item classes (those, for which there is at least another item class representing the same item size $w$), to count the number of items assigned correctly.
**Example 3 (Continuation of Example 1).** At first, we observe that the items are grouped in $\bar{n} = 4$ classes, which can be sorted by non-increasing order of $\bar{t}_i - \bar{w}_i$: $\{1, 4\}, \{2\}, \{5\}, \{3\}$. The associated potential demands are given by $\bar{b}_1 = 3, \bar{b}_2 = 1, \bar{b}_3 = 1$, and $\bar{b}_4 = 1$. A feasible flow, verifying all these artificial demand constraints, is illustrated in Figure 4. Seemingly, no green arc (representing item class $i = 3$) is used in that solution. However, the red arcs caused by item class $i = 1$ in the construction can also represent the allocation of a green item (because the threshold of these items is lower than that of the red items), hence they are used more often than allowed by classic demand to fulfill their potential demand. By that, also the demand of item class $i = 3$ is covered.

<figure><img src="image_6.png" /><figcaption>**Figure 4:** A feasible flow belonging to the network G</figcaption></figure>

To summarize, the combinatorial optimization problem to solve consists in minimizing the value of a $(C, 0)$-flow in the digraph $G$ obtained by Algorithm 2, such that the potential demand of every item class is satisfied. To model this problem as an integer linear program, we introduce the set of variables $x_{uv} \in \mathbb{Z}_+$ representing the units of flow along an arc $(u, v) \in A$.

**Arcflow Model for the BPTP (Model 3)**

$$z^3 = \sum_{(C,v) \in A} x_{Cv} \rightarrow \min$$

$$\text{s.t.} \quad \sum_{(u,v) \in A} x_{uv} - \sum_{(v,w) \in A} x_{vw} = 0, \quad v \in V \setminus \{C, 0\}, \quad (15)$$

$$\sum_{\substack{(u,u-\bar{w}_i) \in A, \\ u \ge \bar{t}_i}} x_{u,u-\bar{w}_i} \ge \bar{b}_i, \quad i \in \bar{I}, \quad (16)$$

$$x_{uv} \in \mathbb{Z}_+, \quad (u, v) \in A. \quad (17)$$

The objective function minimizes the total flow emanating the source node $C$. Constraints (15) ensure the flow conservation at any interior vertex. Moreover, Conditions (16) are responsible for satisfying the potential demands. These constraints ensure that, in a feasible flow, a valid arc can be assigned to each item. Overall, the model involves $O(\bar{n} + C)$ constraints and $O(\bar{n}C)$ variables. The actual solution of the BPTP can then be built from the solution of the flow model with the
greedy procedure detailed in Algorithm 3. This algorithm consists in decomposing the flow into paths, and then assigning the items by non-increasing order of $t_i - w_i$ to the first suitable arc. Eventually, each item is assigned to a path, resulting in a feasible bin.

**Algorithm 3 Reconstruction of the solution**

**Require:** Instance *E* with items sorted according to non-increasing values of $t_i - w_i$, and feasible solution $\boldsymbol{x} = (x_a)_{a \in A}$ of Model 3 with objective value *z*

1: Initialize *z* empty bins $B_1, \dots, B_z$.

2: Transform $\boldsymbol{x}$ into a set of *z* paths $\Gamma = \{\gamma_1, \dots, \gamma_z\}$ using flow decomposition properties.

3: **for** each item $i \in I$ **do**

4: &nbsp;&nbsp;&nbsp;&nbsp;**for** each path $k \in \{1, \dots, z\}$ **do**

5: &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;**if** $\gamma_k$ contains a non-assigned arc $a = (u, u - w_i)$ with $u \ge t_i$ **then**

6: &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;Assign arc $a$ to item $i$ and add $i$ at the corresponding position in bin $B_k$.

7: &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;Move to the next item.

8: &nbsp;&nbsp;&nbsp;&nbsp;**end if**

9: **end for**

10: **end for**

**Ensure:** Items-to-bins assignment via $B_1, \dots, B_z$

**Example 4 (Continuation of Example 1).** Applying the procedure presented in Algorithm 3 to the solution displayed in Figure 4 leads to the following items-to-bins assignment, given as ordered lists: $B_1 = (1, 5)$, $B_2 = (4, 2)$, and $B_3 = (3)$. A corresponding visualization can be found in Figure 5 in *A* where arcs denoted by $\emptyset$ represent unoccupied portions of the bins.

The equivalence between the classical BPP arcflow and set-covering formulation was first observed by Valério de Carvalho in [28], and it is also part of a more general equivalence theorem stated by Martinovic et al. in [21]. Similar arguments could be applied to the pattern-based approach and the flow-based formulation of the BPTP presented in the current section. However, since the flow formulation allows improper pattern, we can only prove that the LP bound of Model 2 is at least as good as the LP bound of Model 3. One of the key ideas of this proof was already collected in Proposition 7. To this end, here we only summarize the corresponding result.

**Proposition 8.** The LP bound of Model 2 is at least as good as the LP bound of Model 3.

Nevertheless, the difference can only come from improper patterns being contained in the LP solution of Model 3. Hence, for many instances, we can expect equality between these two bounds. This is confirmed experimentally in Section 4.6.

Last but not least, we mention that there are other pseudo-polynomial formulations for the classic BPP, such as the one-cut model [13, 21] or alternative views on network-based formulations [3, 4, 11] involving some sort of symmetry-breaking on the packing order. In the case of the BPTP, it is not possible to do so, as the order of the items is critical for a solution to be feasible. As such, these methods cannot be directly adapted for the BPTP.

# 4 Experimental Results

## 4.1 Experimental setup and instance generation

Computational experiments were conducted to evaluate the efficiency of the proposed methods. The implementation was carried out in C++17 using the IBM ILOG CPLEX 22.1.1 API. All experiments were performed on a Dell T640 server equipped with two Intel Xeon Silver 4210R processors (2.4GHz, 10C/20T, 13.75MB cache, 9.6GT/s, 100W, Turbo, HT) and 128GB of RAM,
running Debian GNU/Linux 11 (bullseye).

Due to the novelty of the problem under consideration, the methods were tested on a set of randomly generated benchmark instances covering a wide variety of different aspects. These instances can be fully described using combinations of the following five parameters:

* number of items: $n \in \{50, 100, 200, 500, 750, 1000\}$,

* bin capacity: $C \in \{50, 100, 200, 500, 750, 1000\}$,

* expected item size: $\omega \in \{0.02, 0.05, 0.10, 0.20\}$ referring to the percentage of the capacity,

* threshold generation policy: identical (Idc), inverse (Inv), or random (Rnd),

* instance index: $\kappa \in \{0, \dots, 4\}$.

The purpose of the instance generation protocol is to generate diverse instances with a Gaussian distribution of weights and thresholds and to limit the number of extreme items. For each combination of parameters $n$, $C$, and $\omega$, an item weight vector of size $n$ is generated using a beta distribution with an expected value of $\omega C$ and a variance of $4\omega C$ on the interval $[1, C]$ (values are rounded after generation). A thresholds vector is then generated using a beta distribution with an expected value of $s + \omega C$ and a variance of $4(s + \omega C)$ on the interval $[s, C]$, where $s$ is the value at the same index in the size vector. Note that the beta distribution uses two shape parameters, $\alpha$ and $\beta$, which can be determined by the expected value and variance. Some parameter combinations may be infeasible, such as when $\omega$ and $C$ are small, resulting in insufficient bin capacity to generate items with sufficient weight and threshold dispersion. In such cases, the parameter combination is discarded.

Three instances with $n$ items and bin capacity $C$ are constructed using these vectors. We denote by $s_i$ and $u_i$ the $i$-th elements of the respective weight and threshold vectors:

* identical (Idc): Both size and threshold vectors are sorted in increasing order, and the $i$-th item is assigned the weight $w_i = s_i$ and threshold $t_i = \max\{s_i, u_i\}$.

* inverse (Inv): Both size and threshold vectors are sorted in increasing order, and the $i$-th item is assigned the weight $w_i = s_i$ and threshold $t_i = \max\{s_i, u_{n-i+1}\}$.

* random (Rnd): The $i$-th item is assigned the weight $w_i = s_i$ and threshold $t_i = u_i$. Here, $t_i \ge w_i$, $i \in I$, is ensured by the choice of the distribution parameters from above.

Five instances for each combination of the previous parameters are generated, denoted by the last instance parameter $\kappa$. In total, the benchmark contains 1968 instances and can be found online at https://doi.org/10.5281/zenodo.16962456.

In our computational experiment, the solver was always given a warm start. To this end, we conducted preliminary internal tests (not reported here) to identify the best heuristic from a set of classic candidates. Based on these insights, the following configurations will be investigated:

* the assignment-based model (Reduced Model 1) with the solver being initialized with the solution provided by the first-fit algorithm (with items sorted in decreasing order of $t_i - w_i$),

* the arcflow model with the solver being initialized with the solution provided by the best-fit algorithm (with items sorted in decreasing order of $t_i - w_i$),
* the (heuristic) pattern-based approach, where both the pattern pool (for CG) and the RMH are initialized with the information obtained by the same best-fit algorithm as in the arcflow approach.

Throughout this section, the time limit per instance is established at 5 minutes for the raw models. For the column generation method, the time limit is divided into two segments: 5 minutes allocated for the column generation phase and an additional 5 minutes designated for the RMH.

## 4.2 Preliminary testing on the assignment-based ILP

In the bin packing literature, assignment-based models are well known to be highly inefficient compared to pattern-based or flow-based formulations. For our variant of the problem, we have demonstrated that the relaxation of the assignment-based model can be made arbitrarily weak, exacerbating the inefficiencies observed in classical bin packing problems. Furthermore, we have shown that the relaxations of set-covering and arcflow formulations consistently yield superior results. Consequently, the assignment-based model is anticipated to underperform relative to other methods. To validate this hypothesis and confirm that this model can be excluded from future experiments, preliminary tests were conducted. These experiments were performed on a subset of the generated instances with $\omega = 0.10$ and random threshold policy, encompassing a total of 180 instances of varying sizes. The results, aggregated by item size and bin capacity, are presented in Tables 1 and 2, respectively.

<table><caption>**Table 1:** Optimal and feasible solution rate, and computation time per method when varying n.</caption><thead><tr><th rowspan="2">n</th><th colspan="2">Feasibility</th><th colspan="2">Optimality</th><th colspan="2">Computation time (s)</th></tr><tr><th>Assign. ILP</th><th>Arcflow</th><th>Assign. ILP</th><th>Arcflow</th><th>Assign. ILP</th><th>Arcflow</th></tr></thead><tbody><tr><td>50</td><td>100%</td><td>100%</td><td>87%</td><td>100%</td><td>40.34</td><td>2.60</td></tr><tr><td>100</td><td>100%</td><td>100%</td><td>80%</td><td>100%</td><td>60.66</td><td>5.85</td></tr><tr><td>200</td><td>100%</td><td>100%</td><td>41%</td><td>93%</td><td>178.16</td><td>40.83</td></tr><tr><td>500</td><td>100%</td><td>100%</td><td>17%</td><td>67%</td><td>260.60</td><td>110.75</td></tr><tr><td>750</td><td>100%</td><td>100%</td><td>10%</td><td>67%</td><td>290.45</td><td>112.16</td></tr><tr><td>1000</td><td>20%</td><td>100%</td><td>0%</td><td>70%</td><td>319.59</td><td>124.15</td></tr></tbody><tfoot><tr><th>Total</th><td>87%</td><td>100%</td><td>39%</td><td>83%</td><td>191.70</td><td>66.13</td></tr></tfoot></table>

<table><caption>**Table 2:** Optimal and feasible solution rate, and computation time per method when varying c.</caption><thead><tr><th rowspan="2">C</th><th colspan="2">Feasibility</th><th colspan="2">Optimality</th><th colspan="2">Computation time (s)</th></tr><tr><th>Assign. ILP</th><th>Arcflow</th><th>Assign. ILP</th><th>Arcflow</th><th>Assign. ILP</th><th>Arcflow</th></tr></thead><tbody><tr><td>50</td><td>100%</td><td>100%</td><td>73%</td><td>100%</td><td>105.72</td><td>0.10</td></tr><tr><td>100</td><td>87%</td><td>100%</td><td>50%</td><td>100%</td><td>156.39</td><td>0.32</td></tr><tr><td>200</td><td>83%</td><td>100%</td><td>41%</td><td>100%</td><td>180.72</td><td>0.98</td></tr><tr><td>500</td><td>83%</td><td>100%</td><td>27%</td><td>100%</td><td>226.76</td><td>49.95</td></tr><tr><td>750</td><td>83%</td><td>100%</td><td>27%</td><td>50%</td><td>225.23</td><td>175.98</td></tr><tr><td>1000</td><td>83%</td><td>100%</td><td>17%</td><td>47%</td><td>255.60</td><td>169.44</td></tr></tbody><tfoot><tr><th>Total</th><td>87%</td><td>100%</td><td>39%</td><td>83%</td><td>191.70</td><td>66.13</td></tr></tfoot></table>

In these tables, we list the average rate of success in finding a feasible (or proven optimal solution) within the time limit and the average computation times. Our findings indicate that all instances solved to optimality by the assignment-based ILP were also solved to optimality by the arcflow
method, but typically much faster. The sole exception was an instance with $n = 50$ items and a bin capacity of $C = 1000$. The assignment-based ILP exhibited substantial limitations, failing to find optimal solutions for instances with $n \ge 500$ for any bin capacity, or for instances with $n \ge 200$ and $C \ge 200$. In addition, solver failure was consistently observed for instances with $n = 1000$.

In contrast, the arcflow method demonstrated greater efficiency. For bin capacities $C \le 500$, it successfully solved all instances. However, for $C \ge 750$ and $n \ge 500$, the arcflow method did not always identify the optimal solution. Despite this, arcflow consistently provided a tight absolute optimality gap, never exceeding 1, whereas the assignment-based ILP often resulted in larger gaps, up to 5. As expected, the assignment-based ILP is clearly dominated by the arcflow method and should therefore be discarded from further large-scale experiments.

Two supplementary tables comparing the sizes of the ILP models are included in Tables 9 and 10 in B. Although the model sizes do not scale with the same parameters, these tables indicate that, in practice, the arcflow model involves significantly fewer constraints (with a difference of two orders of magnitude in average). In terms of variables, we observe the following:

* For the assignment model, the number increases drastically with $n$, but does almost not depend on $C$. This behavior is justified, since the number of variables is $O(n^2)$ and independent of $C$ for that type of formulation. The small differences (for varying $C$) are caused by the fact that different instances were considered per fixed value of $C$.

* For the flow model, the number of variables is $O(\bar{n}C)$, so raising either of these parameters leads to a larger number of variables. Typically, doubling $C$ leads to more than twice the number of variables, since almost all partial bin fillings (exceeding a certain bound) can be combined by a given set of items. In terms of $n$, we observe that the growth of the number of variables is stagnating at some point, because $C$ is fixed and, hence, almost all nodes have been activated. This fact also explains why arcflow needs approximately the same average solution time for $n \ge 500$.

Hence, averaging over the entire benchmark set, the size of arcflow is smaller. This observation holds in particular for the difficult scenarios containing many items. These insights, combined with the properties of the LP bounds, explains the observed differences between the two models.

## 4.3 Comparison of the methods and assessment of their limitations

After having figured out that the assignment model is rather unsuitable for dealing with BPTP, we now compare the arcflow model with the pattern-based approach, that is, CG plus RMH, but on the complete benchmark set (not only the subset considered in the previous subsection). The comprehensive benchmark results are presented in Tables 3 and 4.

For each method, we again detail the rate of optimal solution attainment and the average computation time. Both methods successfully provide feasible solutions for all benchmark instances, thus the feasibility rate is omitted from the tables. The performance of the arcflow method is sensitive to both bin capacity and the number of items. This behavior was anticipated, as the ILP model involves $O(\bar{n}C)$ variables. Graph reduction techniques significantly reduce the number of these variables and are more effective on instances with fewer items. On the other hand, the performance of column generation combined with RMH appears to be sensitive primarily to the number of items, because that feature leads to more difficult subproblems (and larger pattern sets) from a combinatorial point of view. This point will also be supported by the data appearing later in Table 5.
<table><caption>**Table 3:** Optimal solution rate and computation time per method when varying *n*.</caption><thead><tr><th rowspan="2">n</th><th colspan="2">Optimality</th><th colspan="2">Computation time (s)</th></tr><tr><th>CG+RMH</th><th>Arcflow</th><th>CG+RMH</th><th>Arcflow</th></tr></thead><tbody><tr><td>50</td><td>11.8%</td><td>100.0%</td><td>0.24</td><td>0.81</td></tr><tr><td>100</td><td>0.0%</td><td>100.0%</td><td>2.68</td><td>3.73</td></tr><tr><td>200</td><td>0.0%</td><td>98.8%</td><td>49.07</td><td>19.43</td></tr><tr><td>500</td><td>0.0%</td><td>89.4%</td><td>193.26</td><td>56.17</td></tr><tr><td>750</td><td>0.0%</td><td>86.5%</td><td>324.77</td><td>62.40</td></tr><tr><td>1000</td><td>0.0%</td><td>84.6%</td><td>368.11</td><td>70.01</td></tr></tbody><tfoot><tr><td>Total</td><td>2.0%</td><td>93.2%</td><td>155.16</td><td>35.26</td></tr></tfoot></table>

<table><caption>**Table 4:** Optimal and feasible solution rate, and computation time per method when varying *c*.</caption><thead><tr><th rowspan="2">C</th><th colspan="2">Optimality</th><th colspan="2">Computation time (s)</th></tr><tr><th>CG+RMH</th><th>Arcflow</th><th>CG+RMH</th><th>Arcflow</th></tr></thead><tbody><tr><td>50</td><td>0%</td><td>100%</td><td>151.43</td><td>0.15</td></tr><tr><td>100</td><td>0%</td><td>100%</td><td>157.31</td><td>0.27</td></tr><tr><td>200</td><td>3.3%</td><td>100%</td><td>147.17</td><td>0.80</td></tr><tr><td>500</td><td>2.5%</td><td>100%</td><td>159.15</td><td>18.78</td></tr><tr><td>750</td><td>3.3%</td><td>87.2%</td><td>157.06</td><td>71.98</td></tr><tr><td>1000</td><td>1.7%</td><td>75.8%</td><td>158.34</td><td>100.89</td></tr></tbody><tfoot><tr><td>Total</td><td>1.98%</td><td>93.2%</td><td>155.16</td><td>35.26</td></tr></tfoot></table>

The results confirm the findings from the previous subsection regarding the arcflow method: feasible solutions are consistently found, and optimal solutions are always identified if $c \le 500$, or $n \le 500$ and $c \le 750$. In other cases, the optimality gap remains very small. As a heuristic method, column generation combined with RMH is not expected to consistently yield optimal solutions. Yet, it is not faster than the arcflow model which is an exact method. Although the current performance of CG plus RMH does not surpass that of the arcflow model, there is some potential for improvement. Techniques such as stabilization and the integration of cutting planes could significantly accelerate the convergence of the method.

## 4.4 Limits of column generation

As observed in the previous subsection, the heuristic approach (CG plus RMH) consistently yields feasible solutions. Having a look at the two components of that solution method, CG successfully finds the optimum of the relaxation for 80.2% of the benchmark instances. Then, solving the pattern-based model with the columns found by CG finds the integer optimum of the RMP for 81.6% of the benchmark instances. Both methods converge on 63% of the instances, indicating a very low opposite correlation between CG and integer RMP stopping criteria, with a Pearson coefficient of -0.16.

It is important to note that when the column generation step does not converge, it is impossible to provide a lower bound for the problem, as the value obtained so far could be higher than the optimal value of the pattern-based integer model. However, the patterns generated can still be used to obtain an upper bound by solving the RMH, without guaranteeing the quality of this bound. The average gap, when CG stops, is 11.33, which reduces to 4.20, when both CG and RMH converge. This suggests that the heuristic solutions found are approximately four bins away
from the optimum, which can be considered a reasonable approximation provided the method is computationally efficient.

Table 5 illustrates the average number of patterns generated by CG as a function of the number of items and bin capacity across the entire benchmark.

<table><caption>**Table 5:** Number of patterns generated depending on the number of items and bin capacity</caption><thead><tr><th>Parameter</th><th>C = 50</th><th>C = 100</th><th>C = 200</th><th>C = 500</th><th>C = 750</th><th>C = 1000</th><th>Average</th></tr></thead><tbody><tr><td>n = 50</td><td>133</td><td>155</td><td>176</td><td>256</td><td>248</td><td>234</td><td>205</td></tr><tr><td>n = 100</td><td>296</td><td>302</td><td>378</td><td>399</td><td>430</td><td>402</td><td>374</td></tr><tr><td>n = 200</td><td>690</td><td>663</td><td>740</td><td>738</td><td>743</td><td>747</td><td>724</td></tr><tr><td>n = 500</td><td>2020</td><td>1933</td><td>1736</td><td>1576</td><td>1562</td><td>1545</td><td>1706</td></tr><tr><td>n = 750</td><td>2511</td><td>2385</td><td>2094</td><td>2030</td><td>2030</td><td>2023</td><td>2141</td></tr><tr><td>n = 1000</td><td>2802</td><td>2640</td><td>2366</td><td>2335</td><td>2341</td><td>2331</td><td>2446</td></tr><tr><td>Average</td><td>1357</td><td>1346</td><td>1248</td><td>1222</td><td>1226</td><td>1214</td><td>1261</td></tr></tbody></table>

As precisely one pattern is added per iteration, these number can also be seen as an indicator, how often the pricing problem had to be solved. As already alluded to earlier, raising the number of items typically leads to much more difficult subproblems (due to the combinatorial nature of the pricing problems), so that more iterations are needed to solve the LP relaxation. The reduction in the number of patterns for large bin capacities is partially attributable to the increased difficulty in solving the pricing problem for larger capacities, which exhibits a pseudo-polynomial time and space complexity of $O(nC)$. Consequently, fewer patterns can be generated within the time limit, potentially biasing the results. However, even when restricting the analysis to instances where CG was successful, it appears that large capacity instances sometimes require fewer patterns, because (on average) more items can be packed into a bin. Overall, an almost linear dependence of the number of patterns generated on the number of items is observed, with a minor impact from bin capacity, resulting in a slight reduction in the number of patterns.

Table 6 presents the performance metrics of the components of the pattern-based approach as a function of the number of patterns. For both column generation and the restricted master heuristic, the table includes the optimality rate (i.e., for CG, stopping with the set partitioning lower bound, and for RMH, stopping with the integer optimum of the RMP) and the average runtime.

<table><caption>**Table 6:** Completion rate and runtime of both steps of the column generation procedure depending on the number of patterns</caption><thead><tr><th>#Patterns</th><th>Relaxation Optimum</th><th>RMP Optimum</th><th>CG Runtime (s)</th><th>RMH Runtime (s)</th></tr></thead><tbody><tr><td>0-499</td><td>100.0%</td><td>100.0%</td><td>0.31</td><td>1.15</td></tr><tr><td>500-999</td><td>100.0%</td><td>87.6%</td><td>2.94</td><td>41.12</td></tr><tr><td>1000-1499</td><td>100.0%</td><td>77.3%</td><td>26.53</td><td>91.55</td></tr><tr><td>1500-1999</td><td>46.6%</td><td>79.9%</td><td>200.37</td><td>63.32</td></tr><tr><td>2000-2499</td><td>61.5%</td><td>66.2%</td><td>194.22</td><td>114.05</td></tr><tr><td>2500-2999</td><td>68.7%</td><td>40.4%</td><td>199.03</td><td>197.10</td></tr><tr><td>3000-3499</td><td>47.5%</td><td>37.3%</td><td>240.37</td><td>195.79</td></tr><tr><td>3500-3999</td><td>62.5%</td><td>0.0%</td><td>222.08</td><td>300.11</td></tr><tr><td>Average</td><td>81.4%</td><td>80.2%</td><td>88.41</td><td>66.75</td></tr></tbody></table>

CG converges with the lower bound on all instances with fewer than 1500 patterns, which, according to Table 5, corresponds to instances with between 200 and 500 items. It is noteworthy that there
is limited data available for instances with 1000-1500 patterns. Even with a significantly larger number of items, CG still converges on the majority of instances. The performance of RMH declines more gradually, with some instances within the 500-999 pattern range already failing to converge. No instances with more than 3500 patterns were solved to the RMP optimum. In this iteration of the method, CG and RMH appear to have a comparable impact on performance, with neither step of the solution method emerging as a clear bottleneck.

## 4.5 Impact of the thresholds

In this subsection, we evaluate how different threshold shapes impact the performance of the proposed solution methods. To this end, Table 7 presents the experimental results for the full benchmark, categorized by threshold types, including the rate of optimal solutions, computation time, and the absolute gap to the best-known lower bound.

<table><caption>**Table 7:** Optimal solution rate, computation time and absolute gap to the best known lower bound per method depending on the threshold generation method</caption><thead><tr><th rowspan="2">Threshold type</th><th colspan="2">Optimality</th><th colspan="2">Computation time (s)</th><th colspan="2">Abs. Gap to Best LB</th></tr><tr><th>CG+RMH</th><th>Arcflow</th><th>CG+RMH</th><th>Arcflow</th><th>CG+RMH</th><th>Arcflow</th></tr></thead><tbody><tr><td>Idc</td><td>1.98%</td><td>95.88%</td><td>172.01</td><td>24.26</td><td>10.10</td><td>0.05</td></tr><tr><td>Inv</td><td>1.98%</td><td>94.21%</td><td>140.51</td><td>37.56</td><td>14.17</td><td>0.08</td></tr><tr><td>Rnd</td><td>1.98%</td><td>89.63%</td><td>152.96</td><td>43.96</td><td>7.61</td><td>0.12</td></tr></tbody></table>

Overall, we observed that the threshold generation method only moderately affects the performance of both methods. For the arcflow method, instances with randomly assigned thresholds appear more challenging. To better understand the reasons of the disparities between different instance types, we summarized in Table 8 some characteristics of the arcflow digraph based on the threshold generation method.

<table><caption>**Table 8:** Average number of vertices, arcs and items classes in the arcflow procedure depending on the threshold generation method</caption><thead><tr><th>Parameter</th><th>Idc</th><th>Inv</th><th>Rnd</th></tr></thead><tbody><tr><td>|V|</td><td>442.63</td><td>442.69</td><td>442.69</td></tr><tr><td>|A|</td><td>14472.85</td><td>17178.30</td><td>19684.68</td></tr><tr><td>n&#773;</td><td>47.84</td><td>47.72</td><td>302.86</td></tr></tbody></table>

Notably, the number of item classes is significantly larger for Rnd instances. Randomly pairing item sizes and thresholds reduces the likelihood of having multiple identical items compared to the other two generation protocols, resulting in a larger number of classes and, consequently, variables. Therefore, the threshold generation method strongly impacts the size of the arcflow ILP model. We also highlight that the starting point obtained by best-fit for Inv and Rnd may be of worse quality than for Idc, meaning that the ILP solver typically has to start with less heuristic information. The reason is that our version of best-fit works on an item list sorted with respect to non-increasing values of $t_i - w_i$, to not exclude any feasible pattern. However, it is known from classic bin packing that best-fit is particularly strong, when the item classes are sorted with respect to non-increasing values of $w_i$, leading to the classic *best-fit decreasing* algorithm. This property is not given, if thresholds do not scale proportionally with the item sizes (like for Inv and Rnd). To be more precise, sorting the items with respect to non-increasing values $t_i - w_i$ potentially leads to a list
starting with many small items. As a consequence of that, these small items are combined cleverly, but then the larger items (coming too late) will lead to many unnecessary bins, deteriorating the quality of the heuristic information passed to the solution methods.

Similarly, for the pattern-based method, the inversely proportional assignment seems to present greater difficulties. Again, we partly attribute this to the initialization of the pattern pool (both, for CG and the RMH) with a best-fit approach, possibly leading to heuristic information of worse quality. While the threshold generation method does have an impact, it is not as pronounced as the impact of instance size parameters. The number of patterns generated appears to be higher on average (1318.05 for Idc, 1212.47 for Inv, and 1251.81 for Rnd) for Idc instances, suggesting that column generation convergence may be slower for these instances.

## 4.6 Lower bounds

Finally, we conducted additional experiments to study the LP bounds of the different models and experimentally confirm the results of Sections 2 and 3. For the sake of exposition, we consider the subset of instances with $n \le 500$ and $C \le 500$, ensuring that the arcflow method consistently yields the optimal solution. This optimal solution serves as a baseline for our comparative analysis. In addition, we compute the LP bound of the three ILP models. For the set-partitioning formulation, we employ column generation, while for the arcflow and assignment-based ILPs, we directly solve the associated LP relaxation with CPLEX.

Out of 840 instances, column generation successfully terminated for 814. All 814 instances satisfy the IRUP condition and, consequently, the MIRUP condition. However, only 41 instances exhibit the integer property. From these empirical observations, we conjecture that MIRUP might be valid for BPTP. For only 11 instances, the set-partitioning formulation provides a stronger lower bound than the arcflow formulation. These instances share a common structure characterized by numerous small items and a large bin capacity, with optimal solutions containing no more than two bins. This outcome can be attributed to the fact that the set-partitioning formulation does not permit improper patterns, unlike the arcflow formulation. In all other cases, the two LP bounds were the same, underlining the theoretical observations made in Section 3.

The results obtained from the assignment-based ILP are consistent with Propositions 4 and 5. In other terms, we observed the weakest LP bounds for that type of model, and the value of this bound was always identical to the closed-form term derived in Section 3. Furthermore, Lemma 3 asserts that, unlike classic BPP/CSP, the lower bound of the assignment-based approach can be arbitrarily weak. Our experimental data also support this statement. By way of example, among the 814 instances tested, we identified one instance with an optimal value of 47, for which the assignment-based model provided a lower bound of 6.04. This lower bound is 7.7 times worse than the integer optimal value, that is, significantly weaker than the worst-case performance ratio for classic BPP. Lastly, we highlight that reducing the symmetric solution space in Model 1 (i.e., the introduction of the set $\Delta$) did not affect the lower bound of the model for the considered subset of instances.

# 5 Conclusion and Perspectives

In this article, a new variant of one-dimensional bin packing, originating from a joint production/maintenance scheduling problem, is introduced. It involves an original threshold constraint that indirectly constrains the order of items within the bins. From a theoretical point of view,
we first study the structural properties of this new problem, leading to a tractable mathematical description of the pattern set, and establish an important result showing that common heuristics (known from classic bin packing) possess an unbounded worst-case when applied to BPTP. Based on these observations, three ILP formulations – an assignment model, a pattern-based approach, and an arcflow formulation – are introduced, and the quality of the associated LP bounds is studied. Most notably, we derive a closed-form term for the LP bound of the assignment model and demonstrate that its approximation factor is arbitrarily bad.

In addition, all three methods were assessed and compared numerically on a wide variety of randomly generated instances. The experiments allowed to establish that arcflow is consistently more efficient than both the other methods. In that regard, the column generation based method is the main candidate for future improvements. Techniques such as stabilization or cutting planes could be implemented to accelerate the convergence of CG, while RMH could be replaced by more sophisticated branching techniques such as branch-and-price for the exact solution or diving heuristics. An alternative column generation procedure based on the arcflow formulation, in which the patterns would be paths in the arcflow network, could also be investigated and compared to the set-partitioning based column generation procedure. From a theoretical point of view, a deeper analysis of the MIRUP conjecture and the construction of heuristics with bounded worst-case seem to be very interesting. In addition, an extension to the *d*-dimensional vector bin packing with thresholds, occurring in the context of maintenance scheduling on multi-component production systems, is a further perspective to investigate. At first glance, the arcflow formulation can easily be extended to that scenario, but at the expense of a significant increase in the size of the underlying network.

## Declarations of interest

Declarations of interest: none

## Acknowledgements

This work has been partially supported by the LabEx PERSYVAL-Lab (ANR-11-LABX-0025-01) funded by the French program Investissement d'avenir.

## References

[1] J. Balogh, L. Epstein, and A. Levin. More on ordered open end bin packing. Journal of Scheduling, 24(6):589-614, 2021.

[2] G. Belov and G. Scheithauer. A branch-and-cut-and-price algorithm for one-dimensional stock cutting and two-dimensional two-stage cutting. European Journal of Operational Research, 171(1):85-106, 2006.

[3] F. Brandão and J. P. Pedroso. Bin packing and related problems: General arc-flow formulation with graph compression. Computers & Operations Research, 69:56-67, 2016.

[4] H. Cambazard and B. O'Sullivan. Propagating the Bin Packing Constraint Using Linear Programming. In D. Cohen, editor, Principles and Practice of Constraint Programming – CP 2010, Lecture Notes in Computer Science, pages 129-136. Springer, 2010.

[5] Y. Chen, L.-H. Su, Y.-C. Tsai, S. Huang, and F.-D. Chou. Scheduling jobs on a single machine with dirt cleaning consideration to minimize total completion time. IEEE Access, 7:22290-22300, 2019.

[6] T.-P. Chung, Z. Xue, T. Wu, and S. C. Shih. Minimising total completion time on single-machine scheduling with new integrated maintenance activities. International Journal of Production Research, 57(3):918-930, 2019.
[7] E. G. Coffman, J. Csirik, G. Galambos, S. Martello, and D. Vigo. Bin packing approximation algorithms: Survey and classification. In *Handbook of Combinatorial Optimization*, volume 1-5, pages 455–531. Springer, New York, NY, 2013.

[8] R. F. F. da Silva and R. C. S. Schouery. Solving cutting stock problems via an extended ryan-foster branching scheme and fast column generation, 2025.

[9] V. L. de Lima, C. Alves, F. Clautiaux, M. Iori, and J. M. V. de Carvalho. Arc flow formulations based on dynamic programming: Theoretical foundations and applications. *European Journal of Operational Research*, 296(1):3–21, 2022.

[10] V. L. de Lima, M. Iori, and F. K. Miyazawa. Exact solution of network flow models with strong relaxations. *Mathematical Programming*, 197(2):813–846, 2023.

[11] M. Delorme and M. Iori. Enhanced Pseudo-polynomial Formulations for Bin Packing and Cutting Stock Problems. *INFORMS Journal on Computing*, 32(1):101–119, 2020.

[12] M. Delorme, M. Iori, and S. Martello. Bin packing and cutting stock problems: Mathematical models and exact algorithms. *European Journal of Operational Research*, 255(1):1–20, 2016.

[13] H. Dyckhoff. A New Linear Programming Approach to the Cutting Stock Problem. *Operations Research*, 29(6):1092–1104, 1981.

[14] M. R. Garey, R. L. Graham, D. S. Johnson, and A. C.-C. Yao. Resource constrained scheduling as generalized bin packing. *Journal of Combinatorial Theory, Series A*, 21(3):257–298, 1976.

[15] M. R. Garey and D. S. Johnson. *Computers and Intractability: A Guide to the Theory of NP-Completeness*. W. H. Freeman & Co., USA, 1990.

[16] P. C. Gilmore and R. E. Gomory. A Linear Programming Approach to the Cutting-Stock Problem. *Operations Research*, 9(6):849–859, 1961.

[17] K. Hoffmann and U. Buscher. Valid inequalities for the arc flow formulation of the railway crew scheduling problem with attendance rates. *Computers & Industrial Engineering*, 127:1143–1152, 2019.

[18] D. S. Johnson. *Near-optimal bin packing algorithms*. PhD thesis, Massachusetts Institute of Technology, 1973.

[19] L. V. Kantorovich. Mathematical Methods of Organizing and Planning Production. *Management Science*, 6(4):366–422, 1960.

[20] E. L. Lawler and J. M. Moore. A functional equation and its application to resource allocation and sequencing problems. *Management Science*, 16(1):77–84, 1969.

[21] J. Martinovic, G. Scheithauer, and J. Valério de Carvalho. A comparative study of the arcflow model and the one-cut model for one-dimensional cutting stock problems. *European Journal of Operational Research*, 266(2):458–471, 2018.

[22] L. Penz, S. Dauzère-Pérès, and M. Nattaf. Minimizing the sum of completion times on a single machine with health index and flexible maintenance operations. *Computers & Operations Research*, 151:106092, 2023.

[23] B. Ramos, C. Alves, and J. Valério de Carvalho. An arc flow formulation to the multitrip production, inventory, distribution, and routing problem with time windows. *International Transactions in Operational Research*, 29(1):526–553, 2022.

[24] M. R. Rao. On the cutting stock problem. *Journal of the Computer Society of India*, 7:35–39, 1976.

[25] H. Rostami, J. Blue, A. Chen, and C. Yugma. Equipment deterioration modeling and cause diagnosis in semiconductor manufacturing. *International Journal of Intelligent Systems*, 36(6):2618–2638, 2021.

[26] R. Sadykov, F. Vanderbeck, A. Pessoa, I. Tahiri, and E. Uchoa. Primal Heuristics for Branch and Price: The Assets of Diving Methods. *INFORMS Journal on Computing*, 31(2):251–267, 2019.
[27] G. Scheithauer and J. Terno. The modified integer round-up property of the one-dimensional cutting stock problem. *European Journal of Operational Research*, 84(3):562–571, 1995.

[28] J. M. Valério de Carvalho. Exact solution of bin-packing problems using column generation and branch-and-bound. *Annals of Operations Research*, 86:629–659, 1999.

[29] J. M. Valério de Carvalho. LP models for bin packing and cutting stock problems. *European Journal of Operational Research*, 141(2):253–273, 2002.

[30] P. H. Vance. Branch-and-Price Algorithms for the One-Dimensional Cutting Stock Problem. *Computational Optimization and Applications*, 9(3):211–228, 1998.

[31] T. S. Wee and M. J. Magazine. Assembly line balancing as generalized bin packing. *Operations Research Letters*, 1(2):56–58, 1982.

[32] L. A. Wolsey. Valid inequalities, covering problems and discrete dynamic programs. In P. Hammer, E. Johnson, B. Korte, and G. Nemhauser, editors, *Studies in Integer Programming*, volume 1 of *Annals of Discrete Mathematics*, pages 527–538. Elsevier, 1977.

## A Additional figure for Example 4

<figure><img src="image_7.png" /><figcaption>**Figure 5:** Path decomposition and bin assignment for the instance from Example 1</figcaption></figure>

## B Complementary tables for Section 4.2
<table><caption>Table 9: Average model size per method for varying n.</caption><thead><tr><th colspan="2">Variables</th><th colspan="2">Constraints</th></tr><tr><th>n</th><th>Assign. ILP</th><th>Arcflow</th><th>Assign. ILP</th><th>Arcflow</th></tr></thead><tbody><tr><td>50</td><td>278.7</td><td>9462.4</td><td>323.0</td><td>448.3</td></tr><tr><td>100</td><td>1052.8</td><td>15310.2</td><td>1141.8</td><td>499.1</td></tr><tr><td>200</td><td>4109.7</td><td>20816.1</td><td>4288.2</td><td>582.9</td></tr><tr><td>500</td><td>25427.1</td><td>26359.3</td><td>25873.5</td><td>799.4</td></tr><tr><td>750</td><td>56759.8</td><td>28591.4</td><td>57429.4</td><td>955.4</td></tr><tr><td>1000</td><td>100952.7</td><td>29894.1</td><td>101845.2</td><td>1096.0</td></tr></tbody><tfoot><tr><td>Total</td><td>31582.5</td><td>21738.9</td><td>31970.7</td><td>730.2</td></tr></tfoot></table>

<table><caption>Table 10: Average model size per method for varying C.</caption><thead><tr><th colspan="2">Variables</th><th colspan="2">Constraints</th></tr><tr><th>C</th><th>Assign. ILP</th><th>Arcflow</th><th>Assign. ILP</th><th>Arcflow</th></tr></thead><tbody><tr><td>50</td><td>29928.1</td><td>819.3</td><td>30317.4</td><td>184.1</td></tr><tr><td>100</td><td>30715.8</td><td>2401.9</td><td>31104.1</td><td>340.6</td></tr><tr><td>200</td><td>32054.0</td><td>6645.2</td><td>32448.6</td><td>527.7</td></tr><tr><td>500</td><td>31919.8</td><td>22932.1</td><td>32306.0</td><td>868.8</td></tr><tr><td>750</td><td>32268.2</td><td>39719.9</td><td>32653.8</td><td>1112.8</td></tr><tr><td>1000</td><td>32624.8</td><td>57915.3</td><td>33010.0</td><td>1347.1</td></tr></tbody><tfoot><tr><td>Total</td><td>31582.5</td><td>21738.9</td><td>31970.7</td><td>730.2</td></tr></tfoot></table>