**American Journal of Industrial and Business Management**

Vol.09 No.01(2019), Article ID:89939,30 pages

10.4236/ajibm.2019.91011

A Mathematical Model: A Flexible Manufacturing System, Prone to Error, Making Two Products Each with Stochastic Demand Schedules

Issa Diop^{1}, Sylvie Nadeau^{2}, Behnam Emami-Mehrgani^{3}^{ }

^{1}Automated Production Engineering Department, École de technologie supérieure, Montréal, Canada

^{2}Mechanical Engineering Department, École de technologie supérieure, Montréal, Canada^{ }

^{3}Synetik Group, Joliette, Canada

Copyright © 2019 by author(s) and Scientific Research Publishing Inc.

This work is licensed under the Creative Commons Attribution International License (CC BY 4.0).

http://creativecommons.org/licenses/by/4.0/

Received: November 12, 2018; Accepted: January 13, 2019; Published: January 16, 2019

ABSTRACT

Lockout/tagout (LOTO) is practiced in manufacturing facilities to ensure safety during machinery maintenance procedures. In flexible manufacturing systems, human error (HE) is a major source of accidents and process deviations. Special measures are needed to minimize occupational risk and increase operational efficiency. In this article, we study a production planning problem involving a failure-prone production system meeting two types of demand and we discuss the associated decision-making process. The aim is to develop an optimal, robust and flexible control strategy that facilitates the integration of LOTO into corrective maintenance (CM) and ultimately into production. The influence of HE on flexible manufacturing systems (FMS) is viewed in terms of production and maintenance planning. The frequency of machine repair depends largely on HE. The intrinsic costs of shortage, inventory and CM are optimized over an unbounded planning horizon. Analytical formalism is combined with discrete event simulation, as well as design of experiments (DOE) and a genetic algorithm (GAs) to identify the optimal planning of production and CM with mandatory LOTO. An illustration and sensitivity analysis are proposed to express, in quantitative terms, the usefulness and efficiency of the proposed approach.

**Keywords:**

Optimal Production Control, Corrective Maintenance, Lockout/Tagout, Human Error

1. Introduction

Quebec regulation on occupational health and safety (RSST) states that maintenance technicians who service hazardous equipment must apply lockout/tagout (LOTO) procedures (cf. art. 188.2). This regulation sets forth safety precautions to follow before undertaking any maintenance activity, repair or unlocking operation involving such equipment. LOTO is defined as a preventive measure taken to avoid injuries caused by premature release of energy during installation, maintenance or repair of industrial machinery, equipment or process devices (D. 1187-2015, a. 3.). In other words, LOTO is a set of procedures and good practices intended to control all of the energy sources of a system during servicing. Only the person authorized to carry out the LOTO procedure on the system is authorized to remove the lock or tag. Inappropriate maintenance activities due to HE increases the total duration and intrinsic cost of servicing [1] and possibly the occupational health and safety risk.

According to Quebec standards, equity, health and work safety commission (CNESST) statistics for the year 2014, a work-related accident costs the life of a worker in Quebec approximately every six days (63 for the year). More than 80,000 accidents occur each year in the workplace. Machinery is involved in 10% - 15% of these. The injuries caused by machines are often severe, causing on average 20 fatalities, 295 amputations and 1125 permanent scars per year. The risk of accident generally increases throughout the machine lifecycle. Particular attention must be paid to risks associated with the production, adjustment and maintenance phases [2].

In September of 2015, the CNESST endorsed a draft proposing stricter occupational health and safety guidelines. The RSST includes now more specific regulations with regard to LOTO and hazardous energy control methods. These new regulations will spell out the duties of each actor and encourage workers and employers to stay mindful of health and safety in the workplace. Nevertheless, the use of LOTO in procedural form alone is an imperfect solution to the problem of accidents involving machinery [3]. A mechanism is needed for monitoring and follow-up of LOTO procedures. Such a mechanism has been proposed [4].

2. Literature Review

Maintenance technicians frequently need to enter hazardous spaces near or even within machinery or processing equipment since these devices break down repeatedly during their lifespan. In recent decades, hundreds of workers have been victims of accidents or incidents while performing routine maintenance on machinery or repairing damaged equipment in manufacturing facilities [5] [6]. The high level of risk inherent in machine maintenance has been confirmed repeatedly [7] [8] [9] [10] [11]. The Quebec health and safety work act (LSST (art. 51)) stipulates that among their obligations, employers are required 1) to ensure that employee working conditions are safe and not unhealthy, 2) to control health and safety risks associated with tasks, 3) to inform workers about these risks, and 4) to train workers to detect and avoid hazards. Furthermore, article 59 provides guidance and recommendations for eliminating risks at the source. On the other hand, workers must help to ensure occupational health and safety (art. 49). One of their principal duties in this respect is to comply with LOTO policies as set forth in occupational health and safety legislation. This is one way to avoid risks associated with equipment maintenance and repair [3]. Researchers point out that the legislation in its current form is of limited effectiveness, its principal shortcomings including inadequate provisions regarding training, inappropriate use of tags or locks, process design weaknesses and others [12]. In addition, many managers continue to regard LOTO as a time-consuming unproductive activity [13].

Recent studies show the numerous conditions that must be satisfied simultaneously for compliance with the recommendations of RSST regarding access to hazardous areas around machinery [14]. This includes alternative access to these areas when LOTO cannot be carried out, for example detection and diagnosis of malfunctions requiring observation of the machine in operation. On the other hand, LOTO is not mandatory when the machine ON/OFF switch is close to and under the exclusive control of the operator. A machine must have only one energy source, and any residual energy remaining after switching the machine off represents a hazard and must be released safely before any examination. It is the duty of the employer to specify these procedures and to ensure compliance with them in accordance with the RSST, in other words to inform and properly educate all workers having access to the hazardous areas with regard to the applicable preventive measures (LSST (1979, c. 63, a. 50)). These measures must indicate procedures to follow when shifts change, when LOTO is forgotten or when more than one employer (e.g. subcontractors) will be present. The LOTO program must guide managers as well as workers, and should cover [15] 1) what hazardous energy is, 2) what types of tools or devices are required for isolating or de-energizing hazardous devices, 3) who is responsible for the tasks, 4) the steps for turning on and off and de-energizing machinery, 5) what steps for locking/tagging the system to put in writing, 6) what training is required, 7) what aspects must be covered by audits. The implementation criteria of these occupational health and safety regulations are described in detail in Canadian standard CSA-Z460 2013―Control of hazardous energy―Lockout and other methods. However, it should be noted that a LOTO system is not an end in itself and must not be viewed as just a padlock [16]. The system must be monitored after its implementation. The documentation and bill of material coding are as crucial as the LOTO mechanism. The documentation must contain the LOTO files and procedures. Software is needed for the audit management, as well as training and reporting. Whereas the success of a LOTO program lies in the quality of daily monitoring (i.e. of reports and audits) [17].

Venkatraman [18] has pondered the connection between safety and maintenance effectiveness. These questions have been of interest to manufacturers for decades and have led to significant progress in understanding of workplace safety [19]. The CSA Z1000 standard indicates that a LOTO system must be integrated into occupational health and safety policies.

We have noted the existence of three categories of study in this subject area: 1) optimizing production planning and maintenance strategies [20] [21] [22] [23] , 2) taking LOTO into account in maintenance activities [24] [25] and 3) taking human factors and human error (HE) into account [26] [27]. The number of studies in this third category is small.

In production scheduling, there is considerable leeway for varying the frequency of machine maintenance. A planner might prefer to stop a machine when its performance has dropped below some threshold or run it to breakdown. In the latter case, the machine may be repaired immediately, or repair may be deferred. Some authors have suggested that manufacturing system performance improves when maintenance task and production activities are scheduled at the same time [28]. A policy based on a stochastic approach to planning production and preventive maintenance activities in an FMS had been developed more than a decade earlier [29]. This strategy incorporated periodic maintenance with the increase in the likelihood of failure depending on machine age. Not long after, the optimality conditions of the hedging point policy were studied [30]. It was thus shown that this policy is optimal when the failure rate does not depend on the production rate; the linearity of the breakdown rate function is sufficient to show that a hedging point policy is optimal. An extension of these models was then proposed, in which maintenance activity scheduling is based on inventory level [31]. This was then expanded to production systems comprising several machines producing several parts [32] [33]. The gains achieved under this policy are measurable in terms of downtime. The problem of planning production and preventive maintenance for a manufacturing operation involving several machines was then treated using analytical formalism combined with simulation and response surface methodology to develop an approximation of optimal control policies [34]. These authors provide an illustration and sensitivity analysis to quantify production and preventive maintenance rates that minimize the production costs and the frequency of maintenance and repair.

Production and maintenance planning have been optimized in conjunction with LOTO policies based on consideration of two types of corrective maintenance (CM): 1) routine failure associated with machine age, in which the mean time between failures decreases with machine deterioration and no LOTO procedure is required, and 2) major failure, with a constant mean time between failures and LOTO is required for repair [25]. In view of the effectiveness of such planning, the study was expanded to include the elimination of possibilities of circumventing LOTO [24]. Considering lockout as a separate and integral machine state while performing maintenance activities is suggested as means of achieving this, and optimal planning for modern manufacturing systems subject to stochastic breakdowns is proposed. This raised the question of the existence of an optimal cost, taking into account reductions in the risk of accidents. Based on numerical illustration and sensitivity analysis, it has been shown that controlling LOTO rather than setting it at a fixed frequency offers more optimal planning [24]. However, the complexity of the problem increases rapidly with the number of machines.

Nowadays, human factors, particularly HE and human reliability, are receiving much attention. In manufacturing settings, HE occurs at a significant rate during equipment maintenance activities [35] and is involved in most of the accidents or incidents that occur during maintenance activities in conjunction with LOTO. Although contributing to more than 25% of the breakdowns that process plant machines undergo during their lifespan, human factors are still not given sufficient consideration [18]. HE or control failures are often behind undesirable events that happen due to unforeseen releases of hazardous energy [36]. For example, most of the deviations (>80%) in quality and productivity in pharmaceutical production appear to be due to HE [37]. Despite technological control of risks inherent in equipment, the probability of HE in machine maintenance remains high, as does production equipment downtime [38]. The findings of more recent studies on the role of HE in risk analysis have been applied to pre-maintenance and post-maintenance procedures. These authors later estimated HE probabilities for several possible failure scenarios associated with the maintenance procedures for a pump [39]. Using HEART (human error assessment and reduction technique, developed in 1988 by J.C. Williams to evaluate the probability of HE occurring while performing a particular task, a reliable technique used in safety analyses), [39] calculated acceptable levels of risk (based on the level of error recovery) beyond which actions should be taken through risk management strategies to raise the level of safety of the maintenance procedure. Other studies show that HE lengthens repair time and increases production costs, inventory costs and shortage as well as the risk of workplace accidents [26]. An optimal policy has been defined to minimize production cost, emphasizing the impact of human on inventory shortage and system capacity [40]. The authors examined occupational health and safety risk acceptability from eight perspectives: economic, personal, cultural, political, social, ethical, psychological, and risk characteristics. These parameters influence the acceptability mechanism and its pervasiveness in workplace safety.

