A study on exponential-size neighborhoods for the bin packing problem with conflicts

TL;DR

This paper introduces an advanced iterated local search method with large neighborhoods and diversification strategies for the bin packing problem with conflicts, demonstrating high-quality solutions and scalability.

Contribution

It presents novel large neighborhood structures, efficient evaluation procedures, and diversification techniques for the bin packing with conflicts, enhancing solution quality and scalability.

Findings

The method produces high-quality solutions on benchmark instances.

0-cost moves and set covering neighborhoods significantly improve search effectiveness.

The approach scales well with increasing problem size.

Abstract

We propose an iterated local search based on several classes of local and large neighborhoods for the bin packing problem with conflicts. This problem, which combines the characteristics of both bin packing and vertex coloring, arises in various application contexts such as logistics and transportation, timetabling, and resource allocation for cloud computing. We introduce $O(1)$ evaluation procedures for classical local-search moves, polynomial variants of ejection chains and assignment neighborhoods, an adaptive set covering-based neighborhood, and finally a controlled use of 0-cost moves to further diversify the search. The overall method produces solutions of good quality on the classical benchmark instances and scales very well with an increase of problem size. Extensive computational experiments are conducted to measure the respective contribution of each proposed neighborhood. In particular, the 0-cost moves and the large neighborhood based on set covering contribute very significantly to the search. Several research perspectives are open in relation to possible hybridizations with other state-of-the-art mathematical programming heuristics for this problem.

Tables

Table 1. Table 1: Comparison of th current state-of-the-art BPPC heuristics for each instance set and size n 𝑛 n .

Instance		HILS-Complete				ILS-Simple				PH			DH-LDS			BKS
Class	$n$	Gap	T(s)	Opt		Gap	T(s)	Opt		Gap	T(s)	Opt	Gap	T(s)	Opt	Bins	Opt
t	60	0.00	1.26	90.0		0.42	0.72	82.1		0.50	41.68	81	0.00	0.20	90	33.40	90
	120	0.01	5.17	89.6		0.63	2.34	66.3		0.66	44.44	66	0.03	1.36	89	66.11	90
	249	0.08	17.27	84.3		0.35	5.09	57.5		0.44	57.61	58	0.00	7.81	90	135.83	90
	501	0.01	23.88	88.9		0.26	13.19	47.9		0.23	65.48	60	0.00	50.43	90	275.69	90
u	120	0.03	3.67	88.8		0.06	1.22	87.1		0.11	24.58	85	0.00	0.76	90	70.38	90
	250	0.01	17.57	88.8		0.13	4.67	74.9		0.19	49.81	71	0.00	3.34	90	143.72	90
	500	0.01	37.56	86.8		0.37	16.47	45.8		0.18	66.44	65	0.00	18.74	90	286.03	90
	1000	0.01	90.07	85.7		0.53	65.07	33.6		0.20	105.34	56	0.00	116.71	90	571.88	90
ta	60	0.00	1.85	90.0		0.09	0.78	88.4		–	–	–	0.00	0.16	90	21.87	90
	120	0.05	23.34	88.1		1.18	3.41	49.4		–	–	–	0.05	1.46	88	41.39	90
	249	0.49	75.17	51.3		0.94	12.88	42.3		–	–	–	0.29	23.67	67	83.80	88
	501	0.40	128.30	44.5		0.57	32.87	44.1		–	–	–	0.15	271.33	48	167.61	70
ua	120	0.00	3.26	90.0		0.52	1.37	66.8		–	–	–	0.00	0.73	90	49.27	90
	250	0.02	28.70	87.1		0.55	6.83	58.9		–	–	–	0.02	6.21	87	100.64	89
	500	0.08	154.09	69.8		0.38	26.62	60.8		–	–	–	0.01	63.10	80	200.77	82
	1000	0.15	399.37	67.8		0.27	55.74	67.6		–	–	–	0.01	516.75	79	400.03	82
d	120	0.00	3.17	90.0		0.04	1.01	88.1		–	–	–	–	–	–	61.82	90
	250	0.00	11.57	89.9		0.08	2.77	80.6		–	–	–	–	–	–	127.93	90
	500	0.00	23.35	89.5		0.10	8.96	69.2		–	–	–	–	–	–	252.79	90
da	120	0.80	12.91	57.7		0.98	1.18	55.0		–	–	–	0.12	6.02	63	23.63	67
	250	1.14	25.50	40.5		1.32	4.88	38.0		–	–	–	0.04	59.60	49	44.66	50
	500	1.83	53.31	37.1		1.85	18.12	37.0		–	–	–	0.10	541.48	47	84.12	49

Table 2. Table 2: Comparison of th current state-of-the-art BPPC heuristics for each instance set and density ρ 𝜌 \rho .

Instance		HILS-Complete				ILS-Simple				PH			DH-LDS			BKS
Class	$\rho$	Gap	T(s)	Opt		Gap	T(s)	Opt		Gap	T(s)	Opt	Gap	T(s)	Opt	Bins	Opt
t	0	0.00	0.93	40.0		1.44	0.97	25.3		0.00	0.00	40	–	–	–	77.50	40
	10	0.00	1.05	40.0		1.26	1.12	26.1		2.08	134.03	2	0.00	11.93	40	77.55	40
	20	0.00	1.77	40.0		0.72	1.67	29.4		1.20	131.98	9	0.00	21.93	40	77.73	40
	30	0.21	19.99	32.8		0.45	7.88	27.2		0.71	115.05	18	0.06	24.15	39	78.33	40
	40	0.00	8.78	40.0		0.23	4.09	30.0		0.14	24.55	36	0.00	19.30	40	94.10	40
	50	0.00	11.35	40.0		0.29	5.00	27.9		0.00	12.73	40	0.00	18.69	40	117.93	40
	60	0.00	13.10	40.0		0.28	5.87	28.3		0.00	21.98	40	0.00	14.18	40	141.53	40
	70	0.00	15.85	40.0		0.27	6.89	23.6		0.00	16.10	40	0.00	10.31	40	164.33	40
	80	0.00	16.94	40.0		0.15	7.36	29.3		0.00	10.18	40	0.00	8.57	40	187.78	40
	90	0.00	18.22	40.0		0.08	8.17	32.0		0.00	4.15	40	0.00	5.51	40	210.58	40
u	0	0.00	0.56	40.0		0.00	0.26	40.0		0.31	75.95	18	–	–	–	188.53	40
	10	0.00	0.60	40.0		0.00	0.28	40.0		0.11	57.33	31	0.00	29.43	40	188.53	40
	20	0.00	0.60	40.0		0.00	0.27	40.0		0.10	57.78	32	0.00	31.63	40	188.53	40
	30	0.05	2.65	38.7		0.12	1.32	36.9		0.27	76.15	25	0.00	37.26	40	188.53	40
	40	0.09	58.14	34.1		0.39	20.12	24.9		0.79	91.95	14	0.00	57.83	40	194.03	40
	50	0.00	41.99	39.1		0.62	22.43	18.2		0.10	46.80	34	0.00	50.78	40	236.80	40
	60	0.00	44.72	39.9		0.36	27.02	21.6		0.06	61.93	34	0.00	34.85	40	284.80	40
	70	0.01	62.13	38.8		0.45	35.54	18.0		0.06	64.40	32	0.00	28.84	40	330.43	40
	80	0.00	68.64	39.5		0.35	42.03	18.0		0.03	59.70	36	0.00	23.78	40	376.73	40
	90	0.00	55.48	40.0		0.16	47.72	23.8		0.00	37.88	39	0.00	19.57	40	423.68	40
