A Finite State Method in the Performance Evaluation of the Bernoulli Serial Production Lines

: Research on the performance measure evaluation of Bernoulli serial production lines is presented in this paper. Important aspects of the modeling and analysis using transition systems within the Markovian framework are addressed, including analytical and approximation methods. The “dimensionality curse” problems of the large scale and dense transition systems in the production system engineering ﬁeld are pointed out as one of the main research and development obstacles. In that respect, a new analytically-based ﬁnite state method is presented based on the proportionality property of the stationary probability distribution across the systems’ state space. Simple and di ﬀ erentiable expressions for the performance measures including the production rate, the work-in-process, and the probabilities of machine blockage and starvation are formulated. A ﬁnite state method’s accuracy and applicability are successfully validated by comparing the obtained results against the rigorous analytical solution.

steady-state probabilities for each state of the system's state space. This property is exploited further in the present paper in order to formulate a finite state method (FSM) that bypasses the system's dimensionality issues and approximates the exact results. The essence of the method reflects the internal relationships between the eigenvector components, allowing further systematic verification of the approximation methods as well as research on system improvability within the Markovian framework of the serial Bernoulli production lines. The method is applied in cases of several serial Bernoulli production lines, providing a possibility of extensive verification of the aggregation method.
The remainder of the paper is structured as follows. A brief literature review is presented in the next chapter. The third chapter outlines a referent analytical solution. It also contains a detailed derivation of the FSM method as the central topic of the present research. The fourth chapter comprises a summary of the FSM validation and discussions on the obtained results. Finally, the fifth chapter presents the main conclusions of the research.

Brief Literature Review
Application of the Markovian modeling in the case of serial production lines has been present in the governing PSE literature body for more than 50 years, starting with the works of Sevast'yanov in 1962. However, it got broader recognition after the development of the approximation methods bridging the exponential growth of the system dimensionality issues associated with the analytical solution formulation. These approximations are usually characterized as the decomposition and the aggregation method.
The decomposition method assumes that a complete serial line can be represented by a set of two machines-one buffer line summarizing the behavior of the upstream and downstream parts of the production flow. Each of the representing lines is composed of machines with the geometric distribution of up and downtimes. The applicability of the method has been demonstrated in cases of several production lines against the performance measures obtained using simulations. In addition, a necessity for more rigorous verification against analytical results is pointed out [15]. The algorithm of the decomposition method was further improved by the evaluation number reduction [16]. An excellent and more detailed review of the decomposition method development is presented in [17]. Further application of the method in cases of the continuous flow lines composed of machines with multiple failure rates is presented in [18,19]. The decomposition approach to the problem was applied in cases of systems with closed loops, two product types, quality failures, automated lines, multiple failure modes, non-linear material flow, assembly systems, etc. [3].
The basic idea of the aggregation method is to reduce a complete production line into a single machine using an algorithm of sequential backward-forward aggregation of two adjacent machines. It was first introduced by Lim et al. [20] in research on simple and analytical modeling of traditional production lines. The method proved to be acceptably accurate and valid for application in the modern mass-production industry. Its governing feature is modeling simplicity and low computational burden. However, the method should be further verified using the analytical results as it is presently validated using extensive simulations in some selected cases [2]. Its further development comprised research in cases like production lines with quality control systems [21], improvability and bottlenecks analysis [22,23], lead time analysis [24], transient problems [25], assembly systems [26], etc.
The analytical solution of the problem is of particular significance as it enables rigorous verification of the approximation and numerical techniques, particularly in the case of decomposition and aggregation methods. Its formulation has been a known issue starting with the work of Sevast'yanov. However, it was solved recently by Hadžić [14] using the concept of the generalized transition matrix and the eigenvalue problem following [27,28]. The solution proved to be exact. However, the evaluation algorithm remained limited due to dimensionality issues. The same issue has also been addressed in [29] using the state ranking transition matrix formulation of the system for arbitrary buffer occupancy and corresponding outcome states.
The "dimensionality curse" problems of the large scale and dense transition systems are present across different scientific disciplines. Consequently, researchers make significant efforts to develop various interpretation algorithms capable of cracking problems into simpler forms while keeping the original state-time framework [30]. Some of the possibilities include external memory storages (storing complete transition matrix), selective matrix criteria, lumping states according to the model properties, and sparse matrix approximation [31]. Additionally, the PSE research community investigated this issue intensively. The development of different approaches and algorithms includes methods for the specification of the system and numerical solution methods. In the first case, the governing research comprises Kronecker and tensor algebra, hierarchical modeling, complex analysis, and other approaches. The main topics in the latter case cover direct, iterative and sparse matrix methods, separable preconditioning, explicit storage algorithms, parallel computing techniques, etc. One of the main drawbacks of these methods is a requirement of compatibility between methods for system specification and the numerical solution of the problem [3].

