A Three-Stage ACO-Based Algorithm for Parallel Batch Loading and Scheduling Problem with Batch Deterioration and Rate-Modifying Activities

Abstract

This paper addresses a batch loading and scheduling problem of minimizing the makespan on parallel batch processing machines. For batch loading, jobs with compatible families can be assigned to the same batch process even if they differ in size; however, batches can only be formed from jobs within the same family, and the batch production time is determined by the family. During the batch scheduling, the deterioration effects are continuously added to batches processed in each parallel machine so that the batch production times become deteriorated. The deteriorated processing time of batches can be recovered to the original processing times of batches by a maintenance or cleaning process of machines. In this problem, we sequentially determine the batching of jobs and the scheduling of batches. Due to the complexity of the problem, we proposed a three-stage ant colony optimization algorithm. The proposed algorithm found an optimal solution for small-sized problems and achieved near-optimal solutions and better performance than a genetic algorithm or a particle swarm optimization for large-sized problems.

1. Introduction

Batch-processing machines (BPMs) have been applied to numerous manufacturing industries such as ceramics, steel, and integrated circuits industries to enhance the productivity of production. Due to this reason, several BPM scheduling problems have been studied in recent years. In general, the BPM sequentially processes batches, a group of jobs processed together in the same machine. The job is the smallest unit of an order requested by a customer. In this paper, we deal with the batch loading and scheduling problem (BLSP) at the diffusion operation in the semiconductor industry [1,2]. In addition to general batching, the concept of job family processing the same operation is adopted. Owing to the chemical nature of the diffusion operation in the semiconductor industry, only jobs with the same family can be assigned to the same batch. Thus, in this manufacturing environment, the batch production time is determined by the family type. In most batch scheduling studies, the batch production time is assumed to be constant. However, for a real-world scheduling problem, the batch production times increase due to the inclusion of activities such as the loading/unloading of jobs and alignment/calibration of tools. The increased batch production times can be recovered to the original production times of each batch by a maintenance or cleaning process of machines. The recovering process is called rate-modifying activity (RMA) [3]; multiple RMAs are considered in the schedule. The RMA time is assumed to be constant, and the RMA can be scheduled between batches. In this paper, the deterioration of batch production time linearly depends on the consecutive batch production runs without the RMA. Thus, it can be formulated as a linear function of the interval between the starting time of the first batch after the previous RMA and the completion time of the last batch before the recent RMA. The interval between RMAs including periods before the first and after the last RMAs is called a bucket, which is defined in Joo and Kim [4].

In this paper, we address BLSP with incompatible job families in parallel BPMs subject to time-dependent batch deterioration and RMAs applied between batches. The BLSP can be decomposed into the two sub-problems of batch loading and batch scheduling [5,6,7]. In the example of a batch loading problem shown in Figure 1a, 11 jobs belonging to three families are shown. All of the jobs must be formed in batches, which can only be assigned to batches within the same family. For example, batch 1 is formed from jobs 1 and 6, which have the same family type. The sum of their job sizes does not exceed the batch capacity. In Figure 1a, job 1–4, 5–8, and 9–11 belong to family 1, 2, and 3, respectively. The size of job is (0.43, 0.49, 0.34, 0.40, 0.49, 0.54, 0.43, 0.66, 0.37, 0.40, 0.46). The machine capacity is 1. The sum of job sizes in a batch does not exceed the batch capacity. For example, the sum of job sizes of job 1 and 4 is 0.83 ($S_1+S_4=0.43+0.40=0.83\le 1)$. After batches have been formed from all of the jobs, the batches must be scheduled to machines, which is referred to as the batch scheduling procedure.

Here, the problem with batch deterioration and RMAs during batch scheduling is considered. An illustrative example of the batch scheduling problem is shown in Figure 1b. As the batch production time including the deterioration increases, the assigning of RMAs between batches should be considered. Since the problem becomes complex, the batch scheduling problem is decomposed into three stages using buckets. The bucket is the set of batches between RMAs. The first stage is to determine the number of buckets. The second stage is to assign batches to each bucket. The last stage is to schedule the buckets to machines. In Figure 1b, the number of buckets is set to 3. The processing time for family is (44, 36, 52). Since batch 1–2, 3–5, and 6–7 consist of jobs for family 1, 2, and 3, respectively, the processing time of batch 1–2, 3–5, and 6–7 is 44, 36, and 52, respectively. The number of machines is 2. The deterioration rate is 0.25 and the processing time of RMA is 30. As a batch is processed, time is added equally to the difference between the last RMA completion time and the batch start time multiplied by the deterioration rate. For example, when batch 4 is processed in bucket 1, since the processing time of batch 2 has elapsed without RMA, 11 ($=44\times 0.25)$ is added to the original processing time of batch 4. On the other hand, batch 6 of bucket 3 is not affected by deterioration because there is an RMA immediately before it. Once the batches have been assigned and sequenced in three buckets, the buckets are scheduled to the corresponding machines. As a result, the makespan can be calculated.

The BLSP with no deterioration was addressed by Jia et al. [2], who proposed ant colony optimization (ACO) and the multi-fit (MF) algorithm for solving the batch loading and batch scheduling problems, respectively. We proposed the same algorithm for the batch loading problem, but the MF algorithm cannot be used in our problem because batch deterioration and RMAs impose complexity on the batch scheduling problem. Therefore, we decompose the batch scheduling problem into the three stages mentioned above. In the first stage, the number of buckets is determined between the lower and upper bounds. The fixed number of buckets determines the positions of the RMAs and helps to reduce the search space. In the second stage, a rule-based ACO is applied to assign batches to the buckets defined in the preceding stage. In the last stage, the buckets are scheduled to machines using a dispatching rule.

The remainder of this paper is organized as follows. A previous related work is reviewed in Section 2. Section 3 presents our mixed-integer programming (MIP) model. In Section 4, Section 5 and Section 6, the overall batch scheduling algorithm is proposed, with the MF algorithm applying the RMA scheduling rule developed in Section 4, the ACO-based three-stage algorithm proposed in Section 5, and the genetic algorithm (GA) and particle swarm optimization (PSO) approaches proposed in Section 6. Section 7 presents our experimental results and, finally, the conclusions are presented in Section 8.

2. Literature Review

To enhance productivity, many BPM studies have attracted the attention of various production areas such as ceramics, steel, and integrated circuit industries. For reviewing the related studies on BPM, we divided BPM studies into seven categories: manufacturing system, production methods, family constraint, job sizes, deterioration constant, recovering process, method, and objective function. The initial studies on BPMs assumed the size of jobs in a batch to be identical. The batch scheduling problem was first introduced by Ikura and Gimple [8]. They proposed a polynomial-time algorithm for minimizing the makespan under an assumption of identical job sizes. In this problem, the total job size in each batch must not exceed the batch capacity. For solving the BLSP with identical job sizes, Lee et al. [9] proposed dynamic programming algorithms to minimize the maximum tardiness and the number of tardy jobs. Cheng et al. [10] proposed a novel ACO for minimizing the makespan. Subsequent studies have extended this approach to single BPMs with non-identical job sizes. Uzsoy [11] proposed the use of heuristics and a B&B algorithm to minimize the total weighted completion time. Two heuristics for solving this problem were proposed by Dupont and Ghazvini [5]. Later studies extended this to the problem of parallel BPMs with non-identical job sizes. Jia et al. [12] considered parallel BPMs with arbitrary capacities carrying out jobs of arbitrary sizes with dynamic arrival times. They proposed the use of two meta-heuristics based on ACO to minimize the makespan. Ozturk et al. [13] considered the problem of parallel BPMs with identical job sizes and release dates. They proposed a branch and bound (B&B) algorithm to minimize the makespan. Zhou et al. [6] considered parallel BPMs with arbitrary release times and non-identical job sizes. They presented different heuristics for solving batch loading and batch scheduling problems. The batch loading problem is tractable under the application of GA, ACO, and greedy randomized adaptive search procedure algorithms. Jia et al. [14] proposed a fuzzy ACO for minimizing the makespans for parallel batch processing machines with jobs having non-identical sizes and fuzzy processing times. Parallel BPM studies have also been extended to unrelated parallel BPMs. Arroyo and Leung [15,16] considered unrelated parallel batch machines with different capacities processing jobs with arbitrary sizes and non-zero ready times. They proposed several heuristics based on rules applying best-fit, first-fit, and meta-heuristic-based iterated greedy algorithms. Arroyo et al. [17] considered the use of unrelated parallel BPMs to minimize the total job flow time. They proposed a simple and effective iterated greedy algorithm and compared them to three meta-heuristic algorithms. Gao and Yuan [18] considered unbounded parallel batch scheduling problems involving jobs with agreeable release dates and processing times. They demonstrated that such problems are binary NP-hard. Zhou et al. [19] considered unrelated parallel BPM processing jobs with non-identical sizes and arbitrary release times. They proposed using a GA with a random key to schedule batches to the machines in such cases.

