SEQUENCING SINGLE MACHINE MULTIPLE-CLASS CUSTOMER ORDER JOBS USING HEURISTICS AND IMPROVED SIMULATED ANNEALING ALGORITHMS

. The multiple job class scheduling problem arises in contexts where a group of jobs belong to multiple classes and in which if all jobs in the same class are operated together, extra setup times would not be needed. On the other hand, the customer order scheduling problem focuses on finishing all jobs from the same order at the same time in order to reduce shipping costs. However, works on customer orders coupled with class setup times do not appear often in the literature. Hence we address here a bicriteria single machine customer order scheduling problem together with multiple job classes. The optimality criterion minimizes a linear combination of the sum of the ranges and sum of tardiness of all customer orders. In light of the high complexity of the concerned problem, we propose a lower bound formula and a property to be used in a branch-and-bound method for optimal solutions. To find approximate solutions, we then propose four heuristics together with a local search method, four cloudy theoretical simulated annealing and a cloudy theoretical simulated annealing hyperheuristic along with five low-level heuristics. The simulation results of the proposed heuristics and algorithms are analyzed.


Introduction
The sequencing (or scheduling) issue involving setup times or setup costs in many manufacturing and service environments is continually receiving more attention in the research community.Allahverdi and Soroush [5] pointed out that the efficacy of reducing setup times or costs contributes to reliable products and services being delivered to customers on time.For more details on applications involving setup times or costs, readers may refer to the five survey studies by Allahverdi et al. [6,7], Allahverdi [3], Cheng et al. [15] and Yang and Liao [62].
Setup times/costs in sequencing models have received wide discussion in the past thirty years because they play a key role in the on-time delivery of customer orders and have crucial impacts on cycle times [30,63,64].In the multiple job class case, the jobs are divided into several classes, and a setup time is needed when a machine switches from one class of jobs to another on account of having to tune the production equipment or change tools.On the other hand, producers process a set of jobs that might come from different orders of distinct clients, and each client's order will contain multiple classes.Our research stems from an actual application involving aluminum extrusion plants.In such a manufacturing environment, if the machine has been set to process a specific set of jobs, it is desirable to run all the jobs of the same class since this will reduce the total setup time to a great extent.In particular, the measurement criterion of holding cost for an order is defined as the difference between the completion times of the first job and the last job of the order.The producer sequences all the jobs of the same order as close as possible to minimize the holding cost (or waiting times) and to reduce the shipping cost [12].
Regarding the optimization of a single objective relating to the setup time studies on a single-machine setting, readers may refer to Psaraftis [46], Monma and Potts [40], Potts and van Wassenhove [45], and Liaee and Emmons [33].For more extensions and improvements such as optimizing the maximum completion time, the makespan, and the mean flow time, we refer readers to van der Veen and Zhang [54], Ahn and Hyun [1], Gupta [24,25], Mason and Anderson [39], and Potts [44].For more literature references on family sequence independent setup times, readers can refer to a comprehensive survey by Allahverdi [3] and to three recent papers by Janiak et al. [29], Muştu and Eren [42], Singh [51], etc.
Regarding the optimization of multiobjective criteria relating to setup time studies on single-machine settings, Dileepan and Sen [16] and Fry et al. [22] were the leaders to discuss three kinds of procedures, including a priori, a posteriori, and interactive procedures.Afterward, Liao [34] developed a branch and bound algorithm to solve a single-machine multiple-objective and multiple-job class model.Gupta et al. [26] proposed two constructive polynomial time algorithms to solve a hierarchical scheduling problem.Erel and Ghosh [18] investigated a customer order scheduling problem with a certain number of products.Lin and Kononov [35] studied another customer order scheduling problem on parallel dedicated machines.Liu [37,38] and Hsu and Liu [28] addressed three different problems of multiple jobs composing a customer order in a job shop environment.More recently, Allahverdi et al. [8] addressed the no-wait flowshop scheduling problem on m machines with separate setup times to minimize total tardiness with an upper bound on makespan.Allahverdi [4] provided a survey of scheduling problems on single machine, parallel machine, flowshop, job shop settings with uncertain processing or sequence independent setup times.As for the sequence-dependent setup times, readers might refer to Alimian et al. [2], Allali et al. [9], Rifai et al. [47], and so on.
Regarding the customer order literature studies with due dates, Blocher et al. [13] examined the performance of order-based dispatching rules in a general job shop to minimize the flow time and tardiness.Erel and Ghosh [18] considered a customer order scheduling model of minimizing the total order lead time in which they assumed that a setup is needed whenever production switches from one family to another.Lee [31] proposed a branch-and-bound method incorporating with some dominance properties and three lower bounds and three heuristics to solve the order scheduling problem of minimizing total tardiness.Following the same problem of Lee [31], but with a position-based learning consideration, Xu et al. [60] applied branch-and-bound method with new properties and lower bound, simulated annealing, and particle swarm optimization (PSO) algorithms to the problem.Adopting the same customer order scheduling model of Lee [31], but with sum-of-processing-time-based learning effect, Wu et al. [56] utilized a branch-and-bound algorithm, four heuristics, and three metaheuristics (the artificial bee colony, a PSO with a varying linear declining inertia weight, and a simulated annealing) to solve it.Framinan and Perez-Gonzalez [20] proposed a new constructive heuristic with a look-ahead mechanism to improve the OMDD heuristic, proposed by Lee [31] and claimed that their proposed method outperformed OMDD.Wu et al. [57] proposed three modified heuristics, a hybrid iterated greedy algorithm, and a particle swarm colony algorithm for an order scheduling problem to minimize the linear sum of the total flowtime and the maximum tardiness.More streamlined scheduling problems of customers ordering multiple products, readers might refer to a review paper covering all concurrent-type scheduling problems by Framinan et al. [21].
Recently, Antonioli et al. [11] considered a customer order scheduling environment in which the setup times are explicit and depend on the production sequence.They proposed a mixed-integer linear programming, several heuristics, and the hybrid matheuristic for minimizing the total tardiness criterion.In view of the fact that the single-machine scheduling issues for orders of multiple work categories have not been addressed, we introduce a bicriteria single-machine scheduling with multiple job classes and customer orders.Our goal is to find a schedule that minimizes a linear combination of the total holding cost and total tardiness cost of given orders.
The remaining part of this study is organized as follows.In Section 2, we introduce the notations and formulate the problem.In Section 3, we propose a dominant property and a lower bound for the branch-and-bound method searching for an optimal solution.In Section 4, we propose four heuristics and a cloudy theoretical simulated annealing (CSA) and a cloudy theoretical simulated annealing hyperheuristic (CSAHH) methods.In Section 5, we experimentally tune the parameters of the proposed CSA and CSAHH.In Section 6, we execute several experiments to evaluate the effectiveness and efficiency of the proposed algorithms.Conclusions and suggestions are summarized in Section 7.

