- Email: [email protected]
A simulation optimization based control policy for failure prone one-machine, two-product manufacturing systems
Computers & Industrial Engineering 46 (2004) 285–292 www.elsevier.com/locate/dsw
A simulation optimization based control policy for failure prone one-machine, two-product manufacturing systems J.P. Kennea, A. Gharbib,* a
Department of Mechanical Engineering, Ecole de Technologie Superieure, University of Quebec, 1100 Notre Dame Ouest, Montreal, Que., Canada H3C 1K3 b Department of Automated Production Enginnering, Production Systems Design and Control Laboratory, Ecole de Technologie Superieure, University of Quebec, 1100 Notre Dame Ouest, Montreal, Que., Canada H3C 1K3
Abstract This paper presents the optimal flow control for a one-machine, two-product manufacturing system subject to random failures and repairs. The machine capacity process is assumed to be a finite state Markov chain. The problem is to choose the production rates so as to minimize the expected discounted cost of inventory/backlog over an infinite horizon. We first show that for constant demand rates and exponential failure and repair time distributions of the machine, the hedging point policy is optimal. Next, the hedging point policy is extended to nonexponential failure and repair time distributions models. The structure of the hedging point policy is parameterized by two factors representing the thresholds of involved products. With such a policy, simulation experiments are coupled with experimental design and response surface methodology to estimate the optimal control policy. Our results reveal that the hedging point policy is also applicable to a wide variety of complex problems (i.e. nonexponential failure and repair time distributions) where analytical solutions may not be easily obtained. q 2004 Elsevier Ltd. All rights reserved. Keywords: Stochastic optimal control; Numerical methods; Simulation; Experimental design; Response surface methodology
1. Introduction We are concerned with the problem of controlling the production rates of a stochastic one-machine, two-product manufacturing system in order to meet the demand rates facing the system at a minimum cost. The stochastic nature of the system is due to machines that are subject to breakdowns and repairs. The decision variables are the production rates of products. Our objective is to choose admissible production rates to minimize the inventory/backlog costs over an infinite horizon. Under some restrictive hypotheses (i.e. exponential distribution of breakdowns and repairs inter-arrival), the problem of * Corresponding author. Tel.: þ1-514-396-8969; fax: þ1-514-396-8595. E-mail addresses: [email protected] (A. Gharbi), [email protected] (J.P. Kenne). 0360-8352/$ - see front matter q 2004 Elsevier Ltd. All rights reserved. doi:10.1016/j.cie.2003.12.004
286
J.P. Kenne, A. Gharbi / Computers & Industrial Engineering 46 (2004) 285–292
controlling a flexible manufacturing system or FMS with unreliable machines was formulated as a stochastic optimal control problem by Older and Suri (1980). Based on the method presented by Rishel (1975) and Kimemia and Gershwin (1983) showed that the optimal control policy of such a problem has a special structure called the hedging point. Investigation in the same direction gives rise to the analytical solution of the one-machine, one-product manufacturing system (Akella & Kumar, 1986). It remains difficult to obtain the optimal control of a large scale manufacturing system (i.e. involving at least two parts). The optimality conditions of the above problems, based on stochastic models, are described by generalized Hamilton Jacobi Bellman (HJB) equations for optimal production policies. However, it is by now well-known that finding the analytical solutions of HJB equations is almost impossible except in a few very simple cases, such as in Akella and Kumar (1986). The crux of the hedging point policy is to maintain surplus levels or thresholds, corresponding to the optimal inventory levels, in order to compensate possible backlogs caused by random failures. We observe that there is still no satisfying method for manufacturing systems modelled by nonMarkov processes (i.e. non-exponential machine up and down times). We will first prove that, with exponential machine up and down times and constant demand rates, the combination of discrete event simulation modelling, experimental design and analytical control theory can provide a control policy that can be compared to the one given by the numerical approach (Kusner & Dupuis, 1992). Following this validation, we will extend our proposal to non-exponential machine up and down time situations. The proposed approach consists of estimating the relationship between the incurred cost and the threshold stock levels considered herein as control factors. The hedging point policy, parameterized by these factors, is then used to conduct simulation experiments. For each configuration of input factors values, the simulation model is used to determine the related output or cost incurred. An input – output data is then generated through the simulation model. The experimental design is used to determine significant factors and/or their interactions and the response surface methodology is applied to the obtained input – output data vectors to estimate the cost function and the related optimum. Details on the combination of analytical approaches and simulation based statistical methods can be found in Abdulnour, Duddek, & Smith (1995), Gharbi and Kenne (2000) and Kenne and Gharbi (1999). The remainder of the paper is organized as follows. Section 2 describes the statement of the optimal control problem. The description of the proposed approach and the simulation model are presented in section 3. Section 4 outlines the experimental design approach and response surface methodology. Section 5 is devoted to extensions of the proposal to non-markovian processes. Section 6 concludes the paper. 2. Problem statement The system under study consists of one machine producing two part types. The operational mode of the machine can be described by a stochastic process zðtÞ: Such a machine is available when it is operational ðzðtÞ ¼ 1Þ and unavailable when it is under repair ðzðtÞ ¼ 2Þ: One can describe zðtÞ statistically with the following state probabilities: ( P½zðt þ dtÞ ¼ blzðtÞ ¼ a ¼
lab dt þ oðdtÞ
if a – b
1 þ lab dt þ oðdtÞ if a ¼ b
ð1Þ
J.P. Kenne, A. Gharbi / Computers & Industrial Engineering 46 (2004) 285–292
287
where lab $ 0 laa ¼ 2 Sb–a lba ða; b [ {1; 2}Þ; limdt!0 oðatÞ=dt ¼ 0 and the transition rate from state a to state b is represented by lab : We assume that lab are given constant transition rates. This process is modeled by a continuous time Markov chain defined on {V; F; P} with a matrix of transition rates Q ¼ {lab }; which is a 2 £ 2 irreducible stochastic matrix. Let uðtÞ be the vector of production rates. The set of feasible control policies depends on the process zðtÞ; such that uðtÞ [ KðaÞ with a ¼ zðtÞ: This set can be given by j ðaÞ; j ¼ 1; 2} KðaÞ ¼ {uðtÞ [ R2 ; 0 # uj ðtÞ # umax
ð2Þ
where ujmax ðaÞ is the maximal production rate of part type Pj at mode a: The system’s behavior is described by a hybrid state comprising both a discrete and a continuous component. The discrete component consists of the above discrete event stochastic process zðtÞ: The continuous component is xðtÞ [ R2 and corresponds to the inventory/backlog vector of products. This state variable is described by the following differential equation: x_ ðtÞ ¼ uðtÞ 2 d
xð0Þ ¼ x
ð3Þ
where x [ R2 and d [ R2 are given initial surplus value and demand rate vectors, respectively. Let ga ðx; uÞ be the cost rate defined as follows: ga ðxðtÞ; uðtÞÞ ¼
2 X
þ 2 2 ðcþ j xj ðtÞ þ cj xj ðtÞÞ
;a [ {1; 2}
ð4Þ
j¼1 2 where cþ j and cj are the incurred cost per unit produced part type Pj for positive inventory and backlog, þ respectively, xj ðtÞ ¼ maxð0; xj ðtÞÞ; x2 j ðtÞ ¼ maxð2xj ðtÞ; 0Þ: Our objective is to control the production rate vector uðtÞ so as to minimize the expected discounted cost given by: ð1 2rt a e g ðx; uÞdtlxðtÞ ¼ x; aðtÞ ¼ a Jða; x; uÞ ¼ E ð5Þ 0
subject to Eqs. (1– 4) for a given discount rate r: The value function of such a problem is: vða; xÞ ¼ min Jða; x; uÞ u[KðaÞ
;a [ {1; 2}
ð6Þ
The continuous dynamic programming set of equations associated with this problem is given by: 8 9 2 < = X › rvðx; aÞ ¼ min ðu 2 dÞ vðx; aÞ þ ga ðx; uÞ þ lab vðx; aÞ ;a; b [ {1; 2} ð7Þ ; u[KðaÞ : ›x b¼1 For the interpretation of Eq. (7), the reader is referred to Sethi and Zhang (1994) where the value function vð·Þ is presented as a viscosity solution of Eq. (7). Therefore, the value function is convex. The optimality conditions of the control problem stated herein are described by Eq. (7), commonly named HJB equations. Due to the complexity of HJB equations, the objective of this paper is not to solve them analytically but to experimentally determine the parameters of a modified hedging point policy, which gives the best approximation of the value vða; xÞ:
288
J.P. Kenne, A. Gharbi / Computers & Industrial Engineering 46 (2004) 285–292
Defining z1 and z2 as the threshold of products of P1 and P2 ; respectively, it is simple to prove that the obtained numerical control policy structure divides the plan ðx1 ; x2 Þ in the following three regions: † No production is required for production of P1 and P2 when x1 $ z1 and x2 $ z2 ; † The machine is allocated to product Pk if xk , zk and xl $ zl ; k – l; † The machine is allocated alternatively to P1 and P2 if x1 , z1 and x2 , z2 : This gives rise to the following control policy: 1. If xk , zk and xl $ zl ; k – l uk ðzÞ ¼ ukmax 2. If xk , zk and xl , zl ; k – l ( j umax if xj , zj and xj , xk uj ðzÞ ¼ 0 otherwise
ð8Þ
ð9Þ
for k ¼ 1; 2 and l ¼ 1; 2: The modified hedging point policy defined by Eqs. (8) and (9) is completely defined for given values of zj ; herein called design factors ðj ¼ 1; 2Þ: 3. Proposed approach and simulation model The proposed approach is summarized by the following steps: † Step 1: Problem statement and analytical approach. The optimal flow control problem for the considered manufacturing system is described and the objective of the study is specified in this step. The objective here is to find the optimal production control policy that minimizes the average discounted surplus cost related to an infinite horizon. The aim of the control problem statement is to develop a mathematical representation of the system based on some simplified hypothesis. By applying an analytical approach, the structure of the production control policy is obtained. Such a policy is parameterized by the factors z1 and z2 and is considered as input of the simulation model. The optimal flow control for the considered manufacturing system is formulated, as in Section 2 with optimality conditions given by HJB Eq. (7). † Step 2: Simulation modeling. For the simulation model, the cost related to each input, given by the values of factors zj ; j ¼ 1; 2; is defined as its output. Recall that the objective of the proposed heuristic approach is to find the best parameters of the control policy uðz1 ; z2 Þ in order to improve the related output (i.e. incurred cost). The proposed simulation model consists of the data input module, the simulation module and the output analysis module. The simulation model describes the dynamic of the system (Pritsker, O’Reilly, & LaVal, 1997) using the control policy parameterized by the factors zj ; j ¼ 1; 2: We refer the reader to Gharbi and Kenne (2000) and Kenne and Gharbi (1999) for details in the application of simulation in the control of FMS. † Step 3: Experimental design and response surface methodology. The experimental design determines, from values of the input factors and related cost values, the input factors and/or their interactions,
J.P. Kenne, A. Gharbi / Computers & Industrial Engineering 46 (2004) 285–292
289
which have significant effects on the output. These factors or interactions are then considered as input of response surface methodology in order to fit the relationship between the cost and input factors. From this estimated relation, the optimal values of the estimated input factors are determined. The related hedging point policy uðzp1 ; zp2 Þ is then an improved hedging point policy to be applied to the manufacturing system. The application of the proposed simulation-based approach to the production rate control problem for the system under study is presented in Section 4, by varying the input factors zj ; j ¼ 1; 2: 4. Experimental design and response surface methodology In this study, we collect and analyze data for a steady state cost which approximates the one defined by the value function, given by Eq. (6) as close as possible. We have considered two independent variables z1 and z2 and one dependent variable (the cost). The levels of independent variables or factors should be chosen carefully so that they represent the domain of interest. The levels of the independent variables were chosen through additional simulation experiments. Due to the convexity of the value function (Section 2), we selected a 32 response surface design since we have two independent variables at three levels each. Such a design leads to the completion of nine experimental trials. Three replications were conducted for each combination; therefore, 27 (i.e., 9 £ 3) simulation runs were made. All possible combinations of different factor levels are provided by the response surface design considered herein. The parameters of the numerical example are presented in Table 1. The experimental design is used to study and understand the effects that some parameters, namely z1 and z2 ; have on the performance measure (i.e. the cost). The resulting standardized Pareto chart for the incurred cost is presented in Fig. 1. From Fig. 1, we see that the main factors z1 ; z2 and their quadratic effects are significant at 0.05 level; however, their interaction is non-significant. The obtained result is related to an R-squared value of 0.95. Such a value states that 95% of the total variability is explained by the model. The second order model is then given by: Cost ¼ b0 þ b1 z1 þ b2 z2 þ b11 z2 þ b22 z22 þ 1:
ð10Þ
where 1 is the experimental error in Cost. The near-optimal control policy to be applied to the considered manufacturing system is defined by zp1 ¼ 10 and zp2 ¼ 10 for b0 ¼ 127:176; b1 ¼ 23:795; b2 ¼ 23:891 b11 ¼ 0:185; b22 ¼ 0:188: A cost value of 88.44 is incurred with such a control policy, as shown in Fig. 2. For the purpose of validation, performances of the obtained control policy have been compared to those given by the numerical approach (zp1 ¼ 9; zp2 ¼ 10 and Cost ¼ 89.76). Table 1 Parameters of the numerical example cþ 1
cþ 2
c2 1
c2 2
d1
d2
u1max
u2max
r
l12
l21
1
1
4
4
0.25
0.25
0.57
0.57
0.09
0.0105
0.4
290
J.P. Kenne, A. Gharbi / Computers & Industrial Engineering 46 (2004) 285–292
Fig. 1. Standardized Pareto chart of the cost.
Through the Student’s t-test at a 95% confidence level using 50 simulation runs, we realized that the incurred cost for the estimated control policy does not differ from the one given by the numerical control policy. Details on the Student test and its application on simulation runs can be found in Kenne and Gharbi (2000). Motivated by our results, we moved on to additional aspects, such as non-exponential machine failure and repair times, in order to control more complex manufacturing systems. 5. Extension to non-exponential failure and repair times distributions In the proposed approach, the simulation model takes into account the nature of the machines’ failures and repair time distributions and gives the corresponding output. These aspects will affect the response surface model. In this section, we explore three additional cases related to commonly used failure and repair time distributions such as gamma, lognormal and Weibull. We refer the reader to Law and Kelton (2000) for details on commonly used failure and repair probability time distributions. The mean
Fig. 2. Estimated response surface.
J.P. Kenne, A. Gharbi / Computers & Industrial Engineering 46 (2004) 285–292
291
Fig. 3. Sensitivity analysis on estimated threshold values.
value and the standard deviation of considered distributions are chosen equal to the basic exponential case (i.e. m ¼ 95 and s ¼ 9:75 for failure times and m ¼ 5 and s ¼ 2:25 for repair time distributions). 2 For product 1, keeping cþ 1 ¼ 1 and varying c1 from 5 to 40 (i.e. varying the instantaneous cost ratio 2 þ c1 =c1 from 5 to 40) the obtained threshold values are depicted in Fig. 3. The same behaviour is obtained for product 2, given that considered products are identical (Table 1). Table 2 illustrates the obtained R-squared values for exponential and non-exponential failure and repair time distributions. For illustrative purpose, we present herein R-squared values of the gamma distribution to represent the class of non-exponential distributions for which obtained R-squared values are greater than 0.99. The relatively small values of R-squared values given by exponential distributions are a result of their high variability. The trend of the curves shown in Fig. 3 confirms the robustness of the proposed approach through sensitivity analysis. This is performed by threshold values analysis versus instantaneous cost. It is wellillustrated in Fig. 3 that large threshold values are obtained for exponential distributions of machine failure and repair time distributions. Such an observation is useful due to the large variability of exponential distributions. One also realizes that the threshold remains constant for large values of c2 1: ; the threshold is also large. With such values of This is due to the fact that, with large values of c2 1 threshold, one cannot reach a backlog situation. Hence, there is no need to increase the threshold to avoid backlog cost. The differences between the obtained four curves (Fig. 3) are mainly based on the variability Table 2 Values of obtained R-squared (%) for different cost ratios Cost ratio
5
10
15
20
25
30
35
40
R2 (exponential) R2 (non-exponential)
91.8 99.8
85.8 99.5
89.3 99.2
90.5 99.4
92.2 99.6
92.7 99.7
93.1 99.7
93.5 99.6
292
J.P. Kenne, A. Gharbi / Computers & Industrial Engineering 46 (2004) 285–292
of involved distributions. Fig. 3 also shows that the proposed approach is not restricted to exponential distribution, as in a few studies based on the control of manufacturing systems. Through the above analysis, it clearly appears that the obtained results make sense and that the proposed approach is robust. 6. Conclusion We have extended the concept of hedging point policy to the control of two-product manufacturing systems with no restriction to exponential machine failure and repair rates. The proposed approach was based on the combination of the analytical model, simulation experiments, experimental design and response surface methodology. A proposed approach was validated by comparing the obtained results with the numerical results through a t-Student test. Motivated by such a validation, we explored the control of manufacturing systems with probability distributions generally used to represent machines failure and repair times, such as gamma, lognormal and Weibull. The Findings of this paper should be very useful for the control of any kind of non-markovian process for which the stochastic optimal control problem remains unsolved. References Abdulnour, G., Duddek, R. A., & Smith, M. L. (1995). Effect of maintenance policies on the just-in-time production system. International Journal of Production Research, 33, 565–583. Akella, R., & Kumar, P. R. (1986). Optimal control of production rate in a failure prone manufacturing system. IEEE Transactions on Automatic Control, 31(2), 116– 126. Gharbi, A., & Kenne, J. P. (2000). Production and preventive maintenance rates control for a manufacturing system: an experimental design approach. International Journal of Production Economics, 65(3), 275– 287. Kenne, J. P., & Gharbi, A. (2000). Production planning problem in manufacturing systems with general failure and repair time distributions. Production Planning and Control, 11, 581– 588. Kenne, J. P., & Gharbi, A. (1999). Experimental design in production and maintenance control of a single machine, single product manufacturing system. International Journal of Production Research, 37(3), 621– 637. Kusner, H. J., & Dupuis, P. G. (1992). Numerical methods for stochastic control problem in continuous time. New York: Springer. Kimemia, J. G., & Gershwin, S. B. (1983). An algorithm for computer control of production in flexible manufacturing systems. IIE Transactions, 15, 353– 362. Law, A. M., & Kelton, W. D. (2000). Simulation modeling and analysis (3rd ed). New York: McGraw-Hill. Older, G. J., & Suri, R. (1980). Time optimal part-routing in a manufacturing system with failure prone machines. Proceedings of 19th Conference on Decision Control, Albuquerque, NM, pp. 722–727. Pritsker, A. A. B., O’Reilly, J. J., & LaVal, D. K. (1997). Simulation with visual SLAM and AweSim. New York: Wiley. Rishel, R. (1975). Dynamic programming and minimal principles for systems with jump Markov disturbances, SIAM. Journal of Control, 13, 338 –371. Sethi, S. P., & Zhang, Q. (1994). Hierarchical control decision making in stochastic manufacturing systems. Boston: Birkhauser.
A simulation optimization based control policy for failure prone one-machine, two-product manufacturing systems J.P. Kennea, A. Gharbib,* a
Department of Mechanical Engineering, Ecole de Technologie Superieure, University of Quebec, 1100 Notre Dame Ouest, Montreal, Que., Canada H3C 1K3 b Department of Automated Production Enginnering, Production Systems Design and Control Laboratory, Ecole de Technologie Superieure, University of Quebec, 1100 Notre Dame Ouest, Montreal, Que., Canada H3C 1K3
Abstract This paper presents the optimal flow control for a one-machine, two-product manufacturing system subject to random failures and repairs. The machine capacity process is assumed to be a finite state Markov chain. The problem is to choose the production rates so as to minimize the expected discounted cost of inventory/backlog over an infinite horizon. We first show that for constant demand rates and exponential failure and repair time distributions of the machine, the hedging point policy is optimal. Next, the hedging point policy is extended to nonexponential failure and repair time distributions models. The structure of the hedging point policy is parameterized by two factors representing the thresholds of involved products. With such a policy, simulation experiments are coupled with experimental design and response surface methodology to estimate the optimal control policy. Our results reveal that the hedging point policy is also applicable to a wide variety of complex problems (i.e. nonexponential failure and repair time distributions) where analytical solutions may not be easily obtained. q 2004 Elsevier Ltd. All rights reserved. Keywords: Stochastic optimal control; Numerical methods; Simulation; Experimental design; Response surface methodology
1. Introduction We are concerned with the problem of controlling the production rates of a stochastic one-machine, two-product manufacturing system in order to meet the demand rates facing the system at a minimum cost. The stochastic nature of the system is due to machines that are subject to breakdowns and repairs. The decision variables are the production rates of products. Our objective is to choose admissible production rates to minimize the inventory/backlog costs over an infinite horizon. Under some restrictive hypotheses (i.e. exponential distribution of breakdowns and repairs inter-arrival), the problem of * Corresponding author. Tel.: þ1-514-396-8969; fax: þ1-514-396-8595. E-mail addresses: [email protected] (A. Gharbi), [email protected] (J.P. Kenne). 0360-8352/$ - see front matter q 2004 Elsevier Ltd. All rights reserved. doi:10.1016/j.cie.2003.12.004
286
J.P. Kenne, A. Gharbi / Computers & Industrial Engineering 46 (2004) 285–292
controlling a flexible manufacturing system or FMS with unreliable machines was formulated as a stochastic optimal control problem by Older and Suri (1980). Based on the method presented by Rishel (1975) and Kimemia and Gershwin (1983) showed that the optimal control policy of such a problem has a special structure called the hedging point. Investigation in the same direction gives rise to the analytical solution of the one-machine, one-product manufacturing system (Akella & Kumar, 1986). It remains difficult to obtain the optimal control of a large scale manufacturing system (i.e. involving at least two parts). The optimality conditions of the above problems, based on stochastic models, are described by generalized Hamilton Jacobi Bellman (HJB) equations for optimal production policies. However, it is by now well-known that finding the analytical solutions of HJB equations is almost impossible except in a few very simple cases, such as in Akella and Kumar (1986). The crux of the hedging point policy is to maintain surplus levels or thresholds, corresponding to the optimal inventory levels, in order to compensate possible backlogs caused by random failures. We observe that there is still no satisfying method for manufacturing systems modelled by nonMarkov processes (i.e. non-exponential machine up and down times). We will first prove that, with exponential machine up and down times and constant demand rates, the combination of discrete event simulation modelling, experimental design and analytical control theory can provide a control policy that can be compared to the one given by the numerical approach (Kusner & Dupuis, 1992). Following this validation, we will extend our proposal to non-exponential machine up and down time situations. The proposed approach consists of estimating the relationship between the incurred cost and the threshold stock levels considered herein as control factors. The hedging point policy, parameterized by these factors, is then used to conduct simulation experiments. For each configuration of input factors values, the simulation model is used to determine the related output or cost incurred. An input – output data is then generated through the simulation model. The experimental design is used to determine significant factors and/or their interactions and the response surface methodology is applied to the obtained input – output data vectors to estimate the cost function and the related optimum. Details on the combination of analytical approaches and simulation based statistical methods can be found in Abdulnour, Duddek, & Smith (1995), Gharbi and Kenne (2000) and Kenne and Gharbi (1999). The remainder of the paper is organized as follows. Section 2 describes the statement of the optimal control problem. The description of the proposed approach and the simulation model are presented in section 3. Section 4 outlines the experimental design approach and response surface methodology. Section 5 is devoted to extensions of the proposal to non-markovian processes. Section 6 concludes the paper. 2. Problem statement The system under study consists of one machine producing two part types. The operational mode of the machine can be described by a stochastic process zðtÞ: Such a machine is available when it is operational ðzðtÞ ¼ 1Þ and unavailable when it is under repair ðzðtÞ ¼ 2Þ: One can describe zðtÞ statistically with the following state probabilities: ( P½zðt þ dtÞ ¼ blzðtÞ ¼ a ¼
lab dt þ oðdtÞ
if a – b
1 þ lab dt þ oðdtÞ if a ¼ b
ð1Þ
J.P. Kenne, A. Gharbi / Computers & Industrial Engineering 46 (2004) 285–292
287
where lab $ 0 laa ¼ 2 Sb–a lba ða; b [ {1; 2}Þ; limdt!0 oðatÞ=dt ¼ 0 and the transition rate from state a to state b is represented by lab : We assume that lab are given constant transition rates. This process is modeled by a continuous time Markov chain defined on {V; F; P} with a matrix of transition rates Q ¼ {lab }; which is a 2 £ 2 irreducible stochastic matrix. Let uðtÞ be the vector of production rates. The set of feasible control policies depends on the process zðtÞ; such that uðtÞ [ KðaÞ with a ¼ zðtÞ: This set can be given by j ðaÞ; j ¼ 1; 2} KðaÞ ¼ {uðtÞ [ R2 ; 0 # uj ðtÞ # umax
ð2Þ
where ujmax ðaÞ is the maximal production rate of part type Pj at mode a: The system’s behavior is described by a hybrid state comprising both a discrete and a continuous component. The discrete component consists of the above discrete event stochastic process zðtÞ: The continuous component is xðtÞ [ R2 and corresponds to the inventory/backlog vector of products. This state variable is described by the following differential equation: x_ ðtÞ ¼ uðtÞ 2 d
xð0Þ ¼ x
ð3Þ
where x [ R2 and d [ R2 are given initial surplus value and demand rate vectors, respectively. Let ga ðx; uÞ be the cost rate defined as follows: ga ðxðtÞ; uðtÞÞ ¼
2 X
þ 2 2 ðcþ j xj ðtÞ þ cj xj ðtÞÞ
;a [ {1; 2}
ð4Þ
j¼1 2 where cþ j and cj are the incurred cost per unit produced part type Pj for positive inventory and backlog, þ respectively, xj ðtÞ ¼ maxð0; xj ðtÞÞ; x2 j ðtÞ ¼ maxð2xj ðtÞ; 0Þ: Our objective is to control the production rate vector uðtÞ so as to minimize the expected discounted cost given by: ð1 2rt a e g ðx; uÞdtlxðtÞ ¼ x; aðtÞ ¼ a Jða; x; uÞ ¼ E ð5Þ 0
subject to Eqs. (1– 4) for a given discount rate r: The value function of such a problem is: vða; xÞ ¼ min Jða; x; uÞ u[KðaÞ
;a [ {1; 2}
ð6Þ
The continuous dynamic programming set of equations associated with this problem is given by: 8 9 2 < = X › rvðx; aÞ ¼ min ðu 2 dÞ vðx; aÞ þ ga ðx; uÞ þ lab vðx; aÞ ;a; b [ {1; 2} ð7Þ ; u[KðaÞ : ›x b¼1 For the interpretation of Eq. (7), the reader is referred to Sethi and Zhang (1994) where the value function vð·Þ is presented as a viscosity solution of Eq. (7). Therefore, the value function is convex. The optimality conditions of the control problem stated herein are described by Eq. (7), commonly named HJB equations. Due to the complexity of HJB equations, the objective of this paper is not to solve them analytically but to experimentally determine the parameters of a modified hedging point policy, which gives the best approximation of the value vða; xÞ:
288
J.P. Kenne, A. Gharbi / Computers & Industrial Engineering 46 (2004) 285–292
Defining z1 and z2 as the threshold of products of P1 and P2 ; respectively, it is simple to prove that the obtained numerical control policy structure divides the plan ðx1 ; x2 Þ in the following three regions: † No production is required for production of P1 and P2 when x1 $ z1 and x2 $ z2 ; † The machine is allocated to product Pk if xk , zk and xl $ zl ; k – l; † The machine is allocated alternatively to P1 and P2 if x1 , z1 and x2 , z2 : This gives rise to the following control policy: 1. If xk , zk and xl $ zl ; k – l uk ðzÞ ¼ ukmax 2. If xk , zk and xl , zl ; k – l ( j umax if xj , zj and xj , xk uj ðzÞ ¼ 0 otherwise
ð8Þ
ð9Þ
for k ¼ 1; 2 and l ¼ 1; 2: The modified hedging point policy defined by Eqs. (8) and (9) is completely defined for given values of zj ; herein called design factors ðj ¼ 1; 2Þ: 3. Proposed approach and simulation model The proposed approach is summarized by the following steps: † Step 1: Problem statement and analytical approach. The optimal flow control problem for the considered manufacturing system is described and the objective of the study is specified in this step. The objective here is to find the optimal production control policy that minimizes the average discounted surplus cost related to an infinite horizon. The aim of the control problem statement is to develop a mathematical representation of the system based on some simplified hypothesis. By applying an analytical approach, the structure of the production control policy is obtained. Such a policy is parameterized by the factors z1 and z2 and is considered as input of the simulation model. The optimal flow control for the considered manufacturing system is formulated, as in Section 2 with optimality conditions given by HJB Eq. (7). † Step 2: Simulation modeling. For the simulation model, the cost related to each input, given by the values of factors zj ; j ¼ 1; 2; is defined as its output. Recall that the objective of the proposed heuristic approach is to find the best parameters of the control policy uðz1 ; z2 Þ in order to improve the related output (i.e. incurred cost). The proposed simulation model consists of the data input module, the simulation module and the output analysis module. The simulation model describes the dynamic of the system (Pritsker, O’Reilly, & LaVal, 1997) using the control policy parameterized by the factors zj ; j ¼ 1; 2: We refer the reader to Gharbi and Kenne (2000) and Kenne and Gharbi (1999) for details in the application of simulation in the control of FMS. † Step 3: Experimental design and response surface methodology. The experimental design determines, from values of the input factors and related cost values, the input factors and/or their interactions,
J.P. Kenne, A. Gharbi / Computers & Industrial Engineering 46 (2004) 285–292
289
which have significant effects on the output. These factors or interactions are then considered as input of response surface methodology in order to fit the relationship between the cost and input factors. From this estimated relation, the optimal values of the estimated input factors are determined. The related hedging point policy uðzp1 ; zp2 Þ is then an improved hedging point policy to be applied to the manufacturing system. The application of the proposed simulation-based approach to the production rate control problem for the system under study is presented in Section 4, by varying the input factors zj ; j ¼ 1; 2: 4. Experimental design and response surface methodology In this study, we collect and analyze data for a steady state cost which approximates the one defined by the value function, given by Eq. (6) as close as possible. We have considered two independent variables z1 and z2 and one dependent variable (the cost). The levels of independent variables or factors should be chosen carefully so that they represent the domain of interest. The levels of the independent variables were chosen through additional simulation experiments. Due to the convexity of the value function (Section 2), we selected a 32 response surface design since we have two independent variables at three levels each. Such a design leads to the completion of nine experimental trials. Three replications were conducted for each combination; therefore, 27 (i.e., 9 £ 3) simulation runs were made. All possible combinations of different factor levels are provided by the response surface design considered herein. The parameters of the numerical example are presented in Table 1. The experimental design is used to study and understand the effects that some parameters, namely z1 and z2 ; have on the performance measure (i.e. the cost). The resulting standardized Pareto chart for the incurred cost is presented in Fig. 1. From Fig. 1, we see that the main factors z1 ; z2 and their quadratic effects are significant at 0.05 level; however, their interaction is non-significant. The obtained result is related to an R-squared value of 0.95. Such a value states that 95% of the total variability is explained by the model. The second order model is then given by: Cost ¼ b0 þ b1 z1 þ b2 z2 þ b11 z2 þ b22 z22 þ 1:
ð10Þ
where 1 is the experimental error in Cost. The near-optimal control policy to be applied to the considered manufacturing system is defined by zp1 ¼ 10 and zp2 ¼ 10 for b0 ¼ 127:176; b1 ¼ 23:795; b2 ¼ 23:891 b11 ¼ 0:185; b22 ¼ 0:188: A cost value of 88.44 is incurred with such a control policy, as shown in Fig. 2. For the purpose of validation, performances of the obtained control policy have been compared to those given by the numerical approach (zp1 ¼ 9; zp2 ¼ 10 and Cost ¼ 89.76). Table 1 Parameters of the numerical example cþ 1
cþ 2
c2 1
c2 2
d1
d2
u1max
u2max
r
l12
l21
1
1
4
4
0.25
0.25
0.57
0.57
0.09
0.0105
0.4
290
J.P. Kenne, A. Gharbi / Computers & Industrial Engineering 46 (2004) 285–292
Fig. 1. Standardized Pareto chart of the cost.
Through the Student’s t-test at a 95% confidence level using 50 simulation runs, we realized that the incurred cost for the estimated control policy does not differ from the one given by the numerical control policy. Details on the Student test and its application on simulation runs can be found in Kenne and Gharbi (2000). Motivated by our results, we moved on to additional aspects, such as non-exponential machine failure and repair times, in order to control more complex manufacturing systems. 5. Extension to non-exponential failure and repair times distributions In the proposed approach, the simulation model takes into account the nature of the machines’ failures and repair time distributions and gives the corresponding output. These aspects will affect the response surface model. In this section, we explore three additional cases related to commonly used failure and repair time distributions such as gamma, lognormal and Weibull. We refer the reader to Law and Kelton (2000) for details on commonly used failure and repair probability time distributions. The mean
Fig. 2. Estimated response surface.
J.P. Kenne, A. Gharbi / Computers & Industrial Engineering 46 (2004) 285–292
291
Fig. 3. Sensitivity analysis on estimated threshold values.
value and the standard deviation of considered distributions are chosen equal to the basic exponential case (i.e. m ¼ 95 and s ¼ 9:75 for failure times and m ¼ 5 and s ¼ 2:25 for repair time distributions). 2 For product 1, keeping cþ 1 ¼ 1 and varying c1 from 5 to 40 (i.e. varying the instantaneous cost ratio 2 þ c1 =c1 from 5 to 40) the obtained threshold values are depicted in Fig. 3. The same behaviour is obtained for product 2, given that considered products are identical (Table 1). Table 2 illustrates the obtained R-squared values for exponential and non-exponential failure and repair time distributions. For illustrative purpose, we present herein R-squared values of the gamma distribution to represent the class of non-exponential distributions for which obtained R-squared values are greater than 0.99. The relatively small values of R-squared values given by exponential distributions are a result of their high variability. The trend of the curves shown in Fig. 3 confirms the robustness of the proposed approach through sensitivity analysis. This is performed by threshold values analysis versus instantaneous cost. It is wellillustrated in Fig. 3 that large threshold values are obtained for exponential distributions of machine failure and repair time distributions. Such an observation is useful due to the large variability of exponential distributions. One also realizes that the threshold remains constant for large values of c2 1: ; the threshold is also large. With such values of This is due to the fact that, with large values of c2 1 threshold, one cannot reach a backlog situation. Hence, there is no need to increase the threshold to avoid backlog cost. The differences between the obtained four curves (Fig. 3) are mainly based on the variability Table 2 Values of obtained R-squared (%) for different cost ratios Cost ratio
5
10
15
20
25
30
35
40
R2 (exponential) R2 (non-exponential)
91.8 99.8
85.8 99.5
89.3 99.2
90.5 99.4
92.2 99.6
92.7 99.7
93.1 99.7
93.5 99.6
292
J.P. Kenne, A. Gharbi / Computers & Industrial Engineering 46 (2004) 285–292
of involved distributions. Fig. 3 also shows that the proposed approach is not restricted to exponential distribution, as in a few studies based on the control of manufacturing systems. Through the above analysis, it clearly appears that the obtained results make sense and that the proposed approach is robust. 6. Conclusion We have extended the concept of hedging point policy to the control of two-product manufacturing systems with no restriction to exponential machine failure and repair rates. The proposed approach was based on the combination of the analytical model, simulation experiments, experimental design and response surface methodology. A proposed approach was validated by comparing the obtained results with the numerical results through a t-Student test. Motivated by such a validation, we explored the control of manufacturing systems with probability distributions generally used to represent machines failure and repair times, such as gamma, lognormal and Weibull. The Findings of this paper should be very useful for the control of any kind of non-markovian process for which the stochastic optimal control problem remains unsolved. References Abdulnour, G., Duddek, R. A., & Smith, M. L. (1995). Effect of maintenance policies on the just-in-time production system. International Journal of Production Research, 33, 565–583. Akella, R., & Kumar, P. R. (1986). Optimal control of production rate in a failure prone manufacturing system. IEEE Transactions on Automatic Control, 31(2), 116– 126. Gharbi, A., & Kenne, J. P. (2000). Production and preventive maintenance rates control for a manufacturing system: an experimental design approach. International Journal of Production Economics, 65(3), 275– 287. Kenne, J. P., & Gharbi, A. (2000). Production planning problem in manufacturing systems with general failure and repair time distributions. Production Planning and Control, 11, 581– 588. Kenne, J. P., & Gharbi, A. (1999). Experimental design in production and maintenance control of a single machine, single product manufacturing system. International Journal of Production Research, 37(3), 621– 637. Kusner, H. J., & Dupuis, P. G. (1992). Numerical methods for stochastic control problem in continuous time. New York: Springer. Kimemia, J. G., & Gershwin, S. B. (1983). An algorithm for computer control of production in flexible manufacturing systems. IIE Transactions, 15, 353– 362. Law, A. M., & Kelton, W. D. (2000). Simulation modeling and analysis (3rd ed). New York: McGraw-Hill. Older, G. J., & Suri, R. (1980). Time optimal part-routing in a manufacturing system with failure prone machines. Proceedings of 19th Conference on Decision Control, Albuquerque, NM, pp. 722–727. Pritsker, A. A. B., O’Reilly, J. J., & LaVal, D. K. (1997). Simulation with visual SLAM and AweSim. New York: Wiley. Rishel, R. (1975). Dynamic programming and minimal principles for systems with jump Markov disturbances, SIAM. Journal of Control, 13, 338 –371. Sethi, S. P., & Zhang, Q. (1994). Hierarchical control decision making in stochastic manufacturing systems. Boston: Birkhauser.