The studies described above focused on problems in which all jobs were compatible. However, the issue of incompatible jobs has occurred in actual situations. For instance, incompatibilities arise from a non-empty intersection on intervals and two-dimensional volumes for individual jobs [20,21]. Uzsoy [1] considered the problem of single BPM with incompatible families. He developed B&B and polynomial-time dynamic programming algorithms to minimize the makespan, lateness, and total weighted completion time. Azizoglu and Webster [22] considered a BPM with incompatible families, identical processing times, arbitrary job weights, and arbitrary job sizes. They proposed a B&B algorithm that can solve this problem with less than 25 jobs in a reasonable time. Dupont and Dhaenens-Flipo [7] considered a BPM with non-identical job sizes for minimizing the makespan. To solve problems with large numbers of jobs, they added two rules to a B&B algorithm to minimize the makespan. Yao et al. [23] considered single BPMs with incompatible job families and dynamic arrivals. They proposed a decomposed B&B algorithm for minimizing the total batch completion time and makespan. Jolai [24] considered BPMs with identical job sizes and incompatible families and applied a polynomial-time dynamic programming algorithm to minimize the number of tardy jobs. Parallel BPMs with incompatible families have also been extensively examined. Balasubramanian et al. [25] considered the problem of identical parallel BPMs with incompatible families. They proposed an approach in which the batch processing time is determined based on the family and developed two GAs for scheduling jobs to batches and machines, respectively. Li et al. [26] further considered the problem of incompatible families and proposed a method for determining the batch processing time based on the length of the longest job. For parallel BPMs with incompatible job families and release times, constraint programming and apparent tardiness cost approaches were presented [27,28].

In previous studies on the BLSP, the processing time was assumed to be constant and known in advance. In several real-world industry situations, however, the processing time increases due to the deterioration phenomenon. The concept of a deteriorating job was introduced by Gupta and Gupta [29], who considered a single BPM governed by a linear deterioration function. The literature on deterioration follows two general tracks: sequence- and time-dependent deterioration. The case of parallel machine scheduling under sequence-dependent deterioration was addressed by Ding et al. [30], who proposed an ejection chain algorithm for minimizing the completion-time-based criteria. The case of time-dependent deterioration was first introduced by Browne and Yechiali [31]. Soleimani et al. [32] proposed cat swarm optimization on unrelated parallel machine scheduling problems with time-dependent deterioration. In general, the degree of deterioration increases with the length of the job process. To recover the inefficient execution of processes under deterioration, it is necessary to carry out the recovering processes called RMAs by Lee and Leon [3]. RMAs return a machine to its initial state. Joo and Kim [33] considered a single BPM with time-dependent deterioration to which RMAs are applied and proposed hybrid meta-heuristic algorithms to minimize the makespan. Woo et al. [34] considered unrelated parallel BLSPs with time-dependent deterioration and RMAs, proposing a GA with a random key to solve the problem. Abdullah and Süer [35] considered the decision-making on the selection of a manufacturing strategy between the classical assembly line and seru according to the skill levels of the operator. Liu et al. [36] proposed exact solutions and heuristic algorithms to minimize the makespan and balance the worker’s workload of divisional and rotating seru in the seru production system. Gai et al. [37] presented an accurate dimensionality reduction algorithm to minimize makespan. It is also compared to the greedy algorithm to verify it.

As a problem becomes more complex, it becomes increasingly important to decompose it. The meta-heuristic approaches developed in previous studies generally decompose a BLSP into sub-problems. Dupont and Dhaenens-Flipo [7] considered unrelated parallel BPMs with non-identically sized jobs. They developed an approach for decomposing the BLSP into two sub-problems of forming and scheduling batches to improve performance. Determining which algorithm to apply to each sub-problem following division is another important challenge. Dupont and Ghazvini [5] assessed GA, ACO, and greedy adaptive search heuristic algorithms as methods for addressing the batch scheduling problem and compared their performances on batch loading situations fixed using the same algorithm. Similarly, Zhou et al. [19] proposed a GA with a random key for solving the batch scheduling problem. In this paper, the BLSP for processes with time-dependent batch deterioration for which RMAs are applied between batches is considered. To the best of our knowledge, this is the first paper to address this problem.

A list of recent studies on BLSPs and scheduling problems under deterioration is provided in

Table 1. The comparisons between recent studies.

	Mfg. System	Production Methods	Family Constraint	Job sizes	Deterioration Constraint	Recovery Process	Method	Objective
Woo et al. [34]	Parallel	Single	No family	N/A	Time-dependent deterioration	RMA	GA	Makespan
Ding et al. [30]	Parallel	Single	No family	N/A	Sequence-dependent deterioration	No recovery	ECA	Completion time
Soleimani et al. [32]	Parallel	Single	No family	N/A	Time-dependent deterioration	No recovery	GA, CSO, IABC	Mean weighted tardiness
Ozturk et al. [13]	Parallel	Batch	No family	Identical	No deterioration	No recovery	B&B	Makespan
Jia et al. [14]	Parallel	Batch	No family	Non-identical	No deterioration	No recovery	ACO	Makespan
Jia et al. [12]	Parallel	Batch	No family	Non-identical	No deterioration	No recovery	ACO	Makespan
Arroyo et al. [17]	Parallel	Batch	No family	Non-identical	No deterioration	No recovery	Iterated greedy	Total flow time
Li et al. [26]	Parallel	Batch	Incompatible family	Non-identical	No deterioration	No recovery	LB	Lateness
Jia et al. [2]	Parallel	Batch	Incompatible family	Non-identical	No deterioration	No recovery	ACO	Makespan
This paper	Parallel	Batch	Incompatible family	Non-identical	Time-dependent deterioration	RMA	ACO	Makespan

. The recent studies are categorized by manufacturing system, production methods, family constraint, job sizes, deterioration constant, recovery process, method, and objective function. To the best of our knowledge, none of the studies dealt with the BLSP simultaneously considering parallel BPMs, incompatible job families, non-identical job sizes, time-dependent deterioration, and RMAs. Among the research, Jia et al. [2] proposed the ACO and MF algorithm for addressing, respectively, the batch loading and batch scheduling problems for parallel BPMs with incompatible families and arbitrary job sizes. For our problem, we adopted the ACO by Jia et al. [2] to solve the batch loading problem. However, the MF algorithm cannot provide the near-optimal makespan during the batch scheduling because time-dependent deterioration and RMAs increase the complexity of the problem. Therefore, in this paper, we decompose the batch scheduling problem into three stages and propose an ACO-based three-stage algorithm for solving the scheduling problem.

3. Mixed-Integer Programming Model

The BLSP is the problem of scheduling jobs to machines. In this problem, we have to make two decisions. We need to determine which jobs to assign in which batches and then schedule the assigned batches and RMA to machines. This is called batch loading problem and batch scheduling problem. In this section, the BLSP for time-dependent deterioration with rate-modifying activities and incompatible families for minimizing the makespan is formulated using MIP. The following parameters and decision variables are used in the mathematical formulation:

Parameters
J	set of jobs
F	set of families
B	set of batches
K	set of buckets
M	set of machines
Jf	set of jobs not belongings to family $f\in F$
Fj	family type of job $j\in J$
Fb	family type of batch $b\in B$
Sj	size of job $j\in J$
Pb	processing time of batch $b\in B$
Q	processing time of RMA
DR	deterioration rate
L	large number
SC	size of machine
Decision variables
Xjb	Equals 1 if job $j\in J$ is assigned in batch $b\in B$
Ybkm	Equals 1 if batch $b\in B$ is assigned to bucket $k\in K$ from machine $m\in M$
Zabkm	Equals 1 if batch $a\in B$ precedes batch $b\in B$ in bucket $k\in K$ from machine $m\in M$
Dependent variables
MS	Makespan
Ck	Completion time of bucket $k\in K$
Cm	Completion time of machine $m\in M$
Tb	Time gap between starting time of batch $b\in B$ and completion time of recent RMA

Based on parameters and decision variables, the MIP is formulated as follows:

$$Minimize MS$$

Subject to

$${\displaystyle \sum }_{b\in B}{X}_{jb}=1 , \forall j\in J$$

(1)

$$X_{jb}=0, \forall b\in B,\forall j\in J^{F_b},$$

(2)

$${\displaystyle \sum }_{j\in J}{S}_j·X_{jb}\le SC, \forall b\in B$$

