Multi-factory Job Shop Scheduling With Due Date Objective

The existing literature on distributed scheduling mainly focuses on the performance measures makespan, total completion time, or costs. Due date related objectives that are gaining in importance in the industry have not been considered to a similar extent. In this contribution, we present a model formulation of the distributed job shop scheduling problem with due date consideration and present adapted greedy heuristics as well as a genetic algorithm to solve large problem instances. Computational experiments are carried out to assess the performance of the model and the algorithms. The heuristics and metaheuristics show promising results in reasonable computation times, with the genetic algorithm outperforming the other heuristics. The results indicate that managers should consider incorporating due date related objectives in the decision-making process of production planning and scheduling in distributed manufacturing.


I. INTRODUCTION
Manufacturing firms are increasingly distributed in geographical terms and operate multiple factories, as a result of globalization and increased customer demands for fast order delivery [1], [2]. This decentralization of production facilities offers advantages such as lower costs, proximity to customers, and flexibility but also poses new challenges regarding planning and scheduling of production in the factories. Scheduling orders in multifactory networks involve more steps than the already complex task of scheduling jobs in a single factory [2], [3]. Before sequencing jobs on machines locally, the orders (as jobs) need to be assigned to a factory. The problems vary in complexity, depending on the manufacturing principle or configuration, similar to the single-factory case. This article focuses on the job shop configuration, which is typical for various industries and also present in distributed networks [4], [5]. In a job shop, a set of n jobs have to be manufactured, each with its own routing on the m machines. Single-factory job shop scheduling has been studied extensively by scholars ever since the heuristic by Jackson in 1956 [6]. Distributed scheduling has recently gained the attention of scholars, for recent reviews refer to [2], [7].
We provide a short review of relevant studies for this article that have been published on the multi-factory or distributed job shop scheduling problem (DJSP): Jia et al. [8] developed genetic algorithms (GA) to minimize makespan and costs. Chan et al. examined the DJSP with flexible routings and makespan objective with a GA [9]. Then, [10] focused on a GA to minimize the total completion time of all jobs in the DJSP. In [11], a tabu search and simulated annealing algorithm were hybridized. Two models and three greedy heuristics were proposed for the DJSP with a makespan objective in [12]. Wu et al. [13] examined different chromosome representations in GAs for the DJSP with a makespan objective.
As this short and selective review indicates, previous studies have primarily focused on the performance measures makespan (Cmax) or total completion time (∑ ). However, the manufacturing industry is evolving quickly and trends like Industry 4.0 led to a situation where for many manufacturing firms, timely deliveries without delays are at least as important as optimal utilization of machine capacities and swift completion of production orders (Cmax) [14]. Earliness and tardiness of orders have different effects on the present and future performance of firms (see Fig. 1). Both objectives should be addressed in order to reach an optimal control of production towards meeting due dates.
To advance in this promising research area, we formulate the problem as a MILP model with a due date objective in this contribution and adapt several heuristics from the literature on classical job shop scheduling for the DJSP Additionally, we propose two adapted greedy heuristics and a genetic algorithm to solve large problem instances efficiently. Experiments are conducted to evaluate the performance of the algorithms using the wellknown instances of Taillard [15] and Fisher and Thompson [16].
The remainder of this paper is structured as follows: In Section 2, the mathematical model is presented and described. Section 3 contains the adapted and developed algorithms and heuristics, whereas Section 4 presents numerical experiments and results. Section 5 concludes the paper with an outlook on further research opportunities. The problem under consideration is modeled as a mixed-integer linear programming model, which adopted the processing constraints (7) and (8) from [12]. Each machine can only process one job at the same time. Although transporting work-in-process orders to another factory in case of disruptions or loss of machine capacity Subject to: , + − = ∀ (3) , ≥ , + , ∀ , , ≠ , , = 1 , ≥ , + , − 3 − , , − , − , ∀ , , < , > , ≥ , + , − 2 + , , − , − , ∀ , , < , > , , ∈ 0, 1 ∀ < , > , Equation (1) is the objective function. The constraint (2) specifies that , is the latest completion time of one of the individual processing steps on one of the machines. The earliness or tardiness of each job is calculated in (3). Then, (4) ensures that each job can only be assigned to one factory. Constraint (5) determines that at least the corresponding processing time must have elapsed before a job can be transferred to the next machine. The following machining step can be terminated at the earliest when the production time for this step has elapsed, in addition to all previous processing times (6). Constraint sets (7) and (8) indicate that a machine can only process one job at a time: If two orders i and j are processed in the same factory and j is processed on machine m after i, then j must be completed at least by its processing time p on machine m later than i. Constraint sets (9), (10) and (11) specify the decision variables of the model.
The model can be extended to focus on minimizing costs of due date deviation (early or late completion). If a specific weighting factor can be assigned to each customer order (indicating the importance of different customer orders), the model (with the remaining parameters and decision variables unchanged) can be adapted to: , Weighting factor for earliness of job i , Weighting factor for tardiness of job i As the determination of the weighting factors is very subjective, we focus on the first model without weighting factors in this contribution. Nevertheless, the adaption highlights the flexibility and simplicity of the developed MILP model.

