Hybrid stochastic and robust optimization model for lot-sizing and scheduling problems under uncertainties

Hybrid stochastic and robust optimization model for lot-sizing and scheduling problems under uncertainties Journal Pre-proof Hybrid stochastic and robust optimization model for lot-sizing and scheduling problems under uncertainties Zhengyang Hu, Guiping Hu PII: DOI: Reference: S0377-2217(19)31070-7 https://doi.org/10.1016/j.ejor.2019.12.030 EOR 16238 To appear in: European Journal of Operational Research Received date: Accepted date: 21 December 2018 20 December 2019 Please cite this article as: Zhengyang Hu, Guiping Hu, Hybrid stochastic and robust optimization model for lot-sizing and scheduling problems under uncertainties, European Journal of Operational Research (2019), doi: https://doi.org/10.1016/j.ejor.2019.12.030 This is a PDF file of an article that has undergone enhancements after acceptance, such as the addition of a cover page and metadata, and formatting for readability, but it is not yet the definitive version of record. This version will undergo additional copyediting, typesetting and review before it is published in its final form, but we are providing this version to give early visibility of the article. Please note that, during the production process, errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain. © 2019 Elsevier B.V. All rights reserved. ฀฀฀฀฀฀฀฀฀ Hybrid stochastic and robust optimization for lot-sizing and scheduling problem under uncertainty Zhengyang Hu and Guiping Hu ∗ Department of Industrial and Manufacturing Systems Engineering, Iowa State University, Ames, IA 50011, United States ∗ Corresponding author: Guiping Hu (gphu@iastate.edu) 1 ฀฀฀฀฀฀฀฀฀ Hybrid stochastic and robust optimization model for lot-sizing and scheduling problems under uncertainties Abstract Uncertainty is among the significant concerns in production scheduling. It has become increasingly important to take uncertainties into consideration for lot-sizing and scheduling. In this paper, we adopt the Hybrid Stochastic and Robust Optimization (HSRO) approach in lot-sizing and scheduling problems in which suppliers have the flexibility of satisfying a fraction of demand based on the market and their policies. Two types of uncertainties have been considered simultaneously: demand and overtime processing cost. Robust optimization is adopted for uncertain demand and Sample Average Approximation (SAA) technique is applied to solve the stochastic program for uncertain overtime processing cost. Numerical results based on a manufacturing company has been conducted to not only validate the proposed hybrid model but also quantitatively demonstrate the merit of our approach. Sample size stability test and sensitivity analyses on various parameters have also been conducted. Keywords: Supply chain management, Stochastic programming, Robust optimization, Lot-sizing and scheduling, Automotive industry 1 ฀฀฀฀฀฀฀฀฀ 1 Introduction Efficient and robust production planning is essential for manufacturing companies to stay competitive. Tomotani pointed out that lot-sizing and scheduling are closely related that both decisions have to be made simultaneously to avoid sub-optimal results arising from considering them separately [40]. As pointed out by Yang et al., the lot-sizing and scheduling decisions are challenging due to the various uncertainties including material shortage, machine breakdowns and demand fluctuation [41]. In addition, replenishing inventory, seeking a new material supplier, purchasing and rearranging machines can be time consuming, so it is almost impossible to make timely response and adjustments to system oscillations [17, 34]. Therefore, the design of lot-sizing and scheduling system must be robust to deal with uncertainties in the production processes. Scenario-based stochastic programming approach has been gaining popularity in studying uncertainties in lot-sizing and scheduling problems. The probability distributions for the uncertain parameters are estimated and then scenarios are generated based on the modeling assumptions. Tas et al. considered stochastic setup time which follows a Gamma distribution and SAA technique was used to solve this stochastic programming problem. In addition, two heuristic algorithms were developed and evaluated in the paper [39]. Hu and Hu studied a two-stage stochastic programming approach under demand uncertainty and later, extended to multi-stage stochastic programing model. The moment matching method was used to generate scenarios [20, 21, 22]. However, this approach is only suitable for situations where the underlying probability distribution of the uncertain variable is known. The major drawbacks of this approach are: First, in some real-world applications, the decision makers may not have enough historical data to fit accurate probability distribution functions for the uncertain parameters. For example, it is almost impossible to predict future demand for a new product due to lack of historical data [25, 37]. One good example is to predict future demand for the fashion industry. A lot of fashion products 2 ฀฀฀฀฀฀฀฀฀ tend to be unique and hence there is not much historical data for fitting an appropriate distribution. Secondly, good approximation for continuous distributions requires a large number of scenarios. In a multi-period problem, overall scenario size grows exponentially with the number of scenarios in each time period. Although techniques like Fast Forward Selection (FFS) and Simultaneous Backward Reduction (SBR) have been used to reduce scenario sample size, important information may be lost during the process [19]. On the other hand, if scenario sample size is limited due to computational complexity, the accuracy of prediction for the future stage can be restricted and solutions may not be feasible for some extreme realizations of uncertainties. Finally, different decision makers may have different attitudes toward uncertainty. Scenario-based stochastic programming approach focuses on the expected value which assumes the decision makers care about the average performance of each scenario. However, in some cases, the decision makers can be more concerned about the worst case scenario than the average scenario. To address these drawbacks, Robust Optimization (RO) has been utilized as an alternative technique to deal with uncertainties. Soyster assumed that all uncertain parameters would take their worst-case values within sets and this approach was perceived as overly conservative for practical implementation [38]. In the mid-1990s, the shortcoming of over-conservatism was addressed by constrainting the uncertain parameters to belong to ellipsoidal uncertainty set. This approach only considers outcomes that are likely to happen but results in a non-linear robust counterpart [3, 4, 5]. More recently, Bertsimas et al. proposed a new robust optimization approach that overcame the issue of high complexity when formulating the robust counterpart. Concretely, the robust counterpart of a linear programming problem is still a linear programming problem [6, 7]. The major advantages of using RO are: First, RO is not based on a probabilistic theory and does not requires extensive amount of historical data to support the parameter estimation. In other words, this approach does not need knowledge of specific probability distribution for the uncertain parameters. Second, this RO approach is computationally tractable because of linear 3 ฀฀฀฀฀฀฀฀฀ robust counterpart. Curcio et al. considered adjustable robust optimization for lot-sizing and scheduling problem under multi-stage demand uncertainty. The results showed that their algorithm is much better than the deterministic model [13]. Diaz et al. considered a production planning problem with uncertain operating and environment conditions. Single and multi-objective formulation for robust optimization were studied. Results show a significant correlation between robustness and sample size in the performance evaluation [15]. There are two types of uncertainties considered in this paper. Historical data is available and reliable to generate scenarios for overtime processing cost and hence we adopt the stochastic programming approach for this type of uncertainty. Due to the unpredictable characteristic, we adopt robust optimization to address demand uncertainty. The major challenge that many companies face is to predict accurate demand volume. If the demand volume is underestimated, customer demand may not be satisfied. On the other hand, if the demand volume is overestimated, then unnecessary inventory cost will incur. Overtime production is highly related to the demand prediction since it helps companies fulfill unpredicted demand in the peak season. The need to study the integration of these two types of uncertainties has been justified by Davis [14]. The author identified three major supply chain uncertainties which are supply uncertainty, process uncertainty, and demand uncertainty. Supply uncertainty depends on suppliers’ reliability. Process uncertainty is related to the production process. Demand uncertainty often arises from inaccurate forecast and market fluctuation. Govindan et al. pointed out that demand quantity, production and transportation costs are the most frequently studied uncertainties in supply chain [18]. Li and Hu studied lot-sizing and scheduling problem under demand and workforce uncertainties [29]. Ramaraj et al. considered the same production problem with demand and yield uncertainties [35, 36]. The difference is that they adopt stochastic programming approach for all uncertain parameters in the model. Alem et al. formulated stochastic and robust models separately and compared the results with Monte 4 ฀฀฀฀฀฀฀฀฀ Carlo simulation. They provided guidelines for decision makers to assess a priori approach based on their preferences [1]. Keyvanshokooh et al. used hybrid approach in the context of a closed loop supply chain problem where uncertain transportation cost was handled by stochastic programming approach and robust optimization was adopted to deal with demand uncertainty [25]. Similar hybrid approach can be found in inventory control and bidding strategy [31, 33]. In this paper, we adopt the HSRO approach introduced by Keyvanshokooh et al. in the context of lot-sizing and scheduling problems [25]. To the best of authors’ knowledge, no existing papers have studied the integration of demand and overtime processing cost uncertainties using HSRO approach in the application of lot-sizing and scheduling. The major contributions of this study are listed as follows: • A lot-sizing and scheduling framework is proposed to study the integration of multiple uncertainties. It provides the flexibility to satisfy only a fraction of customers according to market competition and company’s policies. • We adopt the HSRO approach of Keyvanshokooh et al. in a new context to study two different types of uncertainties simultaneously including stochastic programming approach for overtime processing cost uncertainty and robust optimization for demand uncertainty [25]. • SAA technique is applied to solve the stochastic program which considers overtime processing cost uncertainty. The remainder of this paper is organized as follows: In section 2, we review the basic concept of robust optimization and present the proposed HSRO approach for lotsizing and scheduling formulation. section 3 describes the computational results such as scenario stability test and sensitivity analyses on various parameters that provide managerial insights on the proposed model. Finally, section 4 summarizes the paper and suggests future research directions. 5 ฀฀฀฀฀฀฀฀฀ 2 Problem description In this paper, we consider a multi-product, multi-period, and capacitated lot-sizing and scheduling problem. Raw material are collected from up-stream factories and final products are shipped to customers and down-stream factories. There are multiple steps in the manufacturing process, such as welding, casting and molding. In this study, we assume that the products are perishable and hence the excess amount of inventory in one period cannot be carried over and used to fulfill demand in the future. There are several types of capacities such as maximum available time on the machine, maximum production batch size in both regular time and overtime. The goal is to design a robust production plan to maximize the profits given capacity constraints under demand and overtime processing cost uncertainties. In the HSRO formulation, we introduce the penalty cost for unmet demand and surplus cost for excess production. On one extreme, if there is little market competition and company has significant market share, then losing small fraction of customers can be affordable. On the other extreme, if the market is very competitive then losing any customers may be unacceptable. The proposed hybrid model provides the flexibility to design an optimal plan under any circumstance between these two extreme situations. Penalty cost for unmet demand would be low if market is not competitive. On the other hand, penalty cost can be set to high values if customer satisfaction is critical. Introducing the penalty and surplus costs in the production balance the customer satisfaction and inventory resource requirement [11]. The objective of this study is to design a robust production plan under two different types of uncertainties: One for overtime processing cost and the other for customer demand. We assume that company has complete knowledge of the underlying probability distribution of the uncertain overtime processing cost, so stochastic programming can be used to model this type of uncertainty. On the other hand, predicting the probability 6 ฀฀฀฀฀฀฀฀฀ density function of demand is very challenging due to several reasons. Demand could be influenced by unpredicted situations such as a competitor launches a new product or market fluctuation. Even if market is stable, predicting demand scenarios for new products is very difficult due to insufficient information, therefore, we adopt robust optimization to model demand uncertainty. 2.1 Mathematical notations All notations for the mathematical formulation are listed in Table 1. 2.2 Robust optimization We follow Bertsimas and Thiele to formulate the robust optimization component [9]. Consider a linear programming problem where C ∈ Rn , b ∈ Rm and A is a m ∗ n matrix. Min Cx s.t. Ax ≤ b, x≥0 (1) Without loss of generality, uncertainty is assumed to affect only the elements in the matrix A. We consider a particular row vector i of A and define Ji as the set of uncertain coefficients in that row. To simplify the exposition, every coefficient aij , j ∈ Ji is subject to uncertainty and modeled as independent random variable which belongs to a symmetric interval [aˆij − ∆aij , aˆij + ∆aij ] where ∆aij is the maximum deviation of the uncertain element aij and aˆij is the nominal value. This is reflected in the following formulation of the robust counterpart of Equation 1. Min Cx s.t. max ∀aij ∈Ji A scaled deviation parameter zij = X j aij −âij ∆aij aij xj ! ≤ bi ∀i, x ≥ 0 (2) is introduced. It should be noted that ∆aij , aˆij , aij represent the maximum deviation, nominal value, uncertain parameter for the 7 ฀฀฀฀฀฀฀฀฀ Table 1: Notations used in the HSRO formulation Sets: i j t s Parameters: Di,t s Di,t s Dˆi,t s ∆Di,t + s ηDi,t s ηDi,t − capt ci pti pri 1, 2 · · · 1, 2 · · · 1, 2 · · · 1, 2 · · · ,N ,N ,T ,S Set of raw material comes from up-stream suppliers i and j are alias Set of time periods in the production horizon Set of scenarios for overtime processing cost Stochastic demand for product i at time period t Uncertain demand for product i at time period t in scenario s Nominal demand for product i at time period t in scenario s Maximum demand deviation from nominal value for product i at time period t in scenario s Positive deviation percentage from nominal demand for product i at time period t in scenario s Negative deviation percentage from nominal demand for product i at time period t in scenario s Overall time availability on the machine at time period t Selling price for product i Manufacturing time for product i Regular time manufacturing cost for product i poi posi qi,t sci,j Stochastic overtime manufacturing cost for product i Overtime manufacturing cost for product i in scenario s The maximum regular time batch size for product i at time period t Overall preparation cost when a setup changeover from two different products i, j is taken place sti,j Overall preparation time when a setup changeover from two different products i, j is taken place α Ratio of regular and overtime production batch size N Total amount of products that come from different families probs The probability of scenario s pen Penalty cost per unit of unmet demand sur Surplus cost per unit of unnecessary production ΓsD Overall demand deviation budget for scenario s Decision Variables: Xi,t Production batch size in the regular time for product i at time period t s Oi,t Production batch size in the overtime production for product i at time period t in scenario s Yi,j,t Binary variable which takes value 1 if setup changeover from two different products i, j is taken place at time period t Zi,t Binary variable which takes value 1 if setup for product i is carried over from time period t − 1 to t Vi,t Production sequence at time period t. It takes value from 1 to N 8 ฀฀฀฀฀฀฀฀฀ random variable, respectively. A budget parameter Γi is used to setup a boundary for aggregated scaled deviation. This insight can be incorporated in mathematical terms as follows: X |zij | ≤ Γi , ∀i (3) j∈Ji Based on the description above, set Ji is defined as Ji = {aij aij = aˆij +∆aij zij , ∀i, j, z ∈ P Ψ} where Ψ = {z |zij | ≤ 1, j∈Ji |zij | ≤ Γi , ∀i}. Reformulating each constraint i as P P P P j ∆aij zij xj , the bilevel robust Equation 2 j aˆij xj + j (aˆij + ∆aij zij )xj = j aij xj = can be transformed into: Min Cx s.t. X aˆij xj + Max zi ∈Ψi j∈Jj The lower level problem Maxzi ∈Ψi i is equivalent to Equation 5: Max X ∆aij zij x∗j j∈Ji P s.t. j∈Ji X ∆aij zij xj ≤ bi ∀i, x ≥ 0 (4) j∈Ji ∆aij zij xj for a given x∗ and constraint index X zij ≤ Γi 0 ≤ zij ≤ 1 ∀j ∈ Ji (5) j∈Ji By introducing the dual variables λi and µij , the dual of Equation 5 can be expressed as: Min Γi λi + X µij s.t. λi + µij ≥ ∆aij x∗j , µij , λi ≥ 0 ∀i, Ji (6) j∈Ji By applying the dual Equation 6 to Equation 4, the robust counterpart of Equation 1 9 ฀฀฀฀฀฀฀฀฀ is obtained: Min Cx X X s.t. aˆij xj − Γi λi − µij ≤ bi ∀i j∈Ji (7) j∈Ji λi + µij ≥ ∆aij xj ∀i, j ∈ Ji µij , λi ≥ 0 ∀i, j ∈ Ji For a given solution x∗j , the probability of constraint violation can be calculated by [8]: Pr X j∈Ji aij x∗j > bi  ≤1−Φ   p Γi − 1 / |Ji | ! (8) Reversely, we can setup a maximum violation probability ǫi and Equation 9 gives us the minimum budget Γi to maintain that particular level of violation probability. Γi ≥ 1 − Φ−1 (1 − ǫi ) 2.3 p |Ji | (9) Hybrid stochastic and robust optimization model To introduce our hybrid model, we first introduce a two-stage stochastic programming model which serves as a baseline and the robust optimization is then built on top of the two-stage stochastic programming model to construct the full hybrid model. 2.3.1 Baseline model with uncertainty in overtime processing costs In the two-stage stochastic programming model, demands are set to their nominal values while uncertainties in the overtime processing costs are incorporated with the stochastic programming approach. Sample Average Approximate approach was adopted to generate scenarios. In the lot-sizing and scheduling process, first-stage decision variables identify 10 ฀฀฀฀฀฀฀฀฀ the baseline production i,e. raw material purchase and regular time production planning while second-stage decision variables define possible recourse like overtime production and compensatory actions [2]. The first-stage decisions have to be made in the presence of uncertainties meaning the regular time production decisions are made before we observe the realization of uncertain overtime processing cost [12]. In the objective function, the goal is to maximize profit based on revenue and production costs. Production costs include regular time production cost, overtime production cost, and setup cost. Max ζ = N X T X i=1 t=1 N X − ci ∗ Xi,t + S N X T X X s probs ∗ Oi,t ∗ ci − sci,j ∗ Yi,j,t − pri ∗ Xi,t i=1 t=1 i=1 t=1 s=1 N X T X N X T X S N X T X X (10) s probs ∗ posi ∗ Oi,t i=1 t=1 s=1 i=1 i6=j j=1 t=1 s s Xi,t + Oi,t = Dˆi,t ∀i, t, s (11) Constraint (11) ensures that overall production equals to the nominal demand. Regular time production Xi,t is first-stage decision since it has to be determined before we realize the uncertainty. Overtime production is scenario-based second stage decision since they are served as compensatory actions. In addition, we assume that the products are perishable meaning excess inventory cannot be carried over to satisfy future demand. Backlog is not allowed and hence unmet demand will be lost. In the hybrid model, des mands are random, which can vary in an predetermined interval with mean Dˆi,t and s standard deviation ∆Di,t . Therefore, the overall production can be either greater than or less than the actual realization of demands. This will be discussed in detail in the hybrid model section. Xi,t ≤ qi,t ∗ (Zi,t + N X Yj,i,t ) ∀i, t (12) j6=i Constraint (12) ensures that regular production quantity Xi,t cannot exceed regular 11 ฀฀฀฀฀฀฀฀฀ time production batch size capacity qi,t . Since the setup of a particular product can take P place at most once in each time period meaning Zi,t + N j6=i Yj,i,t ≤ 1. For example, if the setup of product i1 is carried over from t1 to t2 , then Zi1 ,t2 = 1. Furthermore, we know P that N j6=i1 Yj,i1 ,t2 = 0 based on the setup assumption. On the other hand, for any product i2 6= i1 , Zi2 ,t2 = 0 since the carried over setup in time period t2 is not i2 . It provides the opportunity to setup product i2 during this time period. If product i2 is scheduled to be manufactured after product i1 , then we know that Yi1 ,i2 ,t2 = 1. N X pti ∗ Xi,t + i=1 N N X X sti,j ∗ Yi,j,t ≤ capt ∀t (13) i=1 i6=j j=1 Constraint (13) states the overall regular time capacity. The left side represents the total production time includes regular time processing time and setup changeover time. It is assumed that setup is only required when products change from different families. s Oi,t ≤ α ∗ Xi,t ∀i, t, s (14) In constraint (14), overtime production quantity is bounded by the regular time production quantities. Government policy states that the overtime production batch size must be limited to 20% of the regular production batch size. Although overtime exceeds that threshold can be considered under extreme circumstances, but it could lead to injury, fatigue and reduce efficiency [42]. N X Zi,t = 1 ∀t (15) i=1 Xi,t = 0 ∀i, t = T (16) Constraint (15) indicates that only one setup can be carried over to the following time period. Constraint (16) states that no regular time production activity is allowed in the 12 ฀฀฀฀฀฀฀฀฀ dummy period (t = T ). Zi,t + N X j6=i Yj,i,t = Zi,t+1 + N X Yi,j,t ∀i, t = 1 · · · T − 1 (17) j6=i Constraint (17) indicates that setups have to happen in an equilibrium state. One of assumptions is that setups are preservable meaning the last setup in one time period can be carried over to the next period. For example, consider a situation where only products i1 are manufactured in time period t1 and the setup carried over from t0 is i1 , then Zi1 ,t1 = Zi1 ,t2 = 1, since the setup for product i1 is carried over from t1 to t2 . In addition, there is no setup changeover either from other products to i1 or from i1 to any PN PN other products. Therefore, j6=i Yj,i,t = j6=i Yi,j,t = 0 and equality condition is met. Consider a more complicated case where three different products are manufactured in t1 in sequence of i1 − i2 − i3 and the setup carried over from previous time period is i1 . We know Zi1 ,t1 = Zi3 ,t2 = 1 since i3 is the last product on the production line in time period t1 . In addition, Yi1 ,i2 ,t1 = Yi2 ,i3 ,t1 = 1 since two setup changeovers take place in time period t1 . In conclusion, Setups carried over from previous time period and changeover from other products represent flows going into a product node while setups carried over to the next time period and changeover to other products represent flows leaving a product node. Vj,t ≥ Vi,t + 1 − N ∗ (1 − Yi,j,t ) ∀i, j 6= i, t (18) One important assumption is that the setup for a particular product can take place at most once in each time period. Constraint (18) is a sub-tour elimination constraint. N is the number of products and Vi,t is the production sequence of product i in time period t. We explain how this type of constraints avoid sub-tour with an example in Figure 1. Assuming there are five different products to be manufactured in time period t 13 ฀฀฀฀฀฀฀฀฀ and production sequence is 1 − 2 − 3 − 4 − 5. Then N = 5, V1,t = 1, V2,t = 2, V3,t = 3, V4,t = 4 and V5,t = 5. Let’s introduce an extra node 0, and without loss of generality, we assume all feasible paths start from node 0 and end at node 0. Consider there is a sub-tour skipping node 0 and going directly from node 5 back to node 1, then there is a sub-tour. This situation is represented by the dash line in the figure. Since this is a 5-step sub-tour, if we sum up all Yi,j,t = 1, we will have 5 ≤ 4 which is impossible. If all feasible paths have to go through node 0, then there is no sub-tour in the graph. If product i is followed by product j, then we know Yi,j,t = 1 and Vj,t = Vi,t + 1. Otherwise, we have Vj,t − Vi,t ≥ 1 − N . Both inequalities are valid since, for decision variables V , the minimum value is 1 and maximum value is N . 2 3 1 4 5 0 Figure 1: Sub-tour elimination example 2.3.2 Hybrid model with uncertainties in overtime processing costs and demand In this paper, we adopt the HSRO approach proposed by Keyvanshokooh et al. in the lot-sizing and scheduling problem setting [25]. SAA technique is applied to stochastic overtime processing cost. SAA is a Monte Carlo simulation based approach and scenarios are generated by taking random IID samples from the given distribution [26]. Within 14 ฀฀฀฀฀฀฀฀฀ each scenario, we define polyhedral uncertainty sets for demand in each time period and for each product. The details are shown in Figure 2. There are |S| different realizations of stochastic overtime processing cost. Inside each scenario, we define a symmetric ins s terval for robust demand. The nominal demand Dˆi,t and maximum deviation ∆Di,t are s s predetermined. Uncertain demand Di,t is allowed to deviate from the Dˆi,t toward the worst-case value within that interval. Note that strips in Figure 2 are symmetric with respect to the nominal value, but the width of the strips do not need to be the same for different time periods or products. 'Dis,t Dis,t s=1 Overtime Processing Cost Dˆ is,t s=S t=1 t=2 t=3 t=1 t=2 t=3 Product N t=T Product 1 t=T Figure 2: Characterization of uncertain demand and overtime processing cost To develop the uncertainty sets for demand, we introduce the positive and negative deviation percentages from their nominal demands as follows: + s ηDi,t = s s Di,t − Dˆi,t s ∆Di,t s s if Di,t ≥ Dˆi,t s ηDi,t = − 15 s s Dˆi,t − Di,t s ∆Di,t s s if Di,t ≤ Dˆi,t (19) ฀฀฀฀฀฀฀฀฀ + s s Notice that at most one of ηDi,t and ηDi,t will be positive and the other one has to − be zero in (19). The polyhedral uncertainty sets of demand for each scenario of overtime processing costs can be represented by: JDS   s s− s+ s s− s s s+ s ˆ = Di,t Di,t = Di,t + ∆Di,t ∗ ηDi,t − ∆Di,t ∗ ηDi,t ∀i, t ∀ηDi,t , ηDi,t ∈ KD (20) where   X s+ s− s+ s− s+ s− s KD = ηDi,t , ηDi,t 0 ≤ ηDi,t ≤ 1, 0 ≤ ηDi,t ≤ 1, (ηDi,t + ηDi,t ) ≤ ΓD (21) i,t s+ In particular, if a realization of demand is above the mean value Dˆij , then ηDi,t is s strictly positive and ηDi,t equals to 0. Corresponding uncertainty sets of demand becomes   S s s s s+ s ˆ JD = Di,t Di,t = Di,t + ∆Di,t ∗ ηDi,t . Reversely, if a realization of demand is below the − s− s+ mean demand Dˆij , then ηDi,t is strictly positive and ηDi,t equals to 0. Corresponding   − S s s s s s uncertainty sets of demand becomes JD = Di,t Di,t = Dˆi,t − ∆Di,t ∗ ηDi,t . For a given overtime processing cost scenario, the dimension of those sets is |i| ∗ |t|, which means the budget for a given scenario (ΓsD ) can take any value between 0 to |i| ∗ |t|. If ΓsD equals to 0, then there is no protection against uncertain demand and our hybrid model is identical to the two-stage model. If ΓsD equals to |i| ∗ |t|, then corresponding constraint is fully protected. In the HSRO formulation, demands belong to some predefined intervals and the exact values are not known. If we include constraint (11) in the model formulation, then there is no guarantee that these constraints will be satisfied because the exact realization of uncertain parameter is unknown before measuring the optimal regular and overtime production quantities. On the other hand, solving the problem with constraint (11) will 16 ฀฀฀฀฀฀฀฀฀ end up with two results: (1) optimal production quantity is larger than the demand or (2) optimal production quantity is smaller than the demand. In both cases, constraint (11) is not satisfied and the decision making becomes infeasible. In order to avoid this situation, we decide to move this constraint to the objective function with new cost parameters. Penalty cost pen is introduced for one unit of unsatisfied demand and surplus cost sur is introduced for one unit of excess production over demand. These two new cost parameters provide the flexibility for the decision makers based on the market environment and policy. For example, if the products have large market share and there is not much competition, then the company can afford unmet demand and satisfy only proportion of orders. In this case, the decision maker can decrease penalty cost and increase surplus cost. Conversely, if the company is in a very competitive market and losing a customer incurs a high cost, then we should satisfy as much demand as possible since losing customers become expensive. In this case, we can increase penalty cost and decrease surplus cost. Our goal is to minimize the maximum amount of violation cost due to unbalanced production flow in constraint (11). In the HSRO formulation, this set of constraint is incorporated into the objective function with penalty cost pen and surplus cost sur. For a given overtime processing cost scenario s, the violation cost can be formulated as follows: V Cs (X, O) = Max s ∈J S Di,t D X i,t s s Di,t − Xi,t − Oi,t  ∗ pen, X i,t s s Xi,t + Oi,t − Di,t  ∗ sur  (22) s where Di,t is the random demand which belongs to the sets (20) and (21). In equation s (22), Xi,t is the regular time production decision and Oi,t is the overtime production   P s s decisions. The first term in the violation cost i,t Di,t − Xi,t − Oi,t is unmet demand   P s s and the second term in the violation cost i,t Xi,t + Oi,t − Di,t is excess production 17 ฀฀฀฀฀฀฀฀฀ over demand. For each overtime processing cost scenario, equation (22) is the worst-case result for violation cost and we want to minimize this violation cost. By introducing an auxiliary variable Ws for each overtime processing cost scenario, we transform the non-linear programming formulation to a linear programming formulation: Min V Cs (X, O) = Ws  X s s s ∈ JDS s.t. Di,t − Xi,t − Oi,t ∗ pen ≤ Ws , ∀Di,t i,t X Xi,t + s Oi,t − s Di,t i,t  ∗ sur ≤ Ws , s ∀Di,t ∈ (23) JDS Ws ≥ 0 For a given overtime processing cost scenario s, equation (23) should always be feasible for any realizations of uncertain demand within their polyhedral sets. Then we find robust counterparts for each constraints in (23): Max s ∈J S Di,t D X s Di,t − Xi,t − s Oi,t i,t  ∗ pen  ≤ Ws , (24) which can be rewritten as: Max− + s ,ηD s ∈K d ηDi,t i,t X s ∆Di,t s+ ηDi,t ∗ − s ∆Di,t ∗ s− ηDi,t i,t + X s Dˆi,t − Xi,t − i,t + s Oi,t   ∗ pen  (25) ∗ pen ≤ Ws s s ) and negative (ηDi,t ) percentage of deviation We optimize (25) over the positive (ηDi,t − from the nominal value. Let’s focus on the first term in the (25). It should be noted that the differences between the first term in (25) and (26) are: (1) maximization objective function has been changed to a minimization objective function; (2) we add a negative 18 ฀฀฀฀฀฀฀฀฀ sign in front of the coefficients in the objective function. + Min− s ,ηD s ∈K d ηDi,t i,t s.t. X − s ∆Di,t ∗ s+ ηDi,t + s ∆Di,t ∗ s− ηDi,t i,t  ∗ pen  + s 0 ≤ ηDi,t ≤1 s− ηDi,t 0≤ X (26) ≤1 s+ ηDi,t + s− ηDi,t i,t  + ≤ ΓsD s s Notice that the domains of ηDi,t and ηDi,t stay the same, therefore, the objective − value of (26) equals to the negative objective value of the first term in (25). We add an additional negative sign in front of the minimization function in (27) such that it is equivalent to the first term in (25). α1si,t , α2si,t and β1s are the dual variables for constraints that have been transformed into the standard form in (27). − + Min− s ,ηD s ∈K d ηDi,t i,t s.t. X − s ∆Di,t ∗ s+ ηDi,t + s ∆Di,t ∗ s− ηDi,t i,t +  ∗ pen s − ηDi,t ≥ −1 ∀i, t α1si,t s − ηDi,t ≥ −1 ∀i, t   X s+ s− − ηDi,t − ηDi,t ≥ −ΓsD α2si,t −  (27) β1s i,t + s s ≥0 , ηDi,t ηDi,t − Then we take the dual of formulation (27): Max α1si,t ,α2si,t ,β1s s.t.  − ΓsD s ∗ β1 − X α1si,t i,t + α2si,t  s − α1si,t − β1s ≤ −∆Di,t ∀i, t s − α2si,t − β1s ≤ ∆Di,t ∀i, t α1si,t , α2si,t , β1s ≥ 0 19 ∀i, t (28) ฀฀฀฀฀฀฀฀฀ In the linear programming formulation (28), all decision variables are non-negative s and ∆Di,t is the maximum deviation from the nominal value which is also a positive s constant. Therefore, −α2si,t − β1s ≤ ∆Di,t becomes redundant since it can be satisfied in any condition. Due to the strong duality theory, we substitute the objective function (28) without α2si,t in the constraint (25). The robust counterpart of first constraint in (23) can be expressed as follows: X s Dˆi,t − Xi,t − s Oi,t + α1si,t i,t s α1si,t + β1s ≥ ∆Di,t ∀i, t α1si,t , β1s ≥ 0 ∀i, t  + ΓsD ∗ β1 s  ∗ pen ≤ Ws (29) Applying the same procedure to the second constraint in (23), the other robust counterpart is obtained as follows: X Xi,t + s Oi,t s − Dˆi,t + α2si,t i,t s α2si,t + β2s ≥ ∆Di,t ∀i, t α2si,t , β2s ≥ 0 ∀i, t  + ΓsD ∗ β2 s  ∗ sur ≤ Ws (30) Finally, our HSRO formulation of the lot-sizing and scheduling design problem becomes: Max ζ = N X T X i=1 t=1 N X − ci ∗ Xi,t + S N X T X X i=1 t=1 s=1 N X T X sci,j ∗ Yi,j,t − − N X T X pri ∗ Xi,t i=1 t=1 S N X T X X s probs ∗ posi ∗ Oi,t (31) i=1 t=1 s=1 i=1 i6=j j=1 t=1 S X s probs ∗ Oi,t ∗ ci − probs ∗ Ws s=1 subject to constraints (12)-(18), (29) and (30). For each overtime processing cost 20 ฀฀฀฀฀฀฀฀฀ scenario in the hybrid model, the budget parameter ΓsD balances the trade-off between the level of robustness and the degree of conservativeness of the solution. As a consequence, larger ΓsD provides more protection and increases the level of robustness while smaller ΓsD results in higher probability of constraint violation. On the other hand, bigger ΓsD leads to lower expected profit since model allows more deviations toward the worst-case in their uncertainty sets. 3 Case study In order to demonstrate and validate the proposed hybrid model, we conduct a case study on braking equipment production. A manufacturing plant receives two different types of raw material from upstream and produces three different types of braking actuators P1 , P2 , and P3 . These products are directly supplied to customers and the goal is to identify the optimal production strategy such that the total system cost is minimized. 3.1 Data sources This case study considers three different final products (N = 3). According to government policy, ratio of overtime and regular time production batch sizes is set to 20% (α = 20%) meaning overtime production batch size can not bigger than 20% of regular time production batch size [32]. Setup changeovers are product-dependent and hence it is important to identify the optimal production sequence. Setup changeovers times are included in Table 2 and corresponding setup changeover costs are proportional to their setup changeover times [23]. Notice that setup time between products in the same families is 0 and setup changeover matrix is not symmetric due to the fact that setups are productdependent. Table 3 provides the nominal demand, processing time, and regular time processing cost. Nominal demands do not change over the planning horizon and the maximum 21 ฀฀฀฀฀฀฀฀฀ Table 2: Setup changeover times sti,j (mins) P1 P2 P3 P1 0 180 90 P2 270 0 180 P3 90 270 0 s Table 3: Parameters setup for Dˆi,t , pti , and pri s Nominal demand Dˆi,t (unit) Processing time pti (mins/unit) Regular time processing cost pri ($/unit) P1 467.25 6 254.08 P2 33.82 6.6 254.08 P3 149.7 7.2 254.08 demand deviations can vary from 5% to 30% of their nominal values. Overtime processing cost follows a normal distribution with mean equals to 1.5 times regular time processing cost and standard deviation equals to 10% of mean value [16]. Unit selling price, penalty cost, and surplus cost are proportional to regular time processing cost and more details on overall time capacity capt and maximum regular time batch size qi,t can be found in the Hu and Hu [21]. 3.2 Computational results In this section, we first present scenario stability test to validate that the scenario sample size is sufficient to generate stable objective function. The results are illustrated in Figure 3. Five different scenario sample sizes with 20 replications were analyzed. The key idea is that when several scenario samples with the same sample size are generated, the optimal value of objective function should be close if scenario sample size is sufficient enough to represent the distribution. Concretely, we generate 20 different scenario trees ηi ∀i = 1, · · · , 20 and sample sizes are identical across all scenario trees. Suppose we solve the model with each scenario trees generated, the optimal solutions are x∗i ∀i = 1, · · · , 20 and optimal objective values are f (x∗i , ηi ) ∀i = 1, · · · , 20. Stability indicates that f (x∗i , ηi ) ≈ f (x∗j , ηj ) ∀i, j ∈ {1, · · · , 20} and i 6= j. If this type of stability 22 ฀฀฀฀฀฀฀฀฀ Table 4: Model statistics Number of scenarios 20 50 80 100 150 Number of equations 1399 3349 5299 6599 9849 Number of variables 1426 3406 5386 6706 10006 is obtained, then the performance of optimal solution x∗ and f (x∗ , η) are independent of which scenario tree gets selected [24]. Scenario sample sizes = 20, 50, 80, 100, 150 have been tested and illustrated in Figure 3. Intuitively, when the scenario sample size is small, not enough information was collected to generate stable objective values. As scenario sample size increases, more information about the distribution becomes available and objective values become more stable. It can be observed from the figure that when scenario sample size is 20, the optimal objective values are highly variable and hence not stable. When scenario sample size is 150, the minimum objective value is 1,940,176 while the maximum value is 1,945,427. The lack of significant difference between different trials shows sample stability. The model statistics as a function of number of scenarios are shown in Table 4. The computational times are less than 1 minute and the size of problem increases linearly as a function of number of scenarios. From the CPU times, we can see that the number of scenarios does not increase the computational complexity dramatically. In order to investigate the factors that drive the major complexity of the model, the number of time periods and the number of products have been analyzed. The computational time increases linearly as a function of number of time periods. However, if we increase the number of products, the computational time increases much faster than other factors. The solid line in Figure 3 represents the objective value of the deterministic model in which demand and overtime processing cost are fixed and known. It should be noted that roughly half of the objective values of the stochastic models are above that red line 23 ฀฀฀฀฀฀฀฀฀ 1.955 10 6 S = 20 S = 50 S = 80 S = 100 S = 150 Deterministic Bounds Objective value 1.95 1.945 1.94 1.935 1.93 20 40 60 80 100 120 140 160 Scenario size Figure 3: Scenario stability test when the rest are below. The main reason is that we consider overtime processing cost uncertainty in the two-stage stochastic programming model and demands are assumed to be fixed at their nominal values. Therefore, uncertainty only affects the objective function of the two-stage stochastic programming model and feasible region stays the same. Concretely, if the scenarios have overtime processing cost higher than the nominal values, then the two-stage stochastic programming model will have higher objective value than the deterministic model. Reversely, if the scenarios have lower overtime processing cost than the nominal values, then the two-stage stochastic programming model will have lower objective value than the deterministic model. From the 20 random samples in Figure 3, roughly half of the objective values of the two-stage stochastic programming model are higher than the objective value of the deterministic model. The comparison between the deterministic model, two-stage stochastic programming model and hybrid model is shown in Figure 4. The solid line represents the objective value of the deterministic model. When the budget ΓsD = 0, the objective value of the hybrid 24 ฀฀฀฀฀฀฀฀฀ model equals to the objective value of the two-stage stochastic model since no demand uncertainty is allowed. Concretely, the top left point represents the objective value of the two-stage stochastic programming model. However, when we increases the budget, the objective value of the hybrid model decreases dramatically as shown in Figure 4. The key take-away are: First, the objective value of the deterministic model can be either higher or lower than the objective value of the two-stage stochastic programming model depending on the overtime production cost scenarios. Second, as budget increases, the objective value of the hybrid model will be significantly lower than the objective value of the deterministic model which concludes that considering uncertainties are really important. We first perform the sensitivity analysis on parameters related to budget and quantity. The effect of budget uncertainty is studied by varying ΓsD for uncertain demands. Let’s define ρ as the maximum variability of the uncertain demands. Higher ρ results in larger deviation. Note that budget can take any values between 0 to |N |∗|T | = 18. If the budget happens to be integer, then it is the maximum amount of parameters that can deviate from their mean values. ΓsD = 0 indicates that there is no protection against uncertainty and ΓsD = 18 provides fully protection at expense of getting conservative solutions. For a particular trajectory in Figure 4, we fix ρ and vary ΓsD . As budget increases, we allow more variability in the uncertain demand and hence optimal profit becomes lower. One special case is ΓsD = 0, four different ρ values provide the same optimal solution since no deviation is allowed and the level of variability no longer matters. If we fix budget and only focus on the ρ, it is obviously that smaller ρ results in higher optimal objective value. That is, the worst-case value for small ρ is actually better than the one for large ρ since the maximum deviation is smaller. The trajectories in Figure 4 appear piecewise linear because only integer budgets are calculated. Details about relationships between probability of constraint violation and budget are included in the Figure 5. In Figure 5a, we show the probability of violation ǫ with respect to ΓsD . As budget increases, there are more protections toward constraint and the 25 ฀฀฀฀฀฀฀฀฀ 10 20 5 = 5% = 10% = 20% = 30% deterministic Objective value 15 10 5 0 -5 0 2 4 6 8 10 12 14 16 18 Budget Figure 4: Optimal objective value with respect to different ΓsD and ρ probability of violation decreases. When budget reaches the maximum value |N | ∗ |T |, the constraint is completely protected and the probability of violation reduces to 0. Reversely, we can measure the minimum ΓsD with respect to ǫ by taking the inverse function. It is shown in Figure 5b, If the maximum violation probability is set to 0.17, then the minimum budget in order to maintain that violation probability is around 5. In the Figure 6, we investigate how capacity of regular time production batch size (q) and ratio of regular vs overtime production (α) affect the optimal objective value. Clearly, bigger q value means more regular time production resource. In addition, overtime production quantity is bounded by the regular time production batch size and hence bigger q indicates more overtime production resource as well. It is shown in Figure 6 that if extra money can be invested, then production batch size q is much more promising than overtime production α. When we fix q and only focus on the α, bigger α provides more overtime production resource and hence results in higher profit. When we fix α and vary q, it does not only provide more overtime production resource but also regular time 26 ฀฀฀฀฀฀฀฀฀ 0.6 11 0.5 10 Minimum budget Probability of violation 9 0.4 0.3 8 7 0.2 6 0.1 5 0 4 0 2 4 6 8 10 12 14 16 18 0 0.02 0.04 Budget 0.06 0.08 0.1 0.12 0.14 0.16 0.18 0.2 Maximum violation (a) ǫ as a function of ΓsD (b) ΓsD as a function of ǫ Figure 5: Relationship between ΓsD and ǫ production resource. It is why there is big improvement in the optimal profit when we vary q. 1.97 10 6 q = 85% q = 90% q = 95% 1.96 objective value 1.95 1.94 1.93 1.92 1.91 1.9 0.18 0.2 0.22 0.24 0.26 0.28 0.3 Maximum overtime production ratio Figure 6: Optimal objective value with respect to different q and α Next, we conduct the sensitivity analyses on parameters that associate with time re27 ฀฀฀฀฀฀฀฀฀ Table 5: Optimal profit and time utilization (Ut ) for different overall time capacities (capt ) Factor (capt ) 0.7 0.8 0.9 1 Optimal profit 1906807 1931663 1941988 1941988 U1 0.93 0.814 0.724 0.651 U2 1 0.944 0.839 0.755 U3 0.951 0.814 0.724 0.651 U4 0.93 0.814 0.74 0.666 U5 1 1 0.999 0.9 U6 1 1 0.999 0.9 Ut : Utilization of overall time capacity in period t, t = 1, · · · , 6 sources such as total time availability capt , processing time pti and setup changeover time sti,j . In Table 5, the sensitivity analysis for overall time capacity has been illustrated on the machine as well as utilization of time resources. Ui represents the time utilization of time period i. As we increases overall time availability, there is more regular time production resources and hence optimal profit increases. Note that there is no improvement when capt changes from 0.9 to 1 indicating we simply waste 10% of time resource. In addition, we observe that, as we reduce capt , U5 and U6 first approach 1 followed by U2 , U3 , then U1 and U4 . It provides insights on how to invest time resources. Clearly, time period 5 and 6 have the highest priority since their utilizations approach to 1 first. Then time period 2, 3 followed by 1 and 4 should be invested if there is enough budget. Processing time pti and setup changeover time sti,j also play important roles in the decision making process. Different pti and sti,j are tested in Figure 7. Since overall time availability is fixed, both pti and sti,j can affect regular time production plan. From constraint 13, we can see that increasing pti or sti,j reduces regular time production capacity. That is, when we increase factor, optimal profit decreases. Note that there is a point where optimal profit becomes negative meaning current production resources are not able to make profit due to large penalty cost. In addition, we can see that optimal profit is more sensitive to the processing time than setup cost. The main reason is that we assume setup takes place between products from different families. Intuitively, when we produce only one type of product, there is no setup changeover cost but processing cost always exists as long as there is production activity. 28 ฀฀฀฀฀฀฀฀฀ 10 Optimal objective values 2 6 0 -2 -4 -6 -8 1 1.1 1.2 1.3 1.4 1.5 1.6 1.7 Processing time factors 10 6 Optimal objective values 2 1 0 -1 -2 -3 0 5 10 15 Setup time factors Figure 7: Optimal profit for different processing times (pti ) and setup times (sti,j ) Finally, we study the sensitivity for parameters that associate with budget such as regular time processing cost, selling price, penalty and surplus cost. Selling price is originally 3 times more expensive than the regular production cost. That is, when factor is set to 5, pri becomes higher than the profit and optimal profit becomes negative. Smaller pri provides larger revenue and thus profits. Conversely, if we decrease ci , profit becomes smaller. When factor is set to 0.3 meaning ci < pri , optimal profit becomes negative. Since uncertain demand can vary in a predetermined interval where the nominal value and maximum deviation are specified. From Table 6, we can see that production decisions becomes riskier since extra costs for unmet demand and excess amount of production increase as penalty and surplus costs increase. It should be noted that the matrix in Table 6 is symmetric with respect to the main diagonal since uncertain demand vary in a symmetric interval and hence the penalty and surplus costs have identical influence on the optimal profit. 29 ฀฀฀฀฀฀฀฀฀ 10 14 6 Selling price Regular processing cost 12 Objective value 10 8 6 4 2 0 -2 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5 Budget Figure 8: Optimal profit for different regular processing costs (pri ) and selling prices (ci ) Table 6: Optimal objective value as a function of penalty (pen) and surplus (sur) costs Penalty cost (pen) 4 0.5 1 1.5 2 0.5 1496793 1354316 1283077 1240334 Surplus cost (sur) 1 1.5 1354316 1283077 1051597 877072 877072 606401 760721 418141 2 1240344 760721 418141 161205 Conclusion In this paper, we study a multi-product, multi-period and capacitated lot-sizing and scheduling problem under demand and overtime production cost uncertainties. We adopt the HSRO approach of Keyvanshokooh et al. to maximize overall profit under demand and overtime processing cost uncertainties [25]. We assume that overtime processing cost is predictable and hence historical data can be used to generate scenarios. Inside each overtime processing cost, polyhedral sets are introduced for uncertain demands due to 30 ฀฀฀฀฀฀฀฀฀ their unpredictable characteristic. In the case study, braking production from an automotive company is used to illustrate and validate the proposed model and solution. Scenario stability tests have been conducted to identify proper scenario size. As scenario size increases, the objective value becomes more stable since scenarios have better representation of original distribution. When scenario sample size is 150, the relative difference between the maximum and minimum objective value is less than 0.3%. Then we carry out sensitivity analyses for parameters like budget ΓsD , constraint violation probability ǫ and time availability on the machine capt etc. Budget provides the level of robustness at the expense of profit. Higher ΓsD results in better protection against uncertainty, but the corresponding profit becomes lower. Constraint violation probability ǫ can be setup by decision makers as the maximum violation probability of constraints or calculated for given ΓsD . Concretely, ǫ can be written as a strictly decreasing function with respect to ΓsD . Time availability and utilization are also conducted to provide valuable insights. It can be shown that there are at least 10% waste of time resources. In addition, time periods 5 and 6 are bottleneck that can be improved upon if there is an increase in the production load. Impacts of other parameters like setup time, processing cost and selling price have also been studied in the paper. In conclusion, we introduce a framework to investigate multiple source of uncertainties of varying characteristics in the scope of lot-sizing and scheduling production. Future research can be explored in the following directions. First, we assume that a particular setup can be taken place at most once in each time period while in reality, multiple setup times maybe allowed. It increases the complexity of a problem exponentially and hence heuristics can be considered. Second, a lot-sizing and scheduling problem is currently at least as difficult as solving one traveling salesman problem in each time period. Some valid inequalities can be applied for computational performance comparisons when number of products increases. Last, demand is assumed to be independent from previous observations and this may not be true in other context where there are strong seasonalities and 31 ฀฀฀฀฀฀฀฀฀ trends. 32 ฀฀฀฀฀฀฀฀฀ References [1] Alem, D., Curcio, E., Amorim, P., and Almada-Lobo, B. (2018). A computational study of the general lot-sizing and scheduling model under demand uncertainty via robust and stochastic approaches. Computers & Operations Research, 90:125–141. [2] Alfieri, A., Tolio, T., and Urgo, M. (2012). A two-stage stochastic programming project scheduling approach to production planning. The International Journal of Advanced Manufacturing Technology, 62(1-4):279–290. [3] Ben-Tal, A. and Nemirovski, A. (1998). Robust convex optimization. Mathematics of Operations Research, 23(4):769–805. [4] Ben-Tal, A. and Nemirovski, A. (1999). Robust solutions of uncertain linear programs. Operations Research Letters, 25(1):1–13. [5] Ben-Tal, A. and Nemirovski, A. (2000). Robust solutions of linear programming problems contaminated with uncertain data. Mathematical Programming, 88(3):411– 424. [6] Bertsimas, D., Pachamanova, D., and Sim, M. (2004). Robust linear optimization under general norms. Operations Research Letters, 32(6):510–516. [7] Bertsimas, D. and Sim, M. (2003). Robust discrete optimization and network flows. Mathematical Programming, 98(1):49–71. [8] Bertsimas, D. and Sim, M. (2004). The price of robustness. Operations Research, 52(1):35–53. [9] Bertsimas, D. and Thiele, A. (2006). Robust and data-driven optimization: Modern decision-making under uncertainty. INFORMS Tutorials in Operations Research: Models, Methods, and Applications for Innovative Decision Making, 3. 33 ฀฀฀฀฀฀฀฀฀ [10] Biller, S., Chan, L. M. A., Simchi-Levi, D., and Swann, J. (2005). Dynamic pricing and the direct-to-customer model in the automotive industry. Electronic Commerce Research, 5(2):309–334. [11] Birge, Z. and Louveaux, S. (1997). Principles on stochastic programming. [12] Chaharsooghi, S. K., Honarvar, M., Modarres, M., and Kamalabadi, I. N. (2011). Developing a two stage stochastic programming model of the price and lead-time decision problem in the multi-class make-to-order firm. Computers & Industrial Engineering, 61(4):1086–1097. [13] Curcio, E., Amorim, P., Zhang, Q., and Almada-Lobo, B. (2018). Adaptation and approximate strategies for solving the lot-sizing and scheduling problem under multistage demand uncertainty. International Journal of Production Economics, 202:81–96. [14] Davis, T. (1993). Effective supply chain management. Sloan Management Review, 34:35–35. [15] Diaz, J. E., Handl, J., and Xu, D.-L. (2017). Evolutionary robust optimization in production planning–interactions between number of objectives, sample size and choice of robustness measure. Computers & Operations Research, 79:266–278. [16] do Rego, J. R. and de Mesquita, M. A. (2015). Demand forecasting and inventory control: A simulation study on automotive spare parts. International Journal of Production Economics, 161:1–16. [17] Fattahi, M., Mahootchi, M., Moattar Husseini, S., Keyvanshokooh, E., and Alborzi, F. (2015). Investigating replenishment policies for centralised and decentralised supply chains using stochastic programming approach. International Journal of Production Research, 53(1):41–69. 34 ฀฀฀฀฀฀฀฀฀ [18] Govindan, K., Fattahi, M., and Keyvanshokooh, E. (2017). Supply chain network design under uncertainty: A comprehensive review and future research directions. European Journal of Operational Research, 263(1):108–141. [19] Heitsch, H. and Römisch, W. (2003). Scenario reduction algorithms in stochastic programming. Computational Optimization and Applications, 24(2-3):187–206. [20] Hu, Z., Balaskovits, A., and Hu, G. (2017). Multi-stage lot-sizing and production scheduling under demand uncertainty. In IIE annual conference. Proceedings, pages 1054–1059. Institute of Industrial and Systems Engineers (IISE). [21] Hu, Z. and Hu, G. (2016). A two-stage stochastic programming model for lot-sizing and scheduling under uncertainty. International Journal of Production Economics, 180:198–207. [22] Hu, Z. and Hu, G. (2018). A multi-stage stochastic programming for lot-sizing and scheduling under demand uncertainty. Computers & Industrial Engineering, 119:157– 166. [23] James, R. J. and Almada-Lobo, B. (2011). Single and parallel machine capacitated lotsizing and scheduling: New iterative mip-based neighborhood search heuristics. Computers & Operations Research, 38(12):1816–1825. [24] Kaut, M., Vladimirou, H., Wallace, S. W., and Zenios, S. A. (2007). Stability analysis of portfolio management with conditional value-at-risk. Quantitative Finance, 7(4):397– 409. [25] Keyvanshokooh, E., Ryan, S. M., and Kabir, E. (2016). Hybrid robust and stochastic optimization for closed-loop supply chain network design using accelerated benders decomposition. European Journal of Operational Research, 249(1):76–92. 35 ฀฀฀฀฀฀฀฀฀ [26] Kleywegt, A. J., Shapiro, A., and Homem-de Mello, T. (2002). The sample average approximation method for stochastic discrete optimization. SIAM Journal on Optimization, 12(2):479–502. [27] Lei, X. and MacKenzie, C. A. (2019a). Assessing risk in different types of supply chains with a dynamic fault tree. Computers & Industrial Engineering, 137:106061. [28] Lei, X. and MacKenzie, C. A. (2019b). Distinguishing between common cause variation and special cause variation in a manufacturing system: A simulation of decision making for different types of variation. International Journal of Production Economics. [29] Li, Y. and Hu, G. (2017). Shop floor lot-sizing and scheduling with a two-stage stochastic programming model considering uncertain demand and workforce efficiency. Computers & Industrial Engineering, 111:263–271. [30] Lieberman, M. B. and Demeester, L. (1999). Inventory reduction and productivity growth: linkages in the japanese automotive industry. Management Science, 45(4):466– 485. [31] Liu, G., Xu, Y., and Tomsovic, K. (2016). Bidding strategy for microgrid in dayahead market based on hybrid stochastic/robust optimization. IEEE Transactions on Smart Grid, 7(1):227–237. [32] Menezes, A. A., Clark, A., and Almada-Lobo, B. (2011). Capacitated lot-sizing and scheduling with sequence-dependent, period-overlapping and non-triangular setups. Journal of Scheduling, 14(2):209–219. [33] Minoux, M. (2018). Robust and stochastic multistage optimisation under markovian uncertainty with applications to production/inventory problems. International Journal of Production Research, 56(1-2):565–583. 36 ฀฀฀฀฀฀฀฀฀ [34] Pishvaee, M. S., Rabbani, M., and Torabi, S. A. (2011). A robust optimization approach to closed-loop supply chain network design under uncertainty. Applied Mathematical Modeling, 35(2):637–649. [35] Ramaraj, G., Hu, Z., and Hu, G. (2017). A two-stage stochastic programming model for production planning in a kitting facility under demand and yield uncertainties. In IIE Annual Conference. Proceedings, pages 1520–1525. Institute of Industrial and Systems Engineers (IISE). [36] Ramaraj, G., Hu, Z., and Hu, G. (2019). A two-stage stochastic programming model for production lot-sizing and scheduling under demand and raw material quality uncertainties. International Journal of Planning and Scheduling, 3(1):1–27. [37] Scarf, H. E. (1957). A min-max solution of an inventory problem. Rand Corporation. [38] Soyster, A. L. (1973). Convex programming with set-inclusive constraints and applications to inexact linear programming. Operations Research, 21(5):1154–1157. [39] Taş, D., Gendreau, M., Jabali, O., and Jans, R. (2018). A capacitated lot sizing problem with stochastic setup times and overtime. European Journal of Operational Research. [40] Tomotani, J. V. and de Mesquita, M. A. (2018). Lot sizing and scheduling: a survey of practices in brazilian companies. Production Planning & Control, 29(3):236–246. [41] Yang, H., Low, V., Zhang, C., Zheng, L., and Miao, L. (2017). Behaviour perceptionbased disruption models for the parallel machine capacitated lot-sizing and scheduling problem. International Journal of Production Research, 55(11):3058–3072. [42] Zhang, X., Prajapati, M., and Peden, E. (2011). A stochastic production planning model under uncertain seasonal demand and market growth. International Journal of Production Research, 49(7):1957–1975. 37