(3)

$${\displaystyle \sum }_{j\in J}{X}_{jb}\le L·{\displaystyle \sum }_{k\in K}{\displaystyle \sum }_{m\in M}{Y}_{bkm }, \forall b\in B$$

(4)

$${\displaystyle \sum }_{k\in K}{\displaystyle \sum }_{m\in M}{Y}_{bkm}\le 1, \forall b\in B$$

(5)

$${\displaystyle \sum }_{k\in K}{\displaystyle \sum }_{m\in M}{Z}_{bbkm}\le 1, \forall b\in B$$

(6)

$${\displaystyle \sum }_{a\in B}{Z}_{abkm}=Y_{bkm}, \forall b\in B,\forall k\in K,\forall m\in M$$

(7)

$${\displaystyle \sum }_{\begin{array}{c}b\in B\\ b\ne a\end{array}}{Z}_{abkm}\le Y_{akm}, \forall b\in B, \forall a\in B,\forall k\in K,\forall m\in M$$

(8)

$${\displaystyle \sum }_{b\in B}{Z}_{bbkm}\le 1, \forall k\in K,\forall m\in M$$

(9)

$$T_a·\left(1+D\right)+P_a\le T_b+L·\left(1-{\displaystyle \sum }_{k\in K}{\displaystyle \sum }_{m\in M}{Z}_{abkm}\right), \forall a,b\in B,a\ne b$$

(10)

$$T_b·\left(1+D\right)+P_b\le C_k+L·\left(1-Y_{bkm}\right), \forall b\in B,\forall k\in K,\forall m\in M$$

(11)

$${\displaystyle \sum }_{k\in K}{C}_k+Q·({\displaystyle \sum }_{k\in K}{\displaystyle \sum }_{b\in B}{Z}_{bbkm})-Q\le MS, \forall m\in M$$

(12)

The first decision we have to make is which job should be assigned to which batch by Constraint (1), (2) and (3). One job can be assigned to only one batch, and to be assigned to the same batch, the type of family must be the same and the size must not exceed the capacity. Constraints (1) and (2) confirm that each job is assigned to only one batch and that each batch comprises jobs within the same family. Constraint (3) ensures that the total job sizes assigned to each batch do not exceed the capacity of the batch. After all, once batches have been assigned, the assigned batches should be scheduled to the machine considering the RMA by constraint (4) to constraint (9). There are two ways to consider RMAs, i.e., considering scheduling between all batches, and assuming that RMAs are scheduled only between buckets. The bucket means a set of batches between buckets; hence the RMA does not exist between batches in the same bucket. In the batch scheduling problem, we determine which batch should be assigned to which bucket and the order of the batches in the bucket by constraint (4) to constraint (9). Constraints (4) and (5) ensure that each assigned batch is sequenced in one bucket. Constraints (6) and (8) ensure that there cannot be two first-sequence batches in a given bucket. Constraint (7) ensures that, except for the first batch in each bucket, any batches assigned to a bucket must be immediately preceded. Similarly, constraint (8) ensures that if a batch is assigned to a bucket, at most one batch can be performed immediately afterward. When all batches are scheduled, the makespan can be calculated by constraint (10) to constraint (12). Constraint (10) ensures the precedence relation among assigned batches in a bucket and calculates the starting time for each batch. Constraint (11) calculates the completion time of each bucket by computing the maximum time needed to complete each batch in the bucket. Constraint (12) determines the completion times of the machines needed to calculate the makespan.

4. Multi-Fit Batch Scheduling Algorithm with the RMA Scheduling Rule

The MF algorithm has demonstrated good performance in solving the batch scheduling problem with no deterioration [2]. Unlike previous approaches, however, the scheduling of RMAs must be considered in our problem. One of the ways to do so is to consider the scheduling between all batches. Whenever a batch is scheduled one by one, scheduling the RMA should be considered. To determine whether an RMA should be scheduled, an RMA scheduling rule is used. The RMA scheduling rule is that an RMA is scheduled if the deterioration is longer than the RMA processing time. The MF algorithm using this RMA scheduling rule is proposed in Algorithm 1 and the initial lower and upper bounds on completion time, $L{B}_{ct}$ and $U{B}_{ct}$, can be calculated as

$$L{B}_{ct}=max\{\underset{b\in B}{\max }{P}_b,\lfloor {\displaystyle \sum }_{b\in B}\frac{P_b}{M}\rfloor \},$$

(13)

$$U{B}_{ct}=max\{\underset{b\in B}{\max }{P}_b,\lceil {\displaystyle \sum }_{b\in B}\frac{P_b+Q}{M}\rceil \}.$$

(14)

Algorithm 1. Multi-fit algorithm with the RMA scheduling rule for the batch scheduling problem
Input:	The set of assigned batches $B$.
Output:	The makespan.
Begin:	Sort $B$ in non-increasing order of processing times and obtain a batch set $B^{\prime }.$ Compute $L{B}_{ct}$ and $U{B}_{ct}$, respectively. Let iteration $h=1$. While $\left(h<8\right)$
		$A=\left(U{B}_{ct}+L{B}_{ct}\right)/2$  Let batch index $b=1$.
		While $\left(b\le \left|B^{\prime }\right|\right)$
			Select a batch with a sequence $b$ in $B^{\prime }$. Schedule the batch to the machine according to the first-fit rule if the completion time does not exceed $A$. If (Completion time of machine when the batch is scheduled > Completion time of machine when additionally scheduled RMA precedes the batch) then
				Schedule the RMA to precede the batch.
			$b=b+1.$
		End While  If (All batches in $B^{\prime }$ are scheduled) then
			$U{B}_{ct}=A.$
		Else
			$L{B}_{ct}=A.$
		$h=h+1.$
	End While Output the makespan of the global best solution.

Because the positioning of the RMAs and batches has a significant influence on the makespan, the MF algorithm does not perform well in solving this problem; therefore, we divide the batch scheduling problem into three stages and propose an ACO-based three-stage batch scheduling algorithm for solving it.

5. Ant Colony Optimization-Based Three-Stage Algorithm for Batch Scheduling Problem

The batch loading part of the BLSP referred to the solution of Jia et al. [2], and this section deals with the batch scheduling part. Determining the positioning of RMAs between all batches requires lots of calculations. If the number of buckets is pre-determined, however, the RMA positioning will be fixed between buckets. To reduce the computational complexity, we divide the batch scheduling problem into three stages, in which (1) the number of buckets is determined, (2) batches are assigned to the respective buckets, and (3) the buckets are dispatched to the machines.

5.1. Determining the Number of Buckets

In stage 1, the number of buckets is determined. Determining it fixes the positions of the RMAs between the buckets and it can range from zero to $\left|B\right|-1$ in each machine. To reduce the range of this number, we calculate $L{B}_k$ and $U{B}_k$, the lower and upper bounds, respectively, for the number of buckets during batch scheduling with the dispatching rule. The dispatching rule is the rule that selects the bucket with the shortest tentative completion time. Tentative completion time in each machine means the completion time if the current batch is assigned to the machine. The dispatching rule helps to search for the solution with good quality. A detailed explanation of the dispatching rule is given in Algorithm 2.

Algorithm 2. The dispatching rule
Input:	The processing time of batch $P_b.$ The current completion time of buckets $C_k.$
Output:	Selected bucket $k^{\prime },$ which must be assigned by batch $b$
Begin:	$\mathrm{Let} \mathrm{bucket} \mathrm{index} k=1$
	While $(k \le \left|K\right|$)
		Calculate the tentative completion time using the RMA rule
		$k=k+1$
	End While
	Select the bucket $k^{\prime }$ with the smallest tentative completion time

To calculate the $L{B}_k$ and $U{B}_k$ in subsequent algorithm, we must calculate the minimum and the maximum numbers of batches that can be scheduled on one machine. We assume that the batch processing time can be determined between $P_{min}$ and $P_{max}$ during batch scheduling with the dispatching rule, then the minimum number of batches on a machine is calculated in Algorithm 3. $P_{min}$ and $P_{max}$ are the longest and shortest processing time of the given batches, respectively. In determining the minimum number of batches that can be scheduled to the first machine, the first machine should be assigned a minimum number of batches and the rest of the machines should be assigned as many batches as possible. To assign a small number of batches to the machine under the dispatching rule, the completion time of the machine should be tentatively set to a large value. Since the dispatching rule schedules a batch on a machine with a small tentatively completion time, the allocated batches in the first machine must be given the longest processing time as $P_{max}$. The lower bound on the number of batches input to the first machine is calculated using the dispatching rule in Algorithm 3.