Notations and problem formulation
First, some notation used in this study are defined as follows.The considered problem can be stated as follows.Suppose there is a set of  jobs which are grouped into a set of  orders and each order includes several jobs.The jobs can be classified into  different classes and each order contains only one job from each class.Jobs are ready at time zero and will be operated on a single machine.No preemption is allowed during the processing a job.Suppose that each job   has a processing time and belongs to a job class.Furthermore, an order consists of at least one job from each job class.During the processing period, if a job belonging to class  is scheduled immediately following a job belonging to class  ( ̸ = ), then a setup time   is needed.Otherwise, it does not need a setup time.The holding cost HC  () of order  is defined by the range between the completion time of the first job in order  and the completion time of the last job from the same order.In this study, we address the problem of minimizing a linear combination of the total holding cost ∑︀  =1 HC  () and total tardiness cost The proposed problem is NP-hard because for fixed  = 1,  = ,  = 0, and all setup times 0, the resulting problem is an NP-hard problem (see [43]).

A lower bound and property
For the branch-and-bound method, we wish to determine a lower bound of a node  = (,   ), where  is the already scheduled part with   jobs while   is the set of    unscheduled jobs from  orders ( ≤ ).Let   be the completion time of the last job   in , and let ∑︀ ∈ HC  be the total holding cost of those completed orders in .
From the above analysis, we have the following inequality: Therefore, a lower bound on node  = (,   ) can be obtained as follows.
Apart that, let   be a subset of the job classes for jobs in   exclude the last job class in .Another lower bound for the holding cost by sorting only the jobs in   and by considering the setup times is also developed as follows: Therefore, another lower bound on node  = (,   ) can be obtained as follows: )︁ }︃ .
In order to find a better lower bound, we have It is noted that when  = ∅, the lower bound on node  = (,   ) can be obtained as Additionally, a property is proposed to improve the power of the searching ability of the branch-and-bound method.Let  = (,   ,   ,   ) and  ′ = (,   ,   ,   ) be two schedules in which job  is scheduled before job  in , while job  is scheduled before job  in  ′ .Moreover, let job   be the last job in .We have the following property.
Property 1. Assume that jobs   ,   , and   belong to the same job class, and   ∈   ,   ∈   , and  ̸ = .If both jobs   and   are the last assigned jobs of   and   , and the first assigned job of   is scheduled after the first assigned job of   ,   <   ,   <   , then  dominates  ′ .
(1) Because,    ( ′ ) =    (), and The time order   <   +   <   +   <   +   +   and the given condition   <   determine the TD values of ten situations.Three of these ten situations are proven as follows, and the remaining proofs of are similar.