Two approaches appear predominant in system failure research [41]. The first focuses on humans (inattention, forgetfulness and so on) whereas the second focuses on the system (work conditions). Accidents happening to individuals are more recurrent and catastrophic compared to organizational accidents and have adverse side effects in terms of harm, loss of quality of life or of life itself [42]. Studies show the typical limitations on human perception, cognition and physical performance [43]. Factors such as fatigue and stress influence human performance to the point that no amount of effort can eradicate HE [44] [45] [46] [47]. To deal adequately with this reality, Harris et al. [43] recommend keeping two specific goals in mind: 1) reducing HE by implementing a systematic approach and a design methodology focused on human operators, 2) quick identification and correction of sources of malfunction, on the assumption that human error will occur inevitably.

Many authors have militated in favor of flexible manufacturing systems and emphasized their importance in productivity and quality in many manufacturing industries [48] [49]. Indeed, human intervention is less important in such systems compared to those based on conventional material handling. In the latter, the human operator is constantly involved in the transportation of materials between sites, whereas in a flexible manufacturing system, this function is supported fully, for example by automated guided vehicles and material handling systems. This subject has been reviewed [50]. Human intervention is nevertheless required in flexible systems during start-up of operations, as well as in programming, repair activities or loading and unloading of materials and parts. There is a need to take into consideration the system/human mutual impact in terms of HE during these phases. These states lead to the concept of industry 4.0. Elsewhere, Bengtsson [51] warns about this concept very famous in the maintenance community nowadays. These authors show through a case study the necessity of taking basic maintenance concepts and management together. In the same idea, Salonen [52] noticed that on average, 40% of equipment failures are related to improper maintenance practices by human operators.

The purpose of this paper is to develop a production and CM plan with LOTO to improve the safety of a flexible manufacturing system consisting of a failure-prone machine meeting two types of demand. The main contribution of this study is a control strategy that facilitates the integration of LOTO into production and considers the impact of HE during maintenance. Optimal costs of shortage, inventory build-up and CM are identified over an unbounded planning horizon.

Notations and assumptions are presented in Section 3. The FMS under study is presented in Section 4. A numerical illustration is developed in Section 5. A discussion is provided in Section 6. Section 7 is devoted to the conclusion and proposal of future studies.

3. Notations and Assumptions

3.1. Notations

The notations used in this study are as follows:

3.2. Assumptions

The following assumptions are considered:

CM is performed with LOTO;

CM may be completed with HE;

Repairing time with HE is longer than it is without HE;

CM activities restore the system performance to a new one (as good as new).

4. Problem Statement

The subject of this paper is a decision-making problem in an FMS consisting of one failure-prone machine meeting two types of demand. Machine availability is enhanced through CM activities. The decision variables are the frequencies of machine repair with human error (WHE) or without human error (WOHE), which influence the system capacity and the inventory level. HE is considered as an inappropriate action that results in increased CM duration and costs [53]. We are investigating the impact of HE that occurs during CM activities in conjunction with LOTO. We determine an optimal total cost associated with shortage, inventory and CM over an unbound planning horizon. Analytical modeling and numerical resolution approaches were used with discrete events simulation. DOE and a GA were then used to define the optimal plan. Figure 1 summarizes the methodological approach:

Many researchers [13] [24] [25] [54] use in their studies the same approach based on analytical modeling and numerical resolution. To achieve this, they used assumptions that do not take into account phenomena such as the impact of human factors. In response to this situation, a sequential optimization approach combining a numerical method, simulation, experimental designs, genetic algorithm is promising for dealing with optimal control problems. GAs are simple and fast optimization methods for a global exploration and allow to obtain a solution in an interesting time.

Optimization based on simulation makes the resolution of optimal control problems simpler compared to conventional direct search techniques which are too demanding in computing time. Indeed, direct research requires many hypotheses for numerical resolution of control problems compared to simulation models. In addition, simulation techniques provide a detailed description of the dynamic behavior of the manufacturing system [55]. Azadivar [56] recommends simulation as an optimization approach for choosing maintenance strategies. Indeed, in their work, they noticed that the results from the simulation experiment indicate the same forms of response surfaces as those derived from direct search optimization methods. However, the problem of these classical methods (direct search) are not often easy to solve and often involve many assumptions, because the performance of the system depends on the combination of quantitative variables and qualitative variables (choice of maintenance strategy). For more details on these aspects of performance, the reader is referred to reference [56]. A combined approach, which is based on a combination of analytical formalism, simulation modelling, design of experiments and a genetic algorithm presents a more realistic model for industries.

Treating this task as a stochastic dynamic programming problem, we developed

Figure 1. Methodological approach.

a continuous-time Markov chain model integrating possibilities of HE during maintenance activities with mandatory LOTO by human operators over an unbounded planning horizon. Three variables characterize the system in any given state: the stock level $x(.)$ , the stochastic process of machine states $\xi (.)$ and the demand type $d\left(\alpha \right)={\left({d}_{1},{d}_{2}\right)}^{\prime}$ . Raw material is always available (Figure 2).

The dynamics of such a system include continuous and discrete variables. The continuous variable $x(.)$ represents the cumulative inventory/shortage vector, which can be positive (stock build-up) or negative (shortage). The dynamics of $x(.)$ are presented by the differential equation below:

$\frac{\text{d}x(\cdot )}{\text{d}t}=u(\cdot )-d\left(\alpha \right)$ (1)

where: $x\left(0\right)=x$ initial stock level); $u(\cdot )={\left({u}_{1},{u}_{2}\right)}^{\prime}$ ; $d\left(\alpha \right)={\left({d}_{1},{d}_{2}\right)}^{\prime}$ ; $\alpha (\cdot )\in \Omega =\left\{1,2,3,4,5,6\right\}.$

The discrete variable $\xi (.)$ defines machine status (six possible status):

$\Omega \in \left\{\begin{array}{l}\xi \left(t\right)=\text{1:The FMS is working on demand}{d}_{1}\\ \xi \left(t\right)=\text{2:The FMS is shut down for CM without Human Error}\left(\text{WOHE}\right)\text{with demand}{d}_{1}\\ \xi \left(t\right)=\text{3:The FMS is shut down for CM with Human Error}\left(\text{WHE}\right)\text{with demand}{d}_{1}\\ \xi \left(t\right)=\text{4:The FMS is working on demand}{d}_{2}\\ \xi \left(t\right)=\text{5:The FMS is shut down for CM without Human Error}\left(\text{WOHE}\right)\text{with demand}{d}_{2}\\ \xi \left(t\right)=\text{6:The FMS is shut down for CM with Human Error}\left(\text{WHE}\right)\text{with demand}{d}_{2}\end{array}\right\}$ (2)

Its representation as a continuous-time Markov chain (CTMC) is shown in Figure 3.

The mathematical model of the continuous-time stochastic process takes values in the finite state space $\Omega \left(1,2,3,4,5,6\right)$ .

The transition probabilities from state $\alpha \in \Omega \left(1,2,3,4,5,6\right)$ to state $\beta \in \Omega \left(1,2,3,4,5,6\right)$ at time t are given in Equation (3):

Figure 2. The flexible manufacturing system under study.

Figure 3. Transition diagram of machine states.