A. Heuristics (Two-stage iterative approaches)
One way to reduce the complexity of the problem is by taking an iterative solution approach. The two decisions in the DJSP are split into two phases to find near-optimal solutions in reasonable computation times. First, jobs are assigned to factories, followed by the sequencing of jobs on machines. Here, the assignment rule significantly determines the quality of the solution. Allocating the jobs according to descending due dates does not lead to a promising solution. Alternatively, the jobs can be sorted according to their sum of processing times. The difference between this makespan and the due date ("slack") is then used to determine the assignment sequence. In this way, the flexibility in planning is enhanced by preferring urgent jobs -at the expense of longer idle times for jobs with a lot of slack -when creating the plan. However, as the machines are still unscheduled when jobs are assigned, measuring the specific earliness of jobs is not trivial. Therefore, we solely focus on the tardiness of jobs for the used assignment rule, which is conceptually similar to the NEH2 rule [17] that performed well for the distributed permutation flow shop scheduling problem. The main difference is the order in which the jobs are selected, which is minimum slack, i.e., the buffer time. Then, the first s orders are assigned to the s plants one by one in order to separate the most urgent jobs and avoid interference. Each factory receives one urgent job first to ensure that the urgent jobs are started immediately. Subsequently, for the next job in the sorted list, the assignment that leads to the earliest completion time is computed. Various priority/dispatching rules may be used for this calculation. However, rules like FIFO/FCFS or SPT do not yield good results in the context of due date consideration. We focus on rules that include slack and can be used intuitively to meet due dates: • Minimum slack (MSLACK): MSLACK iteratively checks the slack of all jobs and selects the job with the least slack whenever a machine is free. This includes observing the other restrictions like the machine sequence of each job. • Slack per remaining processing time (S/RPT): S/RPT divides the remaining slack by the sum of outstanding raw processing times, to avoid the tendency of MSLACK to prefer long processing steps over shorter ones. The jobs are then sorted in ascending order of this quotient and assigned accordingly. • Slack per remaining operations (S/OPN): Another way to modify MSLACK is by integrating the number of remaining operations in the calculation. Jobs are then sorted according to the remaining slack for each remaining process step.
The following steps are similar for all heuristics. Each job is assigned to a factory and the resulting completion times determined. Once all jobs are assigned, the sequence for each factory is determined, either by again using the three presented rules or by using more sophisticated approaches like the GH1 heuristic [12] presented next.

B. Two-stage greedy heuristic (GH1)
Greedy heuristics have been successfully applied to many combinatorial problems and usually derive a performant solution in a reasonable amount of time. For the problem at hand, we use an encoding scheme for GH1 that consists of as many permutations as factories are available. Each value in the permutation represents an operation and each job appears as often as it needs operations processing. As the quality of the algorithm depends on a performant starting solution, the start permutation is randomly mixed. Subsequently, the operations are iteratively inserted into the sequence of operations in each factory. Naderi and Azab [12] use an additional rule for job-factory assignments based on the smoothing of workload and consider each machine separately. We use GH1 with the MSLACK rule in this article to focus on the minimization of due date deviations.