ta	10	0.00	1.52	40.0		0.74	1.42	29.2		–	–	–	0.00	51.24	40	77.73	40
	20	0.00	2.46	40.0		0.57	1.84	30.9		–	–	–	0.00	42.81	40	77.78	40
	30	0.33	16.99	30.9		0.76	6.73	23.4		–	–	–	0.11	72.27	36	77.80	40
	40	0.28	34.60	30.5		0.28	13.58	30.1		–	–	–	0.27	112.75	28	78.03	40
	50	0.42	74.85	21.0		0.42	25.38	21.0		–	–	–	0.42	120.84	21	78.03	40
	60	0.06	93.88	32.0		0.06	20.04	32.0		–	–	–	0.06	80.08	32	78.40	36
	70	0.00	107.48	32.0		0.10	16.37	30.9		–	–	–	0.00	86.29	32	78.50	32
	80	0.44	82.58	20.3		1.66	16.72	10.5		–	–	–	0.12	9.40	38	79.00	40
	90	0.59	100.13	27.2		1.66	10.27	16.2		–	–	–	0.13	91.72	26	82.75	30
ua	10	0.01	2.46	39.3		0.02	1.14	39.0		–	–	–	0.00	74.71	40	187.58	40
	20	0.01	0.81	39.7		0.01	0.48	39.6		–	–	–	0.00	74.44	40	188.35	40
	30	0.00	0.71	40.0		0.00	0.40	40.0		–	–	–	0.00	69.08	40	187.90	40
	40	0.00	3.69	39.7		0.06	1.79	38.7		–	–	–	0.00	69.72	40	186.78	40
	50	0.00	4.07	40.0		0.05	1.72	39.0		–	–	–	0.00	66.68	40	186.93	40
	60	0.02	13.52	38.6		0.18	8.54	34.1		–	–	–	0.00	81.51	40	187.03	40
	70	0.10	75.84	30.7		0.54	39.83	18.0		–	–	–	0.03	106.85	38	187.48	39
	80	0.18	306.37	21.8		1.08	86.13	2.0		–	–	–	0.00	171.81	34	187.75	34
	90	0.26	909.73	24.9		1.91	63.73	3.7		–	–	–	0.07	605.47	24	189.33	30
d	10	0.00	0.00	30.0		0.00	0.00	30.0		–	–	–	–	–	–	43.90	30
	20	0.00	5.22	30.0		0.06	2.52	28.7		–	–	–	–	–	–	58.67	30
	30	0.01	9.31	29.5		0.13	3.13	26.1		–	–	–	–	–	–	87.67	30
	40	0.00	10.70	30.0		0.17	3.57	25.8		–	–	–	–	–	–	114.30	30
	50	0.00	13.59	30.0		0.02	4.05	28.6		–	–	–	–	–	–	147.63	30
	60	0.00	16.14	30.0		0.11	4.82	24.0		–	–	–	–	–	–	174.07	30
	70	0.00	17.90	29.9		0.11	5.81	22.6		–	–	–	–	–	–	204.40	30
	80	0.00	19.58	30.0		0.03	6.50	26.8		–	–	–	–	–	–	235.07	30
	90	0.00	21.82	30.0		0.04	7.83	25.3		–	–	–	–	–	–	261.93	30
da	10	0.00	0.00	30.0		0.00	0.00	30.0		–	–	–	0.04	92.49	29	43.90	30
	20	0.00	0.00	30.0		0.00	0.00	30.0		–	–	–	0.09	84.50	29	43.80	30
	30	0.38	3.29	28.0		0.38	0.22	28.0		–	–	–	0.00	36.71	30	44.17	30
	40	0.55	18.54	23.6		0.69	4.09	22.4		–	–	–	0.00	79.81	28	44.03	28
	50	1.97	67.05	4.6		2.40	10.05	1.6		–	–	–	0.00	81.00	24	44.43	24
	60	2.30	74.63	0.0		2.65	12.14	0.0		–	–	–	0.04	384.41	2	45.50	3
	70	3.12	47.54	2.1		3.27	15.25	2.0		–	–	–	0.04	543.40	3	49.77	3
	80	2.43	35.15	7.0		2.52	16.01	6.0		–	–	–	0.21	348.24	7	60.30	8
	90	0.55	28.96	10.0		0.57	14.78	10.0		–	–	–	0.34	170.71	7	81.33	10

Table 3. Table 3: Impact of some variations of the shaking and termination parameters.

	$N_{\text{shak}}$	$N_{\text{ls}}$	$S_{\text{shak}}$	(t)		(u)		(ta)		(ua)		(d)		(da)
				Gap	Time	Gap	Time	Gap	Time	Gap	Time	Gap	Time	Gap	Time
Standard	50	100	3	0.02	9.97	0.01	34.44	0.23	52.22	0.05	190.40	0.00	11.78	1.16	31.45
A. More Frequent	100	50	3	0.02	13.76	0.01	42.16	0.22	80.40	0.04	250.58	0.00	13.81	1.16	57.97
B. Less Frequent	25	200	3	0.04	10.96	0.02	32.58	0.26	46.08	0.10	73.39	0.00	13.05	1.40	22.70
C. Longer Run	100	100	3	0.02	21.03	0.00	67.83	0.21	97.67	0.03	382.53	0.00	25.35	1.09	59.95
D. Shorter Run	25	100	3	0.04	5.65	0.02	20.18	0.28	34.69	0.10	50.05	0.01	6.45	1.47	15.60
E. Stronger Shaking	50	100	5	0.03	10.93	0.01	37.71	0.24	56.82	0.09	111.00	0.00	12.87	1.17	31.18
F. Weaker Shaking 1	50	100	1	0.02	10.56	0.01	31.57	0.18	76.25	0.06	166.64	0.00	12.47	1.42	25.78
G. Weaker Shaking 2	50	100	2	0.03	11.03	0.01	32.81	0.24	62.01	0.07	132.53	0.00	12.73	1.29	30.31

Table 4. Table 4: Sensitivity analysis when deactivating some neighborhoods.

	(t)		(u)		(ta)		(ua)		(d)		(da)
	Gap	Time	Gap	Time	Gap	Time	Gap	Time	Gap	Time	Gap	Time
Standard	0.02	9.97	0.01	34.44	0.23	52.22	0.05	190.40	0.00	11.78	1.16	31.45
H. No Assign	0.03	5.21	0.01	24.38	0.27	43.48	0.10	61.89	0.00	4.73	1.29	28.37
I. No Ejection Chains	0.02	9.99	0.01	30.92	0.24	56.67	0.07	108.31	0.00	11.58	1.23	28.87
J. No Grenade	0.02	9.84	0.02	34.99	0.25	46.49	0.07	70.97	0.01	11.58	1.21	29.51
K. No Set Covering	0.37	9.28	0.03	33.16	0.54	33.76	0.30	46.34	0.00	11.88	1.40	13.75
L. No Large Neighborhoods	0.49	4.50	0.24	18.66	0.72	11.35	0.43	18.98	0.07	3.97	1.36	7.39
M. No Large Neighborhoods ++	0.39	50.21	0.15	65.62	0.65	27.61	0.39	47.52	0.05	15.56	1.22	16.33
N. No 0-cost Moves	0.10	9.52	0.02	32.63	0.31	57.97	0.07	172.05	0.00	10.15	3.16	17.06

Table 5. Table 5: Sensitivity analysis on the parameters of the set covering neighborhood.

	Adapt		$S_{\text{pool}}$	$T_{\text{Limit}}$	(t)		(u)		(ta)		(ua)		(d)		(da)
					Gap	Time	Gap	Time	Gap	Time	Gap	Time	Gap	Time	Gap	Time