$P\left(\xi \left(t+\delta t\right)=\beta |\xi \left(t\right)=\alpha \right)=\{\begin{array}{l}{q}_{\alpha \beta}\delta t+o\left(\delta t\right)\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.05em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{if}\text{\hspace{0.17em}}\text{\hspace{0.05em}}\alpha \ne \beta \\ 1+{q}_{\alpha \beta}\delta t+o\left(\delta t\right)\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{if}\text{\hspace{0.17em}}\text{\hspace{0.05em}}\alpha =\beta \end{array}$ (3)

The transition rate is given by the following equations:

$\begin{array}{l}{q}_{\alpha \beta}=\underset{\delta t\to 0}{\mathrm{lim}}\left[P\left(\xi \left(t+\delta t\right)=\beta |\xi \left(t\right)=\alpha \right)/\delta t\right]\\ \text{knowing that:}\alpha \ne \beta \text{,}\underset{\delta t\to 0}{\mathrm{lim}}\left(\frac{o\left(\delta t\right)}{\delta t}\right)=0\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\forall \alpha ,\beta \in \Omega \end{array}$ (4)

The transition rates matrix $Q(\cdot )$ of the continuous-time Markov process (ξ(t)) is defined such that it meets conditions (5) and (6) below:

$Q={q}_{\alpha \beta}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{where}\text{\hspace{0.17em}}{q}_{\alpha \beta}\ge 0\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{if}\text{\hspace{0.17em}}\alpha \ne \beta $ (5)

${q}_{\alpha \alpha}=-\underset{\beta \ne \alpha}{\Sigma}{q}_{\alpha \beta}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{where}\text{\hspace{0.17em}}\alpha ,\beta \in \Omega =\left\{1,2,3,4,5,6\right\}$ (6)

The transitions rates matrix $Q={q}_{\alpha \beta}$ describes the stochastic process (Equation (7)):

$Q(\cdot )=\left|\begin{array}{cccccc}-\left({q}_{12}+{q}_{14}\right)& {q}_{12}& 0& {q}_{14}& 0& 0\\ {r}_{21}& -\left({r}_{21}+{q}_{23}+{q}_{25}\right)& {q}_{23}& 0& {q}_{25}& 0\\ {r}_{31}& 0& -\left({r}_{31}+{q}_{36}\right)& 0& 0& {q}_{36}\\ {q}_{41}& 0& 0& -\left({q}_{41}+{q}_{45}\right)& {q}_{45}& 0\\ 0& {q}_{52}& 0& {r}_{54}& -\left({q}_{52}+{r}_{54}+{q}_{56}\right)& {q}_{56}\\ 0& 0& {q}_{63}& {r}_{64}& 0& -\left({q}_{63}+{r}_{64}\right)\end{array}\right|$ (7)

The admissible production decisions set $B(\cdot )$ and control variables at state $\alpha \in \Omega =\left\{1,2,3,4,5,6\right\}$ are defined as follows:

$B\left(\alpha \right)=\left\{\left[u(\cdot ),{r}_{21}(\cdot ),{r}_{31}(\cdot ),{r}_{54}(\cdot ),{r}_{64}(\cdot )\right]\in {R}^{5},0\le {\Vert u(\cdot )\Vert}_{\infty}\le {u}^{\mathrm{max}}\right\}$ (8)

The total cost function $J(\cdot )$ is given in Equation (9). The total cost function $J(\cdot )$ is convex and non-negative. Our goal is to control the production rate $u(\cdot )$ as demand varies, to minimize the total production cost.

$J(\cdot )=E\left\{{\displaystyle \underset{0}{\overset{\infty}{\int}}{\text{e}}^{-\rho t}}g(.)\text{d}t|x\left(0\right)=x,\xi \left(0\right)=\alpha \right\}$ (9)

$E\{.\}$ is the mathematical expectation $g(\cdot )$ is the instantaneous cost over an unbound planning horizon. The instantaneous cost is given as follows:

$g(\cdot )={c}_{1}^{+}{x}_{1}^{+}+{c}_{1}^{-}{x}_{1}^{-}+{c}_{2}^{+}{x}_{2}^{+}+{c}_{2}^{-}{x}_{2}^{-}+{c}_{\alpha}$

$\begin{array}{l}\text{knowingthat}{x}^{+}=\mathrm{max}\left(0,x\right)\text{,}{x}^{-}=\mathrm{max}\left(-x,0\right)\\ {c}_{\alpha}=\left({c}_{r}+{c}_{LOTO}\right)Ind\left(\alpha =2\right)+\left({c}_{r}+{c}_{LOTO}\right)Ind\left(\alpha =5\right)\\ \text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}+\gamma \left({c}_{r}+{c}_{LOTO}\right)Ind\left(\alpha =3\right)+\gamma \left({c}_{r}+{c}_{LOTO}\right)Ind\left(\alpha =6\right)\\ Ind\langle \Theta (.)\rangle =\{\begin{array}{l}1\text{if}\Theta (.)\text{istrue}\\ 0\text{otherwise}\end{array}\end{array}$ (10)

${c}^{+}$ is the inventory cost per unit; ${x}^{+}$ is the inventory level; ${c}^{-}$ is the shortage cost per unit; ${x}^{-}$ is the shortage level; ${c}^{\alpha}$ is the operation cost under CM; $\gamma $ is the cost index $\left(\gamma \succ 1\right)$

$\begin{array}{l}\Theta \text{isthesystemstate}\left(\alpha \in \Omega =\left\{1,2,3,4,5,6\right\}\right)\text{:}\\ \left\{\begin{array}{l}\alpha =2\to Ind\left(\alpha =2\right)=1\to {c}^{\alpha}={c}_{r}+{c}_{LOTO}\\ \alpha =5\to Ind\left(\alpha =5\right)=1\to {c}^{\alpha}={c}_{r}+{c}_{LOTO}\\ \alpha =3\to Ind\left(\alpha =3\right)=1\to {c}^{\alpha}=\gamma \left({c}_{r}+{c}_{LOTO}\right)\\ \alpha =6\to Ind\left(\alpha =6\right)=1\to {c}^{\alpha}=\gamma \left({c}_{r}+{c}_{LOTO}\right)\end{array}\right\}\end{array}$

The mathematical value function ${\vartheta}_{x}(\cdot )$ is given in Equation (11). The value function should satisfy a set of specific properties represented by the Hamilton-Jacobi-Bellman (HJB) partial differential equations.

$\vartheta (\cdot )=\underset{\left(u,r\right)\in B\left(\alpha \right)}{\mathrm{min}}J\left(x,\alpha ,u,r\right)\text{\hspace{0.17em}}\text{\hspace{0.17em}}\forall \alpha \in \Omega =\left\{1,2,3,4,5,6\right\}$ (11)

The HJB equations characterize the optimality conditions (cf. Appendix 1) for the FMS:

$\begin{array}{l}\rho \vartheta \left(x,\alpha \right)=\underset{\left(u,r\right)\in B\left(\alpha \right)}{\mathrm{min}}\left\{\left[u\left(\alpha \right)-d\left(\alpha \right)\right]\frac{\partial}{\partial x}\vartheta \left(x,\alpha \right)+\left(g\left(x,\alpha \right)+{\displaystyle \underset{\alpha \ne \beta}{\sum}{q}_{\alpha \beta}\vartheta \left(x,\beta \right)}\right)\right\}\\ \forall \alpha ,\beta \in =\left\{1,2,3,4,5,6\right\}\\ \text{where},\{\begin{array}{l}u\left(\alpha \right)=\{\begin{array}{l}{u}_{1}\text{if}\alpha =\text{1}\\ {u}_{2}\text{if}\alpha =\text{4}\end{array}\\ d\left(\alpha \right)=\{\begin{array}{l}{d}_{1}\text{if}\alpha =\text{1,2,3}\\ {d}_{2}\text{if}\alpha =\text{4,5,6}\end{array}\end{array}\end{array}$ (12)

The optimal control policy stands for minimizing the HJB partial differential equations over the set of admissible production capacity $B(\cdot )$ for each $\alpha \in \Omega =\left\{1,2,3,4,5,6\right\}$ . A numerical method for stochastic control problems based on the Kushner and Dupuis method [33] is used to solve the optimality conditions of the HJB equations. Let ${h}_{j},j=1,\cdots ,n$ designate the increment of the finite difference interval of the variable ${x}_{j}$ . Emami-Mehrgani [13] provides a demonstration showing that the value function ${\vartheta}_{i}(\cdot )$ can be approximated by

${\vartheta}_{i}^{h}(\cdot )$ and the partial derivative of the value function ${\left({\vartheta}_{i}\right)}_{{x}_{j}}$ for a small increment of ${h}_{x}$ :

${\left({\vartheta}_{i}\right)}_{{x}_{j}}(.)=\{\begin{array}{l}\frac{1}{{h}_{j}}\left[{\vartheta}_{i}^{h}\left({x}_{1},\cdots ,{x}_{j}+{h}_{j},\cdots ,{x}_{n}\right)-{\vartheta}_{i}^{h}\left({x}_{1},\cdots ,{x}_{j},\cdots ,{x}_{n}\right)\right]\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{if}\text{\hspace{0.17em}}{u}_{j}-{d}_{j}\ge 0\\ \frac{1}{{h}_{j}}\left[{\vartheta}_{i}^{h}\left({x}_{1},\cdots ,{x}_{j},\cdots ,{x}_{n}\right)-{\vartheta}_{i}^{h}\left({x}_{1},\cdots ,{x}_{j}-{h}_{j},\cdots ,{x}_{n}\right)\right]\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{if}\text{\hspace{0.17em}}{u}_{j}-{d}_{j}\prec 0\end{array}$ (13)

From Equation (13) we derive Equation (14):

$\left[u\left(\alpha \right)-d\left(\alpha \right)\right]\vartheta \left(x,\alpha \right)=\{\begin{array}{l}\frac{\left|u\left(\alpha \right)-d\left(\alpha \right)\right|}{{h}_{x}}\left[{\vartheta}^{h}\left(x+{h}_{x},\alpha \right)-{\vartheta}^{h}\left(x,\alpha \right)\right]\text{\hspace{0.17em}}\text{if}\text{\hspace{0.17em}}{u}_{j}-{d}_{j}\ge 0\\ \frac{u\left(\alpha \right)-d\left(\alpha \right)}{{h}_{x}}\left[{\vartheta}^{h}\left(x,\alpha \right)-{\vartheta}_{i}^{h}\left(x-{h}_{x},\alpha \right)\right]\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{if}\text{\hspace{0.17em}}{u}_{j}-{d}_{j}\prec 0\end{array}$ (14)

The dynamic programming equation resulting from the HJB equations therefore can be presented as shown below [29] [32] :

$\begin{array}{l}{\vartheta}^{h}\left(x,\alpha \right)=\underset{\left(u,r\right)\in {B}^{h}\left(\alpha \right)}{\mathrm{min}}\{{\left(\rho +\left|{q}_{\alpha \alpha}\right|+\frac{\left|u\left(\alpha \right)-d\left(\alpha \right)\right|}{{h}_{x}}\right)}^{-1}\\ \text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\times \langle \frac{\left|u\left(\alpha \right)-d\left(\alpha \right)\right|}{{h}_{x}}\left(\Xi \right)+\left(g\left(x,\alpha \right)+Q\cdot {\vartheta}^{h}\left(x,\alpha \right)\right)\rangle \}\\ \text{with},\text{\hspace{0.17em}}\text{\hspace{0.17em}}\Xi =\left({\vartheta}^{h}\left(x+{h}_{x},\alpha \right){K}^{+}+{\vartheta}^{h}\left(x-{h}_{x},\alpha \right){K}^{-}\right)\\ \{\begin{array}{l}{K}^{+}=\{\begin{array}{l}1\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{if}\text{\hspace{0.17em}}\left[u\left(\alpha \right)-d\left(\alpha \right)\right]\ge 0\\ 0\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{otherwise}\end{array}\\ {K}^{-}=\{\begin{array}{l}1\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{if}\text{\hspace{0.17em}}\left[u\left(\alpha \right)-d\left(\alpha \right)\right]\prec 0\\ 0\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{otherwise}\end{array}\end{array}\end{array}$ (15)

5. Numerical Example

A numerical approach [58] is used to solve the HJB equations for the system under study. It is described in details in Appendix 2. The system is feasible if inequality 16 is satisfied:

$\begin{array}{l}\Sigma {\pi}_{i}\left(\alpha \right){U}_{i}^{\mathrm{max}}\left(\alpha \right)\ge d\left(\alpha \right)\\ \text{where,}\\ \{\begin{array}{l}{\pi}_{i}\text{is the probality limi at the state}\alpha \text{,}\\ {U}_{i}^{\mathrm{max}}\text{is the maximum production rate at state}\alpha \text{}\end{array}\\ \left\{\left[{\pi}_{1}{U}_{1}^{\mathrm{max}}Ind\left(\alpha =1\right)+{\pi}_{4}{U}_{4}^{\mathrm{max}}Ind\left(\alpha =4\right)\right]\ge \left[{d}_{1}Ind\left(\alpha =1\right)+{d}_{2}Ind\left(\alpha =4\right)\right]\right\}\end{array}$ (16)

The feasibility in equation is structured so that it satisfies a set of specific properties called probability limits, the details of which are provided below:

$\begin{array}{l}{\pi}_{i}\left(\alpha \right)Q\left(\alpha \right)=0\\ \text{where}Q\text{is the matrix of transition rates}\\ \{\begin{array}{l}{\pi}_{i}\left(\alpha \right)={\pi}_{1}Ind\left(\alpha =1\right)+{\pi}_{{}_{2}}Ind\left(\alpha =2\right)+{\pi}_{3}Ind\left(\alpha =3\right)\\ \text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}+{\pi}_{4}Ind\left(\alpha =4\right)+{\pi}_{5}Ind\left(\alpha =5\right)+{\pi}_{6}Ind\left(\alpha =6\right)\\ {\displaystyle \underset{i=1}{\overset{6}{\sum}}{\pi}_{i}\left(\alpha \right)=1}\end{array}\end{array}$ (17)

$\{\begin{array}{l}\left|{\pi}_{1}\text{\hspace{0.17em}}\text{\hspace{0.17em}}{\pi}_{2}\text{\hspace{0.17em}}\text{\hspace{0.17em}}{\pi}_{3}\text{\hspace{0.17em}}\text{\hspace{0.17em}}{\pi}_{4}\text{\hspace{0.17em}}\text{\hspace{0.17em}}{\pi}_{5}\text{\hspace{0.17em}}\text{\hspace{0.17em}}{\pi}_{6}\right|\left|\begin{array}{cccccc}-\left({q}_{12}+{q}_{14}\right)& {q}_{12}& 0& {q}_{14}& 0& 0\\ {r}_{21}& -\left({r}_{21}+{q}_{23}+{q}_{25}\right)& {q}_{23}& 0& {q}_{25}& 0\\ {r}_{31}& 0& -\left({r}_{31}+{q}_{36}\right)& 0& 0& {q}_{36}\\ {q}_{41}& 0& 0& -\left({q}_{41}+{q}_{45}\right)& {q}_{45}& 0\\ 0& {q}_{52}& 0& {r}_{54}& -\left({q}_{52}+{r}_{54}+{q}_{56}\right)& {q}_{56}\\ 0& 0& {q}_{63}& {r}_{64}& 0& -\left({q}_{63}+{r}_{64}\right)\end{array}\right|\\ =\left|0\text{\hspace{0.17em}}\text{\hspace{0.17em}}0\text{\hspace{0.17em}}\text{\hspace{0.17em}}0\text{\hspace{0.17em}}\text{\hspace{0.17em}}0\text{\hspace{0.17em}}\text{\hspace{0.17em}}0\text{\hspace{0.17em}}\text{\hspace{0.17em}}0\right|\\ {\pi}_{1}+{\pi}_{2}+{\pi}_{3}+{\pi}_{4}+{\pi}_{5}+{\pi}_{6}=1\end{array}$ (18)

A policy improvement approach is used to determine an approximate solution. The search algorithm is presented below:

Step 1:

For a given production control policy ${u}_{n}(\cdot )$ , where n is the number of iterations and h is a given finite difference interval, let ${\varpi}_{{u}_{n}(\cdot )}$ and ${\varpi}_{{u}_{n}(\cdot )}^{\ast}$ be defined by

${\varpi}_{{u}_{n}(\cdot )}\left({\vartheta}^{h}(.)\right)={\vartheta}^{h}(.)$ ; ${\varpi}_{{u}_{n}(\cdot )}^{*}\left({\vartheta}^{h}(.)\right)=\underset{{u}_{n}(\cdot )\in \Gamma (.)}{\mathrm{min}}\left\{{\varpi}_{{u}_{n}(\cdot )}\left({\vartheta}^{h}(.)\right)\right\}$

$\delta \in R$ a given accuracy and ${\vartheta}^{h}(\cdot )$ the value function.

$\begin{array}{l}n=1\\ \delta \in {R}^{+}\\ {\left({\vartheta}^{h}(\cdot )\right)}^{n}:=0\\ \forall \alpha \in \left(t\right)\in \Omega =\left\{1,2,3,4,5,6\right\}\end{array}$

Step 2:

Compute the new value of the value function at rank n and compare it with the old value in memory.

$\begin{array}{l}{\left({\vartheta}^{h}(\cdot )\right)}^{n-1}:={\left({\vartheta}^{h}(\cdot )\right)}^{n}\\ \forall \alpha \in \left(t\right)\in \Omega =\left\{1,2,3,4,5,6\right\}\end{array}$

Step 3:

Determine the production control policy ${u}_{n(.)}$

$\begin{array}{l}{\varpi}_{{u}_{n}(\cdot )}{\left({\vartheta}^{h}(.)\right)}^{n-1}:={\left({\vartheta}^{h}(.)\right)}^{n}:={\varpi}_{{u}_{n}(\cdot )}^{*}{\left({\vartheta}^{h}(.)\right)}^{n-1}\\ \forall \alpha \in \left(t\right)\in \Omega =\left\{1,2,3,4,5,6\right\}\end{array}$

Step 4:

Compute the convergence tests and repeat the above steps until a fixed point $\delta $ precision is found.

$\begin{array}{l}{\epsilon}_{\mathrm{inf}}:=\underset{\begin{array}{l}\alpha \in M\\ {z}_{i}\end{array}}{\mathrm{min}}\left\{{\left({\vartheta}^{h}(.)\right)}^{n}-{\left({\vartheta}^{h}(.)\right)}^{n-1}\right\}\left[\right]\\ {\epsilon}_{\mathrm{sup}}:=\underset{\begin{array}{l}\alpha \in M\\ {z}_{i}\end{array}}{\mathrm{max}}\left\{{\left({\vartheta}^{h}(.)\right)}^{n}-{\left({\vartheta}^{h}(.)\right)}^{n-1}\right\}\\ {\epsilon}_{\mathrm{min}}:=\frac{\rho}{1-\rho}{\epsilon}_{\mathrm{inf}}\end{array}$

${\epsilon}_{\mathrm{max}}:=\frac{\rho}{1-\rho}{\epsilon}_{\mathrm{sup}}$

$\begin{array}{l}\text{if}\left|{\epsilon}_{\mathrm{max}}-{\epsilon}_{\mathrm{min}}\right|\le \delta \to \{\begin{array}{c}\text{stop}\\ {u}^{*}(.):={u}_{n}(.)\\ {u}_{1}={u}_{1}^{*}\text{,}{u}_{2}={u}_{2}^{*}\text{}\end{array}\\ \text{else,}\{\begin{array}{l}{\left({\vartheta}^{h}(.)\right)}^{n}:={\varpi}_{{u}_{n}(\cdot )}{\left({\vartheta}^{h}(.)\right)}^{n}\\ n:=n+1\end{array}\end{array}$

The system capacity is described as a Markov process with six state $\alpha \in \Omega =\left\{1,2,3,4,5,6\right\}$ The system is subject to random failure requiring repair. Machine availability is enhanced through CM activities. The frequencies of machine repair with or without HE influences system capacity and output, which determine the inventory levels and are the decision variables. Table 1 indicates the computational parameters used to run the numerical model.

The following mesh grid is defined to perform the computational domain $D=\left\{x\in Z/-5\le x\le 50\right\}$ .

5.1. Production Policy Structure

The policy outcomes for the manufacturing control variables ${u}_{i}(.)$ are shown in Figure 4 and Figure 5.

Since HJB differential equations cannot be solved analytically, a finite-difference (discretization) method [58] was used to obtain approximate numerical solutions. The resulting value functions are the solutions of the optimality conditions. They represent the cost of operating the production system. Figure 4 and Figure 5 respectively show the outputs of products 1 and 2. These outcomes show that there is no manufacturing urgency for either product (inventory level is adequate). The control policy advises keeping the machine idle in an energy-conserving mode to save power. However, when the inventory drops below a specific threshold, the production control policy recommends setting the machine at its maximal output. Otherwise, its operation is adjusted to meet customer demand.

Based on the results above, the policy illustrated can be summarized as follows:

$\left\{\begin{array}{l}{u}_{1}\left({x}_{1},{x}_{2},1\right)\{\begin{array}{l}{u}_{1}^{\mathrm{max}}\text{if}{x}_{1}\prec {z}_{1}\\ {d}_{1}\text{if}{x}_{1}={z}_{1}\\ 0\text{otherwise}\end{array}\\ {u}_{2}\left({x}_{1},{x}_{2},4\right)\{\begin{array}{l}{u}_{2}^{\mathrm{max}}\text{if}{x}_{2}\prec {z}_{2}\\ {d}_{2}\text{if}{x}_{2}={z}_{2}\\ 0\text{otherwise}\end{array}\end{array}\right\}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{where:}\{\begin{array}{l}{z}_{1}=51.2\\ {d}_{1}=0.775\end{array}\text{}\{\begin{array}{l}{z}_{2}=50.2\\ {d}_{1}=0.650\end{array}$ (19)

Table 1. Values of the numerical parameters.

Figure 4. FMS production rate for product 1.

Figure 5. FMS production rate for product 2.

The above policy is based on the optimal rates of production ${u}_{i}(.)$ that depend on ${z}_{i}$ . It is known as the hedging-point policy (HPP), as described previously [32] [57] [58]. The following sections describe the optimization of ${z}_{i}$ using a simulation approach and DOE combined with a GA.

5.2. Simulation Model

For solving optimal control problems, DOE, statistical analysis, and response surface methodology (RSM) are often used, in combination with simulation, to model the behavior of the studied system and the interactions that have an impact on its performance [20] [55] [58]. At the end, system performance is expressed through linear regression. In doing so, the optimum conditions of the system can be approximated by minimizing the cost function (regression equation) thus obtained. Discrete/continuous event simulation models are used to obtain the real-time production cost.

The simulation software (Arena) and SIMAN language were chosen to perform the discrete simulation. Figure 6 shows the simulation model. The control policy characterized above is used as the simulation model input.

Figure 6. Simulation model.

The inventory trajectories and the time elapsed in the system for products 1 and 2 are presented below in Figure 7 and Figure 8.

The graph of average inventory has a repeating trapezoidal profile because of the limit on the machine output. The stock level increases to this maximum, remains there and then decreases during corrective maintenance (CM) activities. When the machine breaks down, a stock shortage occurs (x drops below zero). After the machine is repaired, the stock begins to increase again and the cycle repeats. This strategy recommends building up inventory while the machine is available at its full capacity to meet demand while production capacity is nil due to CM activities with LOTO. However, this requires accelerating production or paying the workers for overtime, which carry a risk of HE and accidents. The curves of the time spent in the system show that the products are produced at the same pace. We note an initial phase where the time varies considerably, representing the transitional regime, and then varies very little around an average value representing the steady state. The difference in values comes from the fact that each part type has its own distinctive processing time.

5.3. Design of Experiment and Genetic Algorithm Methodologies

The numerical approach enabled us to solve the HJB equations and thereby obtain a first draft of the optimal production control policy, which is a threshold

Figure 7. Average inventory level.

Figure 8. Average time in system.

policy (HPP). The simulation software (Arena) and SIMAN language were then used to perform the discrete event simulation.

Simulation, on the other hand, as a decision support tool alone, is not able to directly perform the solving calculations of production system control problems [59]. However, simulation helps to facilitate the analysis, interpretation and understanding of the behavior of the production system. So, it is necessary to combine this approach with other optimization techniques, in this case meta-heuristics and evolutionary optimization, for example the genetic algorithms. For more details on the simulation-based optimization approach, the reader can refer to the work of Kenné [33] ; Azadivar [56] [60] and their bibliographic references.

The advantage of heuristic resolution approaches, including GA, lies in their ability to provide an acceptable solution in a record time of non-prohibitive calculation. It is worth mentioning that the solution obtained with a heuristic approach is not necessarily the optimal value of the problem to be solved. Nevertheless, the heuristic resolution techniques allow to obtain global optima, allowing the resolution of the problem with a realistic computation time [61].

In this section, the DOE approach and a GA were combined. The DOE model characterized the variation of the control factors to identify the main factor effects (and those of their interactions) on the average cost. The independent variables ranged from a low of 1 to a high of 100 (Table 2).

One response variable and two experimental factors were specified, allowing an experimental design with 27 runs, one sample taken per run. The default model is quadratic with 6 coefficients. P values below 0.05 indicate a probability of less than 5% that the fit with the model is due to chance, while the R2 statistic indicates the percentage of the variation in the response that is explained by the fitted model. In this case, five effects have P values less than 0.05, indicating that it is at least 95% certain that their effects are real (not due to chance), while the model as fitted explains 94.47% (adjusted R2) of the variability in cost. The average cost is a function of z_{1} and z_{2}. In fact, there is a close correlation between the use of DOE and the regression analysis [62]. The regression equation that was fitted to the data is shown below (Equation (20)). The regression coefficients for the cost function are given Equation (20).

$\text{cost}={\beta}_{0}+{\displaystyle \underset{i=1}{\overset{2}{\sum}}{\beta}_{i}{z}_{i}}+{\displaystyle \underset{i=1}{\overset{2}{\sum}}{\displaystyle \underset{\begin{array}{l}j=1\\ j\ge i\end{array}}{\overset{2}{\sum}}{\beta}_{ij}{z}_{i}}}{z}_{j}$ (20)

The regression coefficient for the cost function are given below.

Polynomial coefficients:

${\beta}_{0}=201.046$ ${\beta}_{1}=-0.965064$ ${\beta}_{2}=-0.43045$

${\beta}_{11}=0.00753563$ ${\beta}_{12}=0.00186971$ ${\beta}_{22}=0.00496752$

$\text{cost}={\beta}_{0}+{\beta}_{1}{z}_{1}+{\beta}_{2}{z}_{2}+{\beta}_{11}{z}_{1}^{2}+{\beta}_{12}{z}_{1}{z}_{2}+{\beta}_{22}{z}_{2}^{2}$

The MATLAB Genetic Algorithm “Optimtool” was used to estimate the minimum of the cost function. Table 3 shows the parameters that the GA uses for the computation [61].

The optimal cost (166.87) was found at ${z}_{1}^{*}=51.7375$ and ${z}_{2}^{*}=49.015$ where ${z}_{1}^{*}$ and ${z}_{2}^{*}$ are the optima of independent variables ${z}_{1}$ and ${z}_{2}$ . This policy constitutes an expansion of the hedging-point policy approach in a flexible manufacturing system where the cost function is minimized using GA. Controlled in this manner, the system is able to respond efficiently to customer demand with an infinite planning horizon.

5.4. Sensitivity Analysis

We now provide a sensitivity analysis to verify the usefulness and efficiency of the manufacturing control policy.

Table 4 reveals that the production threshold ${z}_{i}$ increases for product 1 as well as for product 2 while the cost of their shortages increases. However, while

Table 2. Levels on input factors.

Table 3. MATLAB parameters for the GA.

Table 4. Sensitivity analysis outcomes.

the inventory costs for both products increase, the production thresholds ${z}_{i}$ drop to adjust to the new condition.

The effect of the cost of shortages on production hedging point values is shown in Figure 9. Increases in this cost increase the threshold level while all other costs related to production (inventory, repair, and lockout) remain constant. The production policy therefore suggests keeping a large inventory to meet demand in case of machine breakdown.

Figure 10 shows the production threshold as a function of the inventory cost. Increases and decreases in inventory costs both have a direct and significant impact on the threshold. This shows the sensitivity of the inventory cost to the hedging point. Increases in this cost (while other production costs remain constant) accompany decreases in the threshold. At threshold production, the cost of maintaining inventory is closely related to the threshold. Smaller inventory costs less to maintain. However, the production manager must seek the optimal compromise, bearing in mind that a small inventory can easily lead to a shortage in the event of a sudden increase in demand, of machine breakdown or other randomness in production [63]. The outcomes of this analysis suggest increasing inventory to allow for CM and LOTO activities, in spite of the costs associated with accelerated production and the increased risk of HE and accidents.

Figure 9. Changes in the production threshold as a function of the cost of product shortage.

Figure 10. Change in the production threshold as a function of the inventory cost.

Figure 11 shows the production threshold as a function of the duration of CM with LOTO. When the duration of CM and LOTO is very high, machine breakdown cuts into the profit margin of the company.

6. Discussion

The aim of this study is to improve occupational health and safety in a flexible manufacturing system (FMS) by developing a policy of production control that withstands human error (HE) in the practice of lockout/tagout. Analytical formalism was combined with a numerical approach to solve Hamilton-Jacobi-Bellman equations for a system comprising a machine used to manufacture two different products. The outcomes suggest a hedging-point policy (a threshold policy based on stock levels) to maintain a level of inventory that meets demand during machine downtime. The machine is thus set on standby as soon as this threshold

Figure 11. Change in the production threshold as a function of the CM and LOTO duration.

level is reached and operated at its maximal output when the level decreases to (or drops below) the threshold for either product. Hedging-point policies have been described in detail elsewhere [30] [31]. To obtain a model more realistic for flexible systems, the SIMAN processor and simulation language was combined with DOE and a GA to perform a discrete simulation. This yielded a cost optimization policy. A sensitivity analysis was then performed to verify the usefulness and the efficiency of the control policy over an infinite production horizon. It was thus shown that it is possible to integrate CM with LOTO into production in an FMS while optimizing the costs of LOTO, shortage and inventory.

The control method constitutes an extension of the HPP structure. Its outcomes suggest ways for managers to facilitate the integration of occupational health and safety into FMS. However, the proposed model is based on certain simplifying assumptions. It could be improved to allow investigation of systems comprising two or more machines and meeting the demands of several customers. In addition, the likelihood of human error is a constraint on the effectiveness of the model and needs to be addressed properly. The impact of factors such as fatigue and stress also need to be investigated, as do the implications of the limited capacities of humans to perform repetitive tasks without error. Most machinery breakdowns are a result of operator performance and machinery age [64].

As noted in previous studies [24] [27] , LOTO should be included in system planning as a whole rather than treated as an add-on to maintenance. This would allow managers of operations and maintenance to plan their respective activities such that workers have no legitimate excuse for non-compliance with LOTO procedures. Although LOTO does increase the cost of CM, managers should be able to appreciate it as an investment in the reduction of risk (of accidents) rather than as non-value-adding activity. In any case, LOTO is now required by law according to article 188.2 of Quebec regulation respecting occupational health and safety, (RSST - D. 1187-2015, a. 3).

According to CNESST reports for the years 1999 to 2003, machinery was involved in more than 63,000 accidents and more than 100 deaths [65]. HE during maintenance activities must be considered to maintain both productivity and occupational health and safety in a flexible manufacturing system. It increases downtime and hence cost as well as the risk of accidents [26]. HE probabilities in maintenance procedures can be estimated for each possible failure case using the HEART approach [39]. Although a systematic approach focused on human factors does help reduce the likelihood of HE during maintenance activities, error is inevitable and there will be no substitute for awareness and quick adjustment. The fundamental human nature of workers cannot be modified no matter how their working conditions might be [41].

This study shows clearly the impact of HE on a flexible manufacturing system in terms of cost and occupational health and safety when production and CM are planned conjointly. HE during CM activities certainly compromises the safety of maintenance technicians while raising production costs. A seemingly minor flaw in a flexible manufacturing system can lead to economic losses and threaten lives.

7. Conclusions

The aim of this study was to develop a theoretical way of planning production and CM with LOTO procedures that increases safety in a flexible manufacturing system. A system comprising a single machine meeting two types of customer demand and prone to failure was examined. Numerical methods based on finite-difference methods were used to solve Hamilton-Jacobi-Bellman differential equations. An iterative approach was used to improve the manufacturing control policy. A hedging-point policy was thus obtained, and the SIMAN processor and simulation language (Arena software) were then used to carry out the discrete simulation. Finally, DOE and a GA were combined to obtain a more realistic model for industries and to optimize the cost structure.

Our principal contribution is a control strategy that facilitates the integration of LOTO into production and CM and takes into account the impact of HE. The costs of shortage, inventory build-up and CM can thus be optimized over an unbounded planning horizon. This study clearly illustrates the impact of HE on a flexible manufacturing system in terms of production cost and occupational health and safety. The implications for corrective maintenance planning are clear: HE in CM activities compromises the safety of maintenance technicians and increases production cost. A seemingly minor flaw in a flexible manufacturing man-machine system can lead to economic losses and threaten lives.

Inventory must be sufficient to satisfy customer demand during CM. Integrating preventive measures such as LOTO into a flexible manufacturing system can result in a remarkable improvement in occupational health and safety. LOTO helps prevent contact with hazardous spaces during maintenance operations, inadvertent release of hazardous energy and unforeseen start-up of machinery. All activities at risk of HE in a flexible manufacturing system must be identified. This risk can be reduced substantially by paying particular attention to critical maintenance tasks. Accident prevention is a crucial aspect of modern manufacturing systems.

The model is based on certain assumptions and could be improved through study of flexible manufacturing systems comprising a plurality of machines and customers with different demand types. Factors such as fatigue and stress need to be investigated to determine their influence on human performance in such systems.

Acknowledgements

We gratefully acknowledge Natural Sciences and Engineering Research Council of Canada and Association Québécoise pour l’Hygiène, la Santé et la Sécurité au travail. The authors also thank the anonymous reviewers, whose critique and suggestions increased the quality of this manuscript.

Conflicts of Interest

The authors declare no conflicts of interest regarding the publication of this paper.

Cite this paper

Diop, I., Nadeau, S. and Emami-Mehrgani, B. (2019) A Mathematical Model: A Flexible Manufacturing System, Prone to Error, Making Two Products Each with Stochastic Demand Schedules. American Journal of Industrial and Business Management, 9, 139-168. https://doi.org/10.4236/ajibm.2019.91011

References

- 1. Dhillon, B.-S. (2014) Basic Human Factors, Reliability, and Error Concepts. In: Human Reliability, Error, and Human Factors in Power Generation, Springer Series in Reliability Engineering, Springer International Publishing, 27-47. https://doi.org/10.1007/978-3-319-04019-6_3
- 2. Poyard, J.-L. (2016) De la conception au recyclage d’une machine, la sécurité avant tout. Revue Hygiène et sécurité du travail, INRS-DO 15, 5.
- 3. Matsuoka, S. and Muraki, M. (2001) Computer-Aided Planning for Lockout/Tagout Program. Process Safety Progress, 20, 130-135. https://doi.org/10.1002/prs.680200210
- 4. Burlet-Vienney, D., Chinniah, Y. and Pizarro-Chong, A. (2014) Design of an Intelligent Tool for the Observation and Follow-Up of Lockout Procedures during Maintenance Activities on Industrial Machines. Open Journal of Safety Science and Technology, 4, 106-118. https://doi.org/10.4236/ojsst.2014.42012
- 5. Bulzacchelli, M.-T., Vernick, J.-S., Sorock, G.-S., Webster, D.-W. and Lees, P. (2008) Circumstances of Fatal Lockout/Tagout-Related Injuries in Manufacturing. American Journal of Industrial Medicine, 51, 728-734. https://doi.org/10.1002/ajim.20630
- 6. Bulzacchelli, M.-T., Vernick, J.-S., Webster, D.-W. and Lees, P. (2007) Effects of the Occupational Safety and Health Administration’s Control of Hazardous Energy (Lockout/Tagout) Standard on Rates of Machinery-Related Fatal Occupational Injury. Injury Prevention, 13, 334-338. https://doi.org/10.1136/ip.2007.015677
- 7. Grusenmeyer, C. (2000) Interactions maintenance-exploitation et sécurité. Etude bibliographique. 1. Les taches de maintenance: définitions et caractéristiques contribuant à leur criticité. INRS_188.
- 8. Grusenmeyer, C. (2005) Les accidents liés à la maintenance. Etude bibliographique. INRS-NS_248: INRS.
- 9. Lind, S. (2008) Types and Sources of Fatal and Severe Non-Fatal Accidents in Industrial Maintenance. International Journal of Industrial Ergonomics, 38, 927-933. https://doi.org/10.1016/j.ergon.2008.03.002
- 10. Lind, S. (2009) Accident Sources in Industrial Maintenance Operations Proposals for Identification, Modelling and Management of Accident Risks. VTT Publications, No. 710, 1-123.
- 11. Väyrynen, S. (1984) Safety and Ergonomics in the Maintenance of Heavy Forest Machinery. Accident Analysis & Prevention, 16, 115-122. https://doi.org/10.1016/0001-4575(84)90036-8
- 12. Chinniah, Y. and Burlet-Vienney, D. (2013) Study on Lockout Procedures for the Safety of Workers Intervening on Equipment in the Municipal Sector in Quebe. International Journal of Occupational Safety and Ergonomics, 19, 495-511. https://doi.org/10.1080/10803548.2013.11077007
- 13. Emami-Mehrgani, B., Nadeau, S. and Kenné, J.-P. (2011) Lockout/Tagout and Operational Risks in the Production Control of Manufacturing Systems with Passive Redundancy. International Journal of Production Economics, 132, 165-173. https://doi.org/10.1016/j.ijpe.2011.03.028
- 14. Chinniah, Y., Aucourt, B. and Bourbonnière, R. (2015) Prévention des risques mécaniques et physiques: étude sur la sécurité des machines lors des interventions en mode de vitesse ou d’efforts réduits. RAPPORT R-888. Institut de recherche Robert-Sauvé en santé et en sécurité du travail (IRSST), Montreal.
- 15. Burlet-Vienney, D., Jocelyn, S., Chinniah, Y., Daigle, R. and Massé, S. (2009) Vérification du contenu d’un programme de cadenassage. Institut de recherche Robert-Sauvé en santé et en sécurité du travail (IRSST), Montréal, Québec, Canada, Guide RF-617.
- 16. Mario, S. (2012) Le suivi après l’implémentation d’un système de cadenassage: Partie 1: Quelles sont mes obligations? Message envoyé à Blog du Groupe ID.
- 17. Mario, S. (2012) L’après-implantation d’un système de cadenassage Partie 2 L’audit et le suivi de performance. Message envoyé à Blogue du Groupe ID.
- 18. Venkatraman, N. (2012) Business Performance and Maintenance: How Are Safety, Quality, Reliability, Productivity and Maintenance Related? Journal of Quality in Maintenance Engineering, 18, 183-195. https://doi.org/10.1108/13552511211244210
- 19. Bulzacchelli, M.-T., Vernick, J.-S., Webster, D.-W. and Lees, P. (2007) Effects of the Occupational Safety and Health Administration’s Control of Hazardous Energy (Lockout/Tagout) Standard on Rates of Machinery-Related Fatal Occupational Injury. Injury Prevention, 13, 334-338. https://doi.org/10.1136/ip.2007.015677
- 20. Assid, M., Gharbi, A. and Hajji, A. (2015) Production Planning and Opportunistic Preventive Maintenance for Unreliable One-Machine Two-Products Manufacturing Systems. IFAC-PapersOnLine, 48, 478-483. https://doi.org/10.1016/j.ifacol.2015.06.127
- 21. Ben-Salem, A., Gharbi, A. and Hajji, A. (2014) Production Planning and Emission Control for an Unreliable Manufacturing System with Subcontracting Strategy to Achieve Environmental Objectives. MScA. école de technologies supérieure, éTS, Montreal, 201 p.
- 22. Ouaret, S., Kenné, J.-P., Gharbi, A. and Polotski, V. (2015) Age-Dependent Production and Replacement Strategies in Failure-Prone Manufacturing Systems. Proceedings of the Institution of Mechanical Engineers, Part B: Journal of Engineering Manufacture, 231, 540-554.
- 23. Rivera-Gómez, H., Gharbi, A., Kenné, J.-P., Monta-o-Arango, O. and Hernandez-Gress, E. (2016) Production Control Problem Integrating Overhaul and Subcontracting Strategies for a Quality Deteriorating Manufacturing System. International Journal of Production Economics, 171, 134-150. https://doi.org/10.1016/j.ijpe.2015.10.008
- 24. Badiane, A., Nadeau, S., Kenné, J.-P. and Polotski, V. (2016) Optimizing Production While Reducing Machinery Lockout/Tagout Circumvention Possibilities. Journal of Quality in Maintenance Engineering, 22, 188-201. https://doi.org/10.1108/JQME-04-2014-0015
- 25. Charlot, E., Kenné, J.-P. and Nadeau, S. (2007) Optimal Production, Maintenance and Lockout/Tagout Control Policies in Manufacturing Systems. International Journal of Production Economics, 107, 435-450. https://doi.org/10.1016/j.ijpe.2006.09.017
- 26. Emami-Mehrgani, B., Nadeau, S. and Kenné, J.-P. (2014) Optimal Lockout/Tagout, Preventive Maintenance, Human Error and Production Policies of Manufacturing Systems with Passive Redundancy. Journal of Quality in Maintenance Engineering, 20, 453-470. https://doi.org/10.1108/JQME-10-2012-0035
- 27. Emami-Mehrgani, B., Neumann, W.-P., Nadeau, S. and Bazrafshan, M. (2015) Considering Human Error in Optimizing Production and Corrective and Preventive Maintenance Policies for Manufacturing Systems. Applied Mathematical Modelling, 40, 2056-2074.
- 28. Lee, C.-Y. and Lin, C.-S. (2001) Single-Machine Scheduling with Maintenance and Repair Rate-Modifying Activities. European Journal of Operational Research, 135, 493-513. https://doi.org/10.1016/S0377-2217(00)00322-2
- 29. Boukas, E.-K. and Haurie, A. (1990) Manufacturing Flow Control and Preventing Maintenance: A Stochastic Control Approach. IEEE Transactions on Automatic Control, 35, 1024-1031. https://doi.org/10.1109/9.58530
- 30. Hu, J.-Q., Vakili, P. and Yu, G.-X. (1994) Optimality of Hedging Point Policies in the Production Control of Failure Prone Manufacturing Systems. IEEE Transactions on Automatic Control, 39, 1875-1880. https://doi.org/10.1109/9.317116
- 31. Kenné, J.-P. and Boukas, E.-K. (1997) Maintenance and Production Control of Manufacturing Systems with Setups. Lectures in Applied Mathematics—American Mathematical Society, 33, 55-70.
- 32. Kenné, J.-P. and Boukas, E.-K. (2003) Hierarchical Control of Production and Maintenance Rates in Manufacturing Systems. Journal of Quality in Maintenance Engineering, 9, 66-82. https://doi.org/10.1108/13552510310466927
- 33. Kenné, J.-P. and Gharbi, A. (2004) Stohastic Optimal Production Control Problem with Corrective Maintenance. Computers & Industrial Engineering, 46, 865-875. https://doi.org/10.1016/j.cie.2004.05.024
- 34. Gharbi, A. and Kenné, J.-P. (2005) Maintenance Scheduling and Production Control of Multiple-Machine Manufacturing Systems. Computers & Industrial Engineering, 48, 693-707. https://doi.org/10.1016/j.cie.2004.12.007
- 35. Dhillon, B.-S. and Liu, Y. (2006) Human Error in Maintenance: A Review. Journal of Quality in Maintenance Engineering, 12, 21-36. https://doi.org/10.1108/13552510610654510
- 36. Kelley, S.-M. (2001) Lockout/Tagout: A Practical Approach. American Society of Safety Engineers, 208 p.
- 37. Collazo, G.-M. (2010) Reducing Human Error on the Manufacturing Floor. http://www.mastercontrol.com/newsletter/pharmaceutical/reducing_human_error_manufacturing_floor_0310.html
- 38. Whittingham, R.-B. (2004) The Blame Machine—Why Human Error Causes Accidents. Routledge, Elsevier Butterworth-Heinemann, Abingdon-on-Thames, 288 p.
- 39. Noroozi, A., Khan, F., MacKinnon, S., Amyotte, P. and Deacon, T. (2014) Determination of Human Error Probabilities in Maintenance Procedures of a Pump. Process Safety and Environmental Protection, 92, 131-141. https://doi.org/10.1016/j.psep.2012.11.003
- 40. Tchiehe, D.-N. and Gauthier, F. (2017) Classification of Risk Acceptability and Risk Tolerability Factors in Occupational Health and Safety. Safety Science, 92, 138-147. https://doi.org/10.1016/j.ssci.2016.10.003
- 41. Reason, J. (2000) Human Error: Models and Management. Western Journal of Medicine, 172, 393. https://doi.org/10.1136/ewjm.172.6.393
- 42. Reason, J. (2016) Managing the Risks of Organizational Accidents. Routledge, Abingdon-on-Thames. https://doi.org/10.4324/9781315543543
- 43. Harris, S. and Simpson, B. (2016) Human Error and the International Space Station: Challenges and Triumphs in Science Operations. 14th International Conference on Space Operations, Daejeon, 16-20 May 2016, 2406. https://doi.org/10.2514/6.2016-2406
- 44. Givi, Z.-S., Jaber, M.-Y. and Neumann, W.-P. (2015) Modelling Worker Reliability with Learning and Fatigue. Applied Mathematical Modelling, 39, 5186-5199. https://doi.org/10.1016/j.apm.2015.03.038
- 45. Kushner, H. and Dupuis, P. (1992) Numerical Methods for Stochastic Control Problems in Continuous Time. Vol. 24, Springer, Berlin.
- 46. Jaber, M.-Y., Givi, Z.-S. and Neumann, W.-P. (2013) Incorporating Human Fatigue and Recovery into the Learning-Forgetting Process. Applied Mathematical Modelling, 37, 7287-7299. https://doi.org/10.1016/j.apm.2013.02.028
- 47. Neumann, W.-P. and Dul, J. (2010) Human Factors: Spanning the Gap between OM and HRM. International Journal of Operations & Production Management, 30, 923-950. https://doi.org/10.1108/01443571011075056
- 48. Chan, F., Bhagwat, R. and Wadhwa, S. (2007) Taguchi’s Method Analysis of an FMS under Review-Period-Based Operational Controls: Identification of Control Periodicity. IEEE Transactions on Systems, Man and Cybernetics, Part A: Systems and Humans, 37, 212-225. https://doi.org/10.1109/TSMCA.2006.886355
- 49. Chan, F., Bhagwat, R. and Wadhwa, S. (2008) Comparative Performance Analysis of a Flexible Manufacturing System (FMS): A Review-Period-Based Control. International Journal of Production Research, 46, 1-24. https://doi.org/10.1080/00207540500521188
- 50. Raouf, A. and Ben-Daya, M. (1995) Flexible Manufacturing Systems: Recent Developments. Vol. 23, Elsevier, New York, 316 p.
- 51. Bengtsson, M. and Lundström, G. (2018) On the Importance of Combining “the New” with “the Old”—One Important Prerequisite for Maintenance in Industry 4.0. Procedia Manufacturing, 25, 118-125. https://doi.org/10.1016/j.promfg.2018.06.065
- 52. Salonen, A. (2018) The Need for a Holistic View on Dependable Production Systems. Procedia Manufacturing, 25, 17-22. https://doi.org/10.1016/j.promfg.2018.06.052
- 53. Dhillon, B.-S. (2014) Human Reliability, Error, and Human Factors in Power Generation. USA: Springer Series in Reliability Engineering, Springer, Cham, Heidelberg, New York, Dordrecht, London, 188 p. https://doi.org/10.1007/978-3-319-04019-6
- 54. Emami-Mehrgani, B., Kenné, J.-P. and Nadeau, S. (2012) Lockout/Tagout and Optimal Production Control Policies in Failure-Prone Non-Homogenous Transfer Lines with Passive Redundancy. International Journal of Production Research, 51, 1006-1023. https://doi.org/10.1080/00207543.2012.662305
- 55. Hlioui, R., Gharbi, A. and Hajji, A. (2015) Replenishment, Production and Quality Control Strategies in Three-Stage Supply Chain. International Journal of Production Economics, 166, 90-102. https://doi.org/10.1016/j.ijpe.2015.04.015
- 56. Azadivar, F. and Shu, J.-V. (1998) Use of Simulation in Optimization of Maintenance Policies. Proceedings of the 30th Conference on Winter Simulation, Washington DC, 13-16 December 1998, 1061-1068. https://doi.org/10.1109/WSC.1998.745854
- 57. Akella, R. and Kumar, P. (1986) Optimal Control of Production Rate in a Failure Prone Manufacturing System. IEEE Transactions on Automatic Control, 31, 116-126. https://doi.org/10.1109/TAC.1986.1104206
- 58. Assid, M., Gharbi, A. and Hajji, A. (2015) Joint Production, Setup and Preventive Maintenance Policies of Unreliable Two-Product Manufacturing Systems. International Journal of Production Research, 53, 4668-4683. https://doi.org/10.1080/00207543.2015.1030468
- 59. Hajji, A.G. and Kenne, J.P. (2004) Production and Set-Up Control of a Failure-Prone Manufacturing System. International Journal of Production Research, 42, 1107-1130. https://doi.org/10.1080/00207540310001631575
- 60. Azadivar, F. and Tompkins, G. (1999) Simulation Optimization with Qualitative Variables and Structural Model Changes: A Genetic Algorithm Approach. European Journal of Operational Research, 113, 169-182. https://doi.org/10.1016/S0377-2217(97)00430-X
- 61. Chipperfield, A.J. and Fleming, P.J. (1995) The Matlab Genetic Algorithm Toolbox. IEE Colloquium on Applied Control Techniques Using MATLAB, 26 January 1995, 10/1-10/4. https://doi.org/10.1049/ic:19950061
- 62. Robinson, T. (2008) Design and Analysis of Industrial Experiments. Statistical Practice in Business and Industry, 137-161. https://doi.org/10.1002/9780470997482.ch7
- 63. Nadeau, S., Kenné, J.P., Emami-Mehrgani, B. and Badri, A. (2016) Advances in Integration of Equipment Lockout/Tagout, Determination of Actual Production Capacity and Production/Maintenance Planning. Safety Science Monitor, 19, 1-8.
- 64. Njike, A.-N., Pellerin, R. and Kenné, J.-P. (2011) Maintenance/Production Planning with Interactive Feedback of Product Quality. Journal of Quality in Maintenance Engineering, 17, 281-298. https://doi.org/10.1108/13552511111157399
- 65. Chinniah, Y., Champoux, M. and Burlet-Vienney, D. (2008) Comparative Analysis of Lockout Programs and Procedures Applied to Industrial Machines. REPORT R-575. Institut de recherche Robert-Sauvé en santé et en sécurité du travail (IRSST), Montreal.

Appendices

Appendix 1―Optimality Conditions

The value functions and the HJB partial differential equations have properties called optimality conditions. Let $\chi \left(x\right)$ be a function where $x\in R$ , $\chi \left(x\right)$ has

a derivative ${{\chi}^{\prime}}_{\mathcal{l}}\left(x\right)=\underset{x\to 0}{\mathrm{lim}}\frac{\chi \left(x+{\delta}_{\mathcal{l}}\right)-\chi \left(x\right)}{\delta}$ along the direction $\mathcal{l}\in R$ . ${{\chi}^{\prime}}_{\mathcal{l}}(x)$

for every $\mathcal{l}$ and ${{\chi}^{\prime}}_{\mathcal{l}}\left(x\right)=\left(\nabla \chi \left(x\right),\mathcal{l}\right)$ if $\chi \left(x\right)$ is differentiable at x, where $\nabla \chi \left(x\right)$ is the gradient of $\chi \left(x\right)$ . Moreover, if $\chi \left(x\right)$ is convex and continuous on a domain $\psi $ , then $\chi \left(x\right)$ is differentiable and admit a derivative at any point and along direction $\mathcal{l}$ such that $x+\delta \mathcal{l}\in \psi $ for some $\delta \succ 0$ . See Emami-Mehrgani et al. [26] for more in depth details regarding the optimality conditions. Let $\left\{A\stackrel{\u02dc}{u}:\left({u}_{1}{u}_{2}\right)\in B(\cdot )\right\}$ be the admissible directions at x regarding the optimality condition. The HJB equations can be written in terms of the derivative as follows:

$\begin{array}{l}\rho \vartheta \left(x,\alpha \right)=\underset{u\in B\left(x,\alpha \right)}{\mathrm{min}}\left\{{{\vartheta}^{\prime}}_{{A}_{\stackrel{\u02dc}{u}}}\left(x,\alpha \right)+g\left(x,\alpha ,\cdot \right)+{\displaystyle \underset{\alpha \ne \beta}{\sum}{q}_{\alpha \beta}\left[\vartheta \left(x,\beta \right)-\vartheta \left(x,\alpha \right)\right]}\right\},\\ \forall \alpha ,\beta \in M\\ \text{where},\{\begin{array}{l}\rho \vartheta \left(x,\alpha \right)=\underset{u\in B\left(x,\alpha \right)}{\mathrm{min}}\left\{g\left(x,\alpha ,\cdot \right)+{\displaystyle \underset{\alpha \ne \beta}{\sum}{q}_{\alpha \beta}\left[\vartheta \left(x,\beta \right)-\vartheta \left(x,\alpha \right)\right]}\right\}\\ {{\vartheta}^{\prime}}_{{A}_{\stackrel{\u02dc}{u}}}\left(x,\alpha \right)=\left({\stackrel{\u02dc}{u}}_{1}-{\stackrel{\u02dc}{u}}_{4}\right){\vartheta}_{{x}_{1}}\left(x,\alpha \right)+\left({\stackrel{\u02dc}{u}}_{1}-{\stackrel{\u02dc}{u}}_{4}\right)\cdot {\vartheta}_{{x}_{2}}\left(x,\alpha \right)\end{array}\\ \text{with},{\stackrel{\u02dc}{u}}_{1}={u}_{1}{\iota}_{1}^{\alpha}+{u}_{4}{\iota}_{4}^{\alpha},\text{\hspace{0.17em}}\text{\hspace{0.17em}}{\iota}_{1}^{\alpha}=\{\begin{array}{l}1\\ 0\end{array}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}{\iota}_{4}^{\alpha}=\{\begin{array}{l}1\\ 0\end{array}\end{array}$ (21)

The system is considered feasible if:

$\begin{array}{l}\Sigma {\pi}_{i}\left(\alpha \right){U}_{i}^{\mathrm{max}}\left(\alpha \right)\ge d\left(\alpha \right)\\ \text{where,}\\ \{\begin{array}{l}{\pi}_{i}\text{is the probality limi at the state}\alpha \text{,}\\ {U}_{i}^{\mathrm{max}}\text{is the maximum production rate at state}\alpha \end{array}\\ \left\{\left[{\pi}_{1}{U}_{1}^{\mathrm{max}}Ind\left(\alpha =1\right)+{\pi}_{4}{U}_{4}^{\mathrm{max}}Ind\left(\alpha =4\right)\right]\ge \left[{d}_{1}Ind\left(\alpha =1\right)+{d}_{2}Ind\left(\alpha =4\right)\right]\right\}\end{array}$ (22)

The feasibility in Equation (22) is formulated to satisfy a set of specific properties called probability limits, defined in more detail in Equation (23):

$\begin{array}{l}{\pi}_{i}\left(\alpha \right)Q\left(\alpha \right)=0,\text{where}Q\text{is the matrix of transition rates}\\ \{\begin{array}{l}{\pi}_{i}\left(\alpha \right)={\pi}_{1}Ind\left(\alpha =1\right)+{\pi}_{2}Ind\left(\alpha =2\right)+{\pi}_{3}Ind\left(\alpha =3\right)\\ \text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}+{\pi}_{4}Ind\left(\alpha =4\right)+{\pi}_{5}Ind\left(\alpha =5\right)+{\pi}_{6}Ind\left(\alpha =6\right)\\ {\displaystyle \underset{i=1}{\overset{6}{\sum}}{\pi}_{i}\left(\alpha \right)=1}\end{array}\end{array}$ (23)

Let
$\delta \phi $ designate the horizon of
$\phi $ ,
$x\in \delta \phi $ if
${x}_{i}=0$ exists. Let
$\vartheta \left(x,\alpha \right)$ of
$\delta \phi $ be differentiable at an x_{0}. There is a vector
$\stackrel{\u02dc}{\nabla}\vartheta \left({x}_{0},\alpha \right)$ such that
${{\vartheta}^{\prime}}_{\mathcal{l}}\left({x}_{0},\alpha \right)=\langle \stackrel{\u02dc}{\nabla}\vartheta \left({x}_{0},\alpha \right),\mathcal{l}\rangle $ for any acceptable direction at x_{0}. The horizon condition from the continuity on the value function
$\vartheta \left(x,\alpha \right)$ can be presented as below:

$\begin{array}{l}\underset{u\in B\left({x}_{0},\alpha \right)}{\mathrm{min}}\left\{\left(\tilde{\nabla}\vartheta \left({x}_{0},\alpha \right),A\tilde{u}\right)+g\left({x}_{0},\alpha ,\cdot \right)+{\displaystyle \sum _{\alpha \ne \beta}{q}_{\alpha \beta}\left[\vartheta \left({x}_{0},\beta \right)-\vartheta \left({x}_{0},\alpha \right)\right]}\right\}\\ =\underset{u\in B\left(\alpha \right)}{\mathrm{min}}\left\{\left(\tilde{\nabla}\vartheta \left({x}_{0},\alpha \right),A\tilde{u}\right)+g\left({x}_{0},\alpha ,\cdot \right)+{\displaystyle \sum _{\alpha \ne \beta}{q}_{\alpha \beta}\left[\vartheta \left({x}_{0},\beta \right)-\vartheta \left({x}_{0},\alpha \right)\right]}\right\}\\ \forall \alpha ,\beta \in M\end{array}$ (24)

The flexible manufacturing system optimal policy designates a reduction of the value function over the set of admissible production capacities. It is practically impossible to solve the HJB equations analytically. Finite difference methods, based on a numerical approach to stochastic control problems [45] , are used to solve the optimality conditions of the HJB equations.

Appendix 2―Numerical Method

• Value function at mode 1

$\begin{array}{c}{\vartheta}^{h}\left(x,\alpha ,1\right)=\underset{\left(u,r\right)\in {B}^{h}\left(1\right)}{\mathrm{min}}\{{\left(\rho +{q}_{12}+{q}_{14}+\frac{\left|{u}_{1}-{d}_{1}\right|}{{h}_{{x}_{1}}}\right)}^{-1}\\ \text{\hspace{0.17em}}\text{\hspace{0.17em}}\times \langle \frac{\left|{u}_{1}-{d}_{1}\right|}{{h}_{{x}_{1}}}\left({\vartheta}^{h}\left(x+{h}_{{x}_{1}},\alpha ,1\right){K}^{+}+{\vartheta}^{h}\left(x-{h}_{{x}_{1}},\alpha ,1\right){K}^{-}\right)\\ \underset{\underset{}{}}{\overset{}{}}+\left(g\left(x,\alpha ,1\right)+{q}_{12}\cdot {\vartheta}^{h}\left(x,\alpha ,2\right)+{q}_{14}\cdot {\vartheta}^{h}\left(x,\alpha ,4\right)\right)\rangle \}\end{array}$ $\begin{array}{l}\text{where,}\\ g\left(x,\alpha ,1\right)={c}_{j}^{+}{x}_{j}^{+}+{c}_{j}^{-}{x}_{j}^{-}+{c}^{\alpha}\\ {x}_{j}^{+}=\mathrm{max}\left(0,{x}_{j}\right)\\ {x}_{j}^{-}=\mathrm{max}\left(-{x}_{j},0\right)\\ {c}^{\alpha}=0\left(\text{state1}\right)\end{array}$ (25)

• Value function at mode 2

$\begin{array}{c}{\vartheta}^{h}\left(x,\alpha ,2\right)=\underset{\left(u,r\right)\in {B}^{h}\left(2\right)}{\mathrm{min}}\{{\left(\rho +{r}_{21}+{q}_{25}+{q}_{23}+\frac{\left|0-{d}_{1}\right|}{{h}_{{x}_{1}}}\right)}^{-1}\\ \text{\hspace{0.17em}}\text{\hspace{0.17em}}\times \langle \frac{\left|0-{d}_{1}\right|}{{h}_{{x}_{1}}}\left({\vartheta}^{h}\left(x+{h}_{{x}_{1}},\alpha ,2\right){K}^{+}+{\vartheta}^{h}\left(x-{h}_{{x}_{1}},\alpha ,2\right){K}^{-}\right)\\ \underset{\underset{}{}}{\overset{}{}}+\left(g\left(x,\alpha ,2\right)+{r}_{21}\cdot {\vartheta}^{h}\left(x,\alpha ,1\right)+{q}_{25}\cdot {\vartheta}^{h}\left(x,\alpha ,5\right)+{q}_{23}\cdot {\vartheta}^{h}\left(x,\alpha ,3\right)\right)\rangle \}\end{array}$ $\begin{array}{l}\text{where,}\\ g\left(x,\alpha ,2\right)={c}_{j}^{+}{x}_{j}^{+}+{c}_{j}^{-}{x}_{j}^{-}+{c}^{\alpha}\\ {x}_{j}^{+}=\mathrm{max}\left(0,{x}_{j}\right)\\ {x}_{j}^{-}=\mathrm{max}\left(-{x}_{j},0\right)\\ {c}^{\alpha}=\left({c}_{r}+{c}_{LOTO}\right),\left(\text{CM-WOHE}\right)\end{array}$ (26)

• Value function at mode 3

$\begin{array}{c}{\vartheta}^{h}\left(x,\alpha ,3\right)=\underset{\left(u,r\right)\in {B}^{h}\left(3\right)}{\mathrm{min}}\{{\left(\rho +{r}_{31}+{q}_{36}+\frac{\left|{u}_{1}-{d}_{1}\right|}{{h}_{{x}_{1}}}\right)}^{-1}\\ \text{\hspace{0.17em}}\text{\hspace{0.17em}}\times \langle \frac{\left|0-{d}_{1}\right|}{{h}_{{x}_{1}}}\left({\vartheta}^{h}\left(x+{h}_{{x}_{1}},\alpha ,3\right){K}^{+}+{\vartheta}^{h}\left(x-{h}_{{x}_{1}},\alpha ,3\right){K}^{-}\right)\\ \underset{\underset{}{}}{\overset{}{}}+\left(g\left(x,\alpha ,3\right)+{r}_{31}.{\vartheta}^{h}\left(x,\alpha ,1\right)+{q}_{36}.{\vartheta}^{h}\left(x,\alpha ,6\right)\right)\rangle \}\end{array}$ $\begin{array}{l}\text{where,}\\ g\left(x,\alpha ,3\right)={c}_{j}^{+}{x}_{j}^{+}+{c}_{j}^{-}{x}_{j}^{-}+{c}^{\alpha}\\ {x}_{j}^{+}=\mathrm{max}\left(0,{x}_{j}\right)\\ {x}_{j}^{-}=\mathrm{max}\left(-{x}_{j},0\right)\\ {c}^{\alpha}=\gamma \left({c}_{r}+{c}_{LOTO}\right),\left(\text{CM-WHE}\right)\end{array}$ (27)

• Value function at mode 4

${\vartheta}^{h}\left(x,\alpha ,4\right)=\underset{\left(u,\text{}r\right)\in {B}^{h}\left(4\right)}{\mathrm{min}}\{{\left(\rho +{r}_{41}+{q}_{45}+\frac{\left|{u}_{2}-{d}_{2}\right|}{{h}_{{x}_{2}}}\right)}^{-1}$ $\begin{array}{l}\text{where,}\\ g\left(x,\alpha ,4\right)={c}_{j}^{+}{x}_{j}^{+}+{c}_{j}^{-}{x}_{j}^{-}+{c}^{\alpha}\end{array}$

$\begin{array}{c}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\times \langle \frac{\left|{u}_{2}-{d}_{2}\right|}{{h}_{{x}_{2}}}\left({\vartheta}^{h}\left(x+{h}_{{x}_{2}},\alpha ,4\right){K}^{+}+{\vartheta}^{h}\left(x-{h}_{{x}_{2}},\alpha ,4\right){K}^{-}\right)\\ \underset{\underset{}{}}{\overset{}{}}+\left(g\left(x,\alpha ,4\right)+{r}_{41}\cdot {\vartheta}^{h}\left(x,\alpha ,1\right)+{q}_{45}\cdot {\vartheta}^{h}\left(x,\alpha ,5\right)\right)\rangle \}\end{array}$ $\begin{array}{l}{x}_{j}^{+}=\mathrm{max}\left(0,{x}_{j}\right)\\ {x}_{j}^{-}=\mathrm{max}\left(-{x}_{j},0\right)\\ {c}^{\alpha}=0\end{array}$ (28)

• Value function at mode 5

$\begin{array}{c}{\vartheta}^{h}\left(x,\alpha ,5\right)=\underset{\left(u,r\right)\in {B}^{h}\left(5\right)}{\mathrm{min}}\{{\left(\rho +{r}_{54}+{q}_{52}+{q}_{56}+\frac{\left|0-{d}_{2}\right|}{{h}_{x}{}_{2}}\right)}^{-1}\\ \text{\hspace{0.17em}}\text{\hspace{0.17em}}\times \langle \frac{\left|0-{d}_{2}\right|}{{h}_{{x}_{2}}}\left({\vartheta}^{h}\left(x+{h}_{{x}_{2}},\alpha ,5\right){K}^{+}+{\vartheta}^{h}\left(x-{h}_{{x}_{2}},\alpha ,5\right){K}^{-}\right)\\ \underset{\underset{}{}}{\overset{}{}}+\left(g\left(x,\alpha ,5\right)+{r}_{54}\cdot {\vartheta}^{h}\left(x,\alpha ,4\right)+{q}_{52}\cdot {\vartheta}^{h}\left(x,\alpha ,2\right)+{q}_{56}\cdot {\vartheta}^{h}\left(x,\alpha ,6\right)\right)\rangle \}\end{array}$ $\begin{array}{l}\text{where,}\\ g\left(x,\alpha ,5\right)={c}_{j}^{+}{x}_{j}^{+}+{c}_{j}^{-}{x}_{j}^{-}+{c}^{\alpha}\\ {x}_{j}^{+}=\mathrm{max}\left(0,{x}_{j}\right)\\ {x}_{j}^{-}=\mathrm{max}\left(-{x}_{j},0\right)\\ {c}^{\alpha}=\left({c}_{r}+{c}_{LOTO}\right),\left(\text{CM-WOHE}\right)\end{array}$ (29)

• Value function at mode 6

$\begin{array}{c}{\vartheta}^{h}\left(x,\alpha ,6\right)=\underset{\left(u,r\right)\in {B}^{h}\left(6\right)}{\mathrm{min}}\{{\left(\rho +{q}_{63}+{r}_{64}+\frac{\left|0-{d}_{2}\right|}{{h}_{{x}_{2}}}\right)}^{-1}\\ \text{\hspace{0.17em}}\text{\hspace{0.17em}}\times \langle \frac{\left|0-{d}_{2}\right|}{{h}_{{x}_{2}}}\left({\vartheta}^{h}\left(x+{h}_{{x}_{2}},\alpha ,6\right){K}^{+}+{\vartheta}^{h}\left(x-{h}_{{x}_{2}},\alpha ,6\right){K}^{-}\right)\\ \underset{\underset{}{}}{\overset{}{}}+\left(g\left(x,\alpha ,6\right)+{q}_{63}\cdot {\vartheta}^{h}\left(x,\alpha ,3\right)+{r}_{64}\cdot {\vartheta}^{h}\left(x,\alpha ,4\right)\right)\rangle \}\end{array}$ $\begin{array}{l}\text{where,}\\ g\left(x,\alpha ,6\right)={c}_{j}^{+}{x}_{j}^{+}+{c}_{j}^{-}{x}_{j}^{-}+{c}^{\alpha}\\ {x}_{j}^{+}=\mathrm{max}\left(0,{x}_{j}\right)\\ {x}_{j}^{-}=\mathrm{max}\left(-{x}_{j},0\right)\\ {c}^{\alpha}=\gamma \left({c}_{r}+{c}_{LOTO}\right),\left(\text{CM-WHE}\right)\end{array}$ (30)