Analytical Solution
Consider a serial production line composed of M machines m i , i = 1, 2, . . . , M of the Bernoulli(p i ) reliability model, where p i is the probability of machine being in state {up}, Figure 1. In addition, consider a series of buffers b i , i = 1, 2, . . . , M-1 of capacity N i ∈ N 0 , where N 0 denotes the set of natural numbers. Each buffer is placed in between two adjacent machines. Assume that (a) the machine status does not depend on the status of other machines, (b) the time axis is slotted and machines begin operations at the beginning of each slot, (c) machines' and buffers' statuses are determined at the beginning of each slot, (d) the machines have an identical cycle time, and (e) the first machine is never starved and the last one is never blocked [2].
Appl. Sci. 2020, 10, x FOR PEER REVIEW 4 of 22 original state-time framework [30]. Some of the possibilities include external memory storages (storing complete transition matrix), selective matrix criteria, lumping states according to the model properties, and sparse matrix approximation [31]. Additionally, the PSE research community investigated this issue intensively. The development of different approaches and algorithms includes methods for the specification of the system and numerical solution methods. In the first case, the governing research comprises Kronecker and tensor algebra, hierarchical modeling, complex analysis, and other approaches. The main topics in the latter case cover direct, iterative and sparse matrix methods, separable preconditioning, explicit storage algorithms, parallel computing techniques, etc. One of the main drawbacks of these methods is a requirement of compatibility between methods for system specification and the numerical solution of the problem [3].

Analytical Solution
Consider a serial production line composed of M machines mi, i = 1, 2, …, M of the Bernoulli(pi) reliability model, where pi is the probability of machine being in state {up}, Figure 1. In addition, consider a series of buffers bi, i = 1, 2, …, M-1 of capacity Ni ∈ 0 , where 0 denotes the set of natural numbers. Each buffer is placed in between two adjacent machines. Assume that (a) the machine status does not depend on the status of other machines, (b) the time axis is slotted and machines begin operations at the beginning of each slot, (c) machines' and buffers' statuses are determined at the beginning of each slot, (d) the machines have an identical cycle time, and (e) the first machine is never starved and the last one is never blocked [2]. The behavior of such a system can be modeled using stochastic processes. More specifically, it can be modeled as the discrete-time Markov chain {Xn, n = 0, 1, 2, …} over the system's state space S defined according to the buffers' capacities, Figure 1. More specifically, each state in the system's state space represents one possible combination of buffer occupancies. The transition of the system from one state to another is described using transition probabilities, [P(p1,p2,…pM)], also known as the transition matrix, which can be expressed as: The behavior of such a system can be modeled using stochastic processes. More specifically, it can be modeled as the discrete-time Markov chain {X n , n = 0, 1, 2, . . . } over the system's state space S defined according to the buffers' capacities, Figure 1. More specifically, each state in the system's state space represents one possible combination of buffer occupancies. The transition of the system from one state to another is described using transition probabilities, [P(p 1 ,p 2 , . . . p M )], also known as the transition matrix, which can be expressed as: A detailed formulation of the transition probabilities using the generalized transition matrix concept and constitutive matrices [P i (p i )] are available in [14]. The system in hand is a finite, irreducible, and aperiodic Markov chain. Hence, a solution to the Chapman-Kolmogorov equations is the eigenvector associated with the largest eigenvalue of the transition matrix. Such an eigenvector represents a limiting distribution that is identical to the stationary distribution of the Markov chain [32,33]. The stationary distribution of the Markov chain is of great significance for the performance evaluation. It can, therefore, be considered as the DNA of the production line, or more generally of the production system. It is also fundamental for the FSM derivation and validation of the existing approximation techniques.