Heuristics, cloudy theoretical simulated annealing algorithms and branch-and-bound
To find near-optimal solutions, in this section, we first present four heuristics, each of which was iteratively improved by a local pairwise interchange method (pi-method) several times.The pi-method means that taking a schedule with five jobs as  = (1, 2, 3, 4, 5), we generate two random integers  = 2 and  = 4 by randomness and then swap the 2nd job and 4th job in schedule  to generate a new schedule as  ′ = (1, 4, 3, 2, 5).Afterwards, we propose a cloudy theoretical simulated annealing algorithm adopting found job sequences from four heuristics (without improvement by the pi-method) as their initial seeds.Finally, we propose a cloudy theoretical simulated annealing-based hyper heuristic algorithm, along with five low-level heuristics, to solve the problem.Based on the characteristics of the studied problem, multiple-class jobs and customer order of jobs, we develop four heuristics as follows.
Heuristic 1.It is modified from the idea of Gupta et al. [26].
Step 1.We assign  orders based on the smallest value of the sum of processing times of each order plus its setup time to yield "an order schedule".
Step 2. For each order in the yielded order schedule in Step 1, we assign job sequence based on the smallest processing time first to yield job sequence within each order.For simplification, Steps 1 and 2 are termed "OSPTLPT".
Step 3. The solution of the OSPTLPT is improved by using a pi-method iteratively.It is termed "OSPTLPT pi".
Step 1.We assign  orders based on the smallest value of the sum of job processing times in each class to yield "a class schedule".
Step 2. In the class schedule in Step 1, we assign the job sequence based on the largest processing time first in the first class and assign the job sequence based on the smallest processing time first in other classes.For simplification, Steps 1 and 2 are termed "CSLPT".
Step 3. The solution of the CSLPT is improved by using a pi-method iteratively.It is termed "CSLPT pi".Heuristic 3. Step 1.We assign  orders based on the earliest order due date first rule to yield "an order schedule".
Step 2. Given the order schedule in Step 1, we assign the jobs within each order according to the non-increasing order of the sum of processing times of each job plus its setup time to yield a job sequence.For simplification, Steps 1 and 2 are termed "OEDDLPT".
Step 3. The solution of the OEDDLPT is improved by using a pi-method iteratively.It is termed "OED-DLPT pi".