Algorithm 3. Calculating the lower bound on the number of batches input to a machine
Input:	$P_{min};P_{max};$ The number of batches; The number of machines.
Output:	Lower bound on the number of batches input to a machine.
Begin:	Let batch index $b=1$ While$\left(b\le \left|B\right|\right)$
		Select a machine $m^{\prime }$ using the dispatching rule.
		If$\left(m^{\prime }=1\right)$ then
			Assign a batch assuming the processing time is $P_{max}$.
		Else
			Assign a batch assuming the processing time is $P_{min}$.
		$b=b+1$.
	End While Output the number of batches in machine 1

Before calculating the $L{B}_k$ using Algorithm 4, the batch deteriorations must be represented in order. A time-dependent batch deterioration can be represented as the product of $DR$ and $P_i$. In this section, $P_k$ and $D_k$ represent the processing time and deterioration of a batch in sequence $k$, respectively. $D_k$ and $D_{k-1}$ can be represented by the processing time and the deterioration of previous batches as

$$D_k=DR\times \left({\displaystyle \sum }_{i=1}^{k-1}{P}_i+{\displaystyle \sum }_{i=1}^{k-1}{D}_i\right),$$

(15)

$$D_{k-1}=DR\times \left({\displaystyle \sum }_{i=1}^{k-2}{P}_i+{\displaystyle \sum }_{i=1}^{k-2}{D}_i\right).$$

(16)

By subtracting $D_{k-1}$ from $D_k$, $D_k$ can be given in terms of the processing time and the deterioration of the preceding batch as

$$D_k-D_{k-1}=DR\times \left(P_{k-1}+D_{k-1}\right),$$

(17)

and then $D_k$ can be represented as

$$D_k=\left(1+DR\right)\times D_{k-1}+DR\times P_{k-1}.$$

(18)

Using the above equations, the situation in which an RMA must be scheduled is defined through Proposition 1.

Proposition 1.

If the deterioration portion of a batch exceeds the processing time of an RMA, scheduling an RMA before the batch will reduce the completion time of the bucket.

Proof. 

Let $\tilde{B}$ be a set of scheduled batches with no RMAs, which comprises $B_1,B_2\dots ,B_n$. Let $\widetilde{B^{\prime }}$ be a set of the same scheduled batches with one RMA, comprising ${B^{\prime }}_1,{B^{\prime }}_2\dots ,{B^{\prime }}_{l-1}, \mathrm{RMA},{B^{\prime }}_l,\dots$,${B^{\prime }}_n.$ Deterioration of ${B^{\prime }}_i$ can be represented as ${D^{\prime }}_i$. Except for the RMA in $\widetilde{B^{\prime }},$ the processing time of batches in the same sequence as $\tilde{B}$ are equal. Thus, the sums of the deteriorations up to $l-1$ will be equal in both $\tilde{B}$ and $\widetilde{B^{\prime }}$ as follows:

$${\displaystyle \sum }_{i=1}^{l-1}{D}_i={\displaystyle \sum }_{i=1}^{l-1}{D^{\prime }}_i.$$

(19)

As we assume that $D_l$ is greater than $Q$, the sum of ${D^{\prime }}_l$ and $Q$ is less than $D_l$:

$$D_l>{D^{\prime }}_l+Q.$$

(20)

From Equation (19), the deteriorations of $D_{l+1}$ and ${D^{\prime }}_{l+1}$ can be obtained as, respectively,

$$D_{l+1}=\left(1+DR\right)\times D_l+DR\times P_l,$$

(21)

$${D^{\prime }}_{l+1}=\left(1+DR\right)\times {D^{\prime }}_l+DR\times P_l$$

(22)

Because the processing time of each batch is the same and ${D^{\prime }}_l$ = 0, ${D^{\prime }}_{l+1}$ is larger than $D_{l+1},$ repeating the above reasoning, ${D^{\prime }}_{l+2}$ is larger than $D_{l+2}$, ${D^{\prime }}_{l+3}$ is larger than $D_{l+3}$,…, and ${D^{\prime }}_n$ is larger than $D_n$. Adding these deteriorations, we obtain

$${\displaystyle \sum }_{i=l+1}^n{D}_i>{\displaystyle \sum }_{i=l+1}^n{D^{\prime }}_i.$$

(23)

From Equations (20), (21) and (24), $D_1+D_2+\dots +D_l+\dots +D_n$ is larger than ${D^{\prime }}_1+{D^{\prime }}_2+\dots +Q+{D^{\prime }}_l+\dots +{D^{\prime }}_n$. Hence, the deterioration period of $\tilde{B}$ is longer than that of $\widetilde{B^{\prime }}$. □

The above proof gives the condition under which scheduling an RMA reduces the makespan. This can be used to find the minimum number of RMAs. If we set the processing time of each scheduled batch as $P_{min}$, then we can use Proposition 1 to schedule the RMAs starting from the first batch. The lower bound on the number of buckets can then be calculated based on the number of RMAs scheduled per machine using Algorithm 4.

Algorithm 4. Calculating the lower bound on the number of buckets $\left(L{B}_k\right)$
Input:	$P_{max}$; the lower bound on the number of batches in a machine;
Output:	the lower bound on the number of RMAs required by a machine
Begin:	Let batch index $b=1$ Let RMA index $r=0$ While$\left(b\le The lower bound on the number of batches in a machine\right)$
		Calculate a tentative deterioration proportional to the gap between the preceding RMA and $C_{b-1}$  If $(\mathrm{Tentative} \mathrm{deterioration}>\mathrm{Processing} \mathrm{time} \mathrm{of} \mathrm{RMA})$ then
			Schedule RMA after $C_{b-1}$.   $r=r+1$.
		Else
			$C_b=\left(DR+1\right)\times C_{b-1}+P_{max}$
			$b=b+1$.
	End While Output $\left(r+1\right)\times \left|M\right|$

Using these algorithms and propositions, we can obtain the $L{B}_k$ for a machine. In a similar manner, $U{B}_k$ can be calculated. Assuming a batch processing time of $P_{min}$ for the first machine and $P_{max}$ for the remaining machines in Algorithm 3, the upper bound on the number of batches per machine can be applied in Algorithm 4.

5.2. Assigning Batches to Buckets

In stage 2, an ACO algorithm is used to assign batches to each bucket and then a rule is applied to determine the sequencing of the batches. Unlike the problem of assigning batches, the individual buckets do not have size capacities, and the objective function is the makespan. In this case, the batch deterioration has a significant influence on the makespan. Based on the number of buckets determined from stage 1, the batches must be scheduled using a load balancing among buckets.

To assign batches, an ACO algorithm called the min-max ant system (MMAS) is used. MMAS is a constructive meta-heuristic algorithm that builds solutions sequentially, i.e., to solve the batch assignment problem, the batches are assigned one by one. During the batch assignment, we make a probabilistic choice based on pheromone trails and heuristic information until no further batches are available. After assigning the batches, the minimum completion time of each bucket can be calculated using Proposition 2 (Section 5.2.5).

5.2.1. Pheromone Trails

The desirability of unscheduled batches can be calculated using a pheromone trail that gives the relationship between assigned batches. There is a pheromone trail value between all batches and a higher pheromone trail value indicates a higher probability that two batches will be allocated to the same bucket. Defining ${\tau }_{il}\left(t\right)$ as the pheromone trail between $B_i$ and $B_l$ in iteration $t$, ${\tau }_{il}\left(1\right)$ is initialized as

$${\tau }_{il}\left(1\right)={\left(\left(1-\rho \right)\times L{B}_{\tau }\right)}^{-1},$$

(24)

where $\rho$ means evaporation rate and $L{B}_{\tau }$ is the lower bound given by

$$L{B}_{\tau }={\displaystyle \sum }_{bϵB}\frac{P_b}{\left|M\right|}$$

(25)

The desirability of a given batch is calculated using the pheromone trails between the given batch and batches, $B_k\left(t\right)$, that have already been assigned to the selected bucket in iteration $t$. The desirability of assigning batch $i$into the current bucket $k$ during iteration $t$ can be represented as

$${\theta }_{ik}\left(t\right)={\displaystyle \sum }_{B_lϵB_k\left(t\right)}\frac{{\tau }_{il}\left(t\right)}{\left|B_k\left(t\right)\right|}.$$

(26)

5.2.2. Heuristic Information

The load balance of each bucket is important because the deterioration increases rapidly when many batches are allocated per bucket. The heuristic information between batch $i$ and bucket $k$ during iteration $t$ can be represented by

$${\eta }_{ik}\left(t\right)=\frac{1}{\left|\max \left\{c_k\left(t\right)\right\}-c_{k+i}\left(t\right)\right|+p_{max}},$$

(27)