A Finite State Method
Let us first take a look at the properties of a stationary distribution (eigenvector associated with the eigenvalue equal to 1) in the case of a three-machine line with arbitrary properties. Limitation to the case of a three-machine line of specific properties will be used here without loss of generality to illustrate the basic property of the stationary distribution. A similar conclusion can be drawn if a different arrangement of the line or a line with another number of machines is considered. Therefore, assume for simplicity that p 1 = 0.6, p 2 = 0.7, and p 3 = 0.8, and that N 1 = N 2 = 10, Figure 2a. The corresponding stationary distribution can be calculated without significant efforts using the analytical approach [14]. The obtained distribution of the stationary probabilities across the state space of the Bernoulli serial line is presented in Figure 2b including (N 1 + 1) × (N 2 + 1) = 121 possible system states. Proportional relationships between stationary probabilities of the considered system can easily be detected. Moreover, it turns out that such proportionality is a common property of Markov chains describing the behavior of the serial Bernoulli production lines. This property is used as the basis to develop FSM as an analytically-based method capable of breaking the dimensionality issues.
The eigenvector of the considered three-machine line (see Figure 2) can be decomposed into two classes of states. State 0 and states which are a multiple of N 1 + 1 = 11 stand for combinations of buffer occupancy N 1 = N 2 = 0, N 1 = 1 and N 2 = 0, N 1 = 2 and N 2 = 0, ..., N 1 = 10 and N 2 = 0. These combinations represent the external class of the states. States between 0 and 11, as well as states between multiples of N 1 + 1 = 11 cover combinations N 1 = 0 and N 2 = 1, N 1 = 0 and N 2 = 2, . . . , N 1 = 0 and N 2 = 10, N 1 = 1 and N 2 =1, . . . , N 1 = 10 and N 2 = 10. These combinations stand for the internal class of the states. Further analysis of the considered problem reveals that the proportional relationship between the states in the external class is governed by the stationary probability distribution determined analytically in the case of a two-machine line: p 1 -N 1 -p 2 . Similarly, the proportional relationship between the states in the internal class is governed by the stationary probability distribution determined analytically in the case of a two-machine line: p 1 -N 2 -p 3 . Therefore, the corresponding eigenvector can be determined using the well-known and simple analytical solution of the two-machine line problem [2] at considerably lower CPU and memory storage costs using two finite state elements: p 1 -N 1 -p 2 and p 1 -N 2 -p 3 . Notice that p 1 is a common property of both elements and also that p 1 =min{ p 1 , p 2 , p 3 }. Appl. Sci. 2020, 10 Let us now consider a general case of the serial Bernoulli production line (see Figure 3) with a machine m, pm = min{pi}, i = 1, 2, …, M. It is clear that there is a total number of m-1 upstream and Mm downstream machines. Define M-1 finite state elements, (e), such that a complete line is replaced by a set of pe-Ne-pm upstream and pm-Ne-pe+1 downstream elements. Let us now consider a general case of the serial Bernoulli production line (see Figure 3) with a machine m, pm = min{pi}, i = 1, 2, …, M. It is clear that there is a total number of m-1 upstream and Mm downstream machines. Define M-1 finite state elements, (e), such that a complete line is replaced by a set of pe-Ne-pm upstream and pm-Ne-pe+1 downstream elements. More specifically, the stationary distributions of the upstream elements are equal to: where i e = 0, 1, 2, . . . , N e stands for state on the buffer of the element e = 1, 2, . . . , m − 1. The pertaining probabilities can be determined using the analytical solution in the case of a line composed of two machines and one buffer [2]: Similarly, the stationary distributions of the downstream elements can be determined using: where e = m, m + 1, . . . , M − 1. The eigenvector representation of the stationary distribution is the same as in the upstream case, Equation (2). The application of Equations (2)-(4) results in a total number of M − 1 stationary distributions (elementary eigenvectors) that can be used to formulate the stationary distribution (eigenvector) of the complete serial line covering all of the considered system states without loss of generality. The application of the proportional relationships within the state space of the considered transition system yields a general component of the stationary distribution (i.e., the probability that buffers of the system are in state i 1 , i 2 , . . . , i M − 1 ), P i 1 i 2 i 3 ...i M−1 , expressed as: The unknown P 000...0 state of the system (all buffers are empty) from Equation (5) can be determined using the well-known condition . . .
, it follows that: Appl. Sci. 2020, 10, 6602 8 of 23 Equation (6) can be modified by extracting 1/P (e) 0 in the denominator yielding: Since N e i e P (e) i e = 1, the Equation (7) is simplified into the form:

Performance Measures
Performance measures are of particular interest for the PSE and Markov chain framework. It is, therefore, of great significance to develop simple expressions for the production rate (PR), work-in-process (WIP) for each buffer, as well as the probabilities of blockage (BL) and starvation (ST) for each machine [2]. The production rate, PR, is the average number of parts produced by the last machine of the production line. Therefore, it is an intersection of events m M is in state up and b M−1 is not empty , i.e.: Equation (9) can be expressed using the stationary distribution as: . . .
By applying Equation (5), the production rate can be formulated in terms of the finite state elements as: Substitution of Equation (8) into Equation (11) yields: which can be simplified using a relation Similarly, the work-in-process as the average number of parts contained at buffer b i is equal to: Again, using Equation (5), the work-in-process can be expressed as: . . .
Substitution of Equation (8) into Equation (15) after algebraic manipulation and by using identity Consequently, the work-in-process associated with the i th buffer depends solely on the elementary eigenvector P (e) obtained using the corresponding finite state element.
A machine blockage is commonly formulated using union and intersections of several events, i.e.: The probability of blockage, BL i , of the machine m i , Equation (17), can be expressed using the finite state elements as: which after inserting Equations (5) and (8), yields: Since the last machine can never be blocked, i.e., BL M = 0, the probability of the blockage of the penultimate machine equals: The probability of starvation, ST i , of the machine m i can be determined using the intersection of events: Using the finite state elements, Equation (20) can be formulated as: . . .
Again, using Equations (5) and (8), the probability of starvation of the i th machine is equal to: As the first machine cannot be starved, the probability of starvation of the machine m 1 takes a zero value.
Finally, it is important to point out that a complete reconstruction of the stationary distribution of the considered line is not necessary. It is sufficient to evaluate only the stationary distributions at the finite state element level and to collect the first and the last components to determine the production rate as well as the probabilities of blockage and starvation. Evaluation of the work-in-process, however, requires more detailed data on all components of the considered finite state element. The FSM approach is, therefore, very convenient for practical and daily factory-floor applications.

