Scheduling in a two-stage supply chain with uncertain parameters

Abstract

This article considers a scheduling problem of products and distribution in a two-stage supply chain under uncertain environment. The first stage addresses scheduling a set of jobs on parallel machines with uncertain processing times, while the second stage concerns on delivering finished jobs to customers in batches where the transportation times are assumed to be uncertain. The objective of the problem is to minimize the holding cost, tardiness penalty cost, earliness penalty cost and transportation cost under some constraints. In order to obtain a good approximate solution to the problem, a heuristic algorithm is used. A numerical experiment is carried out to demonstrate the efficiency of the approach and make a comparison between three methods.

Keywords

Two-stage supply chain parallel machine uncertain processing time uncertain transportation time heuristic algorithm

1 Introduction

Production and transportation operations are two key segments in supply chain. To satisfy the customers’ demands, manufacturers endeavor to integrate these two segments and scheduling jointly in a coordinated manner. In the past decade, a lot of research has been done on various supply chain scheduling. Many different intelligent algorithms were introduced to solve supply chain scheduling problems. The first study of supply chain scheduling by Hall and Potts [11] considered a variety of scheduling, batching and delivery problems in an arborescent supply chain. They proposed models which minimize the total scheduling and batch delivery cost. Finally, they demonstrated that cooperation between a supplier and a manufacturer may reduce the total system cost depending upon the scheduling objective. Chang and Lee [2] studied a class of two-stage scheduling problems in which the first stage is job production and the second stage is job delivery. The worst-case performance ratio for their heuristic is proven to be 5/3 and the bound is tight. Dawande et al. [7] studied conflict and cooperation issues arising in a supply chain where a manufacturer made products which are shipped to customers by a distributor. They considered several practical scenarios about the level of cooperation between the manufacturer and the distributor and provided algorithms for each of these scheduling problems. At the last, they demonstrated that the cost saving provided by cooperation between the decision makers is usually significant. Chen [3] investigated the integrating production and transportation scheduling problem in a make-to-order environment. Results from simulation showed that their heuristics outperformed the human expert schedule. You and Hsieh [28] proposed a mixed integer nonlinear programming model to deal with an assembly and transportation scheduling problem. They developed a heuristic algorithm to minimize overall costs. The presented heuristic algorithm is shown to perform well compared with the well-known GAMS/BARON and Lingo commercial software on testing randomly-generated problems. Zegordi et al. [30] considered the scheduling of products and vehicles in two-stage supply chain environment. They set a model and proposed a gendered genetic algorithm to solve the model. Their experimental results showed that the proposed gendered genetic algorithm outperformed standard genetic algorithms. Ullrich [26] studied machine scheduling and vehicle routing with time windows integrated problem. His numerical experimental results show that the genetic algorithm outperforms the classic decomposition approaches in some cases. Cakici et al. [1] studied the problem of minimizing the total weighted tardiness and total distribution costs in an integrated problem and distribution environment. They developed different genetic algorithms and computational tests had shown that both the mathematical modeling and heuristics approaches were competitive and can produce a significant number of non-dominated solutions. Gordon et al. [9] reviewed a lot of results on scheduling with due date assignment under such conditions as given precedence constraints or various scenarios of processing time deterioration and learning. Mazdeh et al. [21] described a model for minimizing total weighted flow times plus delivery cost for a set of jobs that are to be processed on a single machine and delivered in batches to a single customer. They proved that when the maximum number of batches is fixed, the branch-and-bound scheme outperforms the dynamic programming algorithm. Pei et al. [22] modeled a supply chain scheduling problem considering multiple manufacturers and non-identical job sizes in a two-stage supply chain to minimize the makespan. Numerical results showed that the proposed modified gravitational search algorithm better than PSO and GA. Hamidinia et al. [12] proposed an effectivegenetic algorithm for minimizing total tardiness and earliness of weighted jobs in a batched delivery system. The proposed genetic algorithm outperformed thetraditional genetic algorithm under non-batched systems. Lei et al. [14] proposed a two-phase solution approach to deal with an integrated production, inventory and distribution routing problem. They provided a two-phase methodology to effectively solve this problem. Phase I solved a mixed integer programming model and phase II solved an associated consolidation problem.

Almost all of the literature mentioned above assume that parameters are known and precisely defined. In other words, they are deterministic. However, in real supply chain scheduling systems, we will encounter various uncertain factors. So many scholars made great efforts on stochastic case. But in many cases, using probability theory to deal with an indeterminate variable is not appropriate [19]. Thus uncertainty theory, based on the experts’ judgement, was proposed by Liu [15]. From then on, uncertainty theory has become a powerful mathematical tool to deal with various problems under incomplete information. For example, Chen [4] derived an American option pricing formula from uncertain financial market [16, 20] and some mathematical properties of these formulas were studied. Ding and Zhu [8] proposed two types of uncertain models for project scheduling problems. One type of models can be transformed to their crisp forms, which may be solved by classical optimization methods. The other type of models could not be transformed to their crisp forms, an uncertain simulation is employed to approximate uncertain functions. Ke [13] studied an uncertain random time-cost trade-off model. He showed that the proposed model can be solved by uncertain random simulation and GA. Sheng and Gao [24] investigated chance distribution of the maximum flow of an uncertain random network. A new method is derived to calculate chance distribution of the maximum flow for an uncertain random network. Yao and Chen [27] proposed a method for solving uncertain differential equations. They designed a numerical method for solving uncertain differential equations, which essentially solves each α-path and produces an inverse uncertainty distribution of the solution. Zhu [31] studied an optimal control model for an uncertain differential equation. The principle of optimality for uncertain optimal control is obtained by applying Bellman’s Principle of Optimality. And a fundamental result called the equation of optimality in uncertain optimal control is given. Finally, a portfolio selection model can be solved by the equation of optimality effectively.