where $c_{k+i}\left(h\right)$ is the completion time of bucket $k$ calculated after assigning batch $i$ to it during iteration $t$.

5.2.3. Forming the Buckets

The pheromone trails and heuristic information are used together to calculate the probability that each batch can be formed. The probability that batch $i$ can be assigned to bucket $k$ during iteration $t$ is expressed as $P_{ik}\left(t\right)$, which represents the desirability of the batches assigned to bucket $k$ as follows:

$$P_{ik}\left(t\right)=\{\begin{array}{l}\begin{array}{cc}\frac{{\theta }_{ik}\left(t\right)\times {\eta }_{ik}{\left(t\right)}^{\beta }}{{{\displaystyle \sum }}_{B_lϵU_k}{\theta }_{ik}\left(t\right)\times {\eta }_{ik}{\left(t\right)}^{\beta }}, if B_iϵU_k& \end{array}\\ \begin{array}{cc} 0, & otherwise\end{array}\end{array} ,$$

(28)

where $\beta$ represents the relative importance of heuristic information and $U_k$ represents the set of unscheduled batches.

5.2.4. Updating the Pheromone Trails

In MMAS, the use of pheromone trails leads to solutions based on experience. Thus, updating the pheromone trails has a significant effect on the performance of the algorithm. Global- and iteration-best solutions are used to update the pheromone trails. Representing the frequency with which jobs i and l are placed in the same batch as $m_{il}$, the pheromone trail updating process can be defined as follows:

$${\tau }_{il}\left(t+1\right)=\left(1-\rho \right)\times {\tau }_{il}\left(t\right)+m_{il}\times \Delta {\tau }_{il}\left(t\right),$$

(29)

$$\Delta {\tau }_{il}\left(t\right)=\frac{Q}{Makespa{n}_{iteration best}\left(t\right)}.$$

(30)

The solutions corresponding to pheromone trails that are too extreme cannot be found easily. Under MMAS, pheromone trail lower and upper bounds are defined as ${\tau }_{min} \mathrm{and} {\tau }_{max}$, respectively, and when the trails are updated, they are modulated as follows:

$${\tau }_{il}\left(t+1\right)=\{\begin{array}{c}\begin{array}{cc}{\tau }_{min},& {\tau }_{il}\left(t+1\right)<{\tau }_{min}\end{array}\\ \begin{array}{cc}{\tau }_{il}\left(t+1\right),& {\tau }_{min}\le {\tau }_{il}\left(t+1\right)\le {\tau }_{max}\end{array}\\ \begin{array}{cc}{\tau }_{max},& {\tau }_{il}\left(t+1\right)>{\tau }_{max}\end{array}\end{array},$$

(31)

$${\tau }_{max}={\left(\left(1-\rho \right)\times Makespa{n}_{global best}\left(t\right)\right)}^{-1},$$

(32)

$${\tau }_{min}=\frac{{\tau }_{max}\times \left(1-\sqrt[\left|B\right|]{0.05}\right)}{(\frac{\left|B\right|}{2}-1)\times \sqrt[\left|B\right|]{0.05}}.$$

(33)

5.2.5. Rule-Based Batch Sequencing

After assigning the batches to the buckets, the completion time of each bucket can be minimized using Proposition 2.

Proposition 2.

Let B be the set of batches in one bucket. Then, scheduling B in ascending order of processing time minimizes the completion time of each bucket.

Proof. 

As the starting time of the first batch is zero, we have

$$D_1=0.$$

(34)

By applying the previous deterioration to Equation (16), we obtain

$$\begin{gathered} D_2=DR\times P_1+\left(1+DR\right)\times D_1 \\ \hspace{1em}=DR\times P_1, \end{gathered}$$

(35)

$$\begin{gathered} D_3=DR\times P_2+\left(1+DR\right)\times D_2 \\ \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}=DR\times P_2+\left(1+DR\right)\times DR\times P_1, \end{gathered}$$

(36)

$$\begin{gathered} D_4=DR\times P_3+\left(1+DR\right)\times D_3 \\ =DR\times P_3+\left(1+DR\right)\times DR\times P_2+{\left(1+DR\right)}^2\times DR\times P_1, \end{gathered}$$

(37)

$$\begin{gathered} D_n=DR\times P_{n-1}+\left(1+DR\right)\times DR\times P_{n-2}+\dots +{\left(1+DR\right)}^{n-2}\times DR\times P_1 \\ ={\displaystyle \sum }_{i=1}^{n-1}{\left(1+DR\right)}^{n-1-i}\times DR\times P_i. \end{gathered}$$

(38)

By using the generalized formula above,$D_n$, the sum of the deteriorations, can be represented as

$$\begin{gathered} {\displaystyle \sum }_{i=1}^n{D}_i=D_2+D_3+D_4+\dots +D_n \\ =\left[DR\times P_1\right]+\left[DR\times P_2+\left(1+DR\right)\times DR\times P_1\right]+\dots \\ +{\displaystyle \sum }_{i=1}^{n-1}{\left(1+DR\right)}^{n-1-i}\times DR\times P_i \\ =\left[1+\left(1+DR\right)+{\left(1+DR\right)}^2+\dots +{\left(1+DR\right)}^{n-2}\right]\times DR\times P_1\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em} \\ +\left[1+\left(1+DR\right)+{\left(1+DR\right)}^2+\dots +{\left(1+DR\right)}^{n-3}\right]\times DR\times P_2 \\ +\dots +\left[1+\left(1+DR\right)\right]\times DR\times P_{n-2}+DR\times P_{n-1} \\ ={\displaystyle \sum }_{i=0}^{n-2}{\left(1+DR\right)}^i\times DR\times P_1+{\displaystyle \sum }_{i=0}^{n-3}{\left(1+DR\right)}^i\times DR\times P_2+\dots \\ +{\displaystyle \sum }_{i=0}^1{\left(1+DR\right)}^i\times DR\times P_{n-2}+DR\times P_{n-1} \end{gathered}$$

(39)

In Equation (39), $P_k,$ the processing time of a batch in sequence k, is multiplied by $\sum _{i=0}^{n-1-k}{\left(1+DR\right)}^i$. $\sum _{i=0}^{n-1-k}{\left(1+DR\right)}^i$ increases as the sequence become slower. As the batch we have to schedule is fixed, the total deterioration can be minimized by scheduling the batches with longer processing times in faster sequences. □

Using the above proof, we can obtain the shortest completion for each bucket. After all, once the shortest completion time of buckets is calculated, buckets must be scheduled to the machines.

5.3. Scheduling Buckets to Machines

In stage 3, buckets are scheduled to machines using the dispatching rule. Tentative timings are used to determine the completion time at which the current batch is scheduled using the RMA rule as Algorithm 2 (Section 5.1). To balance the load among the machines, therefore, the current batch should be scheduled to the machine with the shortest tentative completion time.

5.4. Overall Algorithm

Using the algorithms defined for the three stages in the preceding sections, the ACO-based three-stage algorithm can be defined as Algorithm 5.

Algorithm 5. ACO-based three-stage algorithm for batch scheduling problem
Input:	The set of assigned batches.
Output:	The makespan
Begin:	Compute the $L{B}_K$ and $U{B}_K$, respectively. Compute the $L{B}_{\tau }$ and ${\tau }_{il}\left(1\right)$, respectively. Initialize pheromone trails. Let index of iteration $t=1$ While$\left(t\le t_{max}\right)$
		Let index of iteration $a=1$  While $\left(a\le a_{max}\right)$
			Select the number of buckets between $L{B}_K$ and $U{B}_K$.   While$\left(U_B\ne \varnothing \right)$
				Select bucket $k$ using the dispatching rule.    If (Bucket $k$ has no assigned batches) then
					Select a batch $b$ randomly.     Assigned the batch $b$ to bucket $k$.
				Else
					Select a batch $b$ using pheromone trail and heuristic information.     Assign the batch $b$ to bucket $k$.
			End While   Minimize the completion time of the buckets by sequencing the batches.
			Dispatch the buckets to machines.   Calculate the makespan.   Update the iteration and global best solutions   a = a +1
		End While  Compute ${\tau }_{max}$ and ${\tau }_{min}$.  Update $\Delta {\tau }_{il}\left(h\right)$ and ${\tau }_{il}\left(h\right)$.  $t=t+1$
	End While Output the makespan of the global best solution

The algorithm divides the batch scheduling problem into three-stage. In stage 1, $L{B}_K$ and $U{B}_K$ are used for fast convergence. To demonstrate that $L{B}_K$ and $U{B}_K$ affect the performance of the algorithm, an ACO with no bound (ACO_NB) is proposed.

