Expert Systems With Applications: Said Aqil, Karam Allali
Expert Systems With Applications: Said Aqil, Karam Allali
Expert Systems With Applications: Said Aqil, Karam Allali
Local search metaheuristic for solving hybrid flow shop problem in slabs
and beams manufacturing
Said Aqil ⇑, Karam Allali
Laboratory Mathematics and Applications, Hassan II University of Casablanca, FST, PO Box 146, Mohammedia, Morocco
a r t i c l e i n f o a b s t r a c t
Article history: Solving optimization problems in the modern manufacturing industry is one of the most attractive areas
Received 5 June 2019 of operational research. In this paper, we propose an industrial optimization case of slabs and beams pro-
Revised 8 June 2020 duction problem. The suggested approach consists in modeling the process as a hybrid flow shop with
Accepted 1 July 2020
unrelated parallel machines under constraints of sequence dependent setup time. Jobs launched in pro-
Available online 26 July 2020
duction are of different characteristics and require a machine preparation time at each phase during their
production cycle. The goal is to minimize the total tardiness of all jobs when they are completed on the
Keywords:
last stage. The hybrid flow shop problem is considered as one of the NP-hard issues since it consists of at
Hybrid flow shop
Unrelated parallel machines
least two stages. To solve numerically this problem, we propose two methods, the iterative local search
Sequence dependent setup time (ILS) and the iterated greedy (IG) metaheuristics. In addition, we suggest two new improvements, the first
Total tardiness relates to the initial step while the second concerns the steps of the neighborhood exploration. In the ini-
Iterative local search tial phase, we consider an initial solution set generated from the priority rules whose scheduling is
Iterative greedy algorithm ensured by the Nawaz-Enscore-Ham (NEH) algorithm and the greedy randomized adaptive search proce-
dure (GRASP). In the second phase, we suggest a new version of exploration of the neighborhood. In fact,
from a single solution, we will generate a set of neighboring solutions and we will choose the best solu-
tion. In this new resolution approach, we develop a total of twelve algorithms based on neighborhood
exploration. We note that the proposed algorithms are hybrid metaheuristics and are among the most
powerful optimization algorithms in local search. A simulation study is conducted to verify the effective-
ness of the suggested algorithms.This simulation is performed on a set of instances by varying the num-
ber of jobs, the number of machines per stage and the number of stages. We find that the IG algorithm
based on NEH initialization heuristic gives good results in terms of quality and convergence time to the
best solution.
Ó 2020 Elsevier Ltd. All rights reserved.
https://doi.org/10.1016/j.eswa.2020.113716
0957-4174/Ó 2020 Elsevier Ltd. All rights reserved.
2 S. Aqil, K. Allali / Expert Systems with Applications 162 (2020) 113716
the diameters of the steel wires. The two machines are of different machines is done by workstation or stage containing machines that
performances; they are characterized by a set of rollers with differ- can do the same activity. In addition, the management of the trans-
ent diameters that can be adapted to different wire of metal. The port flow between these stations is ensured by trolleys that can cir-
machines supply is ensured by coil of roll-up wire. Subsequently culate easily to bring the jobs to their destinations. The handling
a bending operation is necessary before the steel wires pass by system consisting of lift trucks and suspension cranes provide
diameter reduction rollers according to the nature of the jobs. logistics in this workshop.
The setting of the rollers, the supply of the coils as well as the cut- The slabs and beams in progress are deposited on specific pal-
ting of the metallic wires is considered as a setup time. In the sec- lets according to their characteristics. Subsequently, lifting machi-
ond phase, three parallel machines M 12 ; M 22 and M 32 of different nes and suspension cranes allow the transfer slabs and beams
performance allow to assemble the steel structure by welding in between stations during processing on machines. They can also
the connection nodes and the addition of suspension structure. transfer finished products to the warehouse for possible delivery
The control of the three machines is different, according to the to customers. Fig. 2, shows the main logistic means used in the pro-
used technology in the servo-controls of the axes and the nature duction cycle. These tools allow to illustrate a real case study, we
of the positioning sensors of the structure. The welding of the sus- represent a small batch of six jobs with different sections and dif-
pension structure will be executed by a manipulator robot, with ferent lengths showing the diversification of manufactured prod-
three degrees of freedom traversing the entire structure. The ucts. We have retained a set of hypotheses to simplify the real
assembly and dis-assembly times, the alignment, as well as the modeling of the flow between the production cells of the work-
adjustment of the positioning sensors of the structure constitute shop. Our goal is to provide the company with a numerical simula-
the setup time of the machines and the jobs according to their tion model that solves the scheduling problem in order to assign
characteristics. In the third phase, three machines M 13 ; M 23 and jobs to the machines and establish the scheduling schedule that
M33 form the last stage of the production unit in the workshop. minimizes the total tardiness.
This set of machines is formed by three mixers of different powers In this study, the processing time of the logistics system as well
with a capacity of the variable enclosure and dependent on the as the adjustment time of the machines will be included in the
power of each mixer. Their main roles are to fill the concrete struc- setup time. Periods of preventive and corrective maintenance of
ture in order to complete the production of the slabs and beams. the machines will not be taken into account and their capacity of
The filling of the hoppers with concrete, water, the necessary the queue is assumed to be always available. We note that in this
ingredients as well as the removal and evacuation of the slabs type of industry we adopt great rigor and precision according to a
and beams at the end of the treatment cycle is considered as a set of international standard norms. The goal is to guarantee secu-
preparation of machines and jobs in this last phase. The beams rity in the planned constructions. The need is to show a great rigor
and slabs are designed with internal prestressing to withstand in scientific study based on mathematical models. This encourages
the different stresses of the materials strengths. They are manufac- business managers and researchers to collaborate to ensure a good
tured with different mechanical dimensional characteristics to bet- quality product that will meet the expectations of customers. We
ter meet the imposed specifications by the company customers. draw upon from this workshop model to study of the hybrid flow
The organization of the production processes in this workshop is shop scheduling (HFS) problem with unrelated parallel machines.
very important at the logistical level. Indeed, the arrangement of We propose a model in the general case of several stages under
S. Aqil, K. Allali / Expert Systems with Applications 162 (2020) 113716 3
the constraint of preparation depending on the sequence. The Ruiz, and Alfaro-Fernández (2017), Zohali, Naderi, Mohammadi,
industrial interest of the modeling, is to propose to the company and Roshanaei (2019) propose various implementations of the
managers a platform based on software allowing the simulation ILS algorithm to solve the HFS problem. In the same context, the
of the production process. This extension will make available to IG algorithm is also one of the most implemented algorithms in
managers a possible future decomposition of the process in several many recent works for HFS problems. We note that Öztop,
stages to better rationalize production resources and logistics. Tasgetiren, Eliiyi, and Pan (2019), Ying and Lin (2018) propose var-
In this paper, we propose a new approach to solve the hybrid ious versions of this algorithm implementation solving HFS and
flow shop scheduling problem with sequence dependent setup distributed HFS, respectively. Our methods consist in implement-
time (SDST/HFS). ing the ILS and IG metaheuristics in various versions with two dif-
We recall that this problem brings together two types of issues, ferent variants of initialization algorithms. The proposed
that of the flow shop and that of parallel machines. Our approach is metaheuristics are considered as hybrid optimization algorithms.
based on two basic concepts, the first concerns the initialization Indeed, the incorporation of neighborhood exploration procedures
phase, the second concerns the neighborhood exploration phase as well as the simulated annealing model in its algorithms has
in metaheuristics implementation. allowed us to consolidate our methods. The ILS and IG algorithms
The resolution of scheduling problems by heuristics is one of are among the most powerful algorithms in optimization for solv-
the most common issues for many research projects. In the case ing scheduling problems. The hybrid methods are more and more
of HFS, these heuristics are usually inspired from the flow shop adopted for the scheduling problem. We note that the good combi-
scheduling resolution methods, for example Fernandez-Viagas, nation and the design of the hybrid algorithms play a determining
Molina-Pariente, and Framinan (2018) develop a set of heuristics role in the achievement of the good solution. Recently the IG algo-
to minimize makspan for HFS with identical parallel machines. In rithm is proposed in Ruiz, Pan, and Naderi (2019), Ribas,
our case, we apply a set of heuristics in the initialization phase. Companys, and Tort-Martorell (2019) for the resolution of the
We start a solution with a priority rule and then apply two main scheduling problem distributed flow shop parallel and flow shop
enhancement algorithms, NEH and GRASP. These two algorithms with blocking respectively. As in Shao, Shao, and Pi (2020) the IG
are among the most promising algorithms for finding good solu- algorithm is implemented with a multiple neighborhood system
tions. Besides, Dios, Fernandez-Viagas, and Framinan (2018) pro- for solving distributed HFS problems. In the same context, the ILS
pose a resolution based on the NEH algorithm of the HFS algorithm has demonstrated a set of applications for solving the
problem with missing operations. Similarly González-Neira and optimization problem in scheduling. In Schulz, Neufeld, and
Montoya-Torres (2017) propose the GRASP algorithm for a multi- Buscher (2019) a hybrid version is proposed to solve the HFS prob-
objective optimization of the HFS problem. Several recent research lem in order to minimize energy consumption. Also in Newton,
studies have addressed the problem of HFS scheduling under dif- Riahi, Su, and Sattar (2019) an improved version of the ILS algo-
ferent constraints (Ahonen & de Alvarenga, 2017; Khare & rithm has been used to solve the flow shop problem with blocking.
Agrawal, 2019). Machine stops for adjusting or mounting tools The application concerns a particular problem of beams and slabs
for specific operations are part of the setup time constraints. These manufacturing, modeled in the form of a HFS. We distinguish three
times are considered unproductive in scheduling; this constraint is types of HFS problem, the HFS with identical parallel machines as
considered a priority for the production managers. Our study in Lei and Guo (2016), the HFS with uniform parallel machines
model is part of this category of problem involving the SDST con- studied in Dessouky, Dessouky, and Verma (1998) and finally the
straint, the presence of unrelated parallel machines is a feature HFS with unrelated parallel machines established in Meng,
of our problem. This model will be applied to the manufacturing Zhang, Shao, Ren, and Ren (2018), Yu, Semeraro, and Matta
case of slabs and beams in order to minimize the total tardiness. (2018). To the best of our knowledge, only few works address
Local search algorithms are used as an efficient method to solve the last case problem with unrelated parallel machines. Motivated
several scheduling problems. Indeed, Lei and Guo (2016), Pan, by the previous works, our contribution consists of tackling unre-
4 S. Aqil, K. Allali / Expert Systems with Applications 162 (2020) 113716
lated parallel machines HFS problem with the presence of setup HFS problems has been presented in Ruiz and Vázquez-Rodríguez
time constraint and with the goal total tardiness. Moreover, we (2010) up to the end of the last decade.
suggest twelve algorithms based on new versions of neighborhood A variety of scheduling problems are represented by the three
exploration in order to enhance the performance of seeking the fields ajbjc in the optimization scheduling literature. This notation
optimum. Our approach will be applied to slabs and beams has been described in Graham, Lawler, Lenstra, and Kan (1979) and
scheduling industrial case. reviewed in Vignier, Billaut, and Proust (1999), it helps to identify
Our paper is organized in the following way, in Section 2 we the nature and characteristics of the problem. Where the first field
present a review of the literature of hybrid flow shop scheduling indicates the nature of the problem, the second represents the con-
problem, in Section 3 we give a description of the chosen model. straints and the last one stands for the one or more criteria to be
Section 4 presents the algorithms used to solve our problem, in optimized. We focus in this review on different implemented
Section 5 we interpret the results of the conducted simulation. Sec- metaheuristics for the resolution of the HFS problem. These meta-
tion 6 will complete the work by emphasizing the strengths of our heuristics are widely used in the optimization domain of schedul-
resolution approaches and by giving some future perspectives of ing. They are applied to optimize several problems in different
this paper. configuration cases for HFS scheduling problem. We find the
genetic algorithm (GA) adapted for solving the HFS problem in
2. Literature review Engin, Ceran, and Yilmaz (2011), Wang and Liu (2013), Jun and
Park (2015), the simulated annealing (SA) for solving the HFS
The slab and beam production scheduling problem can be con- scheduling problem (Wang, Chou, & Wu, 2011; Mirsanei,
sidered as a typical hybrid flow shop problem. In this section, we Zandieh, Moayed, & Khabbazi, 2011) and their combination in
will give a brief review on manufacturing processes with hybrid Figielska (2009); the tabu search (TS) algorithm implemented with
flow shop. Indeed, the manufacturing industry is one of the mod- _
two versions (Wang & Tang, 2009; Bozejko, Pempera, & Smutnicki,
ern attracting research fields because of its high scientific value 2013). We detect also the presence of metaheuristics inspired by
in the area of optimization. It can model continuous flows such the biological phenomenon as the artificial immune system (AIS)
as petrochemicals, metallurgy and plastics or discontinuous flows algorithm (Engin & and Döyen, 2004; Zandieh, Ghomi, &
like logistics, assembly lines, electronics and mechanical industry, Husseini, 2006; Komaki, Teymourian, & Kayvanfar, 2016) and the
etc. Its strong presence in the industry pushes leaders by stream- adaptation of the immunoglobulin-based artificial immune system
lining resources to improve productivity. Indeed, this goal is algorithm (Chung & Liao, 2013; Chung, Sun, & Liao, 2017).
achieved by optimizing one or more criteria in flow processes. In the same situation, we also report the strong implementation
The goal is to meet the expectations of their customers and to be of metaheuristics inspired by the natural behavior of certain living
efficient in the face of their competitors in the industrial market. beings such as the ant colony optimization (ACO) algorithm
The resolution approach is usually based on mathematical models applied in Alaykỳran, Engin, and Döyen (2007), Solano-Charris,
and more specifically on computer simulation techniques. There Montoya-Torres, and Paternina-Arboleda (2011), the artificial bee
are several models representing cases of production management colony (ABC) algorithm (Pan, Wang, Mao, Zhao, & Zhang, 2013;
with various industrial characteristics (Pinedo, 2012). This new Li, Pan, & Duan, 2016) and the migratory bird optimization
link of cooperation between industrial and computer science (MBO) algorithm presented in Pan and Dong (2014), Zhang et al.
develops more and more applications of technological tools in (2017). In addition, we also find other metaheuristics for HFS
the scheduling field. In this approach, the developed models for scheduling like the particle swarm optimization (PSO) algorithm
scheduling problems constitute a relevant reference for the appli- in Tang and Wang (2010), Liao, Tjandradjaja, and Chung (2012).
cation of numerical simulation. In addition, industrial application We note that other algorithms are also present in solving HFS
cases are considered as a goal that enhances scheduling academic problems like the iterative greedy algorithm in Kahraman, Engin,
works. Within this context, the hybrid flow shop or flexible flow Kaya, and Öztürk (2010) and the variable neighborhood search
shop encourage researchers and manufacturers to propose models (VNS) algorithm (Eskandari & Hosseinzadeh, 2014; Jiang, Liu,
of resolutions, in order to optimize some criteria imposed by the Hao, & Qian, 2015; Li, Pan, & Wang, 2014). The constraints are
industrial process. We note that in the case of HFS scheduling diversified according to the industrial context, we distinguish the
problems, several models date back to the last two decades of most frequent in the HFS scheduling the unavailability of jobs or
the last century have been established. machines (Besbes, Teghem, & Loukil, 2010; Allaoui & Artiba,
The suggested approaches and models often study cases with at 2014). This type of constraint appears when periods maintenance
least two stages in the presence of several parallel machines. Often are scheduled or because of the unexpected shutdown of some
the goal is to optimize the makespan criterion under a variety of operating machines. Also, we find like encountered type of con-
constraints, we note that for the case of two stages, the problem straints, the limitation of inter-stage buffer as in Luo, Du, Huang,
is NP-hard in the strong sense when at least one stage contains Chen, and Li (2013), Li and Pan (2015) and the no-wait as in the
at least two machines. The proposed models are usually either industrial studied case in Jolai, Rabiee, and Asefi (2012), Rabiee,
exact methods such as mixed integer linear programming, meth- Rad, Mazinani, and Shafaei (2014). Indeed, we can meet in many
ods with branch and bound based algorithm, or dynamic program- types of industries the no-wait of the jobs because of the nature
ming approach. These approaches are limited for small size of the jobs, which can lose their quality at the end of cycle. This
problems that do not represent real industrial cases; they are also is the case, for example, in the food or pharmaceutical industries.
used to check the complexity of resolution algorithms. The second We note that the type of the used machines radically influences
approach is based on the metaheuristics which consist essentially the optimization of the studied criterion. There exists three kinds
on heuristics in the phase of initialization. Thereafter, a perturba- of HFS problem, the first with identical parallel machines as in
tion procedure is applied either in one or more steps to explore Lei and Guo (2016), the second with uniform parallel machines
the neighborhood of the current solution. These metaheuristics studied in Dessouky et al. (1998) and finally the third with unre-
find their efficiency to solve several problems whose size is rela- lated parallel machines established in Meng et al. (2018), Yu
tively important representing real industrial case. We find several et al. (2018).
works deal with problems of HFS scheduling mono or multi- We inform that the case of HFS scheduling problem in the
criteria. We note that a review of the literature of the various presence of unrelated machines represents the general case of
S. Aqil, K. Allali / Expert Systems with Applications 162 (2020) 113716 5
the studied problem, it will be the subject of our study with the parameters. For each equation, we highlight the necessary rela-
constraint of preparation of the machines. Each type of machine tionships involving the number of machines per stage, processing
can induce specific constraints to the studied scheduling problem. times, inter-jobs setup time per stage and the number of stages.
These constraints are often inspired by real industrial problem, The progress of steps make it possible to determine the assignment
highlighting several methods based essentially on different meta- of jobs to the machines, to sequence the jobs on the machines as
heuristics. In this context, we note that the change of job character- well as the calculation of the end date of each job at each stage.
istics or the change of production tools in machines often require a The objective function is calculated by the sum of the tardiness
setup time. This time can represent a time of adjustment or clean- of all jobs when they are completed.
ing of the machines, it is called the sequence-dependent setup time
it is a constraint that causes a change in the scheduling. The works
discussed in Gómez-Gasquet, Andrés, and Lario (2012), Ramezani,
Rabiee, and Jolai (2015), Pan, Gao, Li, and Gao (2017), 3.1. Description of the problem
Marichelvam, Azhagurajan, and Geetha (2018) have treated this
kind of constraints with a variety of criteria to be optimized. The SDST/HFS scheduling is a problem of planning and organi-
The presence of SDST constraints in workshops production zation of a set of jobs J ¼ fJ 1 ; J 2 ; . . . ; J n g started at the beginning
often affects the scheduling processes compared to the case with- of the production cycle through a set of K stations called stages.
out SDST. Indeed, the repetitive adjustment stops of the machines All jobs have the same production cycle on a set of
expecting the arrival of a new treated job push the production K ¼ f1; 2; . . . ; K g stages arranged in series. Each job will follow K
managers to optimize and rationalize the industrial processes. operations during this process, they will visit all stages in the same
These times are considered as a unproductive dead-time causing order. On each stage, a set MðkÞ ¼ M 1k ; M 2k ; . . . ; Mmk k are able to
a great loss of money for companies and may lose the trust of some process all the present jobs in its queue. We note that
customers in the market. The need to optimize scheduling criteria jM ðkÞ j ¼ mk represents the number of machine on the stage k, where
under SDST constraints becomes a primary goal to gain more cus- jXj represents the cardinal of set X. The machines are parallel
tomers. We are interested in the same problem structures, the unrelated, a job J j requires a processing time pijk on the machine
unrelated parallel machines, our model will be that of SDST/HFS Mik at the stage k, each machine requires a setup time to process
scheduling problem. The longevity of the company will pass a new job. To treat the job J j , we consider shjk , the setup time of a
through an effective management taking into account these con-
machine Mik in a stage k, when it completes the job J h . When the
straints, it will allow it to better rationalize its resources.
job J j is the first to run on this machine in the k stage, this time will
Our model is to minimize the total tardiness in the SDST/HFS
be noted sjjk . The constraint related to the delivery time imposed by
scheduling problem with such constraints in order to satisfy the
the customer will be expressed by dj , the latest date of delivery
customer’s request. For most companies today, this criterion has
beyond this date the job will be considered late. The completion
become a primordial objective comparing to other criteria. Indeed,
time of the job J j on the machine M ik of the stage k is C ijk , when
it allows to respect the constraints of the deadline imposed by the
customers upon receipt of his order from the customer relations the machine is not specified this date will be noted C jk . A sequence
K p ¼ ðp1 ; p2 ; . . . ; pn Þ of jobs is started in the production process by
department. Our model will be noted HFSK ðRMÞk k¼1
Pn starting the first stage and visiting the other stages in the same
jdj ; SDSTj j¼1 T j respecting the notion adopted in scheduling, the order. The problem is to assign and schedule each job on its best
description of the flow in our the studied workshop is given by machine, which means the one that will complete it as soon as pos-
Fig. 3. In this configuration, all jobs will be processed in the same sible. In HFS problem, during the scheduling process and at the end
order, and each job will be processed by the machine that will of each stage, the FIFO rule is applied to schedule the current
complete it as soon as possible during the scheduling process. sequence to the next stage. To simplify this study, we retain the
The goal is to minimize the total tardiness and to complete the usual scheduling hypotheses in the HFS scheduling problem. All
whole job batch as soon as possible by respecting the availability jobs and all machines are available at the initial time, each job will
constraints of the machines in real time scheduling. be handled by one machine at each stage, interrupt running is not
allowed. Each machine processes only one job at a time, inter-stage
3. Problem statement inventory and in the queue of machines are infinite capacity. The
SDST/HSF scheduling problem will be modeled by a set of equa-
We present here the model of the manufacturing process of the tions to highlight the latter constraints and assumptions, this
SDST/HFS scheduling problem, described by set of equations and model is described as follows.
Each job pj in the current sequence is assigned and scheduled d6 ¼ 25. This constraint is often imposed under the customer’s
on a single machine i at stage k, if this machine has not processed requirement when establishing the production schedule. The
any job before, the end date is constraints of inter-job setup in both stages will be expressed by
n o the matrices shj1 ; shj2 and shj3 . In this notation, the job J h is
C ipj k ¼ max spj pj k ; C pj ðk1Þ þ pipj k ð1Þ assumed to be the predecessor of the job J j in the scheduled
sequence. In Fig. 4, we present the Gantt chart of the
During the scheduling process, if a machine i has just completed
sequence p ¼ ð5; 2; 1; 6; 4; 3Þ. According to the plot made on
a job ph on stage k, the completion time of the job pj is
the Gantt diagram we can calculate the tardiness of each job
n o
C ipj k ¼ max C iph k þ sph pj k ; C pj ðk1Þ þ pipj k ð2Þ by applying the formula given by the Eq. (3) which allows to
deduce also the total tardiness TT ¼ 20 units of time given by
The job pj is assigned to the machine i which completes it as the Eq. (4).
soon as possible, with i ¼ argmin16i6mk C ipj k . The tardiness of
4. Solving the SDST/HFS problem
the job pj will be calculated according to its end date on last stage
K by the following expression Our model for solving the problem dealing with SDST/HFS is
based on two main types of approaches. The first is based on
T pj ¼ max 0; C pj K dpj ð3Þ
heuristics, the second is based on metaheuristics. For the first,
The objective is the determination of the tardiness sum of the we propose two heuristics to establish initial solutions that are
all jobs during the scheduling in the workshop, called the total tar- based on four priority rules. Subsequently, we apply the two
diness such as: heuristics; the first is based on the NEH algorithm (Nawaz,
Enscore, & Ham, 1983), the second is based on the greedy random-
X
n
ized adaptive search procedure (Davoudpour & Ashrafi, 2009;
TT ¼ T pj ð4Þ
Abyaneh & Zandieh, 2012).
j¼1
We propose the algorithm flowchart (see Fig. 5) describing
K
In the scheduling problem HFSK ðRMÞk k¼1 jdj ; SDSTjTT, the goal the developed approach for the SDST/HFS problem resolution.
is to find an optimal permutation p that minimizes the TT of all Our approach is based on two main phases: the initialization
the permutations generated in the neighborhood of the search phase and the neighborhood exploration phase. In the first
space. phase a set of priority rules are applied, thereafter the obtained
solutions are improved via two heuristics NEH and GRASP. The
3.2. An illustration example second phase the exploration of the neighborhood system, in
which a set of exploration procedure is implemented and rein-
To illustrate the conducted study in this scheduling situation, forced by the simulated annealing model. For each algorithm, a
we present an example consisting of a small batch of six jobs stopping criterion is applied to limit the search for new solu-
ð Jj; j ¼ 1; . . . ; 6Þ to be scheduled on the three stages forming the tions. Our platform therefore constitutes a software comprising
production unit. The shape of the sections is given in Fig. 2, these twelve resolution algorithms based on the ILS and IG
beams admit different lengths. The first stage is composed of two algorithms.
machines, the second and third stages contain three machines,
on each stage they are unrelated parallel. The processing times 4.1. The suggested heuristics
and the setup time of all jobs in three stages are expressed by
In general, one of the most important phase of solving schedul-
the matrices pij1 ; pij2 ; pij3 ; shj1 ; shj2 and shj3 where
ing problems is the search for the initial solution to start meta-
2 3
" # 5 6 7 8 6 4 heuristics. The goal is to accelerate the convergence to the right
6 2 4 8 5 5 6 7 solution as quickly as possible. In our case study, for solving the
pij1 26 ¼ ; pij2 36 ¼ 4 6 8 4 6 10 85 K P
5 4 6 6 7 4 problem HFSK ðRMÞk k¼1 jdj ; SDSTj nj¼1 T j , we present a set of
8 7 5 5 12 3
heuristics highlighting the problem characteristics that we detail
2 3 later. We give here the four priority rules that are the subject of
3 3 2 3 3 1 the determination of the sequence to be implemented by our
6 7
2 3 62 2 2 1 2 27 heuristics.
8 7 7 6 4 5 6 7
6 7
6 7 64 2 4 2 4 37
pij3 36 ¼ 4 4 9 5 4 5 8 5; shj1 66 ¼6
6
7
7 We calculate the minimum short processing time of each job J j
62 3 3 2 2 27
7 8 5 6 6 7 6 7 on each stage by p jk ¼ min16i6mk pijk and its minimum setup
64 2 1 2 5 37
4 5
time by s
jk ¼ min16h6nðh!¼jÞ shjk which J h represents the assumed
2 2 2 4 2 1 successor in the scheduling sequence. We express a combina-
tion of the two preceding terms in the following expression in
2 3 2 3
3 3 3 3 2 1 3 1 3 2 2 4 order to obtain a ranking order
63 2 2 1 2 27 62 2 3 2 3 37
6 7 6 7
6 7 6 7
X K
X K
61 3 4 2 1 37 64 2 3 1 1 47 ð1Þ
shj2 66 ¼6
62
7; shj3 66 ¼6 7 dj ¼ k
pjk dj
þ ð1 kÞ
sjk dj
ð5Þ
6 3 3 2 3 27
7
62
6 3 2 2 2 37
7
k¼1
k¼1
6 7 6 7
42 2 3 2 3 35 42 1 3 3 3 35
2 2 3 4 2 1 2 4 2 1 4 2 We assume that jobs will be assigned to their best machines
and will be preceded by jobs with shorter setup times. In this
We also give the due dates of this batch constituted by six jobs expressions, we take into account the sum of the most total
launched in this workshop, such as the due dates for all job form- processing and the setup times by calculating the absolute
ing this sequence are: d1 ¼ 23; d2 ¼ 15; d3 ¼ 27; d4 ¼ 29; d5 ¼ 19; difference with the deadlines. With a weighting coefficient
S. Aqil, K. Allali / Expert Systems with Applications 162 (2020) 113716 7
k 2 ½0; 1 for each relative term to the processing time and the For a globally average order, in this expression we calculate the
setup time. average processing time of each job J j at each stage expressed
Pmk
For this rule, we compute the longest processing time of the job p
by pjk ¼ m i¼1 ijk
. We keep the same form established by calculat-
J j at each stage by pþ jk ¼ max16i6mk pijk and we consider the k
þ
ing the following expression
longest setup time by sjk ¼ max16h6nðh!¼jÞ shjk . Subsequently,
X K
X K
dj ¼ k
pjk dj
þ ð1 kÞ
sjk dj
ð7Þ
identical combination to the last form.
k¼1
k¼1
X K
X K
Here we take into account the short setup time in the establish-
ð2Þ
þ
þ
dj ¼ k
pjk dj
þ ð1 kÞ
sjk dj
ð6Þ ment order for each job.
k¼1
k¼1
X K
X K
TT pj 2 ½TT min ; TT min þ q ðTT max TT min Þ ð9Þ
ð4 Þ
þ
dj ¼ k
pjk dj
þ ð1 kÞ
sjk dj
ð8Þ
k¼1
k¼1
Algorithm 2 gives a description of the GRASP heuristic with all
of these features that we adopt in order to solve the SDST/HFS
problem. We can see that the two main heuristics are two con-
We note that we have determined these four expressions by
ðuÞ
structive algorithms where the right solution is built step by step
calculating the terms dj ; u ¼ 1; . . . ; 4. The purpose is to obtain by exploring the neighborhood of the good constructed solution.
an initial sequence for each established rule, we will classify the The diversification of the initial rules of priorities makes it possible
ðuÞ
jobs in descending order of dj . We note that we adopt a weighting to better judge the efficiency of each heuristic during their use in
coefficient k to concertize the effect of the processing times and the resolution phase by metaheuristics.
setup times on the established order of the obtained sequence.
Subsequently, we propose two constructive scheduling heuristics
for our workshop. The first is that given by the NEH algorithm ini-
4.2. The used metaheuristics
tially established for the permutation flow shop problem that we
adopt for the SDST/HFS problem. Algorithm 1 gives a description
Metaheuristics are more and more used in scheduling optimiza-
of the NEH algorithm to establish an initial solution to the SDST/
tion. In most cases, they give good solutions in a reasonable time
HFS problem.
for larger problems. In general, they start with initial solutions
obtained by heuristics, after that a disruption procedure of this
The second is based on the greedy randomized adaptive search solution is generated to give one or more neighboring solutions.
procedure, which is widespread in the scheduling domain. In this The principle is to evaluate this solution among the set of all solu-
procedure, a job is selected from a restricted candidate list (RCL) tions; if it will reduce the optimization criterion (here minimiza-
to be scheduled in the current sequence at an interval established tion case) and then the new solution will be retained. The cycle
by a margin of choice given by the objective function. The selection is repeated until it satisfies an imposed stop condition according
interval of the candidate jobs depends on the parameter q which to the nature of the studied problem. In our case, we implemented
we consider q 2 ½0:5; 0:9. We give in the following equation the two metaheuristics, namely the iterative local search and the iter-
condition to define the condition of the selection of jobs in the list ative greedy algorithms. For these two metaheuristics, we propose
RCL. two main improvements, the first relates to the initial phase, the
S. Aqil, K. Allali / Expert Systems with Applications 162 (2020) 113716 9
second concerns the exploration phase of the neighborhood. In the Where the parameter l0 makes it possible to reduce the value
initial phase, we consider a solution resulting from one of the four of the temperature T e during the evolution of the algorithm; we
rules dðuÞ ; u ¼ 1; . . . ; 4 and improved by one of the heuristics NEH choose l0 2 f0:5; . . . ; 0:9g according to the size of the treated prob-
or GRASP mentioned before. In the local search phase, we start lem. For the exploration of all the neighborhoods we have retained
with the good solution obtained in the initial phase, we generate four types of neighborhood that we present later
a set of neighboring solutions by a specific perturbation procedure,
among which we will retain the best solution. Neighborhood N1 : Insertion with left shift, we selects two posi-
tions p and q with p < q, we will note ðq pÞ ¼ c the step of
4.2.1. Iterative local search algorithm insertion on the left, for i ¼ p; . . . ; q one generates all the
The local search iterative algorithm is one of the powerful algo- sequences obtained by the insertion of all jobs with shift to
rithms in optimization by its simplicity of implementation and its the left. If for example the current sequence is
good performance in scheduling. A description of the adopted p ¼ ð7; 6; 3; 5; 4; 2; 1Þ, with p ¼ 3 and q ¼ 6, all the sequences
model in our approach of resolution is given in Algorithm 3. are X1 ¼ fð7; 6; 5; 3; 4; 2; 1Þ; ð7; 6; 5; 4; 3; 2; 1Þ; ð7; 6; 5; 4; 2; 3; 1Þ;
This algorithm has been implemented in several works espe- ð7; 6; 3; 4; 5; 2; 1Þ; ð7; 6; 3; 4; 2; 5; 1Þ; ð7; 6; 3; 5; 2; 4; 1Þg . In this
cially for flow shop scheduling problems with permutation (PFSP) case the iterative local search algorithm is noted ILS1
(Dong, Huang, & Chen, 2009; Burke et al., 2010), its extension for Neighborhood N2 : Insertion with shift on the right, we select
the same problems with SDST has been studied in Xu, Wu, Yin, two positions p and q with p < q, we note ðq pÞ ¼ c, the
and Lin (2017), or for PFSP with lot streaming in Sang, Gao, and insertion step. For i ¼ q; . . . ; p, the procedure consists of
Li (2014) and just recently with an application in solving a HFS generating the set of solutions obtained by inserting all jobs
problem in Pan et al. (2017), Zohali et al. (2019). This algorithm with shift to the right. Starting from the same previous
is based on two main phases, the first phase start from an initial sequence if p ¼ 3 and q ¼ 6 then the set of sequences
solution, it consists in exploring the neighborhood by a perturba- is X2 ¼ fð7; 6; 3; 5; 2; 4; 1Þ; ð7; 6; 3; 2; 5; 4; 1Þ; ð7; 6; 2; 3; 5; 4; 1Þ;
tion procedure, generally based on the permutation or the inser- ð7; 6; 3; 4; 5; 2; 1Þ; ð7; 6; 4; 3; 5; 2; 1Þ; ð7; 6; 5; 3; 4; 2; 1Þg. In this
tion of the jobs in the sequence. The second phase, is that of the case the local search algorithm is noted ILS2
local search allowing to accept or reject the new solution according Neighborhood N3 : Permutation with exchange of position
to a given criterion. We adopt the simulated annealing model without repetition, we chooses two different positions p and q
developed for the resolution of the PFSF (Ruiz & Stützle, 2007) in the case p < q and for i ¼ p; . . . ; q we generates the set of
and we retain all its parameters. We modify this algorithm by add- sequences obtained by exchange of positions of all the jobs in
ing a neighborhood exploration procedure during the local search this interval. For example, if the current sequence is
to generate a set of neighboring solutions. We retain the best solu- p ¼ ð7; 6; 3; 5; 4; 2; 1Þ, and for p ¼ 3 and q ¼ 6, the set of
tion of this generated set. neighboring solutions is X3 ¼ fð7; 6; 5; 3; 4; 2; 1Þ; ð7; 6; 4; 5; 3;
In this model of the iterative local search algorithm we have 2; 1Þ; ð7; 6; 2; 5; 4; 3; 1Þ; ð7; 6; 3; 4; 5; 2; 1Þ; ð7; 6; 3; 2; 4; 5; 1Þ; ð7; 6;
highlighted several versions of neighborhoods. Moreover, the nec- 3; 5; 2; 4; 1Þg in this case the local search algorithm is noted ILS3
essary adjustment in this implementation is that of the annealing
temperature T e parameter that we choose in the following form In the three neighborhood types, we generate cðc21Þ neighbor-
X
K X
n
hood solutions representing the cardinal of Xg ; g ¼ 1; 2; 3. We
min pijk þ sjjk note that we can vary the neighborhood according to the
16i6mk
T e ¼ l0
k¼1 j¼1
ð10Þ size of the problem to better explore the neighborhoods of
10 K n each type of problem in an efficient manner. We choose only
10 S. Aqil, K. Allali / Expert Systems with Applications 162 (2020) 113716
c 2 fbn=4c; bn=3c; bn=2cg for the different simulated instances phase. To distinguish the different algorithms we will write them
whose objective is to vary the entire generated neighborhood. ILSg ; g ¼ 1; 2; 3 according to the type of neighborhood chosen
In Fig. 6, we explore three neighborhood types with the current Ng ; g ¼ 1; 2; 3.
sequence p ¼ ð7; 6; 3; 5; 4; 2; 1Þ. The sets were given for p ¼ 3 and
q ¼ 6 with a c ¼ b7=2c step, for a total of 6 sequences for each 4.2.2. Iterative greedy algorithm
set Xg ; g ¼ 1; 2; 3. The second algorithm implemented for solving our SDST/
In the same way as the local search iterative algorithm, we HFS problem is the iterative greedy algorithm (see
implemented the three neighborhood versions in the local search Algorithm 4).
S. Aqil, K. Allali / Expert Systems with Applications 162 (2020) 113716 11
This algorithm is based on two main steps, the first consists of vergence towards the right solution. Our interest in this concept of
two phases, the destruction phase and the reconstruction phase. simulator is to allow the managers of the company a flexibility in
The second is that of exploration in local search, this second phase the future at the level of scheduling in the HFS workshops. This
is of great importance for search in the neighborhood of the current guarantees an extension of the platform by adding other criteria
solution. We present in the Algorithm 4 our adopted version for the to be optimized under some constraints.
resolution of our problem. In this version of the greedy iterative
algorithm we add a procedure to explore all neighborhoods based 5. Numerical simulations
on the neighborhood system defined before for the ILS algorithm.
After exploring the set Xg , we apply for each generated solution, 5.1. Characteristics of instances
the iterative greedy algorithm based on the destruction and recon-
struction constituting the first phase of IG. Subsequently, we We conducted an experimental study to test different resolu-
choose the best solution from this neighborhood to exploit it in tion approaches to verify the effectiveness of our algorithms. Our
the second phase. In this metaheuristics, we also distinguish three approach consists of simulating a set of instances and measuring
version of IGg ; g ¼ 1; 2; 3, relatively to the three types of neighbor- certain characteristics by characterizing the algorithms in order
hood exploration Xg ; g ¼ 1; 2; 3. In addition, we consider the same to better judge the relevance of the performances of the proposed
simulated annealing model adopted to validate the acceptance of approaches. The study is based on the generation of a set of
the best solution during the overall course of the algorithm. The instances covering different sizes of problems that can be the sub-
used parameters here are the step c corresponding to the explo- ject of a real industrial case study. The simulation is therefore con-
ration of the neighborhood and the temperature T e that is already ducted for instances of small, medium and large sizes. We divided
taken in ILS algorithm. Moreover, for the parameter d which repre- our instances into two large classes, the relatively small and med-
sents the number of extracted jobs, we consider a variable number ium sized class, and the relatively large class. For the first class of
according to the size of the studied problem. instance, we consider that the number of jobs n 2 f10; 20; . . . ; 90g
Our resolution approach is implemented for a much more glo- and the number of machines per stage mk 2 ½2; 6, the number of
bal vision, developing an industrial platform in the form of a stages will be defined by K 2 ½2; 6. In this category, the processing
scheduling simulator. This realization is of great utility for the times, where pijk 2 ½1; 49, the setup times are also randomize gen-
industrial case studied, indeed allows the company managers to erate, such that sjhk 2 ½1; 20. For the second classes of problems,
anticipate deciding on the scheduling in an efficient way by opti- dealing the relatively large-sized problems, we consider that
mizing the desired criteria according to the external constraints n 2 f100; 120; . . . ; 200g; mk 2 ½6; 12 and K 2 ½6; 10. However, for
in particular the fixed deadlines by their customers. We summarize this class we choose that the processing time and the setup times
the problem-solving tools on a platform modeled by a numerical will be generated respectively according to a uniform law, such as
simulator. In Fig. 7, we propose an optimization simulator pijk 2 ½50; 99 and sjhk 2 ½10; 30. The deadlines are defined accord-
designed for the SDST/HFS scheduling problem, highlighting the ing to the processing time and setup time of the jobs in the stages
different resolution phases. In the first phase, a priority rule is according to the following expression
applied, subsequently an initialization heuristic is chosen to obtain
K
X
an improved initial sequence. We adopt k as a setting parameter to
dj ¼ x max pijk þ sjjk ð11Þ
weight the processing time and the setup time for each established 16i6mk
k¼1
rule. In the second phase, two heuristics based the NEH an an
GRASP algorithms are applied to solve the problem. In the third where x is a variable depending on the studied instance size,
phase, a set of metaheuristics is implemented for current solution such that x 2 ½1; 3 for the instances of small and medium size
in order to find a good solution according to the size of the and x 2 ½4; 6 for the relatively large instances. We have taken into
instance. consideration in this expression the most unfavorable case by
Metaheuristics are based on the iterative local search and the assigning the jobs to their wrong machine on each stage assuming
iterative greedy algorithms and reinforced by the simulated that the processing time would be max16i6mk pijk þ sjjk for each
annealing model. For the six implemented algorithms, we apply job J j . To compare the resolution algorithms, we calculate the aver-
the same type of neighborhood according to the insertion step or age relative percentage deviation by the following expression
permutation procedure. In the four phases a comparative study is !
conducted to verify the effectiveness of the proposed algorithms 1 XR
TT ralgo TT best
ARPD ¼ 100% ð12Þ
by evaluating the performance of each in terms of quality and con- R r¼1 TT best
Where TT ralgo represents the value of the total tardiness for used the same sequence for all the algorithms by adopting the same
algorithm, TT best is the best value found in the set of algorithms and parameters as well as the same type of neighborhood.
R is the number of repeated experiments for each problem. Indeed, The first analysis concerns the comparative study between the
the last parameter is calculated for an average of ten instances different algorithms in terms of performances according to the step
whose purpose is to better judge the difference between the differ- c of the chosen neighborhood. We propose a numerical simulation
ent algorithms. It makes it possible to better quantify the disper- based on the comparative study using the ARPD. We present our
sion of the value of the objective function by calculation to the results with ANOVA variation curves described according to the
average value by differentiating it from the best value found. In NEH or GRASP heuristics and by specifying the studied algorithm
order to define the calculation time, this study is carried out based for each priority rules.
on a limit running time by imposing a criterion according to the In Fig. 8, we present the results of the study for the neighbor-
size of the instance in milliseconds where s ¼ f n M K, and hood step c ¼ bn=4c for a relatively average instance with
M ¼ max16k6K ðmk Þ represents the maximum number of machines n ¼ 50; K ¼ 3; m1 ¼ 3; m2 ¼ 3; m3 ¼ 4 and the parameter k ¼ 0:5.
per stage in all the workshop. The coefficient fwill be chosen The worst result is ARPD ¼ 6:98% given by the dð3Þ rule using
according to the size of the instance being studied such as GRASP for the ILS1 algorithm; While the good result is
f 2 f10; 20; 30g for the problem of low and medium size and ARPD ¼ 2:12% achieved by the dð1Þ rule for NEH and recorded by
f 2 f30; 40; 50gfor large seize problems. IG2 algorithm.
The second comparative study concerns the speed of conver- In Fig. 9, we represent the results of the neighborhood step
gence towards their best solutions for each algorithm. We propose c ¼ bn=3c with the same conditions as the previous case. The worst
a comparative study by recording the evolution of the solution result is ARPD ¼ 5:78% given by the case of dð3Þ using GRASP and
according to the number of iterations for each studied problem performed by the ILS1 algorithm. However the good result is
while fixing the necessary parameters. In order to better calibrate
ARPD ¼ 1:12% recorded by dð1Þ initialized by NEH and realized by
our algorithms we must choose the necessary parameters for each
the algorithm IG2 .
algorithm. This type of adjustment is necessary to limit our case
In the last case and in Fig. 10, we present the case of the neigh-
study, in Table 1 we summarize the main settings parameters.
borhood with c ¼ bn=2c for the same conditions as the two previ-
We note that the number of iterations limit is a fundamental
ous cases. The worst result is ARPD ¼ 2:85% given by dð3Þ in the
parameter when studying the convergence of our algorithms
case of GRASP for the ILS3 algorithm. On the other hand, the good
towards their good solutions. We did several tests to better choose
this number depending on the size of each problem. The estab- result is always done by dð1Þ using NEH and IG2 algorithm with
lished relationships are determined empirically according to the value of ARPD ¼ 0:15%. The first conclusion of this study is that
nature and characteristics of the simulated instance. the algorithm IG2 gives a good result compared to all other algo-
rithms and the best heuristic is the NEH algorithm compared to
the other heuristics based on GRASP algorithm. In the initialization
5.2. The simulation results
phase, we found that the priority rule dð1Þ makes it possible to have
best performance than the other priority rules for a given weight-
In this type of problem, the objective of numerical simulation is
ing factor k.
to characterize the relevance of the proposed resolution algo-
The second analysis of the comparative study is made on the
rithms. The instances are diversified, we have targeted the two
convergence of the algorithms in terms of the quality of the
families of instances cited above to better quantify the perfor-
obtained solution. We propose a simulation study of the objective
mance of each algorithm. We took into account the setting param-
function versus the number of iterations for the different algo-
eters and the chosen priority rule of each algorithm, an
rithms. We choose to represent the simulation of an instance in a
improvement of heuristic is considered and finally the main used
workshop consisting of four stages with a number of jobs n ¼ 40
metaheuristic of resolution will be adopted. We note that to effec-
and the maximum number of machines is M ¼ 5 in all stages.
tively judge the difference between the algorithms, we start with
The iteration limit will be expressed as 40 4 5 ¼ 800, respect-
ing the empirical relationship given in the Table 1. This instance
Table 1 is considered a relatively medium size instance reflecting an indus-
Setting parameters of the algorithms. trially real case.
Instance k Iteration limit f We propose the simulation study based on the four
Family 1 f0:2; 0:5; 0:8g nMK f10; 20; 30g dðiÞ ; i ¼ 1; . . . ; 4 priority rules mentioned before to obtain the initial
Family 2 f0:2; 0:5; 0:8g 1:5 n M K f40; 50; 60g sequence. We record the evolution of our algorithms by indicating
Fig. 8. ARPD with 95% confidence intervals for c ¼ bn=4c and NEH (left) and GRASP (right) initialization algorithms.
S. Aqil, K. Allali / Expert Systems with Applications 162 (2020) 113716 13
Fig. 9. ARPD with 95% confidence intervals for c ¼ bn=3c with NEH (left) and GRASP (right) initialization algorithms.
Fig. 10. ARPD with 95% confidence intervals for c ¼ bn=2c with NEH (left) and GRASP (right) initialization algorithms.
Fig. 11. The value of TT vs number of iteration for NEH dð1Þ (top left), NEH dð2Þ (top right), NEH dð3Þ (bottom left) and NEH dð4Þ (bottom right).
14 S. Aqil, K. Allali / Expert Systems with Applications 162 (2020) 113716
Fig. 12. The value of TT vs number of iteration for GRASP dð1Þ (top left), GRASP dð2Þ (top right), GRASP dð3Þ (bottom left)and GRASP dð4Þ (bottom right).
the two heuristics of improvement namely NEH and GRASP. Subse- the two heuristics. The results are therefore given for the six algo-
quently, for each considered case, we represent the evolution of rithms for the different cases for the four chosen rules at the
the three algorithms of ILSg ; g ¼ 1; 2; 3 for the three IGg algorithms, beginning.
g ¼ 1; 2; 3. We limit ourselves to the representation of the simula- The recorded result concludes that IG2 gives the largest differ-
tion in the case of the neighborhood X2 which we found that it ence between the height value of the objective function and the
gives good results in terms of the ARPD recorded quality in the first minimum value recorded, a difference of 307 per unit of time. This
analysis. difference represents the gain in time unit obtained from the 450th
In Fig. 11, we present the four proposed configurations for a iteration, from which all the algorithms stabilize towards their best
chosen weighting coefficient k ¼ 0:5 to calculate the different cri- solutions. This results confirms the conclusion obtained before,
teria dðiÞ ; i ¼ 1; . . . ; 4 of the priority rules. We find that for heuris- which allows us to conclude that the algorithm IG2 remains the
tics NEH dðiÞ ; i ¼ 1; . . . ; 4, the algorithm IG2 converges to the best best algorithm among the six proposed, whatever the initial
value of the objective function for the four processed configura- solution.
tions. Moreover, for the algorithms IGg ; g ¼ 1; 2; 3, they give very We continue our comparative study between the different algo-
good results compared to the algorithms ILSg ; g ¼ 1; 2; 3 for these rithms in terms of the quality of the obtained solutions by the
three proposed versions. Also, we find that the best result is given determination of ARPD for different instances. These instances
by the heuristic derived from the rule obtained by d1 and that the reveal their characteristics in this case, the number of jobs n, the
maximum deviation recorded the initialization value and the best maximum number of machines M ¼ max16k6K ðmk Þ in the set of
value obtained by IG2 . This difference represents a gain of 818 time stage K as well as the calculation limit time s. This study is done
units made from iteration 400. We have limited the number of iter- with f ¼ 20 for the first category of problem and f ¼ 40 for the sec-
ations to 800 for this instance. We estimate that this number of ond category of problem. The analysis varies the number of stages,
iterations is largely sufficient to judge the convergence of our algo- the number of jobs and the number of machines per stage accord-
rithms for this instance. We also note that after five numerical sim- ing to the classification of two categories of problems according to
ulations of the same instance that all the algorithms converge to their sizes. We retain the respective results of the study for the first
their good solutions around the 600th iteration and we retain the priority rule dð1Þ and the neighborhood step c ¼ bn=2c which gave
best obtained result. Our objective of visualizing the steady-state good results in the previous comparative studies.
convergence for this and for the other types of instance we estab- We report that we have tested several instance configurations
lish the limit value of the iteration number by an empirical study in total for the different priority rules and different initialization
given by the Table 1. In this table we give the relation expressing heuristics. We note that we conducted the study with the four pri-
the iteration number according to the characteristics of the simu- ority rules and improved by the two heuristics NEH and GRASP
lated instance. that are tested for the six algorithms. In Table 2, we represent
In Fig. 12, we also propose the same comparative study, for the the ARPD values of the first category of problems for the value of
same instance, but applying the GRASP heuristics. We also propose f ¼ 20.
the same parameters of the simulation that we have chosen for the The variation of ARPD is between 8:33% and 0:11%, this smaller
NEH heuristic in order to better concertize the difference between value of ARPD is given by the algorithm IG2 . We can see that for the
S. Aqil, K. Allali / Expert Systems with Applications 162 (2020) 113716 15
Table 2
ARPD% values for all algorithms, the best values are in bold.
NEH GRASP
K n M s ILS1 IG1 ILS2 IG2 ILS3 IG3 ILS1 IG1 ILS2 IG2 ILS3 IG3
20 2 2.4 3.51 2.36 4.75 1.24 2.75 4.75 3.37 2.75 4.37 2.67 2.69 3.57
30 3 7.2 5.24 3.45 4.67 1.14 1.56 2.25 2.25 2.67 1.88 1.58 2.12 2.47
3 40 4 9.6 7.23 6.28 4.12 4.15 5.21 4.89 8.33 7.44 5.25 5.27 8.13 5.26
50 5 15 3.22 2.25 4.59 2.14 4.36 3.22 4.56 3.25 5.33 3.22 4.78 4.21
60 6 21.6 2.56 1.22 2.58 1.06 2.47 3.22 2.89 2.28 3.86 2.68 3.44 2.45
50 2 8 2.21 1.57 2.54 1.43 2.57 1.89 2.21 1.65 2.88 1.56 2.68 1.98
60 3 14.4 1.32 0.87 1.54 0.55 1.84 0.67 2.22 1.45 2.74 1.94 2.23 1.44
4 70 4 22.4 4.34 3.24 5.29 3.87 4.25 2.46 5.21 3.32 4.78 3.24 5.23 3.76
80 5 32 1.23 0.33 1.87 0.22 1.64 0.54 1.89 1.22 1.87 1.11 2.02 1.23
90 6 43.2 0.79 0.55 0.87 0.11 0.92 0.57 1.27 0.94 1.69 0.64 1.25 0.64
50 2 10 1.67 1.33 1.98 1.02 1.89 1.54 2.22 1.87 2.86 1.11 2.85 1.14
60 3 18 3.53 2.11 3.45 2.13 3.85 2.32 4.23 3.65 5.21 3.22 5.32 3.38
5 70 4 42 1.01 0.55 1.57 0.67 1.33 0.61 1.78 0.79 1.24 0.76 1.53 0.85
80 5 40 0.89 0.33 0.89 0.25 1.01 0.67 1.25 1.01 2.27 1.03 2.04 0.78
90 6 54 2.33 1.51 2.65 1.42 2.21 1.44 3.22 2.34 4.43 2.35 3.94 2.45
Average 22.65 2.73 1.86 2.89 1.42 2.53 2.06 3.12 2.44 3.37 2.15 3.34 2.37
Table 3
ARPD% values for all algorithms, the best values are in bold.
NEH GRASP
K n M s ILS1 IG1 ILS2 IG2 ILS3 IG3 ILS1 IG1 ILS2 IG2 ILS3 IG3
100 6 144 2.31 1.86 1.75 1.94 2.15 1.75 2.37 1.85 2.17 1.67 2.69 1.57
110 7 184.8 1.04 1.45 2.67 0.74 1.56 1.85 1.25 0.87 1.88 1.18 2.12 1.47
6 120 8 230.4 0.98 0.35 0.35 0.13 0.78 0.29 1.35 0.89 1.22 0.78 1.33 1.45
130 9 280.8 0.53 0.09 0.78 0.06 0.33 0.67 1.31 0.89 1.02 0.67 1.11 0.85
140 10 336 0.35 0.11 0.87 0.05 0.33 0.23 0.98 0.55 0.33 0.25 0.76 0.45
100 6 192 2.34 2.11 1.71 1.52 2.31 1.55 2.33 1.84 2.35 1.22 3.23 1.68
110 7 246.4 3.33 1.21 3.44 1.78 3.67 2.21 3.33 2.78 4.34 2.21 3.45 1.46
8 120 8 307.2 0.47 0.29 0.53 0.22 0.89 0.36 0.89 0.78 0.95 0.33 0.78 0.52
130 9 374.4 1.11 0.78 1.13 0.52 1.44 0.45 1.56 1.22 1.56 1.16 1.98 1.32
140 10 448 0.67 0.15 0.54 0.12 0.25 0.18 1.76 0.95 1.33 0.52 1.45 0.54
100 6 240 2.45 1.64 2.56 1.78 2.67 1.54 3.35 2.98 3.45 2.46 3.62 2.76
110 7 308 5.56 4.35 5.67 3.56 3.67 2.35 7.28 5.34 5.32 5.35 7.66 3.68
10 120 8 384 2.33 1.78 2.67 1.45 2.85 1.25 3.44 2.67 4.68 2.86 4.32 2.32
130 9 468 0.65 0.55 0.74 0.33 0.65 0.62 1.22 1.03 1.67 1.87 1.22 1.04
140 10 560 0.86 0.55 0.78 0.35 0.67 0.45 1.64 1.12 2.33 1.32 1.67 1.12
Average 313.6 1.66 1.15 1.74 0.97 1.61 1.06 2.2 1.71 2.3 1.59 2.49 1.48
15 tested instances, the IG2 algorithm applying NEH gives 11 good To complete this simulation study, we give the success rate per-
results, a success rate of 73:33%. In addition we find that the algo- centage (SRP) of the IG2 algorithm compared to the other algo-
rithms from IG record 93:33% success rate compared to the ILS rithms for the different parameters and the treated family of
algorithms for the same heuristic. problems. Its value is calculated by counting the total number
Also, in Table 3, we show the results of ARDP for the second cat- when ARPD of IG2 is better, divided by the total number of
egory of problems called relatively large with f ¼ 40. This last cat- instances tested, in this case 15.
egory is of great industrial utility in the real case where the We tested 15 problems for each family, we varied
problem requires a resolution by the metaheuristics. We can find f 2 f10; 20; 30g for the first family and f 2 f30; 40; 50g for the sec-
that the values of ARPD vary, for relatively large problems, ond family. The metaheuristics are initialized with the heuristics
between 7:66% and 0:05%. Also this last value is realized by the NEH and GRASP algorithms by varying the weighting coefficient
algorithm IG2 . The latter represents an effective success rate of k 2 f0:2; 0:5; 0:8g obtained by the four priority rules in the initial
80%. Similarly in this category the IG algorithm with its three ver- phase. For each family of problems we have 4 3 3 15 ¼ 540
sion is much more efficient compared to the three versions of ILS simulated problems either in total for both families 1080 problems
algorithm . are treated.
16 S. Aqil, K. Allali / Expert Systems with Applications 162 (2020) 113716
Table 4
SRP for the algorithm IG2 .
Family 1 Family 2
f 10 20 30 30 40 50
k ¼ 0:2 60 66.66 73.33 66.66 73.33 80
k ¼ 0:5 66.66 73.33 80 73.33 80 86.66
k ¼ 0:8 73.33 80 86.66 80 86.66 93.33
Average 66.66 73.33 80 73.33 80 86.66
Finally, in Table 4, we can note the high values of the success worthy to consider, as future perspectives, multi-objective opti-
percentages of the IG2 algorithm for the studied parameters. We mization problem such as minimizing the cost of the production
note that all values exceed 60% and reach the best value of and/or energy consumption. These objectives are usually encoun-
SRP ¼ 93:33%. This high performance recorded by IG2 confirms tered in many manufacturing industries. Moreover, it will interest-
the superiority of the iterative greedy algorithm to schedule this ing to apply the resolution method to other industrial cases other
type of problem. than slabs and beams issue. In addition, we note that the sequence
Those numerical simulations reveal that the iterated greedy of jobs may be independent of setup time. All these questions can
algorithm based on NEH heuristic initialization which is applied be considered as potential future perspectives of the current work.
to solve slabs and beams scheduling problem gives good results
for various neighborhood explorations. According to this analysis,
the IG algorithm in its second version is the best algorithm Declaration of Competing Interest
released from this comparative study. In general, the IG algorithm
is more efficient given its perturbation procedure to generate The authors declare that they have no known competing finan-
neighboring solutions. This quality of better disrupting the current cial interests or personal relationships that could have appeared
solution makes it possible to effectively explore the research direc- to influence the work reported in this paper.
tions in such an optimization problem.
References
6. Conclusion and future perspectives
Abyaneh, S. H., & Zandieh, M. (2012). Bi-objective hybrid flow shop scheduling with
This paper presents an industrial optimization study of slabs sequence-dependent setup times and limited buffers. The International Journal
of Advanced Manufacturing Technology, 58(1–4), 309–325.
and beams production problem. The proposed approach consists Ahonen, H., & de Alvarenga, A. G. (2017). Scheduling flexible flow shop with
in modeling the process as a hybrid flow shop with unrelated par- recirculation and machine sequence-dependent processing times: Formulation
allel machines under constraints of sequence dependent setup and solution procedures. The International Journal of Advanced Manufacturing
Technology, 89(1–4), 765–777.
time. We propose several initialization rules enhanced by NEH Alaykỳran, K., Engin, O., & Döyen, A. (2007). Using ant colony optimization to solve
and GRASP heuristics. In the phase of improvement by metaheuris- hybrid flow shop scheduling problems. The International Journal of Advanced
tics, we suggest two principal algorithms of resolutions, the first Manufacturing Technology, 35(5–6), 541–550.
Allaoui, H., & Artiba, A. (2014). Hybrid flow shop scheduling with availability
one is based on the iterative local search, the second is based on constraints. In Essays in production, project planning and scheduling
the iterative greedy algorithm. For the two algorithms, we present (pp. 277–299). Springer.
three versions based on three types of neighborhoods, the insertion Besbes, W., Teghem, J., & Loukil, T. (2010). Scheduling hybrid flow shop problem
with non-fixed availability constraints. European Journal of Industrial
with right shift, the insertion with shift on the left and the total Engineering, 4(4), 413–433.
permutation without repetitions of the current sequence subset. _
Bozejko, W., Pempera, J., & Smutnicki, C. (2013). Parallel tabu search algorithm for
Our approach allows to explore all the sequences obtained by the the hybrid flow shop problem. Computers & Industrial Engineering, 65(3),
466–474.
neighborhoods of the current solution and later we apply our mod-
Burke, E., Curtois, T., Hyde, M., Kendall, G., Ochoa, G., Petrovic, S., Vázquez-
ified algorithms to retain the best solution. We perform a simula- Rodríguez, J. A. & Gendreau, M. (2010). Iterated local search vs. hyper-
tion of an instance set, in total 1080 instances for several types heuristics: Towards general-purpose search algorithms. In: 2010 IEEE
of problems. The simulation are diversified covering, the small, congress on evolutionary computation (CEC) (pp. 1–8). IEEE..
Chung, T.-P., & Liao, C.-J. (2013). An immunoglobulin-based artificial immune
medium and relatively large instance sizes to better concertize system for solving the hybrid flow shop problem. Applied Soft Computing, 13(8),
the industrial context. We find that the iterative greedy algorithm 3729–3736.
IG2 with its new implementation version gives good results in Chung, T.-P., Sun, H., & Liao, C.-J. (2017). Two new approaches for a two-stage
hybrid flowshop problem with a single batch processing machine under waiting
terms of quality and convergence time to the best solution. time constraint. Computers & Industrial Engineering, 113, 859–870.
Our resolution approach is limited to ten-stage SDST/HFS Davoudpour, H., & Ashrafi, M. (2009). Solving multi-objective sdst flexible flow shop
scheduling problem. This study can be extended by considering using grasp algorithm. The International Journal of Advanced Manufacturing
Technology, 44(7–8), 737–747.
several stages. The number of jobs can be also much greater than Dessouky, M. M., Dessouky, M. I., & Verma, S. K. (1998). Flowshop scheduling with
the one presented in this study. We have limited ourselves to the identical jobs and uniform parallel machines. European Journal of Operational
constraint of the setup time, other constraints can be studied like Research, 109(3), 620–631.
Dios, M., Fernandez-Viagas, V., & Framinan, J. M. (2018). Efficient heuristics for the
the unavailability of the machines during maintenance. The sug- hybrid flow shop scheduling problem with missing operations. Computers &
gested resolution approaches consisting of some well-known Industrial Engineering, 115, 88–99.
heuristics and metaheuristics can be useful for the production Dong, X., Huang, H., & Chen, P. (2009). An iterated local search algorithm for the
permutation flowshop problem with total flowtime criterion. Computers &
managers to find a good solution in a reasonable time. However,
Operations Research, 36(5), 1664–1669.
it can be interesting to develop another hybrid metaheuristic to Engin, O., Ceran, G., & Yilmaz, M. K. (2011). An efficient genetic algorithm for hybrid
solve such case of industrial problems. flow shop scheduling with multiprocessor task problems. Applied Soft
In our work, we have studied HFS problem with unrelated par- Computing, 11(3), 3056–3065.
Engin, O., & Döyen, A. (2004). A new approach to solve hybrid flow shop scheduling
allel machines under the constraint sequence dependent setup problems by artificial immune system. Future Generation Computer Systems, 20
time applied to slabs and beams manufacturing case. It will be (6), 1083–1095.
S. Aqil, K. Allali / Expert Systems with Applications 162 (2020) 113716 17
Eskandari, H., & Hosseinzadeh, A. (2014). A variable neighbourhood search for Öztop, H., Tasgetiren, M. F., Eliiyi, D. T., & Pan, Q.-K. (2019). Metaheuristic
hybrid flow-shop scheduling problem with rework and set-up times. Journal of algorithms for the hybrid flowshop scheduling problem. Computers &
the Operational Research Society, 65(8), 1221–1231. Operations Research.
Fernandez-Viagas, V., Molina-Pariente, J. M., & Framinan, J. M. (2018). New efficient Pan, Q.-K., & Dong, Y. (2014). An improved migrating birds optimisation for a hybrid
constructive heuristics for the hybrid flowshop to minimise makespan: A flowshop scheduling with total flowtime minimisation. Information Sciences,
computational evaluation of heuristics. Expert Systems with Applications, 114, 277, 643–655.
345–356. Pan, Q.-K., Gao, L., Li, X.-Y., & Gao, K.-Z. (2017). Effective metaheuristics for
Figielska, E. (2009). A genetic algorithm and a simulated annealing algorithm scheduling a hybrid flowshop with sequence-dependent setup times. Applied
combined with column generation technique for solving the problem of Mathematics and Computation, 303, 89–112.
scheduling in the hybrid flowshop with additional resources. Computers & Pan, Q.-K., Ruiz, R., & Alfaro-Fernández, P. (2017). Iterated search methods for
Industrial Engineering, 56(1), 142–151. earliness and tardiness minimization in hybrid flowshops with due windows.
Gómez-Gasquet, P., Andrés, C., & Lario, F.-C. (2012). An agent-based genetic Computers & Operations Research, 80, 50–60.
algorithm for hybrid flowshops with sequence dependent setup times to Pan, Q.-K., Wang, L., Mao, K., Zhao, J.-H., & Zhang, M. (2013). An effective artificial
minimise makespan. Expert Systems with Applications, 39(9), 8095–8107. bee colony algorithm for a real-world hybrid flowshop problem in steelmaking
González-Neira, E. M., & Montoya-Torres, J. R. (2017). A grasp meta-heuristic for the process. IEEE Transactions on Automation Science and Engineering, 10(2),
hybrid flowshop scheduling problem. Journal of Decision Systems, 26(3), 307–322.
294–306. Pinedo, M. (2012). Scheduling. Springer.
Graham, R. L., Lawler, E. L., Lenstra, J. K. & Kan, A. R. (1979). Optimization and Rabiee, M., Rad, R. S., Mazinani, M., & Shafaei, R. (2014). An intelligent hybrid meta-
approximation in deterministic sequencing and scheduling: a survey. In: Annals heuristic for solving a case of no-wait two-stage flexible flow shop scheduling
of discrete mathematics (Vol. 5, pp. 287–326). Elsevier.. problem with unrelated parallel machines. The International Journal of Advanced
Hong, W.-K., Park, S.-C., Lee, H.-C., Kim, J.-M., Kim, S.-I., Lee, S.-G., & Yoon, K.-J. Manufacturing Technology, 71(5–8), 1229–1245.
(2010). Composite beam composed of steel and precast concrete (modularized Ramezani, P., Rabiee, M., & Jolai, F. (2015). No-wait flexible flowshop with uniform
hybrid system). Part iii: Application for a 19-storey building. The Structural parallel machines and sequence-dependent setup time: A hybrid meta-
Design of Tall and Special Buildings, 19(6), 679–706. heuristic approach. Journal of Intelligent Manufacturing, 26(4), 731–744.
Jiang, S., Liu, M., Hao, J., & Qian, W. (2015). A bi-layer optimization approach for Ribas, I., Companys, R., & Tort-Martorell, X. (2019). An iterated greedy algorithm for
a hybrid flow shop scheduling problem involving controllable processing solving the total tardiness parallel blocking flow shop scheduling problem.
times in the steelmaking industry. Computers & Industrial Engineering, 87, Expert Systems with Applications, 121, 347–361.
518–531. Ruiz, R., Pan, Q.-K., & Naderi, B. (2019). Iterated greedy methods for the distributed
Jolai, F., Rabiee, M., & Asefi, H. (2012). A novel hybrid meta-heuristic algorithm for a permutation flowshop scheduling problem. Omega, 83, 213–222.
no-wait flexible flow shop scheduling problem with sequence dependent setup Ruiz, R., & Stützle, T. (2007). A simple and effective iterated greedy algorithm for the
times. International Journal of Production Research, 50(24), 7447–7466. permutation flowshop scheduling problem. European Journal of Operational
Jun, S., & Park, J. (2015). A hybrid genetic algorithm for the hybrid flow shop Research, 177(3), 2033–2049.
scheduling problem with nighttime work and simultaneous work constraints: A Ruiz, R., & Vázquez-Rodríguez, J. A. (2010). The hybrid flow shop scheduling
case study from the transformer industry. Expert Systems with Applications, 42 problem. European Journal of Operational Research, 205(1), 1–18.
(15–16), 6196–6204. Sang, H., Gao, L., & Li, X. (2014). An iterated local search algorithm for the lot-
Kahraman, C., Engin, O., Kaya, I., _ & Öztürk, R. E. (2010). Multiprocessor task streaming flow shop scheduling problem. Asia-Pacific Journal of Operational
scheduling in multistage hybrid flow-shops: A parallel greedy algorithm Research, 31(06), 1450045.
approach. Applied Soft Computing, 10(4), 1293–1300. Schulz, S., Neufeld, J. S., & Buscher, U. (2019). A multi-objective iterated local search
Khare, A., & Agrawal, S. (2019). Scheduling hybrid flowshop with sequence- algorithm for comprehensive energy-aware hybrid flow shop scheduling.
dependent setup times and due windows to minimize total weighted earliness Journal of Cleaner Production, 224, 421–434.
and tardiness. Computers & Industrial Engineering, 135, 780–792. Shao, W., Shao, Z., & Pi, D. (2020). Modeling and multi-neighborhood iterated
Komaki, G., Teymourian, E., & Kayvanfar, V. (2016). Minimising makespan in the greedy algorithm for distributed hybrid flow shop scheduling problem.
two-stage assembly hybrid flow shop scheduling problem using artificial Knowledge-Based Systems, 105527.
immune systems. International Journal of Production Research, 54(4), 963–983. Solano-Charris, E. L., Montoya-Torres, J. R., & Paternina-Arboleda, C. D. (2011). Ant
Lei, D., & Guo, X. (2016). Hybrid flow shop scheduling with not-all-machines colony optimization algorithm for a bi-criteria 2-stage hybrid flowshop
options via local search with controlled deterioration. Computers & Operations scheduling problem. Journal of Intelligent Manufacturing, 22(5), 815–822.
Research, 65, 76–82. Tang, L., & Wang, X. (2010). An improved particle swarm optimization algorithm for
Li, J.-Q., & Pan, Q.-K. (2015). Solving the large-scale hybrid flow shop scheduling the hybrid flowshop scheduling to minimize total weighted completion time in
problem with limited buffers by a hybrid artificial bee colony algorithm. process industry. IEEE Transactions on Control Systems Technology, 18(6),
Information Sciences, 316, 487–502. 1303–1314.
Li, J.-Q., Pan, Q.-K., & Duan, P.-Y. (2016). An improved artificial bee colony algorithm Vignier, A., Billaut, J.-C., & Proust, C. (1999). Les problèmes d’ordonnancement de
for solving hybrid flexible flowshop with dynamic operation skipping. IEEE type flow-shop hybride: état de l’art. RAIRO-Operations Research-Recherche
Transactions on Cybernetics, 46(6), 1311–1324. Operationnelle, 33(2), 117–183.
Li, J.-Q., Pan, Q.-K., & Wang, F.-T. (2014). A hybrid variable neighborhood search for Wang, H.-M., Chou, F.-D., & Wu, F.-C. (2011). A simulated annealing for hybrid flow
solving the hybrid flow shop scheduling problem. Applied Soft Computing, 24, shop scheduling with multiprocessor tasks to minimize makespan. The
63–77. International Journal of Advanced Manufacturing Technology, 53(5–8), 761–776.
Liao, C.-J., Tjandradjaja, E., & Chung, T.-P. (2012). An approach using particle swarm Wang, S., & Liu, M. (2013). A genetic algorithm for two-stage no-wait hybrid flow
optimization and bottleneck heuristic to solve hybrid flow shop scheduling shop scheduling problem. Computers & Operations Research, 40(4), 1064–1075.
problem. Applied Soft Computing, 12(6), 1755–1764. Wang, X., & Tang, L. (2009). A tabu search heuristic for the hybrid flowshop
Luo, H., Du, B., Huang, G. Q., Chen, H., & Li, X. (2013). Hybrid flow shop scheduling scheduling with finite intermediate buffers. Computers & Operations Research,
considering machine electricity consumption cost. International Journal of 36(3), 907–918.
Production Economics, 146(2), 423–439. Xu, J., Wu, C.-C., Yin, Y., & Lin, W.-C. (2017). An iterated local search for the multi-
Marichelvam, M., Azhagurajan, A., & Geetha, M. (2018). Minimisation of total objective permutation flowshop scheduling problem with sequence-dependent
tardiness in hybrid flowshop scheduling problems with sequence dependent setup times. Applied Soft Computing, 52, 39–47.
setup times using a discrete firefly algorithm. International Journal of Operational Ying, K.-C., & Lin, S.-W. (2018). Minimizing makespan for the distributed hybrid
Research, 32(1), 114–126. flowshop scheduling problem with multiprocessor tasks. Expert Systems with
Meng, L., Zhang, C., Shao, X., Ren, Y., & Ren, C. (2018). Mathematical modelling Applications, 92, 132–141.
and optimisation of energy-conscious hybrid flow shop scheduling problem Yu, C., Semeraro, Q., & Matta, A. (2018). A genetic algorithm for the hybrid flow shop
with unrelated parallel machines. International Journal of Production Research, scheduling with unrelated machines and machine eligibility. Computers &
1–27. Operations Research, 100, 211–229.
Mirsanei, H., Zandieh, M., Moayed, M. J., & Khabbazi, M. R. (2011). A simulated Zandieh, M., Ghomi, S. F., & Husseini, S. M. (2006). An immune algorithm approach
annealing algorithm approach to hybrid flow shop scheduling with sequence- to hybrid flow shops scheduling with sequence-dependent setup times. Applied
dependent setup times. Journal of Intelligent Manufacturing, 22(6), 965–978. Mathematics and Computation, 180(1), 111–127.
Naaman, A. E. (1982). Prestressed concrete analysis and design: Fundamentals. Zhang, B., Pan, Q.-K., Gao, L., Zhang, X.-L., Sang, H.-Y., & Li, J.-Q. (2017). An effective
McGraw-Hill New York. modified migrating birds optimization for hybrid flowshop scheduling problem
Nawaz, M., Enscore, E. E., Jr, & Ham, I. (1983). A heuristic algorithm for the m- with lot streaming. Applied Soft Computing, 52, 14–27.
machine, n-job flow-shop sequencing problem. Omega, 11(1), 91–95. Zohali, H., Naderi, B., Mohammadi, M., & Roshanaei, V. (2019). Reformulation,
Newton, M. H., Riahi, V., Su, K., & Sattar, A. (2019). Scheduling blocking flowshops linearization, and a hybrid iterated local search algorithm for economic lot-
with setup times via constraint guided and accelerated local search. Computers sizing and sequencing in hybrid flow shop problems. Computers & Operations
& Operations Research, 109, 64–76. Research, 104, 127–138.