When no samples are available to estimate a probability distribution, invite some domain experts to evaluate the belief degree that each event will happen is necessary. In the process of production and processing, many procedures do not require high accuracy. For example, rough turning, wood processing, heavy boring and so on. Because there exist many man-made factors during processing, a worker must estimated the processing time of job based on his past experience and does not estimate the processing time according to a probability distribution. Due to the variability of road conditions and weather, a driver can not estimate transportation time according to a probability distribution. But he can estimate transportation time by his past experience. So it is reasonable that the processing time and transportation time are assumed to be uncertain variables. Moreover, it is inappropriate to model belief degrees by probability theory because it may lead to counterintuitive results. Consider a counterexample presented by Liu [19]. Assume there is one truck and 50 bridges in an experiment. Also assume the weight of the truck is 90 tons and the 50 bridge strengths are independent and identically distributed uniform random variables on the interval [95, 110] in tons. For simplicity, suppose a bridge collapses whenever its real strength is less than the weight of the truck. Now let us have the truck cross over the 50 bridges one by one. It is easy to verify that $\Pr {“ Truck can cross over the 50 {bridges}^{″}} = 1 .$ That is to say, we are 100% sure that the truck can cross over the 50 bridges successfully. However, when there do not exist any observed samples for the bridge strength at the moment, we have to invite some bridge engineers to evaluate the belief degrees about it. As we stated before, human beings usually estimate a much wider range of values than the bridge strength actually takes because of the conservatism. Assume the belief degree function is $Φ (x) = {\begin{matrix} 0, & if x \leq 80 \\ \frac{x - 80}{40}, & if 80 \leq x \leq 120 \\ 1, & if x \geq 120 . \end{matrix}$ Let us imagine what will happen if the belief degree function is treated as a probability distribution. At first, we have to regard the 50 bridge strengths as independent and identically distributed uniform random variables on [80, 120] in tons. If we have the truck cross over the 50 bridges one by one, then we immediately have $\begin{matrix} \Pr {“ the truck cross over the 50 {bridges}^{″}} \\ = 0 . 75^{50} \approx 0 . \end{matrix}$ Thus it is almost impossible that the truck crosses over the 50 bridges successfully. Unfortunately, the two results are at opposite poles. This example shows that, by inappropriately using probability theory, a sure event becomes an impossible one. The error seems intolerable for us. Hence the belief degrees cannot be treated as subjective probability.

In this paper, we will study a two-stage supply chain scheduling problem. In the first stage there are some parallel machines for processing jobs, and in the second stage there are some vehicles for carrying completed jobs to customers. We assume that the processing time and transportation time are uncertain variables.

The remaining of this paper is organized as follows. In the next section, some basic definitions and properties about uncertainty theory will be introduced. In the third section, a mathematical model will be established about the two-stage supply chain scheduling problem with uncertain processing time and transportation time. In the fourth section, a solution methodology of the model will be proposed. In the fifth section, an numerical example is provided.

2 Preliminary

This section will introduce some concepts in uncertainty theory. Let Γ be a nonempty set, and $L$ be a σ-algebra over Γ. Each element $A \in L$ is called an event. A set function $M$ from $L$ to [0, 1] is called an uncertain measure if it satisfies three axioms 15: (normality) $M {Γ} = 1$ ; (duality) $M {Λ} + M {Λ^{c}} = 1$ for any event Λ; and (subadditivity) $M {⋃_{i = 1}^{\infty} Λ_{i}} \leq \sum_{i = 1}^{\infty} M {Λ_{i}}$ for every countable sequence of events Λ₁, Λ₂, ⋯. The triplet $(Γ, L, M)$ is called an uncertainty space. To obtain an uncertain measure of a compound event, a product uncertain measure was defined by Liu [16]: (Product Axiom) Let $(Γ_{k}, L_{k}, M_{k})$ be uncertainty spaces for k = 1, 2, ⋯, the product uncertain measure $M$ is an uncertain measure satisfying $M {\overset{\infty}{\underset{i = 1}{Π}} Λ_{k}} = \overset{\infty}{\underset{i = 1}{\land}} M_{k} {Λ_{k}}$ , where Λ_k are arbitrarily chosen events from $L_{k}$ for k = 1, 2, ⋯, respectively.

An uncertain variable is a function ξ from an uncertainty space $(Γ, L, M)$ to the set $R$ of real numbers such that for any Borel set B of real numbers, the set {ξ ∈ B} = {γ ∈ Γ : ξ (γ) ∈ B} is an event. An uncertain vector ξ is (ξ₁, ξ₂, ⋯ , ξ_n) where ξ₁, ξ₂, ⋯ , ξ_n are uncertain variables. The uncertain distribution Φ of an uncertain variable ξ is defined by $Φ (x) = M {ξ \leq x}$ for any real number x.

Definition 1. [16] The uncertain variables ξ₁, ξ₂, ⋯, ξ_n are said to be independent if $M {⋂_{i = 1}^{m} (ξ_{i} \in B_{i})} = min_{1 \leq i \leq m} M {ξ_{i} \in B_{i}}$ for any Borel sets B₁, B₂, ⋯ , B_m of real numbers.

Example 1. Let ξ₁ be an uncertain variable and let ξ₂ be a constant c. For any Borel sets B₁ and B₂, if c ∈ B₂, then $M {ξ_{2} \in B_{2}} = 1$ and $\begin{matrix} M {(ξ_{1} \in B_{1}) \cap (ξ_{2} \in B_{2})} = M {ξ_{1} \in B_{1}} \\ = M {ξ_{1} \in B_{1}} \land M {ξ_{2} \in B_{2}} . \end{matrix}$ If c ∉ B₂, then $M {ξ_{2} \in B_{2}} = 0$ and $\begin{matrix} M {(ξ_{1} \in B_{1}) \cap (ξ_{2} \in B_{2})} = M {\emptyset} = 0 \\ = M {ξ_{1} \in B_{1}} \land M {ξ_{2} \in B_{2}} . \end{matrix}$

Given an increasing function Φ (x), Peng and Iwamura [23] introduced an uncertainty space $(R, B, M)$ as follows. Let $B$ be the Borel algebra over $R$ . Let $C$ be the collection of all intervals of the form (- ∞ , a], (b, + ∞), ∅ and $R$ . The uncertain measure $M$ is provided in such a way: First, set $\begin{matrix} M {(- \infty, a]} = Φ (a), M {(b, + \infty)} = 1 - Φ (b), \\ M {\emptyset} = 0, M {R} = 1 . \end{matrix}$ Second, for any $B \in B$ , there exists a sequence {A_i} in $C$ such that $B \subset ⋃_{i = 1}^{\infty} A_{i} .$ Thus $M {B} = {\begin{matrix} inf_{B \subset \cup A_{i}} \sum_{i = 1}^{\infty} M {A_{i}} if \\ inf_{B \subset \cup A_{i}} \sum_{i = 1}^{\infty} M {A_{i}} < 0.5, \\ 1 - inf_{B^{c} \subset \cup A_{i}} \sum_{i = 1}^{\infty} M {A_{i}} if \\ inf_{B^{c} \subset \cup A_{i}} \sum_{i = 1}^{\infty} M {A_{i}} < 0.5, \\ 0.5 otherwise . \end{matrix}$ (1) The uncertain variable defined by ξ (γ) = γ from the uncertainty space $(R, B, M)$ to $R$ has the uncertainty distribution Φ. This kind of uncertain variables are called common. An uncertain variable with an uncertainty distribution will be regarded as a common one in order to deal with distributions of its functions.

Let us discuss the distribution of f (ξ) for a common uncertain variable ξ or a common uncertain vector. Assume $C$ is the collection of all intervals of the form (- ∞ , a], (b, + ∞), ∅ and $R$ . Each element A_i emerging in the sequel is in $C$ .

Theorem 1. [32] (i) Let ξ be a common uncertain variable with the continuous distribution Φ (x) and f (x) a Borel function. Then the distribution of the uncertain variable f (ξ) is $\begin{matrix} Ψ (x) = M {f (ξ) \leq x} \\ = {\begin{matrix} inf_{{f (ξ) \leq x} \subset \cup A_{i}} \sum_{i = 1}^{\infty} M {A_{i}} if \\ inf_{{f (ξ) \leq x} \subset \cup A_{i}} \sum_{i = 1}^{\infty} M {A_{i}} < 0.5, \\ 1 - inf_{{f (ξ) > x} \subset \cup A_{i}} \sum_{i = 1}^{\infty} M {A_{i}} if \\ inf_{{f (ξ) > x} \subset \cup A_{i}} \sum_{i = 1}^{\infty} M {A_{i}} < 0.5, \\ 0.5 otherwise . \end{matrix} \end{matrix}$

(ii) Let $f : R^{n} \to R$ be a Borel function, and ξ = (ξ₁, ξ₂, ⋯ , ξ_n) be a common uncertain vector. Then the distribution of the uncertain variable f (ξ) is $\begin{matrix} Ψ (x) \\ = M {f (ξ_{1}, ξ_{2}, \dots, ξ_{n}) \leq x} \\ = M {(ξ_{1}, ξ_{2}, \dots, ξ_{n}) \in f^{- 1} (- \infty, x)} \\ = {\begin{matrix} sup_{Λ_{1} \times Λ_{2} \times \dots \times Λ_{n} \subset Λ} min_{1 \leq k \leq n} M_{k} {Λ_{k}} if \\ sup_{Λ_{1} \times Λ_{2} \times \dots \times Λ_{n} \subset Λ} min_{1 \leq k \leq n} M_{k} {Λ_{k}} > 0.5, \\ 1 - sup_{Λ_{1} \times Λ_{2} \times \dots \times Λ_{n} \subset Λ^{c}} min_{1 \leq k \leq n} M_{k} {Λ_{k}} if \\ sup_{Λ_{1} \times Λ_{2} \times \dots \times Λ_{n} \subset Λ^{c}} min_{1 \leq k \leq n} M_{k} {Λ_{k}} > 0.5, \\ 0.5 otherwise \end{matrix} \end{matrix}$ where Λ = f^-1 (- ∞ , x), and each $M_{k} {Λ_{k}}$ is derived from (1).

Definition 2. [15] Let ξ be an uncertain variable, and α ∈ (0, 1]. Then $ξ_{sup} (α) = sup {r ∣ M {ξ \geq r} \geq α}$ is called the α-optimistic value to ξ, and $ξ_{inf} (α) = inf {r ∣ M {ξ \leq r} \geq α}$ is called the α-pessimistic value to ξ.

The uncertain variable ξ will reach upwards of the α-optimistic value ξ_sup with uncertain measure α, and will be below the α-pessimistic value ξ_inf (α) with uncertain measure α.

Founded by Liu [17] in 2009 and refined by Liu [18] in 2010, uncertain programming is a type of mathematical programming involving uncertain variables. Assume that x is a decision vector, ξ is an uncertain vector, f (x, ξ) is an objective function, and g_j (x, ξ) are constraint functions, j = 1, 2, ⋯ , p. In order to obtain a decision with minimum objective value subject to chance constraints, Liu [17] proposed the following uncertain programming model, ${\begin{matrix} min_{x} min_{\bar{f}} \bar{f} \\ subject to : \\ M {f (x, ξ) \leq \bar{f}} \geq β \\ M {g_{j} (x, ξ) \leq 0} \geq α_{j}, j = 1, 2, \dots, p, \end{matrix}$ where α_j and β are specified confidence levels for j = 1, 2, ⋯ , p, and $min \bar{f}$ is the β-pessimistic value. If we want to minimize the pessimistic value of an objective (such as makespan or cost) subject to some chance constraints, we can use the above model.

A key problem in the research area of uncertain programming is how to solve the above model. Under the framework of uncertainty theory, an uncertain programming model may be transformed into a crisp model if the related functions are monotone [18]. If the related functions are not monotone, we may employ uncertain simulation [32] to approximate the uncertain measure of an event and β-pessimistic value of an uncertain variable. Thus, we can solve the model by simplex method, cutting plane method, tabu search, interior point method, and heuristic algorithms.

3 Problem formulation

We will introduce our problem in two stages. The first stage is processing jobs on machines in a manufactory. And the second stage is to delivery the completed jobs from the manufactory to customs in batch.

At first stage, most of the literature considered single machine. The interested readers may refer to [1 , 29]. However, multiple machines are more closer to the reality. So we assume that there are n independent non-preemptive jobs to be sequenced for processing on m parallel identical machines and processing time of each job is assumed to be an uncertain variable. In addition, the transportation times are assumed to be uncertain variables. The following assumptions are also assumed in advance.

All of the jobs are available for processing at time zero and can be processed on any machine. And each machine can process one job at a time.

There is no precedence relation between jobs. Jobs which are completed before the distribution must be stored.

There are finite vehicles and each vehicle can be reused and has a fixed capacity. The loading quantity of a vehicle is related to the size of jobs.

Jobs after completed will be delivered in batches to customers. The batch delivery start time is equal to or greater than completion time of the last job in a batch.

Our objective is to minimize total cost, which consists of inventory holding cost, delivery earliness penalty cost, transportation cost, and delivery tardiness penalty cost. The structure of the problem is shown in Fig. 1.

Our problem considered in this paper is a parallel-machine scheduling problem in which jobs are ultimately delivery to customers in batches. There are n jobs have to be processed by m machines. M vehicles transport jobs in b batches. Each job has a size v. The capacity of each vehicle is Q. Each job has a processing time ξ_i for i = 1, 2, ⋯ , n and each trip has a transportation time η_r. As previously mentioned, we assume that the parameters ξ_i and η_r are uncertain variables. They have uncertainty distribution Φ_{ξ
_i} and Φ_{η
_r}, respectively. The cost delivery each batch to customer is denoted by R. The completion time of job i is denoted by C_i. Completion time of a job in kth location on machine j is denoted by T_jk. C_max denotes the makespan. The actual departure time of job i is denoted by G_i. Due date of batch r is denoted by d_r. The job is on-time and thus no penalty is incurred if the job is delivered to the customer at the due date, while it will incur an earliness penalty if the job is arrived before the due date and a tardiness penalty if the job is arrived after the due date. The job will incur a holding penalty if the job is not delivered to the customer at its completion time. Moreover, associate with each job is a unit earliness penalty cost α_r, a unit tardiness penalty cost β_r and a unit inventory holding cost h_i. In addition, decision variables is denoted by x_ijk.

For the sake of convenience, we list the symbols used in the sequel as follows:

number of jobs

number of machines

number of batches

job index, i = 1, 2, ⋯ , n

machine index, j = 1, 2, ⋯ , m

batch index, r = 1, 2, ⋯ , b

the total number of vehicles

size of a job

capacity of each vehicle

ξ _i

processing time of job i (uncertain variable) Φ

_{ξ
_i}

uncertainty distribution ofξ_i

η _r

transportation time of a trip for jobs in batch r (uncertain variable)

Φ _{η
_r}

uncertainty distribution of η_r

cost of each vehicle taking a trip

C _i

completion time of job i

T _jk

completion time of a job in kth location on machine j

C _max

makespan

G _i

the actual departure time of job i

d _r

due date of batch r

α _r

delivery earliness penalty cost per unit time of batch r

β _r

delivery tardiness penalty cost per unit time of batch r

h _i

inventory holding cost per unit time of job i

Decision variables:

x_ijk:

1 if job i is processed on machine j at & the kth position, 0 otherwise.

Objective function:

$\begin{matrix} f = & \sum_{i = 1}^{n} \sum_{r = 1}^{b} δ_{ir} \cdot max {0, G_{i} - C_{i}} \cdot h_{i} \\ + \sum_{r = 1}^{b} max {0, max_{i \in B_{r}} {G_{i}} + η_{r} - d_{r}} \cdot β_{r} \\ + \sum_{r = 1}^{b} max {0, d_{r} - (max_{i \in B_{r}} {G_{i}} + η_{r})} \cdot α_{r} \\ + R \cdot (\sum_{r = 1}^{b} ⌈ \frac{v \cdot δ_{ir}}{Q} ⌉) \end{matrix}$ (2) where B_r indicates the set of indices of jobs which belong to batch r, δ_ir = 1 if job i is assigned to batch r and 0 otherwise.

In the right hand side of (2), the first term indicates holding cost; the second term indicates delivery tardiness penalty cost; the third term indicates delivery earliness penalty cost; and the fourth term indicates transportation cost.

The completion time T_jk of the kth job on machine j: $T_{j 1} = \sum_{i = 1}^{n} x_{ij 1} ξ_{i}$ when job k is assigned at the first position; otherwise, $T_{jk} = \sum_{i = 1}^{n} x_{ijk} ξ_{i} + T_{j, k - 1}$ , k = 2, 3, ⋯ , n.

The completion time C_i of job i: C_i = T_jk where j and k are unique indices that satisfy x_ijk = 1.

The earliest time of batch r being sent is equal to the completion time of the last job belonging to batch r if there are enough vehicles. Otherwise, jobs must wait for a vehicle back. That is, $G_{i} = max_{i^{'} \in B_{r}} {C_{i^{'}}}$ if there are enough vehicles and $G_{i} = max_{i^{'} \in B_{r}} {C_{i^{'}}} + H \cdot η_{r}$ for some positive integer H otherwise.

Due to the size of jobs and the capacity of vehicles, when there are enough vehicles, we can carry away all jobs at once; otherwise, we can carry away a portion of jobs and wait for the vehicles coming back to carry the remaining portion. A job can be delivered when all other jobs of the batch are completed. Otherwise, the job must wait for other jobs. A job should be stored from accomplishment to delivery. The period will emerge storage costs. Arrival time of the entire batch is equal to arrival time of the last job of the batch. The symbol ⌈· ⌉ indicates rounding up operation.

3.1 First stage

The formulation of the first stage of the problem is presented as follows: Model 1 $\min_{x_{i j k}} \min_{\bar{f}} \bar{f}$ (3) $\begin{array}{l} subject to \\ M {f \leq \bar{f}} \geq τ \end{array}$ (4) $\sum_{i = 1}^{n} x_{i j k} = 1, j = 1, \dots, m; k = 1, \dots, n$ (5) $\sum_{j = 1}^{m} x_{i j k} = 1, i = 1, \dots, n; k = 1, \dots, n$ (6) $\sum_{k = 1}^{n} x_{i j k} = 1, i = 1, \dots, n; j = 1, \dots, m$ (7)

The model is to minimize the pessimistic value $\bar{f}$ of f subject to some constraints. Equation (4) ensures that the total cost will be less than or equal to $\bar{f}$ on τ confidence level. Equation (5) ensures that there is only one job in a location on a machine. Equation (6) ensures that each job can appear on only one machine for a position. Equation (7) ensures that each job can occur in only one position on a machine.

Because objective function f is not a monotone function with respect to ξ_i and η_r for i = 1, 2, ⋯ , n and r = 1, 2, ⋯ , b, so inequality (4) can not be transformed to a crisp form by using the inverse uncertain distribution method [18]. We will employ the uncertain simulation [32] to obtain the result of $min \bar{f}$ . Let ξ = (ξ₁, ξ₂, ⋯ , ξ_n) be a common uncertain vector where ξ_i is a common uncertain variable with continuous uncertainty distribution Φ_i (x) for i = 1, 2, ⋯ , n, and g : Rⁿ → R be a Borel function. The following algorithm is used to simulate the pessimistic value: $g_{inf} = inf {r ∣ M {g (ξ) \leq r} \geq α}$ where α ∈ (0, 1) is a predetermined confidence level.

Algorithm (Uncertain simulation for g_inf) [32]:

Step 1. Randomly generate $u_{k} = (γ_{k}^{(1)}, γ_{k}^{(2)}, \dots,$ $γ_{k}^{(n)})$ with $0 < Φ_{i} (γ_{k}^{(i)}) < 1$ , i = 1, 2, ⋯ , n, k = 1, 2, ⋯ , m .

Step 2. Set a = g (u₁) ∧ g (u₂) ∧ ⋯ ∧ g (u_m), b = g (u₁) ∨ g (u₂) ∨ ⋯ ∨ g (u_m) .

Step 3. Set r = (a + b)/2 .

Step 4. If $M {g (ξ) \leq r} \geq α$ , then b ← r .

Step 5. If $M {g (ξ) \leq r} < α$ , then a ← r .

Step 6. Repeat the third to fifth steps until b - a < ɛ for a sufficiently small number ɛ.

Step 7. g_inf = (a + b)/2 .

In the Steps 4 and 5, the $M {g (ξ) \leq r}$ can be obtained according to the Algorithm (Uncertain simulation for uncertainty distribution) in [32].

3.2 Second stage

After scheduling at the first stage, completion time of the jobs in each batch has been determined. At the second stage, completed jobs in a batch can be delivered according to the weight of each batch. To begin with, it is required to determine the number of vehicles needed and the number of remaining vehicles at the beginning of each batch. Assume that N_vehicle represents the number of vehicles needed and A_batch represents the total amount of jobs of the current batch. Then $N_{vehicle} = ⌈ \frac{A_{batch}}{Q} ⌉$ . If $w = \frac{N_{vehicle}}{M} = 1$ , we carry away all jobs at once; if w < 1, we use M - N_vehicle to represent the number of remaining vehicles; if w > 1, we use $⌈ \frac{A_{batch}}{Q \cdot M} ⌉$ to represent the number of required transportation trip. The number of vehicles have three kinds of cases. Case 1: there are no idle vehicles. Delivery time of jobs is equal to the earliest free time of vehicles; Case 2: there are idle vehicles and w ≤ 1. All of jobs can be transported to the destination at once, and arrival time is equal to start time plus transportation time; Case 3: there is idle vehicles. But vehicles can not take all jobs away one-time. There is a part of jobs must wait for other vehicles. Then we can schedule vehicles at the second stage to ensure the maximum utilization rate. When there are idle vehicles but can not carry away all jobs one time, we transport the job which has higher inventory holding cost at the first, then transport the rest of these jobs.

4 Solution methodology

The complexity of production and transportation problem is NP-complete [10]. Due to the complexity, heuristic approaches are used to generate good approximate solutions. Genetic algorithm (GA) and simulated annealing algorithm (SA) have been shown to be promising technique by many researchers for solving supply chain scheduling problems. GA can reach the optimization solution quickly, but it is liable to be trapped in a local optima. SA can jump out of the local optimization and search for the best solution. The hybrid algorithm (GSA) combining GA and SA may also be an efficient method for optimization problems. In this paper, we use GSA to deal with the first stage of the problem. The second stage of the problem can be handled according to the method described in Section 3.2. Before using GSA, we must use the aforementioned uncertain simulation algorithm to approximate the objection function $min \bar{f}$ . Most researchers are familiar with GA and SA. Now we list the steps of GSA.

Step 1. Solution representation

Solution is expressed as a string of real numbers. The length of a string of real numbers is equal to the number of jobs. The integral part of each real number is generated randomly from 1 to m (the number of machines). The fractional part is generated from the interval (0, 1), which is used to determine the sequence of jobs. An example of solution representation with 2 machines and 4 jobs is shown in Fig. 2. The decoding result of the solution is shown in Fig. 3. The integral part means the index of machine. If some real numbers have identical integral part, sort these real numbers in ascending order. The order of corresponding jobs is determined by sequence of real numbers. If two real numbers are equal, sequence of jobs will be randomly generated.

Step 2. Initialization

Generate chromosomes randomly in a range according to weight of batch.

Step 3. Calculate the target

Calculate objective function value of each chromosome by uncertain simulation.

Step 4. Selection procedure

Use roulette wheel method to select the parent chromosomes.

Step 5. Crossover and mutation

In the crossover operation, we select a position from parents’ genes randomly and then exchange all genes before this position. An example of crossover operation as shown in Fig. 4. In the mutation operation, we randomly select two position of a chromosome and exchange their values. An example of mutation operation as shown in Fig. 5.

Step 6. Update optimal value

We can obtain a neighbor of current optimal solution by choosing one position at random and alter its current value according to upper and lower bounds. Let Z_cur be optimal value of the objective function for the current solution and Z_new be the value of the objective function for the new candidate. If Z_new ≤ Z_cur, then new solution is accepted as the current solution. If Z_new > Z_cur, the new solution is accepted with probability ${Prob}_{accept} = exp (\frac{- (Z_{new} - Z_{cur})}{T_{t}})$ , where t is the iteration number and T_t is the current temperature.

Step 7. Update temperature

GSA starts with an initial temperature T₀ = T_max. And the temperature decreases according to T_t = T_max0 . 95^t.

Step 8. Stopping rule

If the maximal generation is reached, then stop; otherwise, go to Step 3 to perform another iteration.

5 A Numerical example

Because our problem is different from the traditional deterministic problem, we can’t solve our problem with general heuristic method directly. First, we must put it into deterministic model. But the objective function of our problem is not a monotone function with respect to uncertain variables, so we can’t transform our model to a crisp form by the inverse uncertain distribution method in [18]. We may deal with our objective function by uncertain simulation [32]. Then our problem can be solved by heuristic algorithm after being transformed to a crisp form. In order to test the efficiency of heuristic method, we used GSA, HGA [28] and GA [26] deal with our problem, respectively.

Assume that there are 5 machines and 20 jobs (denote job i by J_i for i = 1, 2, ⋯ , 20). Consider six batches as shown in Table 1.

The weights of batches are from batch 1 to batch 6 in descending order. The time is specified in minutes. The cost is specified in yuan. Processing time of job i is assumed to be a linear uncertain variable: $ξ_{i} \sim L (i, i + 1)$ for i = 1, 2, ⋯ , 20, transportation time of a trip for jobs in the rth batch is assumed to be a linear uncertain variable: $η_{r} \sim L (r, (r + 1)^{2})$ , the size of each job is 10, the number of vehicles is 3, the capacity of each vehicle is 10, a transportation cost of a vehicle for a trip is 100.

Due date d = [20, 25, 70, 100, 130, 200], inventory holding cost h = [20, 15, 20, 25, 30, 10, 35, 30, 25, 20,10, 30, 25, 35, 30, 20, 30, 15, 10, 20], delivery earliness penalty cost α_r = [30, 20, 35, 40, 25, 30] , delivery tardiness penalty cost β_r = [40, 30, 20, 25, 35, 15] , confidence levels τ = 0.8.

Population size of GSA is 30. Mutation rate is 0.25. Crossover rate is 0.75. T_max is 100. Maximum iteration number is 200.

The optimal objective value obtained by running GSA is 16562. The processing and transport time of the case are come from linear distribution of uncertain variable:

ξ = [1.06, 2.09, 3.27, 4.40, 5.48, 6.90, 7.60, 8.32,

9.48, 10.58, 11.16, 12.83, 13.94, 14.60, 15.03,

16.81, 17.61, 18.70, 19.10, 20.43],

η = [2.13, 3.16, 13.83, 21.61, 19, 47.13] .

Tables 2–4 list the optimal values of the decision variables x_ijk of Model 1 by using GSA.

First stage:

According to the optimal decision variables, we sort jobs on parallel machines. It is assumed that there is no time gap among these jobs on each machine. In other words, the optimal scheduling plan is

Machine 1: 6 → 8 →11 → 18

Machine 2: 4 → 10 → 13 → 17

Machine 3: 2 → 7 →12 → 20

Machine 4: 5 → 1 →14 → 16

Machine 5: 3 → 9 →15 → 19

Second stage:

Having completed scheduling at the first stage, jobs must be sent in batches. Assume that the earliest beginning processing time at the first stage is 8 a.m. The processing time of a job and transportation time are considered at confidence level τ = 0.8.

Batch 1 will be distributed at the first. Initially there are three idle vehicles. Inventory holding cost of {J₁, J₃, J₄} are higher than J₂, so we transport {J₁, J₃, J₄} preferentially. Because makespan of batch 1 is 6′31″, so the earliest transporting time of {J₁, J₃, J₄} is 8 : 06′31″. After three vehicles come back, use one vehicle to transport J₂. Departure time of J₂ is 8 : 14′15″. The arrival time of whole batch is 8 : 16′21″.

Makespan of batch 2 is 9′42″. That is to say that batch 2 can be delivered after 8 : 09′42″. However, the earliest vehicles are usable at 8 : 14′15″ and only have two. Thus the delivery time of {J₅, J₇} is 8 : 14′15″. The vehicle which transports J₂ in batch 1 will free at 8 : 18′27″. So J₆ can be delivered at 8 : 18′27″. The arrival time of whole batch is 8 : 21′36″.

Makespan of batch 3 is 26′24″. Batch 3 can be distributed after 8 : 26′24″. The earliest two idle vehicles which transport {J₅, J₇} at batch 2 are valid until 8 : 17′21″. The departure time of jobs {J₈, J₁₂} is 8 : 26′24″. The vehicle which transports J₆ in batch 2 will free at 8 : 21′33″. Departure time of job J₉ is 8 : 26′24″. The free time of the two vehicles which transport {J₈, J₁₂} is 8 : 40′02″. So the departure time of jobs J₁₀ and J₁₁ is 8 : 40′02″. The arrival time of whole batch is 8 : 53′50″.

Makespan of batch 4 is 29′. Batch 4 can be distributed after 8 : 29′. The earliest two vehicles are valid until 8 : 44′. Then the departure time of job J₁₄ is 8 : 44′. The vehicle which transports {J₁₀, J₁₁} are usable at 9 : 07′38″. So departure time of job J₁₃ is 9 : 07′38″. A vehicle free after 9 : 07′38″. The arrival time of whole batch is 9 : 29′14″.

Makespan of batch 5 is 46′36″. Batch 5 can be distributed after 8 : 46′36″. The earliest idle vehicle is valid until 9 : 27′12″ at batch 4, so departure time of job J₁₅ is 9 : 27′12″. The vehicle which transports J₁₄ at batch 4 is valid until 9 : 27′12″. Departure time of job J₁₇ is 9 : 27′12″. The vehicle which transports J₁₃ at batch 4 is valid until 9 : 50′12′. Departure time of job J₁₆ is 9 : 50′12″. The arrival time of whole batch is 10 : 09′.

Makespan of batch 6 is 46′54″. Batch 6 can be distributed after 8 : 46′54″. The earliest idle vehicle which transports J₁₅ at batch 5 is valid until 10 : 05′12″. Departure time of job J₂₀ is 10 : 05′12″. Departure time of job J₁₈ is 10 : 28′12″. Departure time of job J₁₉ is 10 : 28′12″. The arrival time of whole batch is 11 : 15′18″. The frequency of total schedule vehicle is 20.

Finally, we use Tables 5–7 to schedule the transportation start time in the second stage.

To compare the proposed method with the existing methods of HGA [28] and GA [26], we assume that the parameters of HGA [28] and GA [26] are the same as our algorithm.

The objective value obtained by running uncertain simulation and HGA is 16852. The scheduling result at the first stage is

Machine 1: 7 → 4 →15 → 17

Machine 2: 3 → 10 → 13 → 18

Machine 3: 5 → 6 →12 → 14

Machine 4: 8 → 1 →16 → 20

Machine 5: 2 → 9 →11 → 19

The objective value obtained by running uncertain simulation and GA is 16963. The scheduling result at the first stage is

Machine 1: 5 → 9 →13 → 14

Machine 2: 6 → 8 →10 → 19

Machine 3: 1 → 12 → 15 → 17

Machine 4: 7 → 3 →11 → 20

Machine 5: 4 → 2 →16 → 18

The results show that the three algorithms mix with uncertain simulation respectively are effective to deal with our problem, but GSA would reach a better result.

6 Conclusions

In this paper, we established a scheduling and distribution model for an uncertain two-stage supply chain problem. The model for the first stage problem was proposed in the form of pessimistic value under the assumptions that processing time and transportation time of the jobs are uncertain variables. The objective function includes holding cost, tardiness penalty cost, earliness penalty cost and transportation cost. It is the non-monotone function of processing time and transportation time. A novel method consisting of an uncertain simulation and GSA algorithm is used for solving the model. At last, we give a numerical example and make a comparison between our method and two methods in literatures. Results show that the three heuristic algorithms are effective and GSA is a better one.

Footnotes

Acknowledgement

This work is supported by the National Natural Science Foundation of China (No. 61273009).

References

Cakici

, Mason

and Kurz

, Multi-objective analysis of an integrated supply chain scheduling problem, International Journal of Production Research50(10) (2012), 2624–2638.

Chang

and Lee

, Machine scheduling with job delivery coordination, European Journal of Operational Research158(2) (2003), 470–487.

Chen

, Integrating production and transportation scheduling in a make-to-order environment, Ph.D. Thesis, Cornell University, New York, 2000.

Chen

, American option pricing formula for uncertain financial market, International Journal of Operations Research8(2) (2011), 32–37.

Cheng

, Kovalyov

and Lin

, Single machine scheduling to minimize batch delivery and job earliness penalties, SIAM Journal on Optimization7(2) (1997), 547–599.

Cheng

, Gordon

and Kovalyov

, Single machine scheduling with batch deliveries, European Journal of Operational Research94(2) (1996), 277–283.

Dawande

, Geismar

, Hall

and Sriskandarajah

, Supply chain scheduling: Distribution systems, Production and Operations Management15(2) (2006), 243–261.

Ding

and Zhu

, Two empirical uncertain models for project scheduling problem, Journal of the Operational Research Society66(9) (2015), 1470–1480.

Gordon

, Due date assignment and scheduling under special conditions on job processing, Information Control Problems in Manufacturing13(1) (2009), 522–527.

10.

Guo

and Cheng

, Parallel machine scheduling with batch delivery costs, International Journal of Production Economics68 (2000), 177–183.

11.

Hall

and Potts

, Supply chain scheduling: Batching and delivery, Operations Reseach51(4) (2003), 566–584.

12.

Hamidinia

, Khakabimamaghani

, Mazdeh

and Jafari

, A genetic algorithm for minimizing total tardiness/earliness of weighted jobs in a batched delivery system, Computers & Industrial Engineering62(1) (2012), 29–38.

13.

, Uncertain random time-cost trade-off problem, Journal of Uncertainty Analysis and Applications2(23) (2014), 1–10.

14.

Lei

, Liu

, Ruszczynski

and Park

, On the integrated production, inventory, and distribution routing problem, IIE Transactions38(11) (2006), 955–970.

15.

Liu

, Uncertainty Theory, 2nd Edition, Springer-Verlag, Berlin, 2007.

16.

Liu

, Some research problems in uncertainth theory, Journal of Uncertain Systems3(1) (2009), 3–10.

17.

Liu

, Theory and Practice of Uncertain Programming, 2nd Edition, Springer-Verlag, Berlin, 2009.

18.

Liu

, Uncertainty Theory: A Branch of Mathematics for Modeling Human Uncertainty, Springer-Verlag, Berlin, 2010.

19.

Liu

, Why is there a need for uncertainty theory?, Journal of Uncertain Systems6(1) (2012), 3–10.

20.

Liu

, Toward uncertain finance theory, Journal of Uncertainty Analysis and Applications1(1) (2013), 1–15.

21.

Mazdeh

, Shashaani

, Ashouri

and Hindi

, Single-machine batch scheduling minimizing weighted flow times and delivery costs, Applied Mathematical Modelling35(1) (2011), 563–570.

22.

Pei

, Pardalos

, Fan

, Yang

and Wang

, Application of an effective modified gravitational search algorithm for the coordinated scheduling problem in a twostage supply chain, The International Journal of Advanced Manufacturing Technology70(1-4) (2014), 335–348.

23.

Peng

and Iwamura

, A sufficient and necessary condition of uncertainty distribution, Journal of Interdisciplinary Mathematics13(3) (2010), 277–285.

24.

Sheng

and Gao

, Chance distribution of the maximum flow of uncertain random network, Journal of Uncertainty Analysis and Applications2(15) (2014), 1–14.

25.

Selvarajah

, Steiner

and Zhang

, Single machine batch scheduling with release times and delivery costs, Journal of Scheduling16(1) (2013), 69–79.

26.

Ullrich

, Integrated machine scheduling and vehicle routing with time windows, European Journal of Operational Research227(1) (2013), 152–165.

27.

Yao

and Chen

, A numerical method for solving uncertain differential equations, Journal Intelligent & Fuzzy Systems25(3) (2013), 825–832.

28.

You

and Hsieh

, A heuristic approach to a single stage assembly problem with transportation allocation, Applied Mathematics and Computation218(22) (2012), 11100–11111.

29.

Yuan

, A note on the complexity of single-machine scheduling with a common due date, earliness-tardiness, and batch delivery costs, European Journal of Operational Research94(1) (1996), 203–205.

30.

Zegordi

, Abadi

and Nia

, A novel genetic algorithm for solving production and transportation scheduling in a two-stage supply chain, Computers & Industrial Engineering58(3) (2010), 373–381.

31.

Zhu

, Uncertain optimal control with application to a portfolio selection model, Cybernetics and Systems41(7) (2010), 535–547.

32.

Zhu

, Functions of uncertain variables and uncertain programming, Journal of Uncertain Systems6(4) (2012), 278–288.