Validation and Application of the Developed Theory
Validation of the developed FSM is performed in cases of the serial Bernoulli lines L 1 -L 9 with probabilities of the state {up} specified in Table 1. Performance measures were calculated using Equations (12), (15), (18), (19) and (22) for each line using the analytical solution (AN), the aggregation procedure (AGG) and the FSM using ProLab, an in-house software developed by the authors, and PSEToolbox (Production System Engineering Toolbox) [2]. Evaluation of the performance measures was performed for the specified lines as well as for their permutations including even and uneven distribution of the buffer occupancy. However, only selected and the most interesting results are presented here. In addition, to enable a simple graphical presentation of the results, all of the buffer occupancies were selected as fractions or multiplies of the first or the last buffer's capacity. Nevertheless, the considered approach was valid in case of the arbitrary buffer occupancies. Buffer occupancy was specified separately for each considered case. As the first step of the FSM validation, a line L 1 was considered to check the fundamental accuracy of the method against the well known analytical solution of the line composed of two machines and one buffer. The performance measured in the case of the line L 1 and its reverse were evaluated using AN, AGG, and FSM approaches and are compared in Figures 4 and 5. Both AGG and FSM approaches agreed very well with the analytical results in this simplest case, except for some negligible discrepancies due to numerical reasons. An excellent agreement could also be noticed in the asymptotic values approached by the performance curves in all three cases. In addition, a reversibility property of the line was excerpted nicely.
The performance measures of the line with three machines were evaluated in several cases. The first case, L 2 , took into account the original arrangement of machines (see Table 1) and even buffer capacity distribution along the line. Figure 6 presents the distribution of the performance measures as a function of buffer capacity. A good agreement of both AGG and FSM with the analytical results was evident, except for slight discrepancies in PR and ST 3 between AN and FSM solution for the lowest level of buffers' capacity.
Uneven distribution of the buffer capacity and perturbed arrangement of the machines were evaluated in the case of the line labeled as L 2 . The third case, L 2 R , considered a reverse of the L 2 line. The performance measures obtained using AN, AGG, and FSM approaches are presented in Figures 7 and 8 as functions of the buffer capacity. It can be seen that generally, all three methods agreed well, particularly in asymptotic values of the performance measures. It can also be noticed that the reversibility property held for lines L 2 and L 2 R .
Performance measures of the line with four machines were evaluated in three cases of different machine arrangements including the uneven distribution of the buffer capacity. A comparison between the obtained results is presented in Figures 9-11 in cases of lines L 3A , L 3B , and L 3C . Again, the discrepancies were related to the smallest state spaces of buffers with small capacity. As the state space increased, the methods approached the same asymptotic values. It is interesting to notice a relationship between the position of the worst machine, WIP, BL, and ST. In a case when the worst machine was in position i = 2, 3, . . . , M, WIP i−1 altered from the asymptotic function to an almost linear curve, while BL i−1 and ST i altered pertaining asymptotes.
The above examples proved a quite good agreement between the aggregation procedure, the finite state method, and the analytical solution. Therefore, the analytical solution could be omitted from further evaluations of lines L 4 -L 9 as it would require considerably more CPU time as compared to the AGG and FSM algorithms. Further, to avoid the presentation of extensive data generated by the evaluation, only asymptotic values of the performance measures of lines L 4 -L 9 are presented in Figures 12 and 13. The machine arrangement is presented in Table 1, while the buffer occupancy is considered to be even along the lines. The obtained asymptotic value of the production rate was, as expected, the same for both lines and was equal to 0.6. The probability of machine blockage was also equal to 0 since the first machine was also the worst one. The asymptotic values of the work-in-process for each buffer of lines L 4 -L 9 are presented in Figure 12. It can be seen that both methods yielded almost the same values. Additionally, the asymptotic values of the probability of starvation for each machine of the considered lines are presented in Figure 13. A nice agreement between the results can be noticed.

Discussion
Three different methods were considered in the research: The analytical solution, the aggregation procedure, and the newly developed finite state method. Both approximation methods demonstrated respectable accuracy as compared to the analytical results. Therefore, both the AGG and the FSM approaches could be recommended for further application for the performance evaluation of the serial Bernoulli production line. This also includes the decomposition method since it was already verified successfully against the aggregation method [34].
Some drawbacks of the FSM as compared to the aggregation algorithm manifest in the range of the small buffer occupancy due to relatively rough system state discretization using the finite elements. This drawback is more pronounced as the probabilities p i approach the same value. However, this issue diminishes with the augmentation of the state space. The presented FSM was a powerful and analytically-based tool that enabled validation of the aggregation and decomposition procedure in complex cases that were previously unreachable to the conventional analytical approach. Its CPU requirements were comparable to other approximation algorithms, while the analytical approach remained an extremely CPU-demanding approach. Additional advantages of the FSM in comparison to other approximation methods are the differentiability and reconstruction of the production system eigenvector. The differentiability is of great significance in the field of improvement of production lines. This feature was enabled by the FSM since the performance measures were analytically related to the governing line properties. Additionally, the reconstruction of the system's stationary probability distribution vector (eigenvector) was of great value to the researcher in the PSE field as it enabled a deeper understanding of the behavior of complex production systems.