Heuristic 4.
Step 1.We assign  orders based on the earliest order due date first rule to yield "an order schedule".
Step 2. Given the order schedule in Step 1, we assign the jobs in each order according to the non-decreasing order of the sum of processing times of each job plus its setup time to yield a job sequence.For simplification, Steps 1 and 2 are termed "OEDDSPT".
Step 3. The solution of the OEDDSPT is improved by using a pi-method iteratively.It is termed "OED-DSPT pi".
There are three stopping rules used presented in used approximate methods such as the maximum number of iterative cycles, specified CPU time limit, and the maximum number of cycles between two improvements of the global best solution (refer to [50,55]).However, the stopping rule of proposed CSA or CSAHH are adopted an initial temperature  0 and a final temperature  f in this study.In what follows, to find the good quality of approximate solutions, we set the best solutions (job schedules) found from the four heuristics (OSPTLPT, CSLPT, OEDDLPT, and OEDDSPT) as the four initial seeds in the cloudy theoretical simulated annealing (CSA) algorithm [53].The CSA algorithm equipped with four different seeds is recorded as algorithms CSA1, CSA2, CSA3, and CSA4.As mentioned in Torabzadeh and Zandieh [53], the CSA algorithm uses several parameters, including an initial temperature  0 , a final temperature  f , an annealing index , and a number of improvement repetitions  r .In addition, let  < 10 −8 (be a very small number) and  be a random number selected from a uniform distribution  (0, 1).Then, the procedure of the CSA algorithm is summarized as follows.
In addition, a cloudy theoretical simulated annealing hyperheuristic algorithm (CSAHH) along with five low-level heuristics is also proposed to solve this problem.A hyperheuristic with the operators of the low-level heuristics and high-level strategy has been widely utilized in searching for near-optimal solutions (Dowsland et al. [17], Burke et al. [14], Gascón-Moreno et al. [23], Anagnostopoulos et al. [10], Wu et al. [58], Zhao et al. [65], and so on.).Instead of solving the problem directly, a low-level heuristic is directly applied for searching the solution space.The high-level policy includes the heuristic selection policy and acceptance criterion.The high-level policy is adopted to find a low-level heuristic to yield a new solution.Wu et al. [58] proposed five local improvement methods and separately used each of them in the cloud theory-based simulated annealing algorithm to solve a single-machine past sequence setup scheduling with two scenarios.For the diversity of neighborhood search, we modify these five local improved methods to be the five low-level heuristics (called LH 1 , LH 2 , LH 3 , LH 4 , and LH 5 ) embedded in the CSA algorithm.Applying a low-level heuristic on a job sequence, it first needs to choose two jobs randomly from the sequence.Then, these five low-level heuristics consist of (LH 1 ) swapping two chosen jobs; (LH 2 ) switching both jobs with their immediately succeeding jobs; (LH 3 ) switching both jobs with their closet-most preceding jobs; (LH 4 ) switching the job in front with its immediately succeeding job and switching another with its closet-most preceding job; and (LH 5 ) switching the job in front with its close-most preceding job and switching another with its immediately succeeding job.
Let iCmax denote the total number of iterations in performing CSAHH.In the first run of performing the CSAHH, the probability of all five low-level heuristics contributing to the CSAHH are equally set to 1/5.We recode the accumulated performance times   of each LH  in the first run and then set   / ∑︀ 5 =1   as the selection probability of LH  for the subsequent run.For diversity, we set   = max{1,   } to keep all five lowlevel heuristics in the pool from the first run to the end of the last run in performing CSAHH.Furthermore, the notations  0 ,  f , and  used in CSAHH are the same as those in CSA, while Nr denotes the total number of times for five low-level heuristics called in each run.The description of the CSAHH is summarized as follows.
The procedures of the CSAHH method: 01: Set initial seed  =  0 , and its objective function ℎ() = ℎ( 0 ); 02: Input parameters Nr, iCmax,  0 ,  f ;  = 1 03: Set  =  0 ; 04: Set   = 1,  = 1, . . ., 5; 05: To solve the proposed problem, we applied the branch-and-bound (B&B) method for finding an exact solution.Performing the B&B method, we begin the root node with no job and branch the root node to obtain new nodes by adding the unscheduled jobs to the root node.We adopt a depth-search first rule and schedule the jobs from the first position to the last.
For each considered node, new jobs (nodes) selected from those unexplored jobs are one by one fabricated into a complete schedule by appending the selected jobs.A lower bound is calculated for each such node.Then, the dominance Property 1 and these lower bounds are used to remove those nodes that satisfy the dominance property or calculated lower bounds are greater than an established upper bound (obtained from the heuristics or CSA algorithm).Further branching is continued with the remaining nodes that are not considered yet.Until all possible nodes are either considered or removed, this procedure is repeated to obtain an optimal solution.On the basis of the above description, the steps of the proposed B&B method are described below.
Input: A set of (= ) jobs with processing times { 1 ,  2 , . . .,   } assigned independently into two sets: Step 2. For each considered active node, compute the lower bounds using the equation LB(), determinate if the lower bounds are equal to or larger than the incumbent upper bound, and remove those nodes and all beneath nodes in the branching tree.
Step 3. Applying the dominance Property 1 to purge the unwanted nodes from the branching tree.
Step 4. If the considered node is a complete schedule (say ), then compute its objective value ℎ() and if this objective value is smaller than the upper bound, then update the upper bound and the current solution by the ; otherwise, branch from the node with the depth search rule on the remaining nodes to create a set of active nodes.
Step 5. Conduct Steps 2 through 4 repeatedly, until all nodes have been considered.Let the final complete schedule be  opt .Output.The  opt is an optimal solution.Besides, an illustrative example and node trees also were provided to explain the proposed model and the steps of the proposed Branch and Bound algorithm, respectively.
Example.The following processing times of four jobs can be classified 2 classes and 2 orders.Furthermore, let  = 0.25, and let  1 = 5 and  2 = 6 denote the setup time for class 1 and class 2 while  1 = 70 and  2 = 62 denote the due date for order 1 and order 2, respectively.
Given  = ( 3 ,  1 ,  2 ,  4 ), the completion times are computed as follows: [1] () =  2 +  The best solution ℎ() among all the proposed methods was designed as the upper bound of B&B. Figure 1 presents the details of the steps of the proposed Branch and Bound algorithm.Thus, the lower bound of node  = ( 1 ,  4 , −, −) is 28.5.