C. Greedy heuristic (integrated approach)
So far, the presented heuristics and rules focused on a two-stage approach. We now turn to an integrated approach of assigning jobs to factories and sequencing jobs on machines in each factory. A greedy heuristic (GH3 by [12]) was adapted to the examined problem with a due date focus. The job-factory assignment is based on minimum slack and the shortest cycle times in the plants. The processing sequences of the jobs are calculated in each iteration using the GH1 procedure, where all plants are considered in each step. The procedure of the algorithm is as follows: First, one job is assigned to each factory. Again, jobs with minimum slack are selected first and the makespan is calculated based on their processing times. The next job is then assigned to the factory with the shortest makespan after the assignment of the respective job. This is done by adding each operation to the previous sequence and inserting it in all possible positions, without, however, changing the previously determined sequence. For each position in the sequence, the makespan is recorded and then compared with the other position to select the best solution. If several positions lead to the same (best) solution, the earliest is chosen. This selection enables the jobs that are scheduled subsequently to be completed as early as possible and potential delays can be avoided. Unlike the two-stage algorithm GH1, GH3 does not depend on a random start sequence.
To specifically focus GH3 on due dates, we created a further adaption, denoted as GH3_SlackMin: The jobfactory assignment is made after calculating the slack of all orders. Each job is iteratively assigned to each factory and scheduled using the exchange procedure. Then the remaining slack of all jobs is cumulated for the sequence. The job is assigned to the factory where the slack is maximized, that is, where flexibility is higher compared to the other locations. In this way, avoiding tardy orders is especially considered.

D. Genetic algorithm adapted for the DJSP with due date objective
Genetic algorithms (GA) are popular meta-heuristics for scheduling problems due to their ability to investigate large solution spaces and simultaneously apply local search methods. We use a simple chromosome encoding for the adapted GA (based on [8]) aimed at the DJSP with a due date objective: The first number indicates the factory, followed by the job number. The first appearance of the job number indicates the first operation, the second appearance the second operation and so on. The fitness value that is calculated for each chromosome refers to the due date deviation of each developed schedule. Rather than randomly assigning a job to a factory during initialization and only changing this assignment rarely using global mutation, we modify the crossover and mutation operators to maintain a more extensive search space for the algorithm.
The initial job assignment is done separately for each chromosome in GA_adapted (see Fig. 2). In the iterations, two chromosomes are selected for the crossover process and the subsequences are exchanged. The feasibility of the new solution is checked and in case of a violation (i.e., job operations now spread among more than one factory), the crossover procedure overwrites the previous factory assignment for all further genes of this job. Thus, each generation explores new options for factory assignment using the crossover operator. Additionally, the best chromosome from the last iteration and the best chromosome overall are transferred to the next generation. The succeeding chromosome is assessed regarding its fitness value after crossover and after mutation. The old chromosome is replaced if the new one outperforms it. This facilitates the movement of chromosomes towards a better solution and the crossover procedure prevents an excessive restriction of the search space.
In mutation, a better chromosome is taken over from the next generation with a certain probability. The selection probability is determined using a ranking of chromosomes based on their fitness value that is then transformed into a vector of selection probabilities (p) using the following equation [18]: The greater is, the greater the probability of selecting the best chromosome. If = 0, the worst chromosome is never selected for mutation or crossover. For = 0.01, the best chromosome is twice as likely to be selected than the median, which does not restrict the search space too much as all chromosomes might be selected. All steps in the iteration procedure occur for a fixed set of iterations and the algorithm is terminated if no improvement is determined for several iterations. The parameters of the GA were determined using several small instances of the MILP model and led to the following specifications:

IV. EXPERIMENTAL EVALUATION AND DISCUSSION
In this section, we use conducted experiments to evaluate the performance of the model and the proposed heuristics. We start with the evaluation of the model on small instances before we move on to larger instances to compare the performance of the heuristics. The model was formulated and solved using CPLEX 12.10, whereas the algorithms were coded in MATLAB. All runs were conducted with a PC with 1.6 GHz i5-8250 CPU and 8 GB of RAM. The performance measure used in this contribution is the relative percentage deviation (RPD), with Alg as an indicator for the objective value of the tested algorithm and Min as the best solution that was obtained by any algorithm for the given problem instance: So far, there are no benchmark sets for the DJSP taking due dates into account. Therefore, we created some smaller instances to test the MILP model. Parameters include jobs (J), machines (M), factories (S) and the processing time , , which was randomly set to U (1,99). Setting realistic due dates has a significant impact on the quality of the solution. We experimented with this parameter and set a value of 120% of the sum of processing times of each job for the tested small instances. For the larger instances, this value was increased to a maximum of 200% to take into account the increasing number of orders and potentially more scheduling conflicts. The time limit for the model and the algorithms was set at 3600 seconds. Table II indicates the results of the MILP model for the small instances. The model is able to solve the smaller instances to optimality, whereas for the larger models, CPLEX fails to determine the optimum solution (at least in the case of 2 factories). Regarding the complexity of the model, it is quadratic in J, but linear in M and S. For example, for the instance 10x4x2, the number of constraints is 870. The variables amount to 200 binary, 70 integers, and 60 continuous. The most influential parameter on the complexity of the model is the number of jobs J.  All the heuristics were able to solve the small instances in ≤1 s, mostly with the optimal solution. Therefore, and due to space restrictions, we only present the values of the MILP model in Table II and do not present the heuristics. To evaluate the heuristics and meta-heuristics for larger instances, we derived a second set of ten instances based on the benchmark instances of Taillard [15], denoted as "tai" in Table III, as well as Fisher and Thompson [16], denoted as "ft". The results of the best performing heuristics and the meta-heuristic are displayed in Table III. The values in bold represent the best performance value found (Min). Table IV indicates the different instances with the parameter values used. Here, DDO represents the due date offset used.  Instance  J  M  S  DDO  10x10  10  10  2  500 time units  ft06  6  6  2  120%  ft06  6  6  3  120%  ft10  10  10  2  120%  ft10  10  10  3  120%  ft10  10  10  4  120%  ft20  20  5  5  150%  tai15  15  15  2  150%  tai15  15  15  3  150%  tai50 50 15 5 200% The genetic algorithm GA_adapted performed best overall and was able to determine the best solution in 7 of 10 instances. The total RPD for GA_adapted is also the lowest at 15%. Due to the high number of generations, the computational times of GA_adapted are also the longest on average. The second best performing heuristic is GH1_MSLACK, which was run for 20 times for each problem instance to generate different random start sequences and assess the solution quality comprehensively. If we specifically compare the two-stage algorithm GH1_MSLACK with the integrated algorithms GH3 and GH3_SlackMin, we find that for medium problem sizes the integrated approach does not perform better than the twostage approach. The integrated approaches only outperform the two-stage approach on larger instances with more than ten orders where the advantage of iterative assignment to factories and sequencing is exploited.
The simple heuristic rules MSLACK, S/RPT, and S/OPN show promising results in very short computation times (< 1 ) . Interestingly, they even outperform the greedy heuristics in some instances. This may be caused by the due date objective, which is explicitly considered in the three rules. The greedy heuristics were originally designed to minimize the makespan and therefore arrange jobs to be completed as early as possible, which might lead to earliness and thus to a deterioration of the objective value.
Further adjustments or new algorithm designs seem promising in this context. Other metaheuristics could also be meaningfully applied to the DJSP. For example, adaptations of simulated annealing or iterated greedy heuristics have performed well on related problem settings and may be considered [19], [20]. The need for standardized benchmark instances for the DJSP with a focus on due date related objectives has also become apparent.

V. CONCLUSION
Research on distributed scheduling in a job shop configuration that also considers due dates is scarce. However, ensuring on-time delivery is becoming increasingly important for several manufacturing industries and there is a need for models and algorithms that focus on due dates. In this contribution, we developed and presented a MILP model for the DJSP with a due date objective. Several heuristics from the single-factory scheduling literature were adapted to the specific problem at hand. Additionally, a genetic algorithm was presented to solve larger instances of the DJSP effectively. The scheduling model and the algorithms can be readily applied in realworld settings, where manufacturing is distributed among several factories. In this way, industries can advance towards modern configurations of manufacturing systems. The experimental results indicate that the heuristics and metaheuristics provide promising solutions in reasonable computation times. They are better suited for large problem instances than the exact model. Future research should deal with heterogeneous factories, where different machine types or processing times are present in the factories as well as stochastic and dynamic demand. Besides, additional characteristics should be considered in the model to make the problem setting more realistic.
This could also include multi-objective problems, as many practical settings in the industry are subject to multiple objectives. For instance, it might be worth striving to find a balance between efficient manufacturing with a high machine utilization, focusing on makespan, and the fulfillment of customer requested delivery dates, e.g. tardiness.