Standard	Yes		1500	20	0.02	9.97	0.01	34.44	0.23	52.22	0.05	190.40	0.00	11.78	1.16	31.45
O. Static Pool –	No		1000	10	0.05	9.80	0.02	34.14	0.33	41.57	0.16	44.97	0.00	11.84	1.30	16.58
P. Static	No		1500	20	0.04	9.85	0.02	35.18	0.26	54.75	0.14	48.92	0.00	11.80	1.20	24.84
Q. Static Pool +	No		2000	60	0.03	11.42	0.02	34.29	0.20	95.66	0.11	71.90	0.00	11.93	1.28	58.21
R. Set Partitioning	Yes		1500	20	0.02	10.19	0.01	33.85	0.22	52.82	0.07	122.83	0.00	12.13	1.20	27.71

Equations

$$\min$$

$$\sum _{i=1}^n{w}_i{x}_{ik}⩽Q{y}_k$$

$$\sum _{k=1}^n{x}_{ik}=1$$

$$x_{ik}+x_{jk}⩽1$$

$$y_k\in \{0,1\}$$

$$x_{ik}\in \{0,1\}$$

$$\Phi (B)={\omega }^{\text{textsc}c}C(B)+{\omega }^{\text{textsc}w}W(B)$$

$$\Delta =\text{textsc}Conf[i][B{}^{\prime }]+\text{textsc}Conf[j][B]-\text{textsc}Conf[i][B]-\text{textsc}Conf[j][B{}^{\prime }]-2\times \text{textsc}isConflict(i,j)$$

$$\Delta$$

$$+2\times (N_{\text{textsc}intra}(S)+N_{\text{textsc}intra}(S{}^{\prime })-N_{\text{textsc}inter}(S,S{}^{\prime }))$$