Exploring the parameters used in CSA algorithms
In this section, we explore the appropriate values of parameters built into the CSA and CSAHH algorithms.These parameters are the annealing index (), the number of improvements (Nr), the initial temperature ( 0 ), and the number of runs (iCmax).To choose values for these parameters, we design a scheme consisting of 9 jobs, which are grouped into three orders and three classes, and generate one hundred tested instances for each parameter combination of , Nr,  0 , and iCmax.The job processing times are generated uniformly from ]︁ × 100%.In addition,   is the solution obtained from the CSAHH and OB  is the solution obtained from the B&B method.
The parameter adjustment process is based on the one-factor at a time experimental method [41] on a generated set of problem instances.To explore the value of , after several trials, we fixed iCmax = 20,  0 = 0.1,  f = 10 −10 and Nr = 950 and tested the value of  first over the interval (0.1, 0.9) by an increment of 0.1 and second over the interval (0.9, 0.99) by an increment of 0.01. Figure 2a shows that the AEP approaches a lower point at 0.9 on the interval (0.1, 0.9) case and then approaches a lower point at 0.97 on the interval (0.9, 0.99) case.Thus, we adopt  = 0.97 in the later simulation experiments.
To explore the value of Nr, we set  = 0.97, fix the other parameters as iCmax = 20,  0 = 0.1,  f = 10 −10 , and test the value of Nr over the interval (50,1500) by an increment of 50 per trial.Figure 2b shows that the trend of AEP declines rapidly as the value of Nr increases and becomes flat when Nr is larger than 1350.Thus, we adopt Nr = 1350 in the later tests.
Following exploring parameters  and Nr, we set  = 0.97 and Nr = 1350 and fix iCmax = 20, and  f = 10 −10 .We test the value of  0 from 0.00001 to 0.1 by a common rate of 0.1 per trial.Figure 2c shows that the AEP declines rapidly to a lower point at  0 = 0.1.
For the value of parameter iCmax, we set  = 0.97, Nr = 1350, and  0 = 0.1 and fix  f = 10 −10 .We then explore the values of iCmax from 10 to 40 (times of runs) by an increment of 2 runs per trial.Figure 2d shows that there are three local lower points in those tests.The lowest point is located at 34.Therefore, iCmax = 34 is used in later simulation experiments.
The number of jobs was set at  = 8, 10 and 12 for small-sized cases.For each , the combination of number of orders and number of classes (mo × mk) were set differently.There were 2, 2, and 4 combinations for  = 8, 10, 12, respectively.For each combination of the aforementioned parameters, a set of 100 problem instances were generated randomly.Thus, in total, 28 800(= × 100) instances were tested.To avoid too much computational time, as the number of nodes searched exceeded 10 9 , the B&B method jumped to the next set of problem instances.
With regard to computational experiments on large-sized job cases, the number of jobs was set at  = 60 and 100.A set of 100 random problem instances were generated for each combination of parameters ,  , , , mo × mk, and "setup".In total, 61 200 = × 100 instances were tested.