6. Genetic Algorithm- and Particle Swarm Optimization-Based Batch Scheduling Algorithms

In this section, two population-based meta-heuristic algorithms for solving the batch scheduling problem—a GA-based three-stage algorithm and a PSO-based three-stage algorithm—are proposed. Unlike constructively generated ACO solutions, these two algorithms search solutions to develop an encoded solution that is represented using a random key. The two proposed population-based meta-heuristic algorithms (Algorithms 6 and 7) are defined as follows.

Algorithm 6. GA-based three-stage algorithm for batch scheduling problem
Input:	The set of assigned batches from the batch loading problem
Output:	The makespan
Begin:	Compute the $L{B}_K$ and $U{B}_K$ respectively. Set the initial decoded solution. Let index of generation $t=1$ While $\left(t\le t_{max}\right)$
		Let index of chromosome $c=1$  While$\left(c\le c_{max}\right)$
			Select the number of buckets between $L{B}_K$ and $U{B}_K$.   Let random value $r_1=Unif\left(0,1\right).$   Let random value $r_2=Unif\left(0,1\right).$   Select two parent chromosomes, c, and c + 1, in sequence.
			If$(r_1$ < Crossover rate) then
				Do one-point crossover operations.
			If$(r_2$ < Mutation rate) then
				Do swap mutation operation.
			Decode the offspring chromosome $c$ and calculate the makespan.
			$c=c+1$
		End While  Reproduce the next generation parents from current-generation offspring  Update the iteration and global best solutions  $t=t+1$
	End While Output the makespan of the global best solution

Algorithm 7. PSO-based three-stage algorithm for batch scheduling problem
Input:	The set of assigned batches from the batch loading problem
Output:	The makespan
Begin:	Compute the $L{B}_K$$\mathrm{and} U{B}_K$, respectively. Set the initial solution. Let index of iteration $t=1$ While $\left(t\le t_{max}\right)$
		Let index of particle $p=1$  While $\left(p\le p_{max}\right)$
			Select the number of buckets between $L{B}_K$ and $U{B}_K$   Update the factor velocity   Update the factor position   Decode the particle $p.$   $p=p+1$
		End While  Decode the particles.  Update the iteration and global best solutions;  $t=t+1$
	End While Output the makespan of global best solution

7. Computational Experiments

To evaluate the performance of the ACO-based three-stage batch scheduling algorithm, extensive computational experiments were conducted using the solutions obtained by an ACO-based batch loading algorithm presented in Jia et al. [2]. In the first experiment, the absolute differences between the optimal solutions produced by CPLEX and the solutions produced by the ACO, ACO_NB, GA, PSO, and MF algorithms were compared for small-sized problems. In the second experiment, ACO and ACO_NB were compared to validate the $L{B}_k$ and $U{B}_k$ used in the proposed algorithm. In the third experiment, the relative differences between the solutions obtained by ACO and other meta-heuristic algorithms, GA and PSO, were compared for large-sized problems. All meta-heuristic algorithms are implemented in C#.

7.1. Problem Parameter Settings

In the experiments, BLSPs corresponding to the batch scheduling problem with varying complexity were randomly generated. The complexity of a batch scheduling problem with batch deterioration and RMAs is highly dependent on $\left|M\right|, \left|K\right|, \left|B\right|$, and $\left|F\right|$ [24,25,26]. Therefore, four problem parameters, $\left|M\right|, \left|K\right|/\left|M\right|,$ $\left|B\right|/\left|K\right|$, and $\left|F\right|$, were used to generate problems with different complexities, where $\left|K\right|/\left|M\right|$ is the average number of buckets per machine and $\left|B\right|/\left|K\right|$ is the average number of batches per bucket. The machines and families were generated as $\left|M\right|$ and $\left|F\right|$, respectively, and the batches were generated using $\left|M\right|\times \left|K\right|/\left|M\right|\times \left|B\right|/\left|K\right|$. To assign $\left|B\right|/\left|K\right|$ batches to each bucket, $Q$ was generated as $\frac{Q_{min}+Q_{max}}{2}$. As $Q$ decreases, RMAs will occur more frequently and the number of batches per bucket will be reduced. To assign $\left|B\right|/\left|K\right|$ batches to each bucket, $Q$ should be adjusted so that it falls between $Q_{min}$ and $Q_{max}$, where those values mean the minimum and maximum sums of deteriorations that occur when there are $\left|B\right|/\left|K\right|$batches, respectively. The sums can be calculated by assuming the processing time of batches as $P_{min}$ and $P_{max}$, respectively. $P_f$ is randomly generated in $\left[P_{min},P_{max}\right]$ with $P_{min}$ and $P_{max}$ fixed at 40 and 60, respectively. $D$ and $\overline{S_j}$ are fixed at 0.2 and 0.4, respectively. Using the four problem parameters defined above, two groups of problems (small-sized problems with fewer than 10 batches and large-sized problems with more than 10 batches) were randomly generated.

7.2. Algorithm Parameter Settings

To find the major parameters that affect the performance of algorithms, the Taguchi method was applied in this section. The method can be used to conduct experiments more rapidly than a full factorial experiment because the method can carry out experiments using fewer scenarios. Using the Taguchi method, five control factors for the ACO-based batch scheduling algorithm were analyzed. The levels of each control factor are listed in

Table 2. Factors levels for MMAS based three-stage batch scheduling algorithm.

Factor 1	Factor 2	Factor 3	Factor 4	Factor 5
Number of Ants	Evaporation Rate	Iteration Limit	Updating Parameter	Relative Importance of Heuristic Information
$\mathit{A}\mathit{n}{\mathit{t}}_{\mathit{m}\mathit{a}\mathit{x}} \left(\mathit{A}\right)$	$\mathit{\rho } \left(\mathit{B}\right)$	$\mathit{G} \left(\mathit{C}\right)$	$\mathit{\delta } \left(\mathit{D}\right)$	$\mathit{\beta } \left(\mathit{E}\right)$
A$\left(1\right)\text{:} 1.0\times \left|\mathrm{B}\right|$	B (1): 0.5	C (1): 10	D$\left(1\right)\text{:} 0.5\times L{B}_{\tau }$	E (1): 4
A$\left(2\right)\text{:} 1.5\times \left|\mathrm{B}\right|$	B (2): 0.6	C (2): 30	D$\left(2\right)\text{:} 1.0\times L{B}_{\tau }$	E (2): 6
A$\left(3\right)\text{:} 2.0\times \left|\mathrm{B}\right|$	B (3): 0.7	C (3): 50	D$\left(3\right)\text{:} 1.5\times L{B}_{\tau }$	E (3): 8
A$\left(4\right)\text{:} 2.5\times \left|\mathrm{B}\right|$	B (4): 0.8	C (4): 70	D$\left(4\right)\text{:} 2.0\times L{B}_{\tau }$	E (4): 10
A$\left(5\right)\text{:} 3.0\times \left|\mathrm{B}\right|$	B (5): 0.9	C (5): 90	D$\left(5\right)\text{:} 2.5\times L{B}_{\tau }$	E (5): 12

. The first factor, A, is the number of ants ($An{t}_{max})$. When a problem becomes more complicated, more ants are needed to solve it. Therefore, we set $An{t}_{max}$ to increase according to the complexity of the problem; as the complexity was proportional to $\left|B\right|$, $An{t}_{max}$ is varied over the range $\left\{1.0\times \left|B\right|,1.5\times \left|B\right|,2.0\times \left|B\right|,2.5\times \left|B\right|,3.0\times \left|B\right|\right\}$. The second factor, B, is the evaporation rate ($\rho )$, which is used to calculate the initial pheromone trails and update later trails. Because the evaporation must be a real value between [0, 1] and values lower than 0.5 are relatively insignificant to the parameter settings, $\rho$ is varied over $\left\{0.5, 0.6, 0.7, 0.8, 0.9\right\}.$ The third factor, C, is the iteration limit ($G$), which indicates how many iterations are needed to reset a pheromone trail in the ACO. If there is no change in the global solution over $G$ iterations, the trail must be reset and, therefore, smaller values of $G$ mean that trails must be more often reset. In this paper, $G$ is varied over $\left\{10, 30, 50, 70, 90\right\}$. The fourth factor, D, is the updating parameter ($\delta$), which is used to update the pheromone trails. In a previous paper [2], $\delta$ performs better when it is based on $L{B}_{\tau }$ rather than a constant parameter. Therefore, $\delta$ is varied over $\left\{0.5\times L{B}_{\tau },1.0\times L{B}_{\tau },1.5\times L{B}_{\tau }, 2.0\times L{B}_{\tau }, 2.5\times L{B}_{\tau }\right\}$. The final factor, E, is the relative importance of heuristic information ($\beta$), which is used in calculating desirability. $\beta$ is varied over {4, 6, 8, 10, 12}.