$$c_{ij}=\{\begin{array}{ll}\text{the cost difference of bin B(j) when replacing item j by item i}& \text{mbox}ifi\in \mathcal{V}\text{ and }j\in \mathcal{V},\\ \text{the cost difference of bin B(j) when removing item j}& \text{mbox}ifi\in {\mathcal{V}}^{\text{textsc}zero}\cup v^{\text{textsc}source}\text{ and }j\in \mathcal{V},\\ \text{the cost difference of bin B(j) when inserting item i}& \text{mbox}ifi\in \mathcal{V}\text{ and }j\in {\mathcal{V}}^{\text{textsc}zero}.\end{array}$$

$$\min$$

$$\sum _{j\in R_i}{u}_j⩾1$$

$$\sum _{j\in R}{u}_j⩽K$$

$$u_j\in \{0,1\}$$

Full text

A study on exponential-size neighborhoods for the bin packing problem with conflicts

Renatha Capua, Yuri Frota, Luiz Satoru Ochi, Thibaut Vidal

A study on exponential-size neighborhoods for

the bin packing problem with conflicts

Renatha Capua

Instituto de Computação – Universidade Federal Fluminense

Rua Passo da Pátria, 156 - São Domingos, Niterói - RJ, 24210-240, Brazil

rcapua@ic.uff.br

Yuri Frota

Instituto de Computação – Universidade Federal Fluminense

Rua Passo da Pátria, 156 - São Domingos, Niterói - RJ, 24210-240, Brazil

yuri@ic.uff.br

Luiz Satoru Ochi

Instituto de Computação – Universidade Federal Fluminense

Rua Passo da Pátria, 156 - São Domingos, Niterói - RJ, 24210-240, Brazil

satoru@ic.uff.br

Thibaut Vidal

Departamento de Informática, Pontifícia Universidade Católica do Rio de Janeiro (PUC-Rio)

Rua Marquês de São Vicente, 225 - Gávea, Rio de Janeiro - RJ, 22451-900, Brazil

vidalt@inf.puc-rio.br

Technical Report – Universidade Federal Fluminense – May 2017

Abstract. We propose an iterated local search based on several classes of local and large neighborhoods for the bin packing problem with conflicts. This problem, which combines the characteristics of both bin packing and vertex coloring, arises in various application contexts such as logistics and transportation, timetabling, and resource allocation for cloud computing. We introduce ${\mathcal{O}}(1)$ evaluation procedures for classical local-search moves, polynomial variants of ejection chains and assignment neighborhoods, an adaptive set covering-based neighborhood, and finally a controlled use of 0-cost moves to further diversify the search. The overall method produces solutions of good quality on the classical benchmark instances and scales very well with an increase of problem size. Extensive computational experiments are conducted to measure the respective contribution of each proposed neighborhood. In particular, the 0-cost moves and the large neighborhood based on set covering contribute very significantly to the search. Several research perspectives are open in relation to possible hybridizations with other state-of-the-art mathematical programming heuristics for this problem.

Keywords. Metaheuristics; Bin packing with conflicts; Large neighborhood search; Ejection chains; Assignment; Set covering

1 Introduction

The Bin Packing Problem with Conflicts (BPPC) is a difficult combinatorial optimization problem which brings together the features of bin packing and vertex coloring problems. We are given a set of bins of identical capacity $Q$, a set of items ${\mathcal{I}}=\{1,\dots,n\}$ with weights $w_{1},\dots,w_{n}$ and a conflict graph $G=(V,E)$. Any edge $(i,j)\in E$ represents a conflict between items $i$ and $j$. The objective of the BPPC is to assign each item to a bin, in such a way that the number of bins in use is minimized and there are no two items in conflict in the same bin. This problem can be formulated as follows (Gendreau et al. 2004):

[TABLE]

In this model, each binary variable $x_{ik}$ is set to $1$ if and only if item $i$ is assigned to bin $k$, and each binary variable $y_{k}$ is set to $1$ if and only if the bin $k$ has been used. Inequality (2) enforces the capacity constraints, Equation (3) states that each item should be assigned to exactly one bin, and Inequality (4) models the conflict constraints.

The BPPC is linked to several important applications. In the field of transportation and logistics, restrictions and compatibility constraints are imposed on several classes of products, including food, medicine or hazardous materials. Combinations of flammable, explosive, or toxic products should never be carried jointly in the same load (Hamdi-Dhaoui et al. 2014, Minh et al. 2013), leading to difficult item-to-vehicle assignment decisions during operational planning. The problem also arises in examination timetabling, where we search for a schedule that respects room capacities and avoids conflicts (Laporte and Desroches 1984). No student, for example, should take more than one exam at the same time. Finally, another application comes from the field of parallel computing, where we aim to assign a set of tasks to a minimum number of processors, subject to restrictions on the collocation of some conflicting tasks in some machines (Jansen 1999, Masson et al. 2013).

From a methodological standpoint, the problem is also of high interest. It includes two types of constraints: capacity restrictions and conflicts. Each constraint class, alone, leads to a NP-hard subproblem: bin packing or vertex coloring. It is also noteworthy that the research on the BPPC has been in majority concentrated on approximation algorithms (Epstein and Levin 2008) and mathematical-programming techniques. Mathematical programming has been successful on the BPPC up to this date, since decomposition methods such as branch-and-price (Fernandes Muritiba et al. 2010, Elhedhli et al. 2011, Sadykov and Vanderbeck 2013) are very efficient for problem instances with few items per bin (a common characteristic of most of the test sets available in the literature). Even in these conditions, the current mathematical programming algorithms remain applicable only to medium-scale problems, and many instances with one thousand items cannot be exactly solved to this date. This is unfortunately insufficient for several applications of BPPC, for cloud computing and large-scale task allocation problems (Masson et al. 2013), which can involve tens or hundreds of thousands of items.

The objective of this article is to progress on the understanding and development of metaheuristics for this class of problems. A significant part of recent research on metaheuristics, in the last five years, has been focused on producing novel search concepts, often from natural analogies. However, as in Sörensen (2015), we believe that part of this research is misled. Beyond natural phenomenons, the success of many metaheuristics come from a careful definition of their neighborhoods, efficient techniques to search them, and an exploration strategy which finds a sensible balance between intensification and diversification (Blum and Roli 2003, Vidal et al. 2013). Metaheuristics are also better understood in terms of their basic operations: neighborhoods, local search, restarts and perturbations, recombinations, short and long-term memories and diversification procedures, among others. Following this general logic, we develop a simple metaheuristic framework, an iterated local search, but enhance this search procedure with a variety of local and large (exponential-sized) neighborhoods, with the aim of carefully investigating their relative contribution to the search. Overall, this leads to a simple and efficient algorithm, plus a more complete version with additional classes of neighborhoods, which both return results of good quality on classic instances for the considered problem. We also observe that set covering-based neighborhoods had the biggest impact on the search performance, as well as the proposed diversification strategies based on 0-cost moves.

Overall, the main contributions of this paper are:

A simple and efficient ILS approach for the BPPC; 2. 2.

New ejection chain, assignment, and grenade neighborhoods tailored for the problem; 3. 3.

Efficient amortized ${\mathcal{O}}(1)$ move evaluations for standard enumerative neighborhoods; 4. 4.

A controlled use of 0-cost moves to better explore the search place; 5. 5.

An investigation of the relative contribution of all proposed neighborhoods.

2 Related literature

Although conflict constraints have a long history in the domain of timetabling (Laporte and Desroches 1984), the first studies dedicated to the bin packing problem with conflicts appeared in Jansen and Öhring (1997) and Jansen (1999). Early works were focused on approximation algorithms, a research line which was later pursued in Epstein et al. (2008) and Epstein and Levin (2008). Unfortunately, no bounded polynomial approximation scheme can be obtained for the BPPC with general conflict graphs (Jansen 1999). As such, most approximation results concern specific graph structures, such as perfect graphs, interval graphs, and bipartite graphs, with approximation ratios of $2.5$, $2.33333$, and $1.75$, respectively (Epstein and Levin 2008).

More recently, Gendreau et al. (2004) investigated some heuristics (without performance guarantee) for the problem. The authors proposed six different constructive heuristics, as well as two lower bounds, and a first set of benchmark instances based on the data set of Falkenauer (1996) for the bin packing problem. Unfortunately, these instances are not available anymore to this date, and thus were re-generated later on (Fernandes Muritiba et al. 2010, Elhedhli et al. 2011).

Subsequently, a significant research effort has been dedicated to mathematical programming algorithms for the problem. Khanafer et al. (2010) proposed new lower bounds based on the concept of data-dependent dual feasible solutions. Furthermore, three recent articles proposed exact branch-and-price procedures (Fernandes Muritiba et al. 2010, Elhedhli et al. 2011, Sadykov and Vanderbeck 2013). These articles use different branching rules, initial columns, and sub-procedures for the pricing problem (a knapsack problem with conflicts). When the conflicts form an interval graph, the pricing problem can be solved efficiently via dynamic programming (Sadykov and Vanderbeck 2013), while branch-and-bound procedures can be used in more general cases. These exact methods can solve many benchmark instances, with up to 1000 items in some cases. Still, instances with a large number of items per bin (e.g., set (da) of Sadykov and Vanderbeck 2013) remain challenging. The computational effort can also vary significantly from one instance to another, and becomes impracticable for large data sets with general graphs.

Metaheuristics for the BPPC have received, to this date, less attention than mathematical-programming algorithms. An advanced population-based metaheuristic was proposed in Fernandes Muritiba et al. (2010), and used to generate good upper bounds and initial columns. The method is a complex combination of a tabu search based on the impasse class neighborhoods of Morgenstern (1996), with a genetic algorithm using the crossover operator of Galinier and Hao (1999). Sadykov and Vanderbeck (2013) also introduce a diving heuristic based on a controlled partial exploration of the branch-and-price search tree. This method is a good alternative between exact solution approaches and heuristics, as it combines the benefits of mathematical programming for problem instances with few items per bin, with a shorter CPU time. Overall, metaheuristics for the BPPC still deserve further investigation.

Finally, as the BPPC is a generalization of both the bin packing problem and vertex coloring, results about metaheuristic searches for these two subproblems can be a good source of inspiration. The literature on these two subproblems is extensive, and we refer to Malaguti and Toth (2010), Lewis et al. (2012), Quiroz-Castellanos et al. (2015), Delorme et al. (2016) and Lewis (2016) for some detailed reviews. Some of the current best metaheuristics for the bin packing problem include the grouping genetic algorithm of Falkenauer (1996), the extension of the Perturbation-MBS method by Fleszar and Charalambous (2011), and the grouping genetic algorithm with controlled gene transmission of Quiroz-Castellanos et al. (2015). Population-based heuristics tend to be very frequently used, due to their capability of exploring very diverse areas of the search space. For the vertex coloring problem, some of the state-of-the-art metaheuristics include the TabuCol approach of Hertz and de Werra (1987), the PartialCol algorithm of Blöchliger and Zufferey (2008) with fluctuation of the objective function, the hybrid evolutionary algorithms of Galinier and Hao (1999), Malaguti et al. (2008), and the ant colony optimization metaheuristic of Dowsland and Thompson (2008).

3 Methodology

We opted for a simple metaheuristic framework, which allows to carefully investigate the impact of several local and large neighborhoods. Our method is thus built on an Iterated Local Search (ILS) metaheuristic, with several neighborhood classes for solution improvement, and one single perturbation procedure. Its pseudo-code is presented in Algorithm 1.

An initial solution is first generated via a constructive heuristic (Line $1$). Subsequently, to decrease the number of used bins, the method iteratively removes one bin (Line $4$), randomly assigning its items to other bins, and attempts to solve the resulting conflicts or capacity excesses by means of iterations of local search with 0–cost moves, ejection chains, grenade moves, resolutions of set covering or assignment formulations, and perturbation procedures (Lines $7$–$17$). Using the terminology of the vertex coloring problem, and as surveyed in Lewis et al. (2012), this type of strategy exploits the search space of complete, but improper colorings. The algorithm terminates when no feasible solution is found after $N_{\textsc{shak}}$ perturbations, or if a trivial lower bound $K_{\textsc{lb}}$ on the number of bins (the sum of items weight divided by bin capacity) is reached. We now describe each component of the algorithm in further details.

3.1 Initial Solution

The initial solution is generated by means of the modified first fit decreasing heuristic (Gendreau et al. 2004). Items are enumerated in non-decreasing order of weight, and inserted in the first bin which has enough residual space and does not contain conflicting items. This leads to a feasible initial solution with $K$ bins.

3.2 Local Search

The local search procedure aims to improve any infeasible solution produced by the removal of one bin or by the shaking operator. To define an improvement, the cost of each solution ${\mathcal{S}}$ is evaluated as the sum of the cost of its bins $B\in\mathcal{S}$ given by Equation (7). In this equation, $C(B)$ represents the current number of conflicts in the bin, $W(B)$ represents the total excess weight, and $\omega^{\textsc{c}}$ and $\omega^{\textsc{w}}$ are the penalty factors associated to each unit conflict and weight violation.

[TABLE]

Neighborhoods and exploration. Three classic neighborhoods are used: Swap, Relocate and Swap2vs1. A Swap move exchanges two items from different bins. As suggested by its name, a Relocate move reassigns an item to a different bin. Finally, a Swap2vs1 move exchanges a pair of consecutive items from one bin with a single item from another bin. The neighborhoods are explored as detailed in Algorithm 2, enumerating the pairs of bins in random order to test the associated moves (Lines 4–5). For each bin pair, the best improving move is applied if such a move exists (Line 9), and the local search terminates when no improving move can be found.

This evaluation policy allows to save computational effort, since it is unnecessary to re-evaluate the moves between a pair of bins when the evaluation has been done in the past without success, and if the bins have not undergone any change. Moreover, at least one bin $B$ involved in a move must satisfy $\Phi(B)>0$ to hope for an improvement. This condition is included as an additional filter (Line 4), allowing to save a significant amount of CPU time. Finally, during the first loop of the local search, the method accepts 0-cost moves (Line 8). These moves do not improve the objective but contribute to diversify the search.

Efficient Move Evaluations. The computational efficiency of the move evaluations is of critical importance, as it is the main bottleneck of most recent local search-based heuristics. Any significant reduction in the CPU time needed to evaluate moves can be translated into a direct gain in terms of number of moves which can be evaluated during one run of the algorithm, and thus into enhanced solution quality. For the BPPC, move evaluations are not trivial since they require to compute the number of conflicts in the new solution. A straightforward implementation of Swap, for example, between bins $B$ and $B\,^{\prime}$ would require ${\mathcal{O}}(|B|+|B\,^{\prime}|)$ elementary operations per move evaluation, where $|B|$ is the current number of items in a bin $B$. This effort is due to the evaluation of conflicts between the items exchanged and those existing in the bins. This effort grows linearly with the average number of items per bin, and can be responsible for the majority of the CPU time for problem instances with a large number of items per bin.

To address this issue, we propose a pre-processing and move evaluation mechanism which supports amortized ${\mathcal{O}}(1)$ move evaluations. We rely on an array which gives ${\mathcal{O}}(1)$ access to the current number of conflicts $\textsc{Conf}[i][B]$ between an item $i$ and all items in bin $B$. Let $d_{i}$ be the degree of $i$ in the conflict graph. This array can be built in ${\mathcal{O}}(Kn+\sum_{i=1}^{n}d_{i})$ elementary operations for the initial solution. Subsequently, any solution modification can be assimilated to a sequence of item relocations. For each relocation of an item $i$ from a bin $B$ to a bin $B\,^{\prime}$, the Conf array is updated in ${\mathcal{O}}(d_{i})$ by means of Algorithm 3.

This data structure can now be used to perform efficient move evaluations. Consider a Swap move between an item $i$ from bin $B$, with item $j$ from bin $B\,^{\prime}$. The difference $\Delta$ in the number of conflicts after and before the move can be computed as follows:

[TABLE]

Remark the corrective term $2\times\textsc{isConflict}(i,j)$, which takes value $2$ if and only if there is a conflict between the items $i$ and $j$. This term comes from the fact that the array $\textsc{Conf}[i][B\,^{\prime}]$ contains the number of conflicts when moving $i$ in $B\,^{\prime}$, but before the removal of item $j$, and vice-versa for the array $\textsc{Conf}[j][B]$.

This approach can be generalized to evaluate in ${\mathcal{O}}(1)$ any move involving an exchange of bounded sets of items $S$ and $S\,^{\prime}$ between a pair of bins $B$ and $B\,^{\prime}$, hence allowing efficient evaluations for all the considered moves. An example of a move, which exchanges three items of each bin, is displayed in Figure 1. Let $N_{\textsc{inter}}(S,S\,^{\prime})$ be the number of conflicts between pairs of items in $S$ and $S\,^{\prime}$, let $N_{\textsc{intra}}(S)$ be the number of conflicts within the set $S$, and let $N_{\textsc{intra}}(S\,^{\prime})$ be the number of conflicts within the set $S\,^{\prime}$. The difference in the number of conflicts after the application of the move can be obtained from Equation (9).

[TABLE]

Capacity excesses can also be simply evaluated in ${\mathcal{O}}(1)$. Overall, this leads to move evaluations in ${\mathcal{O}}(1)$ at the price of a slightly higher computational effort when updating the solution. In all computational experiments, we observed that hundreds of moves are tested, in average, before even applying one of them. The complexity needed to initialize and update the data structures remained largely dominated by the effort dedicated to move evaluations, leading to the announced amortized ${\mathcal{O}}(1)$ evaluations.

3.3 Large Neighborhoods

Our hybrid ILS is built on simple perturbation and local search concepts, but it relies on four large neighborhoods to enhance the search performance. The first two neighborhoods – Assignment and Ejection Chains – exploit restrictions of the problem which have the benefit to be combinatorial but still polynomial. This leads to neighborhoods of exponential size which can be searched in polynomial time. The third neighborhood – Grenade – is based on the enumeration of a larger move. Finally, the fourth neighborhood – Set Covering – uses the fact that effective integer programming solvers can address the set covering formulation of the problem in the presence of a limited pool of columns (item combinations in a bin). These solvers are used to combine high-quality columns from local minima to produce better solutions.

Some of these families of neighborhoods are well known in the operations research literature, and we refer to Deineko and Woeginger (2000) and Ahuja et al. (2002) for detailed surveys and analyses. Not all combinatorial optimization problems allow for an efficient search of such exponential-sized neighborhoods. For the case of the BPPC, such a search is possible. We now describe these large neighborhoods, tailored for the problem, as well as a specific form of ejection chains, which also allows a polynomial exploration of an exponential subset of solutions, under the condition that ejections between bins should follow a pre-specified ordering.

3.3.1 Assignment Neighborhood

The assignment neighborhood selects some items for removal, and reinserts these items in the holes thus formed by solving an assignment problem in a weighted bipartite graph. Such a methodology was successfully used for traveling salesman and vehicle routing problems in Sarvanov and Doroshko (1981), Gutin (1999) and Toth and Tramontani (2008).

We adapt this approach to the BPPC as follows. In our context, only a few problematic items usually present conflicts or belong to a bin with excess weight during the course of the search. To improve the solution, better allocation decisions should be sought for some of these items with simultaneous movements of other items. Hence, the proposed algorithm randomly selects one problematic item, call it $k$, as well as $N_{\textsc{assign}}$ random items, at most one in each bin. Together, these items form a vertex set $V^{\prime}_{1}$. Let $V^{\prime}_{2}$ be a copy of this set of vertices, and define the edge set $E^{\prime}=\{(i,j)\text{ such that }\allowbreak i\in V^{\prime}_{1}\text{ and }j\in V^{\prime}_{2}\}$. The cost of each edge $(i,j)$ is defined as the cost of inserting item $i$ in the bin $B(j)$ of item $j$ (which has been removed beforehand). Finally, to break possible ties and favor a movement of $k$, an $\epsilon$-penalty is added to the cost of the edge linking $k$ with its copy. This defines a bipartite graph $G^{\prime}=(V^{\prime}_{1},V^{\prime}_{2},E^{\prime})$, illustrated in Figure 2.

Then, finding a best re-allocation of items to bins is done by solving a linear-sum assignment problem on this graph. The optimal solution is found in ${\mathcal{O}}(n^{3})$ by means of an efficient implementation of the Hungarian algorithm (Kuhn 1955). Any non-trivial optimal solution of this assignment model (different from the pairing of each node $i$ with its clone) corresponds to a move of the BPPC solution. This move is applied, leading to a new (better or equal) incumbent solution in the ILS.

3.3.2 Ejection Chains

Ejection chains allow to explore a different subset of neighbor solutions obtained via chained relocations of items (see, e.g., Thompson and Psaraftis 1993, Glover 1996 and Deineko and Woeginger 2000). A key difference with the previous assignment neighborhood is that the choice of items to be relocated is not a-priori fixed. On the other hand, to allow for a polynomial exploration procedure, we restrain item relocations so as to comply with a pre-defined bin ordering.

Our ejection chains algorithm works as follows. First, choose a random order $\Pi$ of the bins in the current solution. Then, define an auxiliary graph $G^{\prime\prime}=(V^{\prime\prime},A^{\prime\prime})$ with $V^{\prime\prime}={\mathcal{V}}\cup{\mathcal{V}}^{\textsc{zero}}\cup v^{\textsc{source}}$, illustrated in Figure 3. The set ${\mathcal{V}}$ contains one node for each item, while ${\mathcal{V}}^{\textsc{zero}}$ contains one node $v^{\textsc{zero}}_{k}$ for each bin $k$ in the current solution. Each node in ${\mathcal{V}}$ will model the possible replacement of one item by another, and each node in ${\mathcal{V}}^{\textsc{zero}}$ will model the possible insertion of an item and the end of an ejection chain. The node $v^{\textsc{source}}$ represents a source. We denote by $B(i)$ the bin that contains vertex $i$ in the current solution. The set $A^{\prime\prime}$ contains one arc $(i,j)$ for any pair of nodes $i\in{\mathcal{V}}\cup{\mathcal{V}}^{\textsc{zero}}$ and $j\in{\mathcal{V}}\cup{\mathcal{V}}^{\textsc{zero}}$ such that bin $B(i)$ precedes bin $B(j)$ in the order $\Pi$, and one arc $(v^{\textsc{source}},j)$ for any $j\in{\mathcal{V}}$. The costs of the arcs are defined as follows:

[TABLE]

A shortest path in this graph from the source to any node in ${\mathcal{V}}^{\textsc{zero}}$ is equivalent to a sequence of combined item insertions and removals. Note that, thanks to the use of the nodes in ${\mathcal{V}}^{\textsc{zero}}$, we can possibly obtain a collection of several disjoint ejection chains. The solution is finally updated in case of improvement.

3.3.3 Grenade Neighborhood

The Grenade neighborhood has been previously used for the vertex coloring problem in Avanthay et al. (2003), and consists in the enumeration of some larger moves. In our context, we consider in turn each item $k$ in the set of problematic items $\mathcal{P}$ (as defined in the previous section: with either a conflict or in a bin with excess capacity). For each possible bin $B\neq B(k)$, the method evaluates the possibility of relocating $k$ in $B$ and jointly relocating any conflicting item $i$ of $B$ to another bin by means of a best insertion criterion. The total cost of these combined relocations is evaluated and the best grenade move for $k$ is applied in case of improvement.

The overall scheme is illustrated in Figure 4, and the pseudo-code of the neighborhood exploration is detailed in Algorithm 4. Note that the search consists of a single loop over the items $k\in\mathcal{P}$, as a complete descent towards a local minimum of this neighborhood would be time consuming.

3.3.4 Adaptive Set Covering

Set covering (or set partitioning) formulations have been successfully used to solve a variety of combinatorial optimization problems to optimality and to generate larger neighborhoods (see, e.g., Monaci and Toth 2006, Muter et al. 2010 and Subramanian et al. 2013 for recent applications on a wide range of packing, coloring and routing problems). In our context, we investigate a set covering formulation in which the maximum number of bins is fixed, and one aims to minimize the sum of column costs associated to possible capacity of conflict constraints violations, and use it as a base for an additional large neighborhood in the ILS. To limit the computational effort, this procedure is used only a few times during the search, after each $N_{\textsc{sc}}$ consecutive iterations of the local search procedure without improvement.

The considered set covering formulation is expressed in Equations (10–13). It involves the set $R$ of all possible columns, where each column represents a combination of items, associated with a cost $c_{j}$. The cost $c_{j}$ is based on the same weighted penalties for conflicts and capacity excess as in Equation (7), and let $R_{i}\subseteq R$ be the subset of columns containing an item $i\in{\mathcal{V}}$. The goal of the formulation is to find a set of up to $K$ columns that covers all items with minimum cost. We highlight the fact that the positive-cost columns (with $c_{i}>0$) correspond to sets of items with excess capacity or conflicts, and that any solution with cost [math] would correspond to a feasible packing in $K$ bins.

[TABLE]

This formulation involves an exponential number of decision variables ($2^{n}-1$). To circumvent this issue, we restrict the model to a subset of columns $R^{\prime}\subset R$, obtained from local minima of the ILS. More specifically, before every shaking iteration, the columns of the current solution are added to the pool $R^{\prime}$, filtering possible duplicates (this can easily implemented by hashing), and controlling the size of the pool to remain below a limit $S_{\textsc{pool}}$. The new columns start to replace the oldest columns when this limit is attained. The formulation is solved by means of an integer programming solver subject to a CPU time limit $T_{\textsc{Limit}}$. The pool-size limit is adapted during the search, with the aim to balance the difficulty of the problem and thus obtain good solutions with a formulation that remains tractable. If the model was solved to the optimality during the allowed time, then $S_{\textsc{pool}}$ is increased by $15\%$. On the contrary, if the time limit was attained without an optimal solution, then $S_{\textsc{pool}}$ is reduced by $15\%$. The optimal solution provided by this model is taken as new incumbent solution.

3.4 Perturbation

Shaking is only operated after $N_{\textsc{ls}}$ consecutive iterations of the [math]-cost moves, local and large neighborhood search without improvement of the best solution. Indeed, the 0-cost moves themselves already allow significantly modifications to the solutions, which help to progress towards different areas of the search space.

During shaking, $N_{\textsc{shak}}$ items are randomly relocated, considering $50\%$ items with conflicts or capacity excess, and $50\%$ random items. The size of the shaking operator is a key parameter for the search. Its value has been defined as a result of our preliminary computational experiments, as described in the next section.

4 Computational Experiments

The goal of these experiments is to evaluate the performance of the proposed method on the BPPC, as well as to identify the relative contribution of each neighborhood and the impact of its key design choices to drive future research. The classical benchmark instances for the BPPC are divided into six sets: (t) and (u) from Fernandes Muritiba et al. (2010), and (ta), (ua), (d), (da) from Sadykov and Vanderbeck (2013), for a total of 2060 instances ranging from 60 to 1000 items.

The instances of the group (t – triplet) and (u – uniform) were first derived from the bin packing instances of Falkenauer (1996) by Gendreau et al. (2004). The generation procedure, however, contains some randomized choices and the original instances are no longer available. Thus, these instances were generated again in Fernandes Muritiba et al. (2010) and Elhedhli et al. (2011). In this paper, we rely on the files of Fernandes Muritiba et al. (2010), as this allows a comparison with the only existing population-based metaheuristic in the literature, as well as the current state-of-the-art mathematical programming – heuristic and exact – algorithms of Sadykov and Vanderbeck (2013). The set (t) includes four groups of instances with a number of items ranging from 60 to 501. One remarkable characteristic of this instance group is that, in the optimal bin packing solution, without conflict constrains, the bins are filled by exactly three items each. In a similar manner, the set (u) includes four classes of instances, with 120 to 1000 items. In both groups, the conflict graph is an interval graph, with a density ranging from 0% to 90%. Previous articles did not always report results for the instances with 0% density, leading to 360 instances per set (9 density levels $\times$ 4 size levels $\times$ 10 instances).

The instances of sets (ta) and (ua) have the same characteristics as the instances of sets (t) and (u), but with an arbitrary conflict graph. These instances were generated by Sadykov and Vanderbeck (2013) with edge densities from 10% to 90%, for a total of 360 instances in each set. Finally, the instances of sets (d) and (da) were generated with item weights uniformly distributed in [500, 2500], and a bin capacity of 10000. These instances include either 120, 250 or 500 items, with interval (d) or arbitrary (da) conflict graph. The objective of these instances was to include a higher average number of items per bin and obtain more difficult problems for column generation-based methods. This objective was, however, only partially attained, since the average number of items per bin in the best known solutions remain relatively small for all data sets: 2.18, 1.96, 2.91, 2.48, 2.76, and 5.83 for the classes (t), (u), (ta), (ua), (d), and (da), respectively. For this reason, these instance sets remain ideal for existing branch-and-price methods.

4.1 Performance Evaluations

The first experiment aims to measure the performance of the proposed ILS with large neighborhoods relatively to the performance of the best current metaheuristic, the population heuristic (PH) of Fernandes Muritiba et al. (2010), and the best mathematical-programming based heuristic, the diving approach with limited discrepancy search (DH-LDS) of Sadykov and Vanderbeck (2013). As the total number of instances for the BPPC is large, the results are usually presented in an aggregated form. For a precise comparison, the detailed results of PH and DH-DLS were kindly provided by the authors in a personal communication. The proposed algorithm was coded in C++, uses CPLEX 12.5.1, and was run on a single thread of a Xeon X5675, 3.07 GHz CPU. PH and DH-LDS were run on Pentium IV 3.0GHz and Xeon X5460 3.16GHz CPUs, respectively.

A preliminary calibration of the search parameters – size and frequency of the shaking operator, number of iterations of the method, size of the set covering subproblem – was conducted with the aim of producing good solutions in a CPU time which is comparable with previous research for medium-sized instances. This led to a standard set of parameters: $N_{\textsc{shak}}=50$, $N_{\textsc{ls}}=100$, $N_{\textsc{sc}}=25$, $S_{\textsc{pool}}=1500$, $T_{\textsc{Limit}}=20$ and $S_{\textsc{shak}}=3$. The impact of any deviation from these search parameters will be investigated in deeper details in Section 4.2.

The proposed method was then run 10 times for each instance, with different random seeds, in order to measure its average performance. The complete version of the algorithm will be called HILS-Complete. To examine the impact of the large neighborhood search on solution quality, we will also report the results of a simplified configuration, ILS-Simple, which corresponds to the same algorithm without the large neighborhoods. These experiments are presented in Tables 1–2. In Table 1, the results are aggregated for each instance class $\times$ number of items, while Table 2 aggregates the results per instance class $\times$ density. The following metrics are reported:

• Gap –

The average gap from the best known solutions (BKS) in the literature. For each instance, the gap is calculated as $100(z-z_{BKS})/z_{BKS}$, where $z$ is the value of the solution found by the method and $z_{BKS}$ is the value of the BKS;

• T(s) –

The average CPU time per instance until termination of the search;

• Opt –

The average number of optimal solutions found for this group of instances in one run.

For each group of instances, the best results are highlighted in boldface. The last two columns indicate the average number of bins in the best known solutions of the literature, and the number of known optimal solutions for the group.

Two factors are determining for the performance of an optimization approach: the quality of its solutions, and the scalability of the method, i.e., how fast the CPU time grows for larger problem instances. From these results, the following observations can be made.

In terms of solution quality, HILS-Complete provides solutions of higher quality than the current best metaheuristic (PH), and slightly lower quality than the column generation-based diving algorithm (DH-LDS). The magnitude of the differences depends on the sets of instances. For the sets (t) and (d), which are common to the three algorithms, HILS-Complete finds an average gap of $0.02\%$, compared to $0.29\%$ and $0.00\%$ for PH and DH-LDS, respectively. For the sets (ta) and (ua), the difference of solution quality is also small, and both HILS-Complete and DH-LDS return near-optimal solutions with average gaps below $0.15\%$. For the (da) instances, DH-LDS produces solutions of higher quality ($1.17\%$ gap difference in average) than HILS-Complete, albeit at the expense of a larger CPU time. HILS-Complete attained 1819 out of the 1917 known optimal solutions, and a total of 9 new best upper bounds were also found for some open (da) instances.

DH-DLS performs a truncated search in the branch-and-price tree, which loses its optimality certificate but often keeps finding near-optimal solutions when the problem size allows for it. Yet, as $n$ grows, its CPU grows very quickly as noted in Sadykov and Vanderbeck (2013): “a disadvantage of our primal heuristic is that the running time increases rapidly with the number of items”. To better evaluate the scalability of the methods, Figure 5 reports the CPU time of DH-DLS and HILS-Complete for each instance set and size. Moreover, we fitted the CPU time as a power law $f(n)=\alpha\times n^{\beta}$ of the number of items $n$ for each instance class (least-squares regression of an affine function on the log-log graph). The observed CPU time of HILS-Complete grows as ${\mathcal{O}}(n^{1.49})$. In contrast, the time of DH-DLS grows at a rate which is faster than cubic for some instance sets, with power laws of the form $n^{\beta}$ where $2.45\leq\beta\leq 3.62$.

We also analyzed the impact of some instance characteristics on the performance of HILS-Complete. Figure 6 displays the boxplots of the percentage gap for each combination of class $\times$ number of items. For the data sets (t), (u), and (d) with interval conflict graph, the optimal solution was attained by HILS-Complete on nearly all instances. As a consequence, the boxplots are represented as a single line at a gap of $0.00\%$ and any other gap value appears as an outlier. The other data sets, especially (ta) and (da) with arbitrary conflict graph, are more difficult, such that a gap of $0.00\%$ does not represent the average behavior of the method anymore. The difficulty of the instances tends to grow with their size, although this effect is counterbalanced by the way the “gap” metric is defined: on small instances with few bins in the optimal solution, a mistake of one bin directly translates into a large gap (e.g., $5\%$ if the BKS includes 20 bins) leading to some outliers for small problems.

Instance class and size are not the only factors which impact the performance of the methods. The density of the conflict graph is also important. This effect is illustrated in Figure 7, which displays the boxplots of the percentage gap as a function of density. Generally, the most difficult instances are neither those with low density (equivalent to a bin packing problem), nor those with high density (equivalent to a tightly-constrained coloring problem), but those which combine the effects of both packing and conflict constraints at a density level between $50\%$ and $80\%$. HILS-Complete appears to perform very well for problem instances with low density (up to $20\%$). On these instance, the BKS has been systematically reached and even sometimes improved, as highlighted by the two outliers located below $0.00\%$.

We also evaluated the impact of the large neighborhoods, comparing the results of ILS-Simple and HILS-Complete. We observe a significant difference of solution quality, confirmed by a pairwise Wilcoxon test (p-value $\leq 10^{-64}$). The large neighborhoods have a positive impact on solution quality, but consume some additional CPU time (49.85 seconds in average per instance with the large neighborhoods, and 12.53 seconds without). Note that we also conducted side experiments to investigate whether an increase in the number of iterations of ILS-Simple could help to profit from this additional CPU time to reach better solutions, but the quality improvements were only minor when doubling the number of iterations (see Section 4.2).

Finally, we analyzed the time consumption of each component of the method. The results of this analysis are reported in Figure 8, which displays the percentage CPU time of each main component of the approach. This metric has been aggregated per instance class $\times$ number of items. The printed article contains the figure in black and white, and we refer to the on-line article for a color figure. The search effort appears to be well distributed among the search components, i.e., no specific neighborhood consumes the majority of the time. Overall, the assignment neighborhood, local search and zero-cost moves are the largest time consumers. In contrast, the ejection chains and grenade neighborhoods use less than $6\%$ and $4\%$ of the CPU time, respectively. The computational burden of the 0-cost moves mostly relates to the fact that it constitutes the first loop of the local search, and thus involves many moves evaluations which are later known as non-improving, and thus pruned in subsequent local search iterations. In line with the previous discussions, we also observe that the neighborhood based on set covering requires a larger amount of CPU time for difficult instances with an arbitrary conflict graph, such as (ta) and (da).

Overall, from these experiments, HILS-Complete is a promising alternative method for the BPPC, producing high quality solutions with an overall gap of $0.22\%$ from the BKS, and improving upon the previous best metaheuristic for the problem in a more controlled CPU time than existing mathematical programming-based methods. Having two alternative methods with their own strengths (near optimality for DH-LDS, and good scalability for HILS) will be a good asset for future works, which can explore new hybridizations. This option was kept aside in this work in order to better concentrate on the neighborhood aspects. Since these neighborhoods are central in the design of the method, we now go a step further to identify the most promising search components and measure the impact of its key parameters. This is the subject of a detailed analysis in the next section.

4.2 Sensitivity Analysis – Search Parameters and Neighborhoods

This section analyses the sensibility of the method to a change of its parameters, and evaluates the contribution of its key components and neighborhoods. Starting from the standard HILS-Complete configuration, we varied one parameter/factor at a time (OFAT approach) and tested the resulting configurations. These experiments are summarized in Tables 3 to 5. Each line presents the performance of an alternative configuration of the method (indexed from A to R), in terms of average gap and time, on the six classes of instances. To maintain a reasonable computational effort, a single run was performed for each configuration on the 2060 instances.

Impact of the general search parameters. Our first investigation concerns the shaking operator and the termination criterion, as these parameters are recurrent in all metaheuristics based on iterated local search. In the first two alternative configurations, we modified the periodicity of the shaking operator, i.e., the number $N_{\textsc{ls}}$ of zero-cost moves and local search phases before shaking. This parameter is halved in configuration A, and doubled in configuration B. Since the product $N_{\textsc{shak}}\times N_{\textsc{ls}}$ governs the total number of local search descents before termination, we kept it constant by updating $N_{\textsc{shak}}$ accordingly. The configurations C and D aim to investigate the impact of longer and shorter runs, respectively, by doubling or halving the value $N_{\textsc{shak}}$ without impacting $N_{\textsc{ls}}$. Finally, the configurations E, F and G investigate alternative strengths for the shaking operator, which may involve to perform stronger shaking for a better exploration of the search space ($S_{\textsc{shak}}=5$), or weaker shaking for further intensification ($S_{\textsc{shak}}=1$ or $2$).

The changes of parameter setting for the configurations A, B, E, F and G are generally detrimental for the performance of the method. It appears that strengthening the shaking, or calling upon this operator more frequently has only a small impact on performance, but the inverse process would have a larger negative impact. The configuration D allows a longer run, and thus naturally leads to better solutions. This improvement of solution quality, however, remains minor in comparison to the additional CPU time.

Contribution of the large neighborhoods. The second experiment aims to evaluate the respective role of each large neighborhood as well as the 0-cost moves, which are critical for search diversification. We thus deactivated, in the configurations H to K, each one of the four large neighborhoods. In configuration L, we deactivated all four neighborhoods. Finally, no 0-cost moves are applied in configuration N.

All search components appear to contribute positively to the final solution quality. By order of magnitude, the largest loss occurs when deactivating the 0-cost moves, as they contribute very significantly to diversity the search. In second position, deactivating all large neighborhoods (as in ILS-Simple) also leads to a large loss of solution quality, but also to a gain of computational time. We tested longer runs (configuration M, via a fourfold increase of $N_{\textsc{shak}}$) without the large neighborhoods, and this was not sufficient to regain the solution quality of the complete method.

The large neighborhood based on set covering has the largest impact on the search, especially for the instance classes (t) and (ua). As noted in Figure 8, on these particular instances, the set covering problem is conveniently solved in small CPU time. This allows to use a larger pool of columns, thus increasing the chances of solution improvement. Finally, the last three neighborhoods, ejection chains, assignment, and grenade have a more moderate impact on solution quality. Their contributions are only observed on the difficult instances: (ta), (ua) and (da). Nevertheless, note that the removal of these three neighborhoods together is more largely detrimental (compare configurations K and L). The three neighborhoods provide different forms of solutions improvement and diversification, which are not critical when taken one by one, but much more significant when considered as a whole.

Impact of the set-covering parameters. The last parameters of the method concern the set covering neighborhood. Our preliminary calibration led to a time limit $T_{\textsc{Limit}}=20$ seconds for the solver, as well as an adaptive mechanism for the size of the pool, starting from a pool size of $S_{\textsc{pool}}=1500$ columns. We thus tested the impact of this adaptive mechanism in the configurations O to Q, by evaluating three static pool size values of 1000, 1500 and 2000 columns, along with a time of 10, 20 and 60 seconds, respectively, for the resolution of the set covering problem. Finally, the last configuration R considers a possible replacement of the set covering by a set partitioning formulation, using an equality in Equation (11).

In these experiments, we observe that the size of the pool is a very sensible parameter. A small value leads to simple set covering problems which are not likely to result in solution improvements, even when solved to optimality, while a larger value opens the way to more improvement opportunities, at the risk of not solving the mathematical model in the allowed time. For the instances of class (da), an intermediate pool size of 1500 appears to be desirable, while for other instances, e.g., (ta) and (ua), a larger pool size with 2000 columns leads to solutions of better quality at the cost of additional computational time. Clearly, the best parameter setting is dependent on the instance characteristics. This led us to use the adaptive mechanism, which performs better overall than each of the static configurations. Finally, using either the set covering or the set partitioning model did not lead to a significant difference of solution quality.

5 Conclusions

In this paper, we introduced an ILS based on several classes of local and large (exponential-sized) neighborhoods for the BPPC. We intentionally relied on simple metaheuristic concepts so as to draw the complete focus on the large neighborhoods and accurately measure their contribution to the search. We proposed ${\mathcal{O}}(1)$ move evaluation procedures for the local search, polynomial variants of ejection chains and assignment neighborhoods, as well as an adaptive search mechanism for the set covering-based neighborhood. We also introduced a controlled use of 0-cost moves to further diversify the search.

The resulting ILS produces solution of good quality on the classic BPPC benchmark instances, with an average gap to the BKS of $0.22\%$. This is a better solution quality than the previous best metaheuristic, the genetic algorithm of Fernandes Muritiba et al. (2010), and a slightly lower quality than the diving heuristic with limited discrepancy search (DH-LDS) of Sadykov and Vanderbeck (2013). Still, the current instances are an ideal testing ground for column-generation-based methods, e.g. DH-LDH, as they contain very few items per bin (2.87 in average). On the other hand, the proposed metaheuristic exploits different search concepts, and scales very well with an increase of problem size, with a measured CPU time in ${\mathcal{O}}(n^{1.49})$. In comparison, the CPU time of DH-LDH tends to increase by a factor ten for each twofold increase of problem size. Finally, we conducted extensive experiments to measure the contribution of each neighborhood as well as the impact of the key search parameters. As underlined by these experiments, the 0-cost moves are critical for the search performance, as well as the set covering-based neighborhood. The other large neighborhoods, i.e., ejection chains, assignment, and grenade, have a smaller individual impact but contribute significantly to the search when used together.

Many research perspectives are open on this class of problems. First, future hybridizations should allow to combine the strengths of both neighborhood-centered search and column-generation based-heuristics, with the aim of keeping both near-optimality and tractability. Future research should thus be done along the path of matheuristics, aiming to harness the strength of mathematical programming and neighborhood search. Second, populations of solutions and more advanced memory structures can be a strong asset to help diversifying the search and favoring the discovery of better solutions. Finally, as the current benchmark instances are now solved to near-optimality, new and larger test sets, in terms of overall number of items, and number of items per bin, should also be investigated.