Analysis of small-𝑛 computational results
The average error percentage (AEP) is employed to evaluate the capability of searching optimal job sequences for the proposed 8 heuristics (OSPTLPT, CSLPT, OEDDLPT, OEDDSPT, OSPTLPT pi, CSLPT pi, OED-DLPT pi, OEDDSPT pi) and the five CSA algorithms (CSA1 to CSA4 and CSAHH).The AEP is defined as 100(mean[(  − OB  )/OB  ])[%], where   is the value of ℎ() found from each heuristic/algorithm and OB  from the B&B method.
Tables 1 and 2 convey the performance of the B&B method.Column 3 (Tabs.1 and 2) shows a sharp increase in the mean of the explored nodes as  rises from 8, 10 to 12.The average CPU times (in seconds), presented in Tables 1 and 2, dramatically increase as the number of jobs  increases to 12.In addition, all of the generated problem instances were solvable within 10 9 nodes.The ratio of feasible solution (marked as RF, RF = the number of tested instances solved out by B&B or by heuristics/total number instances tested per case) is presented in Column 5 of Tables 1 and 2. Table 2 summarizes the CPU time and the number of explored nodes for each of the Factors ,  , , , mo × mk, and "setup".From Columns 3 and 4 of Table 2, it is seen that as the numbers of job classes ("mk" in mo × mk) increase, the number of nodes explored and CPU times increase; the greater the number of job classes, the greater the number nodes explored.It can be observed that our proposed B&B was not so powerful when number of jobs increased to  = 12.But it was noted in Table 1 that on average, proposed B&B took less 700 s to solve an instance at cases mo × mk = 4 × 3 and 6 × 2. It was also noted that overall, the average of CPU time B&B took was only 550.32 s.
Table 3 also summarizes the means of AEP for the levels of the Factors ,  , , , mo × mk, and "setup".Table 3 shows that the AEP of all heuristics and algorithms increases slightly in general as  increases from 8 to 12. Furthermore, Figure 3 displays the boxplots (distributions) of AEP for the 8 heuristics and five CSA algorithms.The diamond in the boxplot represents the mean of the AEP for each heuristic/algorithm.Regarding the CPU times, it can be seen in Figure 4 that the group (CSA1, CSA2, CSA3, CSA4 and CSAHH 13) consumed more CPU time than other eight heuristics because metaheuristics are more complex than heuristics.However, all of proposed heuristics/algorithms executed all less than 1 s for small-sized instances.
To verify whether the difference among the 8 heuristics and five CSA algorithms are significant or not in a statistical sense, based on rank-sums of AEP for the 288 groups of problem instances, which are formed by the combinations of ,  , , , mo × mk, and "setup", the Freidman test (a nonparametric method) was executed (in SAS 9.4).The value of the statistic is 3287.8 with degrees of freedom 12 and a  value less than 0.0001.According to the test results, it is concluded that the AEP samples do not follow the same distribution at the