Using the control factors described above, an experiment involving 25 scenarios was conducted. The number of scenarios was chosen so that the best-fit design for five factors with five levels.

Table 3. Factor levels of array $L_{25}$.

Scenario No.	Factor Levels
	A	B	C	D	E
1	A (1)	B (1)	C (1)	D (1)	E (1)
2	A (1)	B (2)	C (2)	D (2)	E (2)
3	A (1)	B (3)	C (3)	D (3)	E (3)
4	A (1)	B (4)	C (4)	D (4)	E (4)
5	A (1)	B (5)	C (5)	D (5)	E (5)
6	A (2)	B (1)	C (2)	D (3)	E (4)
7	A (2)	B (2)	C (3)	D (4)	E (5)
8	A (2)	B (3)	C (4)	D (5)	E (1)
9	A (2)	B (4)	C (5)	D (1)	E (2)
10	A (2)	B (5)	C (1)	D (2)	E (3)
11	A (3)	B (1)	C (3)	D (5)	E (2)
12	A (3)	B (2)	C (4)	D (1)	E (3)
13	A (3)	B (3)	C (5)	D (2)	E (4)
14	A (3)	B (4)	C (1)	D (3)	E (5)
15	A (3)	B (5)	C (2)	D (4)	E (1)
16	A (4)	B (1)	C (4)	D (2)	E (5)
17	A (4)	B (2)	C (5)	D (3)	E (1)
18	A (4)	B (3)	C (1)	D (4)	E (2)
19	A (4)	B (4)	C (2)	D (5)	E (3)
20	A (4)	B (5)	C (3)	D (1)	E (4)
21	A (5)	B (1)	C (5)	D (4)	E (3)
22	A (5)	B (2)	C (1)	D (5)	E (4)
23	A (5)	B (3)	C (2)	D (1)	E (5)
24	A (5)	B (4)	C (3)	D (2)	E (1)
25	A (5)	B (5)	C (4)	D (3)	E (2)

shows $L_{25}\left(5^5\right)$ orthogonal array. In each scenario, 12 test problems were solved three times each for a total of 900 runs for the experiment. Under the Taguchi method, variation in performance is measured using the mean signal-to-noise ($\mathrm{S}/\mathrm{N}$) ratio. Generally, the S/N ratio is expressed as $-10\log {\left(objective function\right)}^2$, but this formula was not directly applicable to this case because the experimental tests had varying objective functions and sizes. Therefore, we adopted the relative percentage deviation ($RPD$) of the objective function as follows:

$$RPD(\%)=\frac{OB{J}_{sol}-OB{J}_{best}}{OB{J}_{best}}\times 100,$$

(40)

where $OB{J}_{sol}$ is the objective function given by an algorithm and $OB{J}_{best}$ is the objective function of the best solution. Using the $RPD$, the S/N ratio can be computed as follows:

$$S/N rati{o}_k=-10\log \left(\frac{1}{12}{\displaystyle \sum }_{i=1}^4{\displaystyle \sum }_{j=1}^{3}RP{D}_{ijk}^2\right) \forall k\in 1,2,\dots ,2,D h the 12525.$$

(41)

To analyze statistically significant factors, analysis of variance (ANOVA) was performed using the S/N ratio data. The factor D, which had the smallest sum squares (SS), was selected as the error term. After pooling this factor, ANOVA was performed, the results of which are listed in

Table 4. Analysis of variance for the S/N ratio.

Factor	SS	df	V	${\mathit{F}}_0$	$\mathit{p}$-Value
A	0.2593	4	0.0648	2.6456	0.1844
B	1.4578	4	0.3644	14.8721	0.0114
C	0.1809	4	0.0452	1.8459	0.2836
D(Error)	0.0980		0.0245
E	0.2144	4	0.5360	21.8756	0.0055
Total	4.1404	20

. Factors B and E exceeded the significance level, indicating that they had a large effect on the algorithm performance. The levels of these factors with the highest S/N ratios—B (5) and E (5)—were identified based on an examination of the S/N results in Figure 2. The factors and levels A (4), C (3), and D (1) had the least effect on the performance, as reflected by the $RPD$ values shown in Figure 3. From these results, we identified an optimal combination of A (4), B (5), C (3), D (1), and E (5).

The GA and PSO parameters were also set using the Taguchi method. For the GA, the parameters $\mathrm{population} \mathrm{size}=3\times \left|B\right|$, $\mathrm{crossover} \mathrm{rate}=0.9$, and $\mathrm{mutation} \mathrm{rate}=0.2$ were set. For the PSO, the parameters $\mathrm{population} \mathrm{size}=3\times \left|B\right|$, $\mathrm{inertia} \mathrm{weight}=0.9$, social cognitive 1$=1.0$, and social cognitive 2 $=0.5$ were set.

7.3. Experimental Results

Three experiments were conducted to validate the ACO-based three-stage batch scheduling algorithm. All of the experiments were carried out on a PC with an Intel Core i7-7700 CPU running at 3.6 GHz with 16 GB of memory. In each experiment, the batch scheduling problem with the same problem instance was solved using the same given solution obtained from the ACO-based batch loading algorithm. All the algorithms were run until iteration reached 1500 or performance did not increase for 100 iterations.

In the first experiment, the ACO, ACO_NB, GA, PSO, and MF algorithms were validated by comparing their solutions with optimal solutions produced by ILOG CPLEX. The experiment was conducted on eight small-sized problems with 30 replications, with the results summarized in

Table 5. The test results of small-sized problems.

$\left|\mathit{M}\right|$	$\frac{\left|\mathit{K}\right|}{\left|\mathit{M}\right|}$	$\frac{\left|\mathit{B}\right|}{\left|\mathit{K}\right|}$	CPLEX		ACO		ACO_NB		GA		PSO		MF
			Opt.	Time	RPD	Time	RPD	Time	RPD	Time	RPD	Time	RPD
1	1	1	65	0.02	0.00	0.01	0.00	0.01	0.00	0.02	0.00	0.01	0.00
1	1	2	217.2	0.24	0.00	0.01	0.00	0.01	0.00	0.02	0.00	0.01	0.73
1	2	1	149	0.14	0.00	0.01	0.00	0.01	0.00	0.02	0.00	0.01	0.00
1	2	2	423	5.93	0.00	0.01	0.00	0.01	0.00	0.02	0.00	0.01	3.30
2	1	1	136	0.12	0.00	0.01	0.00	0.01	0.00	0.03	0.00	0.01	0.73
2	1	2	246	8.85	0.00	0.01	0.00	0.01	0.00	0.03	0.00	0.01	9.75
2	2	1	219	9.412	0.00	0.01	0.00	0.01	0.00	0.03	0.00	0.01	2.73
2	2	2	N/A	7200+	N/A		N/A		N/A		N/A		N/A
Avg.					0.00	0.01	0.00	0.01	0.00	0.02	0.00	0.01	2.16

. To compare the performance of the respective algorithms, the RPD values were calculated using the optimal solution obtained from CPLEX. In the table, “N/A” denotes cases in which CPLEX was unable to obtain the optimal solution before 7200s. The problems with more than eight batches could not be solved by ILOG CPLEX. For problems with fewer than eight batches, the average RPD values of ACO, ACO_NB, GA, and PSO were zero whereas the average for the MF algorithm was 2.16. The result indicates that the meta-heuristic algorithms provided optimal solutions to small-sized problems. However, the MF algorithm performed poorly even on small-sized problems.

In the second experiment, the relative performances of ACO and ACO_NB in solving large-sized problems were compared to validate the $L{B}_k$ and $U{B}_k$ values obtained in stage 1. The experiment was conducted on 36 large-sized problems with 30 replications. The results are summarized in

Table 6. The test results of large-sized problems.

$\left|\mathit{F}\right|$	$\left|\mathit{M}\right|$	$\frac{\left|\mathit{K}\right|}{\left|\mathit{M}\right|}$	$\frac{\left|\mathit{B}\right|}{\left|\mathit{K}\right|}$	Best	ACO		ACO_NB		GA		PSO
					RPD	Time	RPD	Time	RPD	Time	RPD	Time