Discussion
Three different methods were considered in the research: The analytical solution, the aggregation procedure, and the newly developed finite state method. Both approximation methods demonstrated respectable accuracy as compared to the analytical results. Therefore, both the AGG and the FSM approaches could be recommended for further application for the performance evaluation of the serial Bernoulli production line. This also includes the decomposition method since it was already verified successfully against the aggregation method [34].
Some drawbacks of the FSM as compared to the aggregation algorithm manifest in the range of the small buffer occupancy due to relatively rough system state discretization using the finite elements. This drawback is more pronounced as the probabilities pi approach the same value. However, this issue diminishes with the augmentation of the state space. The presented FSM was a powerful and analytically-based tool that enabled validation of the aggregation and decomposition procedure in complex cases that were previously unreachable to the conventional analytical approach. Its CPU requirements were comparable to other approximation algorithms, while the analytical approach remained an extremely CPU-demanding approach. Additional advantages of the

Discussion
Three different methods were considered in the research: The analytical solution, the aggregation procedure, and the newly developed finite state method. Both approximation methods demonstrated respectable accuracy as compared to the analytical results. Therefore, both the AGG and the FSM approaches could be recommended for further application for the performance evaluation of the serial Bernoulli production line. This also includes the decomposition method since it was already verified successfully against the aggregation method [34].
Some drawbacks of the FSM as compared to the aggregation algorithm manifest in the range of the small buffer occupancy due to relatively rough system state discretization using the finite elements. This drawback is more pronounced as the probabilities pi approach the same value. However, this issue diminishes with the augmentation of the state space. The presented FSM was a powerful and analytically-based tool that enabled validation of the aggregation and decomposition procedure in complex cases that were previously unreachable to the conventional analytical approach. Its CPU requirements were comparable to other approximation algorithms, while the analytical approach remained an extremely CPU-demanding approach. Additional advantages of the

Conclusions
Manufacturing is of great importance for the global economy and society. It is, therefore, of great significance to master the analysis and design of various production systems. In that respect, research on the performance measures evaluation of the Bernoulli serial production lines was presented in this paper. Several important aspects of the modeling and analysis using transition systems within the Markovian framework were addressed, including analytical and approximation methods. The "dimensionality curse" problems of the large scale and dense transition systems in the PSE field were pointed out as one of the main research and development obstacles.
Given that, a new analytically-based FSM approach was developed based on the proportionality property of the stationary probability distribution across the systems' state space. An analytical solution of the two machine-one buffer line was exploited to formulate finite state elements used to model a complete Bernoulli serial production line. Simple and differentiable expressions for the performance measures including the production rate, the work-in-process, and the probabilities of machine blockage and starvation were developed. The FSM accuracy and applicability were successfully validated by comparing the obtained results against the rigorous analytical solution. In addition to that, a thorough validation of the aggregation method was provided, proving its accuracy and applicability. Currently, the FSM is limited to the evaluation of the single product lines as the cycle time was assumed to be equal at each machine along the line. Other limitations are related to the assumption of a Bernoulli reliability model of each machine along the line as well as to the assumption of occurrence of the machine breakdowns and repairs at the beginning of each cycle which may not be in complete agreement with the real production system.
Further research in the PSE field, as well as further development and application of the FSM, should focus on the analytical formulation of problems like improvability analysis, design of lean production lines, closed Bernoulli lines, Bernoulli lines with rework, assembly systems, the steady-state and transient behavior of transition systems, etc. Additionally, further effort should be put into FSM modeling of the production lines involving multiple products of different processing times. Such research will make it possible to model complex stochastic relationships in cases of systems like ship production or other job shop production systems. A significant impact on the PSE research body would also be accomplished in case of validation of the evaluation methods against the factory floor data. Last but not the least, it would be interesting to research a possibility to apply the FSM approach in cases of other large scale transition systems typically encountered in fields of physics, biology, chemistry, ecology, etc.

Conflicts of Interest:
The authors declare no conflict of interest.

AGG
AGGregation procedure AN Analytical solution CPU Central processing unit DNA DeoxyriboNucleic acid FSM Finite state method PSE Production system engineering