Analysis of large 𝑛 experimental results
For a large-sized job instance, the optimal solution of the test instance and its optimal objective value were rarely found.Thus, an average relative percentage deviation (RPD) is reported for the performed results.RPD = 100(mean[(  −   )/  ])[%], where   is the value of objective function ℎ() for a job sequence searched by a heuristic or by a CSA algorithm and   is the smallest value among the values found among 13 heuristics/CSA algorithms.Regarding the performances of the proposed 8 heuristics and five CSA algorithms,  and CSAHH algorithms, respectively.The CSAHH algorithm performed the best according to the RPD measure, followed by OSPTLPT pi and OSPTLPT, CSA1, CSA3, CSA4, and CSA2, and CLSPT pi and CLSPT.It is noted that the problem-based heuristic improved with a pairwise interchange method (OSPTLSPT pi) could be as powerful on producing (near-) optimal solutions as a metaheuristic (CSA algorithm) does.Table 5 also summarizes the means of RPD for levels of the factors  , , , mo × mk, and "setup" for  = 60 and 100, respectively.Table 5 shows that the RPDs of all heuristics/algorithms increase slightly in general as  increases from 60 to 100.Furthermore, Figure 6 shows the boxplots (distributions) of RPD for the 8 heuristics and five CSA algorithms.
To verify whether the difference among the 8 heuristics and five CSA algorithms are statistically significant, based on rank-sums of RPD for the 612 groups of problem instances, formed by the combinations of ,  , , , mo × mk, and "setup", the Freidman test was executed (in SAS 9.4).The value of the statistic is 7039.8 with degrees of freedom 12, and the  value is less than 0.0001.Based on the test results, it is concluded that the RPD samples do not follow the same distribution at the level of significance  = 0.05.Table 4 (Column 4) reveals sums of RPD ranks of the 612 sets of problem instances for the 8 heuristics and five CSA algorithms.The rank sums of the CSA1 to CSA4 and CSAHH algorithms are 1758.0,4082.5, 3214.5, 3706.5, and 612.0, respectively.The rank sums of the heuristics for OSPTLPT pi and OSPTLPT are 1452.5 and 2646.5, respectively, and those for the other heuristics range from 4750.5 (OEDDLPT pi) to 7956.0 (CLSPT).It is observed that the heuristics OSPTLPT pi performed better than that of any one of four CSA algorithms but worse than that of CSAHH for the large-sized job cases.To further analyze the performances among these heuristics and CSA algorithms, a multiple comparison procedure, WNMT, was used.There were 78 pairwise difference comparisons among the 8 heuristics and five CSA algorithms.The best-performing group is CSAHH, and OSPTLPT pi and CSA1 are in the second-best group at a level of significance of 0.05.The phenomenon that the performance of OSPTLPT pi performed better than that of CSA1 to CSA4 in this case reflects the advantage of Gupta's method [26], especially for a large number of jobs.Regardless, the CSAHH has the smallest RPD (0.01%) and rank-sum (612.0) for a large number of jobs.
Pertaining to the computational time or CPU times (in sec.)cost on searching near-optimal solutions, Figure 7 exhibits the violin plots of the CPU times of the heuristics and CSA algorithms.The differences between CPU times within groups are no more than 0.01, 0.001, and 0.00001 (s) for CSA algorithms, four heuristics-pi, and four heuristics, respectively.On average, the CSA algorithms, the four heuristics with the pi-method, and the four heuristics spent approximately 0.66, 0.03, and 0.0001 (s), respectively.