4	3	3	3	964.68	0.08	0.50	0.08	0.48	0.59	0.22	2.16	0.15
4	3	3	4	1458.448	0.00	1.09	0.00	1.10	0.96	0.42	1.10	0.28
4	3	4	3	1402.04	0.29	1.19	0.29	1.17	0.07	0.55	1.55	0.39
4	3	4	4	1776.96	1.36	2.62	1.39	3.37	0.15	1.03	1.89	0.56
4	4	3	3	859.28	0.00	0.65	0.00	0.69	0.43	0.38	0.47	0.27
4	4	3	4	1498.52	0.22	3.49	0.30	3.23	0.48	1.21	2.73	0.60
4	4	4	3	1218.04	0.00	1.40	0.00	1.48	1.54	0.53	1.28	0.54
4	4	4	4	1806.64	0.51	7.96	0.53	6.81	0.14	2.34	1.80	1.28
4	5	3	3	953.6	0.00	1.00	0.00	1.14	0.39	0.72	0.39	0.46
4	5	3	4	1383.896	0.00	3.50	0.00	3.19	1.08	1.01	1.07	0.92
4	5	4	3	1324.6	0.64	3.72	0.67	3.74	0.10	1.97	1.96	0.94
4	5	4	4	1924.68	0.00	7.12	0.00	6.75	0.73	2.62	0.82	1.76
5	3	3	3	944	0.00	0.28	0.00	0.31	0.79	0.22	0.91	0.14
5	3	3	4	1347.904	0.08	1.28	0.09	1.30	0.32	0.53	2.29	0.25
5	3	4	3	1168.44	0.12	1.12	0.14	1.08	0.27	0.54	1.58	0.29
5	3	4	4	1922.64	0.00	3.32	0.01	2.91	0.86	1.17	1.03	0.58
5	4	3	3	966	0.00	0.64	0.00	0.78	0.56	0.43	0.59	0.27
5	4	3	4	1583.16	0.61	3.52	0.60	3.26	0.27	1.31	1.27	0.53
5	4	4	3	1206.48	0.00	1.86	0.00	1.94	1.42	0.56	1.13	0.55
5	4	4	4	1733.96	0.00	5.08	0.00	4.52	1.54	1.68	1.67	1.03
5	5	3	3	977.6	0.73	1.80	0.78	2.18	0.22	0.89	3.51	0.52
5	5	3	4	1360.248	0.06	5.22	0.10	5.39	1.05	2.51	2.31	0.97
5	5	4	3	1152	0.00	2.76	0.00	2.98	1.06	0.89	0.73	0.90
5	5	4	4	1827.6	0.11	15.07	0.16	14.06	0.10	4.65	1.64	2.09
6	3	3	3	990.92	0.72	0.55	0.73	0.57	0.36	0.31	2.20	0.16
6	3	3	4	1449.64	0.06	1.47	0.07	1.52	0.31	0.66	0.92	0.27
6	3	4	3	1300.6	0.00	0.62	0.00	0.71	0.45	0.39	0.52	0.28
6	3	4	4	1932.56	0.02	3.25	0.02	2.99	0.80	1.43	0.95	0.53
6	4	3	3	915.04	0.08	1.14	0.11	0.95	0.25	0.54	1.72	0.26
6	4	3	4	1382.528	0.26	3.27	0.27	3.73	1.48	1.35	2.98	0.52
6	4	4	3	1286.6	0.00	1.52	0.00	1.67	0.35	0.93	0.41	0.66
6	4	4	4	1959.84	0.06	8.84	0.08	7.50	1.33	3.27	2.65	1.33
6	5	3	3	962.8	0.00	1.35	0.01	1.66	0.64	0.73	0.67	0.47
6	5	3	4	1437.832	0.32	5.16	0.38	5.81	1.74	1.83	2.33	0.90
6	5	4	3	1306.64	0.06	4.24	0.06	4.35	0.21	1.90	1.19	1.02
6	5	4	4	1919.76	0.08	13.48	0.09	12.55	0.40	4.64	1.88	1.91
Avg.					0.18	3.36	0.19	3.27	0.65	1.29	1.51	0.68

. To validate the results, a paired $t$-test $\left(\alpha =0.05\right)$ was conducted to statistically evaluate the differences in performance between ACO and ACO_NB. The results of the paired $t$-tests are listed in

Table 7. Paired $t$-test between ACO_NB-ACO.

.	$\mathit{N}$	Mean	St. Dev.	St. e. Mean	Lower	$\mathit{t}$-Value	$\mathit{p}$-Value
Paired $t$-test for ACO_NB-ACO					0.16	3.75	0.001
ACO_NB	36	1386.2	337.0	56.2
ACO	36	1385.8	336.8	56.1
Difference	36	0.35	0.56	0.09

. ACO significantly outperformed ACO_NB because it is seen from the table that all of the p-values are less than 0.05. The differences between the results obtained by ACO and ACO_NB are dependent upon the calculated values of $L{B}_k$ and $U{B}_k$, respectively. To demonstrate that the use of $L{B}_k$ and $U{B}_k$ quickly converge $\left|K\right|$ and the makespan, the performance improvements obtained at different levels of $\left|F\right|$ are shown in Figure 4, in which the average makespan and $\left|K\right|$ are presented for five iterations of ACO and ACO_NB, respectively. The results demonstrate that the makespan converges faster under ACO than under ACO_NB because the $\left|K\right|$ from ACO rapidly converges to an optimal value between $L{B}_k$ and $U{B}_k$.

In the third experiment, the results produced by ACO, GA, and PSO in solving large-sized problems were relatively compared. The experiment was conducted on the same large-sized problems used in the second experiment, and the average RPD values obtained by ACO, GA, and PSO were 0.18, 0.65, and1.51 in Table 6, respectively. These results indicate that ACO obtained better RPD values than GA and PSO. However, it took more computing time for it to obtain the best feasible solutions because, despite the fact that the algorithm generally converges quickly to the best solution, ACO requires a considerable amount of computing time to execute the proposed MMAS algorithm within one pheromone search iteration. In our batch scheduling problem, the calculation of desirability under the MMAS algorithm exponentially increases the algorithm running time when $\left|B\right|/\left|K\right|$ increases.

To determine whether the application of GA or PSO could result in an improved solution, an additional experiment in which the computing times for GA and PSO were extended until the ACO solution converged was conducted. Using the input data from the large-sized problem, GA and PSO were retested, with the results summarized in

Table 8. The mean RPD of large-sized problems with an equal running time.

$\left|\mathit{F}\right|$	ACO	GA	PSO
4	0.35	0.51	1.36
5	0.40	0.42	1.66
6	0.23	0.51	1.63
Average	0.33	0.48	1.55

. The average RPD values produced by ACO, GA, and PSO were 0.33, 0.48, 1.55, respectively, indicating that applying GA and PSO could not improve the ACO solution quality even if they were provided with additional computing time. To statistically validate the results in Table 8, $t$-tests were conducted to verify the differences in performance between ACO and the other meta-heuristics. As the same instance of the problem was used in each case, a paired $t$-test $\left(\alpha =0.05\right)$ was conducted. The results of paired $t$-testing between PSO-ACO and GA-ACO are shown in

Table 9. Paired t-test between PSO-ACO and GA-ACO.

	N	Mean	St. Dev.	St. e. Mean	Lower	t-Value	$\mathit{p}$-Value
Paired t-test for PSO-ACO					15.577	18.04	0.000
PSO	180	1403.6	339.0	25.3
ACO	180	1386.1	332.6	24.8
Difference	180	17.490	13.009	0.970
Paired t-test for GA-ACO					0.209	2.21	0.029
GA	180	1388.1	332.5	24.8
ACO	180	1386.1	332.6	24.8
Difference	180	1.987	12.090	0.901

. In both cases, the $p$-values are less than 0.05, suggesting that ACO performed significantly better than either PSO or GA.

8. Conclusions

In this paper, we considered a BLSP for parallel BPMs with batch deterioration and applied RMAs. In the proposed BLSP, the time-dependent deterioration occurring during batch processing and the need to schedule RMAs between all batches increases the complexity of the batch scheduling problem. As an MF algorithm from a previous paper could not find an optimal solution to even small-sized versions of this problem, we solved the batch scheduling problem by dividing it into three stages, determining the number of buckets, assigning the batches to buckets, and scheduling the buckets to machines. Determining the number of buckets fixes the position of the RMAs and reduces the complexity of the batch scheduling problem. The lower and upper bounds of the number of buckets, $L{B}_k$ and $U{B}_k$, respectively, are calculated to improve the solution performance and increase the speed of convergence. To schedule batches into a fixed number of buckets, an ACO is used to assign batches and a derived rule is used to sequence the batches. In solving the batch assignment problem, the ACO outperformed both GA and PSO. Finally, a dispatching rule is used to schedule buckets to the machines. Using this three-stage ACO-based batch scheduling algorithm, the proposed method finds optimal solutions for small-sized problems and provides better-quality solutions for large-size problems than can be obtained using other meta-heuristics.