The performance of lower bound at small-𝑛 and large-𝑛 experimental results
In order further to evaluate the performance of the lower bound at small- and large- experimental results, we tested some additional problem instances at  = 8, 10, 12 for the small- case and at  = 60, 100, 200   Notes.Note that (•) denotes total number of optimal solutions obtained from each approximate method.for large- case, respectively.Other settings are included  = 0.5, six combinations of ( , ): (0.25, 0.25), (0.25, 0.50), (0.25, 0.75), (0.50, 0.25), (0.25, 0.50), and (0.25, 0.75), setup times was generated from uniformly distribution  (1, 10), number of classes (mo × mk) is only one case for each , and a set of 100 problem instances were generated randomly for each case.For the small- case, we reported AEPs for the performance of each algorithm, the number of the optimal solutions (in (•)) obtained from each algorithm obtained, and the mean or maximum of gap ratio of the optimal solution obtained from the B&B with respect to the lower bound (LB1 or LB2) at the root node (100(mean[(OB  − ()  )/OB  ])[%],  = 1, 2. For the large- case, we reported the average gap ratio of proposed approximate methods to the lower bound LB1 (100(mean[(  − (1)  )/1  ])[%]), where OB  is obtained from the B&B method, ()  is the value of the lower bound 1 and lower bound 2 at the root node from each instance, and   is the value of objective function ℎ() for a job sequence searched by a heuristic.To avoid too much computational time, as the number of nodes searched exceeded 10 9 , the B&B method jumped to the next set of problem instances.Thus, in total, 1800 (= × ×××(mo × mk)×setup = 3 × 2 × 3×1 × 1 × 1×100) instances were tested for the small- case and for the large- case, respectively.To avoid too much computational time, as the number of nodes searched exceeded 10 9 , the B&B method jumped to the next set of problem instances.The results were reported in Table 6 for the small- case and in Table 7 for the large- case, respectively.It can be observed in Table 6 that on average, the mean of AEPs or the numbers of the optimal solutions obtained from each algorithm of the group (CSA1, CSA2, CSA3, CSA4, CSAHH) were better than those of the group (OSPTLPT pi, CSLPT pi, OEEDLPT pi, OEEDSPT pi), and while the group (OSPTLPT pi, CSLPT pi, OEEDLPT pi, OEEDSPT pi) were better that those of the group (OSPTLPT, CSLPT, OEEDLPT, OEEDSPT).The last four columns in Table 6 reported that on average, the difference was only 6% between the improved lower bound LB2 and the original lower bound LB1.Regarding to the performances of proposed approximate methods to the lower bound LB1, it can be seen in Table 7 that OSPTLPT was the best than other three simple heuristics (CSLPT, OEEDLPT, OEEDSPT); OSPTLPT pi was the best that other three improved heuristics (CSLPT pi, OEEDLPT pi, OEEDSPT pi); and CSA1 was the best than other three metaheuristics (CSA2, CSA3, CSA4).It can also be observed that OSPTLPT or OSPTLPT pi did better than the metaheuristics (CSA2, CSA3, CSA4).The reason was that the performance of CSA may be affected by the initial heuristic or seed.Table 7 also reported that CSAHH with the smallest gap ratio did the best among all proposed 13 approximate methods.It is also noted that the required CPU time of the above lower bound is less than one second.

Conclusions
In this study, we discuss the customer order scheduling problem with jobs belonging to different classes.To balance the costs incurred by holding the finished products and total tardiness, the optimality goal is to find a job sequence that minimizes a linear combination of the holding cost and total tardiness objective.Embedded with a dominance property, a branch-and-bound method is adopted to produce solutions for a small number of jobs; it performs well up to  = 12 jobs, which consist of different combinations of the number of job classes and the number of customers.For a relatively large number of jobs ( = 60, 100), four problem-based heuristics, each coupled with a local (pairwise interchange) iteratively improved searching method, a cloud simulated annealing (CSA) algorithm with four different seeds (found from four heuristics), and a CSAHH algorithm were proposed for finding near-optimal solutions.It is found that the CSAHH performed best among all 8 heuristics and five CSA algorithms at the number of jobs  = 60 and 100.This study also observed that a well-designed problembased heuristic improved with a simple pi-method (OSPTLSPT pi) could be as powerful on producing (near-) optimal solutions (sequences of jobs) as a metaheuristic (CSA algorithm).
Further research can consider extending the single-machine processing environment to a two-stage production process scenario with multiple job classes and customer orders.The first stage might be with multiple machines in parallel (identical, uniform or unrelated) or in flow-shop, and the second stage, an assembly production line.

Figure 1 .
Figure 1.The steps of B&B for illustrative example.

Figure 4 .
Figure 4. Boxplots of CPU time for heuristics and CSA algorithms (small ).

Figure 5 .
Figure 5. Variation of probabilities called for five low-level heuristics.

Figure 7 .
Figure 7. Violin plots of RPD for heuristics and CSA algorithms (large ).

Table 1 .
Performance of the branch-and-bound method.

Table 2 .
Summary of the performance of the branch-and-bound method.

Table 3 .
AEP of Heuristics and CSA algorithms (small ).

Table 4 .
The rank-sum of heuristics and CSA algorithms.

Table 5 .
RPD of Heuristics and CSA algorithms (large ).

Table 6 .
The performances of approximate methods and two lower bounds for small- case.

Table 7 .
The performances of approximate methods with respective to the lower bound LB1 for big- case.