# A Max-Plus Algebra Approach to Study Time Disturbance Propagation within a Robustness Improvement Context

Authors: César Martínez-Olvera; Jaime Mora-Vargas
Journal: Mathematical Problems in Engineering
 (2018)

---

## Abstract

Industry 4.0 aims to ensure the future competitiveness of the manufacturing industry, where one of the major challenges faced by its implementation is the manufacturing/production system robustness (that is, able to perform in the presence of noise), as they may not be able to absorb input disruptions without bending or breaking. In this paper we propose to use the Max-Plus algebra approach to study the propagation of manufacturing disturbances (i.e., processing time variations), presenting a case study and performing a sensitivity analysis, with the idea of understanding under which conditions disturbance propagation takes place. Findings show that the impact propagation depends on where the variation source is located within the manufacturing system. Two are the main original contributions of this paper: the use of Max-Plus algebra to study the impact propagation of processing time variations and a four-step methodology to derive the equations representing the deterministic manufacturing system.


---

## Body

## 1. Introduction

As e-commerce has changed the distribution process, i.e., delivery distances are decreasing, the product variety is increasing, and more products are being sold in smaller quantities, the current supply chain needs to be rethought [1]. According to [2, 3], Industry 4.0 technologies will contribute to the more robust functioning of supply chains, by allowing a quicker and more reliable recognition of the potential external and internal disruptions and disturbances and the minimization or avoidance of their negative consequences. Different authors, i.e., [2–12], identify these value creation enabler technologies of Industry 4.0 as the Internet of Things, Big Data Analytics, and Cyber-Physical Systems. The real time Internet-networking of a growing number of physical objects with intelligent sensor/actuator technology is triggering the next stage of the industrial revolution known as Industry 4.0 [13]. Industry 4.0 aims to radically combine manufacturing, automation, and information and communications technology (ICT), in order to ensure the future competitiveness of the manufacturing industry [11, 14]. This requires the horizontal integration of the intelligent digital networking from different companies, with vertical integration of the autonomous, rules-based decision-making, and performance management within a single company, based on Big Data analysis [8], to create organizational transparency and seamless production [13]. More specifically, Industry 4.0 involves the networking and computerization of industrial production areas, in order to organize and optimize their operations, via the real time coordination of the production stages within value creation network, so they can boost customer satisfaction by producing customized products at serial production price [12, 15].


## 2. Literature Review

### 2.1. Industry 4.0 Reliability, Robustness, and Resilience

On the other hand, one of the major challenges faced by an Industry 4.0 implementation is the manufacturing/production system reliability, robustness, and resilience, that is, that seldom fails, able to perform in the presence of noise, that remains elastic under changing load, and that autonomously recovers from fault or failure situations [14]. Many manufacturing/production systems are not resilient to disruptions of their normal operations, as plenty of operations management techniques depend on assumption of stable and predictable environment, so their transforming processes are tuned in such way that they may not be robust enough to absorb input disruptions without bending or breaking [16–18]. References [19, 20] state that manufacturing/production systems need to be planned to be robust and resilient enough to (
1
) maintain their basic properties and ensure execution and (
2
) be able to adapt their behavior in the case of disturbances in order to achieve planned performance with the help of actions for recovery. Reference [21] summarizes this last as being resistant to different kinds of disturbances and able to react properly to them, while [22, 23] state that, in order to assure the effectiveness and efficiency necessary to keep a permanent high performance level, a manufacturing/production system must exhibit well behaved dynamics, that is, responsiveness and robustness. Now, robustness and resilience do not necessarily coexist, but both of them provide manufacturing/production systems with a sustainable competitive advantage in the global marketplace [24]. Even though robustness has been used interchangeably with resilience, they present some basic differences [25–27]:(i)
Robustness, defined as the ability of a system to perform its intended function and retaining its original level of performance while sustaining some damage (reasonable variations/disturbances) without failure. Robustness is a crucial feature of a reliable manufacturing/production system [21].
(ii)
Resilience, defined as the ability of a system to adapt to major variations/disturbances and gradually return to its original level of performance (or migrate to desired different one).
(iii)
Reliability, defined as the probability that a system will be able to perform its function without failing for a specific time period, under certain operating conditions.According to [28, 29], robustness is a static concept (the system can resist disruptions and retain its previous stable state), while resilience is a dynamic concept (the system can survive a disruption and return to a new different stable state). Finally, Table 1 shows the relationship between these robustness, resilience, and reliability, based on the work presented by [30–32].Table 1
Robustness, resilience, and reliability relationship.
Disturbance severity
Solution validity over time
No stress
Robustness
Reliability
Survival stress
Transient perturbation
Resilience
Sustainability
Permanent perturbation
New steady-state reference
Lethal stress
New system


### 2.2. Robustness, Resilience, and Production Disturbances Propagation

Robustness of a manufacturing/production system should aim for both a target value of the process outcome and a stable or consistent performance with minimum deviation or variation [33], while its resilience, to cope with disturbances and maintain their original state [28, 29, 34–36]. Now, robustness is related to the propagation by a given disturbance [37], so a realistic manufacturing/production system’s robustness study requires identifying each disturbance’s specific cause (unexpected/unwanted event) and effect (substantial deviation from a planned quality, quantity, time, cost value deviation), in terms of quantity of events and manufacturing downtimes they cause [21]. Internal/external industrial disturbances interfere not only with efficiency but also with the lead time of the manufacturing/production system [38], as their occurrence often leads to undesired perturbations degrading the system’s performance [21, 31, 32, 39–41], so understanding their dynamics is of paramount importance [22, 23]. According to [37], manufacturing disturbances can be classified as upstream, internal, and downstream disturbances, and more likely to propagate in highly networked manufacturing environments [38, 42]. Reference [43] refers to this propagation as ripple effect, that is, the failed nodes transfer the failure to their downstream neighbor nodes through connectivity links, with this process continuing through the entire supply network. This propagation of manufacturing disturbances affect the behavior of the production network [23], and their impact on the shop floor, the most sensitive part of a manufacturing system and where customer value is created [38], depends on the associated consequences derived from their propagation [31, 32, 39]. For this reason, authors like [44–47] have proposed models that try to predict the effects of manufacturing disturbance propagation (in this case, dimensional variation propagation), representing the disturbance propagation in terms of state-time equations (of some sort). If we take into account that (i) robustness and resilience are necessary to respond to disruptions [24], (ii) the topology of the manufacturing system affects its robustness and resilience [48], and (iii) in order to increase its robustness, a manufacturing system must be modelled and the disturbances’ negative effects studied [21], it becomes evident the need of a tool that facilitates the modelling of the propagation of a disturbance, from the manufacturing system structure point of view. For this reason, in this paper we propose to investigate the feasibility of using state-time equations to study, the propagation of manufacturing disturbances (i.e., processing lead time variation) downstream the production network. With this idea on mind, next section presents the literature review of relevant topics. Section three presents a case study and the evolution equations derived from the Max-Plus algebra approach are put to the test in Section four. Section five presents the conclusions derived from the research findings and establishes future research venues


### 2.3. Max-Plus Algebra and Modelling of Manufacturing Systems

Manufacturing/production systems can be modelled from a qualitative or quantitative point of view [49]: qualitative models (i.e., Automata, finite state machines, and Petri Nets) capture logical aspects of a system, while quantitative models (i.e., discrete-event simulation and Max-Plus algebra) highlight the quantitative system performance. Regarding these last, according to [50], discrete-event simulation models are dynamic in nature and are often easier to apply than analytic models, while Max-Plus algebra models are static in nature and express the event timing dynamics in terms of a set of algebraic linear equations analogous to conventional state-space linear equations [51]. However, even though discrete-event computer simulation is a widely used technique to analyze a manufacturing system [52], an important drawback is that it does not supply equations needed to analyze system’s behavior, making the understanding of how parameter changes affect important system properties such as stability, robustness and optimality of system performance very difficult [53, 54]. On the other hand, state-time dynamic models are an effective method to make quantitative analysis and evaluation of manufacturing systems, as it allows analyzing the influence of random disturbance factors on the system’s outputs [22, 23], so regarding this last issue, Max-Plus algebra can provide a much better insight than discrete-event simulation [52]. Max-plus algebra is an effective tool for modelling the event timing dynamics of a deterministic discrete-event system like a manufacturing line or a job shop manufacturing system [55–57]. According to [58], research related to Max-Plus algebra has been done in both the areas of developing the tool and the applications of the tool itself. A complete detailed description and analysis of the mathematics behind Max-Plus algebra can be found in [59, 60], a review of the basic concepts in [51–54, 61, 62]. A review of the Max-Plus algebra tool and its applications can be found in [58] and some of the work done on the modelling of manufacturing systems in [53, 54, 63, 64]. In summary, in Max-Plus algebra, the maximum operator (denoted by ⊕) replaces the addition operator, and the addition operator (denoted by ⊗) replaces the multiplication operator. The ⊕ called operation has a neutral element ϵ called zero; the ⊗ operation has a neutral element e, called unity. For all a, b ∈ real numbers R ∪ −∞ the Max-Plus operators are defined according to the following [49]:(1)
a
⊕
b
=
max
⁡
a
,
b
(2)
a
⊗
b
=
a
+
bIf we take into account the fact that state-space models are an effective method to quantitatively analyze the effect of manufacturing disturbances [23] and that Max-Plus algebra provides linear state-space-like equations to model manufacturing systems [51], then it becomes natural to follow a Max-Plus algebra approach to study the effects of manufacturing disturbances propagating downstream the production network. In fact, [58] mentions that an useful application of Max-Plus algebra equations is in not only evaluating the effect of variation in processing times on the total line idle time, but also possible performance measures as throughput, throughput time, and stations utilization. Next section presents a brief tutorial of the basic concepts of Max-Plus algebra and introduces the manufacturing system used as case study.


## 2.1. Industry 4.0 Reliability, Robustness, and Resilience

On the other hand, one of the major challenges faced by an Industry 4.0 implementation is the manufacturing/production system reliability, robustness, and resilience, that is, that seldom fails, able to perform in the presence of noise, that remains elastic under changing load, and that autonomously recovers from fault or failure situations [14]. Many manufacturing/production systems are not resilient to disruptions of their normal operations, as plenty of operations management techniques depend on assumption of stable and predictable environment, so their transforming processes are tuned in such way that they may not be robust enough to absorb input disruptions without bending or breaking [16–18]. References [19, 20] state that manufacturing/production systems need to be planned to be robust and resilient enough to (
1
) maintain their basic properties and ensure execution and (
2
) be able to adapt their behavior in the case of disturbances in order to achieve planned performance with the help of actions for recovery. Reference [21] summarizes this last as being resistant to different kinds of disturbances and able to react properly to them, while [22, 23] state that, in order to assure the effectiveness and efficiency necessary to keep a permanent high performance level, a manufacturing/production system must exhibit well behaved dynamics, that is, responsiveness and robustness. Now, robustness and resilience do not necessarily coexist, but both of them provide manufacturing/production systems with a sustainable competitive advantage in the global marketplace [24]. Even though robustness has been used interchangeably with resilience, they present some basic differences [25–27]:(i)
Robustness, defined as the ability of a system to perform its intended function and retaining its original level of performance while sustaining some damage (reasonable variations/disturbances) without failure. Robustness is a crucial feature of a reliable manufacturing/production system [21].
(ii)
Resilience, defined as the ability of a system to adapt to major variations/disturbances and gradually return to its original level of performance (or migrate to desired different one).
(iii)
Reliability, defined as the probability that a system will be able to perform its function without failing for a specific time period, under certain operating conditions.According to [28, 29], robustness is a static concept (the system can resist disruptions and retain its previous stable state), while resilience is a dynamic concept (the system can survive a disruption and return to a new different stable state). Finally, Table 1 shows the relationship between these robustness, resilience, and reliability, based on the work presented by [30–32].Table 1
Robustness, resilience, and reliability relationship.
Disturbance severity
Solution validity over time
No stress
Robustness
Reliability
Survival stress
Transient perturbation
Resilience
Sustainability
Permanent perturbation
New steady-state reference
Lethal stress
New system


## 2.2. Robustness, Resilience, and Production Disturbances Propagation

Robustness of a manufacturing/production system should aim for both a target value of the process outcome and a stable or consistent performance with minimum deviation or variation [33], while its resilience, to cope with disturbances and maintain their original state [28, 29, 34–36]. Now, robustness is related to the propagation by a given disturbance [37], so a realistic manufacturing/production system’s robustness study requires identifying each disturbance’s specific cause (unexpected/unwanted event) and effect (substantial deviation from a planned quality, quantity, time, cost value deviation), in terms of quantity of events and manufacturing downtimes they cause [21]. Internal/external industrial disturbances interfere not only with efficiency but also with the lead time of the manufacturing/production system [38], as their occurrence often leads to undesired perturbations degrading the system’s performance [21, 31, 32, 39–41], so understanding their dynamics is of paramount importance [22, 23]. According to [37], manufacturing disturbances can be classified as upstream, internal, and downstream disturbances, and more likely to propagate in highly networked manufacturing environments [38, 42]. Reference [43] refers to this propagation as ripple effect, that is, the failed nodes transfer the failure to their downstream neighbor nodes through connectivity links, with this process continuing through the entire supply network. This propagation of manufacturing disturbances affect the behavior of the production network [23], and their impact on the shop floor, the most sensitive part of a manufacturing system and where customer value is created [38], depends on the associated consequences derived from their propagation [31, 32, 39]. For this reason, authors like [44–47] have proposed models that try to predict the effects of manufacturing disturbance propagation (in this case, dimensional variation propagation), representing the disturbance propagation in terms of state-time equations (of some sort). If we take into account that (i) robustness and resilience are necessary to respond to disruptions [24], (ii) the topology of the manufacturing system affects its robustness and resilience [48], and (iii) in order to increase its robustness, a manufacturing system must be modelled and the disturbances’ negative effects studied [21], it becomes evident the need of a tool that facilitates the modelling of the propagation of a disturbance, from the manufacturing system structure point of view. For this reason, in this paper we propose to investigate the feasibility of using state-time equations to study, the propagation of manufacturing disturbances (i.e., processing lead time variation) downstream the production network. With this idea on mind, next section presents the literature review of relevant topics. Section three presents a case study and the evolution equations derived from the Max-Plus algebra approach are put to the test in Section four. Section five presents the conclusions derived from the research findings and establishes future research venues


## 2.3. Max-Plus Algebra and Modelling of Manufacturing Systems

Manufacturing/production systems can be modelled from a qualitative or quantitative point of view [49]: qualitative models (i.e., Automata, finite state machines, and Petri Nets) capture logical aspects of a system, while quantitative models (i.e., discrete-event simulation and Max-Plus algebra) highlight the quantitative system performance. Regarding these last, according to [50], discrete-event simulation models are dynamic in nature and are often easier to apply than analytic models, while Max-Plus algebra models are static in nature and express the event timing dynamics in terms of a set of algebraic linear equations analogous to conventional state-space linear equations [51]. However, even though discrete-event computer simulation is a widely used technique to analyze a manufacturing system [52], an important drawback is that it does not supply equations needed to analyze system’s behavior, making the understanding of how parameter changes affect important system properties such as stability, robustness and optimality of system performance very difficult [53, 54]. On the other hand, state-time dynamic models are an effective method to make quantitative analysis and evaluation of manufacturing systems, as it allows analyzing the influence of random disturbance factors on the system’s outputs [22, 23], so regarding this last issue, Max-Plus algebra can provide a much better insight than discrete-event simulation [52]. Max-plus algebra is an effective tool for modelling the event timing dynamics of a deterministic discrete-event system like a manufacturing line or a job shop manufacturing system [55–57]. According to [58], research related to Max-Plus algebra has been done in both the areas of developing the tool and the applications of the tool itself. A complete detailed description and analysis of the mathematics behind Max-Plus algebra can be found in [59, 60], a review of the basic concepts in [51–54, 61, 62]. A review of the Max-Plus algebra tool and its applications can be found in [58] and some of the work done on the modelling of manufacturing systems in [53, 54, 63, 64]. In summary, in Max-Plus algebra, the maximum operator (denoted by ⊕) replaces the addition operator, and the addition operator (denoted by ⊗) replaces the multiplication operator. The ⊕ called operation has a neutral element ϵ called zero; the ⊗ operation has a neutral element e, called unity. For all a, b ∈ real numbers R ∪ −∞ the Max-Plus operators are defined according to the following [49]:(1)
a
⊕
b
=
max
⁡
a
,
b
(2)
a
⊗
b
=
a
+
bIf we take into account the fact that state-space models are an effective method to quantitatively analyze the effect of manufacturing disturbances [23] and that Max-Plus algebra provides linear state-space-like equations to model manufacturing systems [51], then it becomes natural to follow a Max-Plus algebra approach to study the effects of manufacturing disturbances propagating downstream the production network. In fact, [58] mentions that an useful application of Max-Plus algebra equations is in not only evaluating the effect of variation in processing times on the total line idle time, but also possible performance measures as throughput, throughput time, and stations utilization. Next section presents a brief tutorial of the basic concepts of Max-Plus algebra and introduces the manufacturing system used as case study.


## 3. Case Study

In order to introduce the use of Max-Plus algebra approach for studying the propagation of disturbances in a manufacturing system environment, we will use the scenario presented in Figure1, based on the one mentioned by [59]. Following is the notation used in further analysis made throughout the rest of the document:(i)
M
i: variable denoting a certain machine type. In example, machine type 3 is denoted as M3.
(ii)
P
i: variable denoting a certain part type. In example, part type 1 is denoted as P1.
(iii)
P
i
M
i; variable denoting a part’s routing, that is P
a
r
t
i using M
a
c
h
i
n
e
i for a certain amount of time. In example, Part P2 routing is denoted as P2M1 (1 min), P2M2 (2 min), P2M3 (3 min).
(iv)
p.a.s.: variable denoting a part’s availability sequence, that is, the sequence of machines M
i a part P
i will make use. In example, thep.a.s. for P1 is M1 – M3 – M2 – M1 – M3.
(v)
m.a.s.: variable denoting a machine’s availability sequence, that is, the sequence of parts P
i a machine M
i will process. In example, them.a.s.for M3 is P2 – P1 – P2 – P1.
(vi)
X
i: variable denoting the position of each P
i
M
i has on its correspondingp.a.s.
(vii)
X
i (k): variable denoting the starting time ofX
i for the first unit of P
i to be processed.
(viii)
X
i (k+1): variable denoting the starting time ofX
i after a previous unit of P
i has been processed.
(ix)
U
i (k+1): variable denoting the initial arriving time of part P
i to the manufacturing system.
(x)
C
t
i: variable denoting the cycle time of part P
i, that is, the amount of time it takes P
i from its arrival until its departure from the manufacturing system.
(xi)
U
i: variable denoting the utilization of machine M
i, that is, the ratio between the actual time machine M
i being busy processing units of part P
i and the actual time machine M
i being available for processing units of part P
i.
(xii)
w.t.: variable denoting the waiting time a unit of partP
i, that is, the amount of time it takes P
i from its arrival to machine M
i until it starts being processed by machine M
i.
(xiii)
m.a.t.: variable denoting the machineM
i availability time, that is, the time at which certain machine M
i becomes available to process the next unit of part P
i of its machine’s availability sequence (m.a.s.).
(xiv)
p.a.t.: variable denoting the partP
i availability time, that is, the time at which a certain part P
i becomes available to be processed by the next machine M
i on its part’s availability sequence (p.a.s.).Figure 1
Part routing and machine sequencing.As stated before, Figure1 shows the sequence of machines each part follows, denoted as part routing P
i
M
i (that is, P
a
r
t
i using M
a
c
h
i
n
e
i), where columns represent P
i and rows represent M
i. The processing times per part type are as follows:(i)
Part P1: M2 (3 min), M3 (4 min).
(ii)
Part P2: M1 (1 min), M2 (2 min), M3 (3 min).
(iii)
Part P3: M1 (5 min), M2 (3 min).where transportation and set-up times are negligible, and all processing times are deterministic in minutes.Now, according to [56] there are four methods to derive Max-Plus algebra models for a manufacturing system:(i)
Timed event graph or a directed graph representation of the system [49, 55].
(ii)
Max-plus linear queuing networks [57].
(iii)
The system specifications.
(iv)
Block diagrams [53, 54, 56].Reference [63] claims that obtaining the evolution equations from a timed event graph can be cumbersome, and, for that reason, he presents the following recursive formulation:(3)
x
k
+
1
=
A
0
x
k
+
1
⊕
A
1
x
k
⊕
B
0
k
+
1where(4)
A
0
=
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
1
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
1
ɛ
3
ɛ
ɛ
ɛ
ɛ
ɛ
5
ɛ
2
ɛ
ɛ
ɛ
ɛ
ɛ
3
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
4
ɛ
6
ɛ
A
1
=
ɛ
5
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
3
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
3
4
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
3
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
B
0
=
e
ɛ
ɛ
ɛ
ɛ
e
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
e
ɛ
e
ɛ
e
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
e
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
ɛNow, in this paper we are more interested in the practical application of this modelling approach, so the industrial practitioner can truly exploit the benefits of it, independently the level of knowledge of the mathematics behind Max-Plus algebra. With this idea on mind, we present an easy-to-follow four-step algorithm to derive the evolution equations of the manufacturing system.


### 3.1. Proposed Algorithm to Derive the Evolution Equations

In order to exemplify the benefits of the proposed algorithm, we will use the following reentrant manufacturing system (next section will highlight the features of such type of system):(i)
Part P1: M1 (t1), M3 (t3), M2 (t2), M1 (t1), M3 (t3).
(ii)
Part P2: M3 (t3), M2 (t2), M1 (t1), M3 (t3), M1 (t1).Step 1.
For eachP
i, establish its M
i sequence (in a column format) and label each M
i
P
i (from left to right, top to bottom) in terms ofX
i. This will be known as the part availability sequence,p.a.s. (Figure 2).Figure 2
Part availability sequence (p.a.s.) of products P1, P2.Step 2.
Based on the previous step, for eachM
i, identify the P
i sequence (from left to right, top to bottom) in terms of its correspondingX
i. This will be known as the machine availability sequence,m.a.s. (Figure 3).Figure 3
Machine availability sequence (m.a.s.) of machines M1, M2, and M3.Step 3.
Based on Steps1 and 2, build them.a.s./p.a.s.sequencing matrix (Table 2). For example,(i)
P1p.a.s.hasX1 as the 1
s
t element in the sequence (corresponding to M1), followed byX3 as the 2
n
d (corresponding to M3),X5as the 3
r
d (corresponding to M2),X7 as the 4
t
h (corresponding to M1), andX9 as the 5
t
h (corresponding to M3)
(ii)
P1m.a.s.hasX1 as the 1
s
t element in the sequence, followed byX6as the 2
n
d,X7as the 3
r
d, andX10as the 4
t
h, all of them corresponding to M1Table 2
Machine availability/part availabilitysequencing matrix.
P
1
P
2
M
1
X1
X7
X6
X10
1
s
t
1
st
3
r
d
4
th
2
n
d
3
rd
4
t
h
5
th
t1
t7
t6
t10
M
2
X5
X4
2
n
d
3
rd
1
s
t
2
nd
t5
t4
M
3
X3
X9
X2
X8
2
n
d
2
nd
4
t
h
5
th
1
s
t
1
st
3
r
d
4
th
t3
t9
t2
t8Step 4.
Put together the evolution equations for eachP
i
M
i correspondingX
i (see (5)–(14)), where the first element of eachX
i evolution equation corresponds to am.a.s. relationship sequence, and the second element to ap.a.s. relationship sequence. It must be highlighted that (
1
) any n+1 => n relationship is denoted as (k), and (
2
) any n => n+1 is denoted as (k+1).(5)
X
1
k
+
1
=
t
10
⊗
X
10
k
⊕
t
9
⊗
X
9
k
(6)
X
2
k
+
1
=
t
9
⊗
X
9
k
⊕
t
10
⊗
X
10
k
(7)
X
3
k
+
1
=
t
2
⊗
X
2
k
+
1
⊕
t
1
⊗
X
1
k
+
1
(8)
X
4
k
+
1
=
t
5
⊗
X
5
k
⊕
t
2
⊗
X
2
k
+
1
(9)
X
5
k
+
1
=
t
4
⊗
X
4
k
+
1
⊕
t
3
⊗
X
3
k
+
1
(10)
X
6
k
+
1
=
t
1
⊗
X
1
k
+
1
⊕
t
4
⊗
X
4
k
+
1
(11)
X
7
k
+
1
=
t
6
⊗
X
6
k
+
1
⊕
t
5
⊗
X
5
k
+
1
(12)
X
8
k
+
1
=
t
3
⊗
X
3
k
+
1
⊕
t
6
⊗
X
6
k
+
1
(13)
X
9
k
+
1
=
t
8
⊗
X
8
k
+
1
⊕
t
7
⊗
X
7
k
+
1
(14)
X
10
k
+
1
=
t
7
⊗
X
7
k
+
1
⊕
t
8
⊗
X
8
k
+
1


### 3.2. Benefits of the Proposed Algorithm

According to [58], the main difficulty in modelling reentrant manufacturing systems using Max-Plus algebra is the main underlying assumption that the system (to be modelled) must be representable by a decision free, Timed Event Graph. This means, in Petri nets terms, that there is only one upstream and one downstream transition for each place. This is something a reentrant manufacturing system clearly does not comply with. Moreover, the same author presents the following reentrant manufacturing system:(i)
Part P1: M1 (operation A), M2 (operation B), M1 (operation C).And the corresponding recursive formulation is(15)
X
k
=
G
⊗
X
k
⊕
H
⊗
X
k
-
1where(16)
G
=
A
1
E
E
B
3
B
1
E
C
3
C
2
C
1
⊕
C
5
H
=
E
E
A
2
E
B
2
E
E
E
E
A
1
=
ɛ
ɛ
t
A
ɛ
A
2
=
t
C
ƞ
1
ɛ
t
C
B
1
=
ɛ
ɛ
t
B
ɛ
B
2
=
ɛ
t
B
ɛ
ɛ
B
3
=
t
A
ɛ
ɛ
t
A
C
1
=
ɛ
ɛ
t
C
ɛ
C
2
=
ɛ
t
A
ɛ
ɛ
C
3
=
t
B
ɛ
ɛ
t
B
C
4
=
ɛ
ɛ
ƞ
2
ɛ
C
5
=
ɛ
ɛ
t
C
ƞ
2
ƞ
1
ƞ
2As we can see, for a very simple reentrant manufacturing system (one product – one machine used twice), the Max-Plus algebra involved can be cumbersome to the common average industrial practitioner, where, on the other hand, the proposed four-step methodology allows the depiction of fairly more complex reentrant manufacturing systems (two products – two machines used twice). For this reason, we consider the proposed methodology as a value-adding original contribution of this paper.


### 3.3. Simulation Study

As we are interested in studying the impact of varying the processing time (at eachP
i
M
i,) has on the rest of the manufacturing system, we apply the recursive formulation (see (3)) or use the four-step algorithm proposed in the previous section to obtain evolution equations of the manufacturing system as follows:(17)
X
1
k
+
1
=
5
X
2
k
⊕
3
X
7
k
⊕
U
1
k
+
1
⊕
U
5
k
+
1
(18)
X
2
k
+
1
=
1
X
1
k
+
1
⊕
3
X
5
k
⊕
U
6
k
+
1
(19)
X
3
k
+
1
=
3
X
5
k
⊕
4
X
6
k
⊕
U
2
k
+
1
⊕
U
4
k
+
1
(20)
X
4
k
+
1
=
1
X
1
k
+
1
⊕
3
X
3
k
+
1
(21)
X
5
k
+
1
=
5
X
2
k
+
1
⊕
2
X
4
k
+
1
(22)
X
6
k
+
1
=
3
X
3
k
+
1
⊕
3
X
7
k
⊕
U
3
k
+
1
(23)
X
7
k
+
1
=
2
X
4
k
+
1
⊕
4
X
6
k
+
1Table 3 shows the starting time of each product P
i at each machine M
i, for a total processing of fifteen units per P
i (that is, fifteen units of P1, P2, and P3, respectively). Table 4 shows the seven different tested scenarios, consisting of varying the processing time of each P
i
M
i, one at the time, while maintaining the original processing times of rest of the P
i
M
i. For example, Scenario 1 consists of varying t1 from 1 to 5, while maintaining unchanged the original processing times of t2 through t7. Moreover, in order to test the validity of (17)–(23), a discrete-event simulation model of the case study was developed using the simulation software ARENA [65], validated, and used to collect sample data, following an approach similar to [66]. The running conditions of the seven different tested scenarios were the following: fifteen units per P
i were simulated; a new unit of P
i can enter the manufacturing system, only after the previous unit of P
i departed already; the running time was long enough to allow all the fifteen units per P
i to be processed and depart from the manufacturing system.Table 3
Starting time of each productP
i at each machine M
i.
1st unit
……..
15th unit
(K+1)
(K+1)
(K)
X1
X2
X3
X4
X5
X6
X7
(K)
X1
X2
X3
X4
X5
X6
X7
X2
ɛ
ɛ
ɛ
ɛ
ɛ
X2
130
131
131
136
133
X5
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
X5
133
133
136
138
136
140
X6
ɛ
ɛ
ɛ
ɛ
ɛ
X6
130
133
135
133
135
X7
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
X7
133
134
134
139
133
136
U1(K+1)
0
1
1
6
3
U1(K+1)
0
1
1
6
3
U2(K+1)
0
3
5
3
5
U2(K+1)
0
3
5
3
5
U3(K+1)
0
4
U3(K+1)
0
4
U4(K+1)
0
3
5
3
5
U4(K+1)
0
3
5
3
5
U5(K+1)
0
1
1
6
3
U5(K+1)
0
1
1
6
3
U6(K+1)
0
5
U6(K+1)
0
5
starting time
0
1
0
3
6
3
5
starting time
133
134
133
136
139
136
140Table 4
Processing time variations of the tested scenarios.
Processing time per P
i
M
i
Scenario
t1
t2
t3
t4
t5
t6
t7
1
1-5
5
3
2
3
4
3
2
1
1-5
3
2
3
4
3
3
1
5
1-5
2
3
4
3
4
1
5
3
1-5
3
4
3
5
1
5
3
2
1-5
4
3
6
1
5
3
2
3
1-5
3
7
1
5
3
2
3
4
1-5The impact of varying the processing time (at eachP
i
M
i) has on the rest of the manufacturing system was evaluated in terms of he following:(i)
Cycle time (C
t
i) of each P
i: Table 5 shows the calculation of the C
t
i (for the case of Scenario 1 with t1 = 1), based on Table 3. For example, Ct1 = starting time of P1M3 + processing time t6 – starting time of P1M2; Ct2 = starting time of P2M3 + processing time t7 – starting time of P2M1; Ct3 = starting time of P3M2 + processing time t5 – starting time of P3M1, and so on.
(ii)
Utilization (U
i) of each M
i: Table 5 shows the calculation of the U
i (for the case of Scenario 1 with t1 = 1), based on Table 3. For all three U
i, the utilization was calculated by adding the processing time of each P
i
M
i using the same machine, multiplied by the fifteen processed units, and divided by the highest departure time of the 15
t
h processed unit. For example, in Scenario 1 with t1 = 1, M1 consumed 6 minutes (t1 = 1 + t2 = 5), for a total of 90 minutes (15 ∗ 6), and was divided by 143 (starting time of X7 = 140 plus t7 = 3), for a utilization of 90/143 = 0.629.
(iii)
Waiting time (w.t.) before eachP
i
M
i: Table 6 shows the calculation of the w.t. (for the case of Scenario 1 with t1 = 1) before each P
i
M
i, according the following rules:
(a)
XI; IF t7 X7 (K+1) ≥ t2 X2(K+1) THEN w.t. = 0 ELSE w.t. = t2 X2(K+1) - t7 X7 (K+1).
(b)
X2; IF t5 X5 (K+1) ≥ t1 X1(K+2) THEN w.t. = 0 ELSE w.t. = t1 X1(K+2) - t5 X5 (K+1).
(c)
X3; IF t6 X6 (K+1) ≥ t5 X5(K+1) THEN w.t. = 0 ELSE w.t. = t5 X5(K+1) - t6 X6 (K+1).
(d)
X4; IF t1 X1 (K+1) ≥ t3 X3(K+1) THEN w.t. = 0 ELSE w.t. = t3 X3(K+1) - t1 X1 (K+1).
(e)
X5; IF t2 X2 (K+1) ≥ t4 X4(K+1) THEN w.t. = 0 ELSE w.t. = t4 X4(K+1) - t2 X2 (K+1).
(f)
X6; IF t3 X3 (K+2) ≥ t7 X7(K+1) THEN w.t. = 0 ELSE w.t. = t7 X7(K+1) - t3 X3 (K+2).
(g)
X7; IF t4 X4 (K+1) ≥ t6 X6(K+1) THEN w.t. = 0 ELSE w.t. = t6 X6(K+1) - t4 X4 (K+1).Table 5
Scenario 1 (with t1 = 1),C
t
i, and U
i calculations.
t1 = 1
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
X1
0
10
19
29
38
48
57
67
76
86
95
105
114
124
133
X2
1
11
20
30
39
49
58
68
77
87
96
106
115
125
134
X3
0
9
19
28
38
47
57
66
76
85
95
104
114
123
133
X4
3
12
22
31
41
50
60
69
79
88
98
107
117
126
136
X5
6
16
25
35
44
54
63
73
82
92
101
111
120
130
139
X6
3
12
22
31
41
50
60
69
79
88
98
107
117
126
136
X7
7
16
26
35
45
54
64
73
83
92
102
111
121
130
140
Ct
P1
7
7
7
7
7
7
7
7
7
7
7
7
7
7
7
7
7
P2
10
9
10
9
10
9
10
9
10
9
10
9
10
9
10
9.4
9.53
P3
8
8
8
8
8
8
8
8
8
8
8
8
8
8
8
8
8
U
M1
0.63
0.62
M2
0.84
0.83
M3
0.73
0.73Table 6
Waiting times before each PiMi, for the case of Scenario 1 (with t1 = 1).
P2M1 (XI), t1 = 1
0
10
19
29
38
48
57
67
76
86
95
105
114
124
133
w.t.
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
P3M1 (X2), t2 = 5
1
11
20
30
39
49
58
68
77
87
96
105
115
125
134
w.t.
2
1
2
1
2
1
2
1
2
1
2
1
2
1
1.500
1.466
P1M2 (X3), t3 = 3
0
9
19
28
38
47
57
66
76
85
95
104
114
123
133
w.t.
2
3
2
3
2
3
2
3
2
3
2
3
2
3
2.500
2.333
P2M2 (X4), t4 = 2
3
12
22
31
41
50
60
69
79
88
98
107
117
126
136
w.t.
2
1
2
1
2
1
2
1
2
1
2
1
2
1
1.500
1.400
P3M2 (X5), t5 = 3
6
16
25
35
44
54
63
73
82
92
101
111
120
130
139
w.t.
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0.000
0.000
P1M3 (X6), t6 = 4
3
12
22
31
41
50
60
69
79
88
98
107
117
126
136
w.t.
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0.000
0.000
P2M3 (X7), t7 = 3
7
16
26
35
45
54
64
73
83
92
102
111
121
130
140
w.t.
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2.000
1.866Regarding Table5, the row with the C
t
i calculation presents two extra values, i.e., for Ct2, the 9.4 and 9.53 values. The first value refers to the average Ct2 value obtained from Table 5 and the second value from the simulation run output. For both cases, the results are practically the same. The same applies to the U
i calculation. Regarding Table 6, the row with the w.t. calculation presents two extra values, i.e., for X2 and t2 = 5, the 1.5 and 1.466 values. The first value refers to the average w.t. value obtained from Table 6 and the second value from the simulation run output. For both cases, the results are practically the same.


### 3.4. Experimental Results

Figure4 shows the impact of varying the processing times of P1M1, from t1 = 1 to t1 = 5, have on the waiting times right in front of P3M1. Figures 5(a)–5(g) show a similar analysis for the tested scenarios presented in Table 4.Figure 4
Impact of varying the processing time of P1M1, on the waiting time of P3M1.Figure 5
Scenarios 1 – 7.
(a)
Scenario 1
(b)
Scenario 2
(c)
Scenario 3
(d)
Scenario 4
(e)
Scenario 5
(f)
Scenario 6
(g)
Scenario 7Three things must be noted regarding Figures5(a)–5(g):(i)
EachP
i
M
i is located according to the layout presented in Figure 1.
(ii)
TheC
t
i is presented at the end of each column, and the U at the end of each row.
(iii)
All the m.a.t. relationships are always presented on the horizontal direction, as each row represents a certain machineM
i.
(iv)
All the p.a.t. relationships are always on the vertical direction, as each column represents a certain partP
i.Figure6 shows an analysis of the impact of varying the processing time of P2M1 has on the rest of the manufacturing system (Scenario 1 of Table 4, Figure 5(a)), where such impact is typified in terms of degree link types and INC/CON effects:(i)
1
s
t, 2
n
d, 3
r
d, and 4
t
h degree link type refers to the degree of closeness between the affecting P
i
M
i (that is, the source of processing time variation) and the affected P
i
M
i. In this sense, a 1
s
t degree link type means the affecting and the affected P
i
M
is are adjacent, where, in a 2
n
d link type, there is a P
i
M
i in between the affecting and the affected P
i
M
is. This idea of using 1
s
t, 2
n
d, 3
r
d, and 4
t
h connectivity degree links is taken from the supply network risk propagation model proposed by [43].
(ii)
INC (increase) and CON (constant) effect refers to the effect the affectingP
i
M
i has on the m.a.t/p.a.t. of the affected P
i
M
i and the t
i value at which these last start increasing, or the value at which these last remain constant.Figure 6
Waiting time values behind Scenario 1.The analysis of Figure6 shows the following:(i)
X1 (P2M1): its m.a.t increases due to a 2
n
d degree link from P3M1; its p.a.t. remains constant (and then starts increasing) due to a 3
r
d degree link from P2M3. This results in the pattern A of Figure 7.
(ii)
X2 (P3M1): both its m.a.t and its p.a.t. increase due to a 1
s
t degree link from P2M1 and a 3
r
d degree link from P3M2. This results in pattern B of Figure 7.
(iii)
X3 (P1M2): its m.a.t. increases due to a 3
r
d degree link from P3M2; its p.a.t. remains always constant due to a 4
t
h degree link from P1M3. This results in the pattern C of Figure 7.
(iv)
X4 (P2M2): its m.a.t. remains constant due to a 4
t
h degree link from P1M2; its p.a.t. increases due to a 1
s
t degree link from P2M1. This results in the pattern D of Figure 7.
(v)
X5 (P3M2): its m.a.t. remains constant (and then starts increasing) due to a 2
n
d degree link from P2M2; its p.a.t. increases due to a 2
n
d degree link from P3M1. This results in the pattern E of Figure 7.
(vi)
X6 (P1M3): its m.a.t. remains constant (and then starts increasing) due to a 3
r
d degree link from P2M3; its p.a.t. increases due to a 4
t
h degree link from P1M2. This results in the pattern F of Figure 7.
(vii)
X7 (P2M3): its m.a.t. remains constant due to a 4
t
h degree link from P1M3: its p.a.t. remains constant due to a 2
n
d degree link from P2M2. This results in the pattern G of Figure 7.Figure 7
Scenario 1 in the context of the m.a.t./p.a.t. analysis.


### 3.5. Results Analysis and Derived Conclusions

The previous analysis seems to suggest that the further a certain affectedP
i
M
i is from the affecting P2M1, the less impact its varying processing time has on it. As the ultimate objective of this research effort is to study the conditions under which the manufacturing disturbance propagation takes place, we performed a similar analysis on the seven tested scenarios of Table 4 (Table 7 presents a summary of the results obtained from this analysis, while the Appendix presents the full details of the results obtained). For this matter, we defined the type of impact the affecting P
i
M
i has on the m.a.t. and p.a.t. of the affected P
i
M
i: primary (TOT), the m.a.t./p.a.t. of the affected P
i
M
i increases; secondary (PAR), the m.a.t./p.a.t. of the affected P
i
M
i remains constant and then increases; or null (NON), the m.a.t./p.a.t. of the affected P
i
M
i remains constant. Also, we classify each affected P
i
M
i in terms of its location: level I, for X1 (P2M1) and X2 (P3M1); level II, for X3 (P1M2), X4 (P2M2), and X5 (P3M2); level III for X6 (P1M3) and X7 (P2M3). An analysis of the results presented in Table 7 reveals the following facts:(1)
Independently at which level the affectingP
i
M
i is located (i.e., I, II, and III), it has a primary impact on the m.a.t. of its immediate affected P
i
M
is (1
s
t degree link), 85.71% of the times. From here, its frequency starts decreasing: 21.42% for a 2
n
d degree link, 11.76% for a 3
r
d degree link, and 9.09% for a 4
t
h degree link. Moreover, 50% of the times there is a primary impact on an affected P
i
M
is, is due to a 1
s
t degree link affecting P
i
M
i, and only 25% of the times from 2
n
d degree link, 16.66% from a 3
r
d degree link, and 8.33% from a 4
t
h degree link. A similar situation presents the p.a.t. case. These results lead us to conclude that the primary impact from the affecting P
i
M
i is mostly contained by the immediate, 1
s
t degree link affected P
i
M
is.
(2)
Independently at which level the affectingP
i
M
i is located (i.e., I, II, and III), it has a secondary impact on the m.a.t. of its immediate affected P
i
M
is (1
s
t degree link), 14.28% of the times. From here, its frequency increases (42.85% for a 2
n
d degree link) and then decreases again (23.52% for a 3
r
d degree link). Moreover, 14.28% of the times there is a secondary impact on an affected P
i
M
is, is due to a 1
s
t degree link affecting P
i
M
i, 42.85% of the times from 2
n
d degree link, and 23.52% from a 3
r
d degree link. A similar situation presents the p.a.t. case. These results lead us to think that there must be some residues of the primary impact (absorbed by the 1
s
t degree link affected P
i
M
is) that keeps propagating further. This issue is reflected on the fact that the secondary impact affecting the 2
n
d degree link affected P
i
M
is is almost three times higher than the frequency of the secondary impact affecting the 1
s
t degree link affected P
i
M
is (42.85% versus 14.28%).
(3)
Independently at which level the affectingP
i
M
i is located (i.e., I, II, and III), it has a null impact on the m.a.t. of its immediate affected P
i
M
is (1
s
t degree link), 0% of the times. From here, its frequency starts increasing: 35.71% for a 2
n
d degree link, 64.70% for a 3
r
d degree link, and 90.9% for a 4
t
h degree link. Moreover, 0% of the times there is a null impact on an affected P
i
M
is, due to a 1
s
t degree link affecting P
i
M
i, 19.23% of the times from 2
n
d degree link, 42.30% from a 3
r
d degree link affecting, and 38.46% for a 4
t
h degree link. A similar situation presents the p.a.t. case. As these results go on the opposite direction of those regarding the frequency of primary impact, and the fact that a null impact is the complement of a primary impact, we arrive to the same conclusion already stated: the primary impact from the affecting P
i
M
i is mostly contained by the immediate 1
s
t degree link affected P
i
M
is.Table 7
m.a.t./p.a.t. analysis, degree link and impact types context.
m.a.t. analysis
p.at. analysis
TOT
PAR
NON
TOT
PAR
NON
TOT
PAR
NON
TOT
PAR
NON
1st
85.71
14.28
100
1st
50
9.09
1st
85.71
14.28
100
1st
46.15
10
2nd
21.42
42.85
35.71
100
2nd
25
54.54
19.23
2nd
28.57
42.85
28.57
100
2nd
30.76
60
15.38
3rd
11.76
23.52
64.70
100
3rd
16.66
36.36
42.30
3rd
12.5
12.5
75
100
3rd
15.38
30
46.15
4th
9.09
90.90
100
4th
8.33
38.46
4th
8.33
8.33
83.33
100
4th
7.69
38.46
100
100
100
100
100
100Figure8 shows the same results of Table 7, where Series1 refer to the TOT impact type, Series2 to NON impact type, an Series3 to the PAR impact type; 1
s
t, 2
n
d, 3
r
d, and 4
t
h refer to the connectivity link type.Figure 8
m.a.t. and p.a.t. values behavior (frequency), according to Table7.


### 3.6. Implications for Robustness Improvement

According to [1, 19, 20, 37], the improvement of robustness in a production system is provided by the presence of inventory (to make up production short-falls) and capacity buffers (slack available at resources). On the other hand, [67] mentions that, in order to improve the robustness of a production/supply chain network, its topology needs to be modified. Author in [68] expands this idea when he states that a dynamic network reconfiguration bring robustness, as new machines can be plugged to provide redundancy. Reference [37] claims that the robustness of the system depends on the resource with the minimum ability to recover. Reference [40] proposes a dynamic network of the manufacturing system, where this critical resource is the one with the highest correlation between performance and system disturbance. For [43], there is the node that provokes two or more downstream nodes to be disrupted. It is the authors belief that, within a context of a dynamic reconfiguration of the network topology, typical of an Industry 4.0 environment, where the production system reconfiguration takes place through the easy integration of plug-and-produce, fully automated, digitized, highly cost-efficient, smart new manufacturing units [69], following an approach like the proposed in this paper would allow(i)
identifying under which conditions a resource will become critical
(ii)
establishing the appropriate stock/inventory and time/capacity buffers countermeasures to minimize it


## 3.1. Proposed Algorithm to Derive the Evolution Equations

In order to exemplify the benefits of the proposed algorithm, we will use the following reentrant manufacturing system (next section will highlight the features of such type of system):(i)
Part P1: M1 (t1), M3 (t3), M2 (t2), M1 (t1), M3 (t3).
(ii)
Part P2: M3 (t3), M2 (t2), M1 (t1), M3 (t3), M1 (t1).Step 1.
For eachP
i, establish its M
i sequence (in a column format) and label each M
i
P
i (from left to right, top to bottom) in terms ofX
i. This will be known as the part availability sequence,p.a.s. (Figure 2).Figure 2
Part availability sequence (p.a.s.) of products P1, P2.Step 2.
Based on the previous step, for eachM
i, identify the P
i sequence (from left to right, top to bottom) in terms of its correspondingX
i. This will be known as the machine availability sequence,m.a.s. (Figure 3).Figure 3
Machine availability sequence (m.a.s.) of machines M1, M2, and M3.Step 3.
Based on Steps1 and 2, build them.a.s./p.a.s.sequencing matrix (Table 2). For example,(i)
P1p.a.s.hasX1 as the 1
s
t element in the sequence (corresponding to M1), followed byX3 as the 2
n
d (corresponding to M3),X5as the 3
r
d (corresponding to M2),X7 as the 4
t
h (corresponding to M1), andX9 as the 5
t
h (corresponding to M3)
(ii)
P1m.a.s.hasX1 as the 1
s
t element in the sequence, followed byX6as the 2
n
d,X7as the 3
r
d, andX10as the 4
t
h, all of them corresponding to M1Table 2
Machine availability/part availabilitysequencing matrix.
P
1
P
2
M
1
X1
X7
X6
X10
1
s
t
1
st
3
r
d
4
th
2
n
d
3
rd
4
t
h
5
th
t1
t7
t6
t10
M
2
X5
X4
2
n
d
3
rd
1
s
t
2
nd
t5
t4
M
3
X3
X9
X2
X8
2
n
d
2
nd
4
t
h
5
th
1
s
t
1
st
3
r
d
4
th
t3
t9
t2
t8Step 4.
Put together the evolution equations for eachP
i
M
i correspondingX
i (see (5)–(14)), where the first element of eachX
i evolution equation corresponds to am.a.s. relationship sequence, and the second element to ap.a.s. relationship sequence. It must be highlighted that (
1
) any n+1 => n relationship is denoted as (k), and (
2
) any n => n+1 is denoted as (k+1).(5)
X
1
k
+
1
=
t
10
⊗
X
10
k
⊕
t
9
⊗
X
9
k
(6)
X
2
k
+
1
=
t
9
⊗
X
9
k
⊕
t
10
⊗
X
10
k
(7)
X
3
k
+
1
=
t
2
⊗
X
2
k
+
1
⊕
t
1
⊗
X
1
k
+
1
(8)
X
4
k
+
1
=
t
5
⊗
X
5
k
⊕
t
2
⊗
X
2
k
+
1
(9)
X
5
k
+
1
=
t
4
⊗
X
4
k
+
1
⊕
t
3
⊗
X
3
k
+
1
(10)
X
6
k
+
1
=
t
1
⊗
X
1
k
+
1
⊕
t
4
⊗
X
4
k
+
1
(11)
X
7
k
+
1
=
t
6
⊗
X
6
k
+
1
⊕
t
5
⊗
X
5
k
+
1
(12)
X
8
k
+
1
=
t
3
⊗
X
3
k
+
1
⊕
t
6
⊗
X
6
k
+
1
(13)
X
9
k
+
1
=
t
8
⊗
X
8
k
+
1
⊕
t
7
⊗
X
7
k
+
1
(14)
X
10
k
+
1
=
t
7
⊗
X
7
k
+
1
⊕
t
8
⊗
X
8
k
+
1


## 3.2. Benefits of the Proposed Algorithm

According to [58], the main difficulty in modelling reentrant manufacturing systems using Max-Plus algebra is the main underlying assumption that the system (to be modelled) must be representable by a decision free, Timed Event Graph. This means, in Petri nets terms, that there is only one upstream and one downstream transition for each place. This is something a reentrant manufacturing system clearly does not comply with. Moreover, the same author presents the following reentrant manufacturing system:(i)
Part P1: M1 (operation A), M2 (operation B), M1 (operation C).And the corresponding recursive formulation is(15)
X
k
=
G
⊗
X
k
⊕
H
⊗
X
k
-
1where(16)
G
=
A
1
E
E
B
3
B
1
E
C
3
C
2
C
1
⊕
C
5
H
=
E
E
A
2
E
B
2
E
E
E
E
A
1
=
ɛ
ɛ
t
A
ɛ
A
2
=
t
C
ƞ
1
ɛ
t
C
B
1
=
ɛ
ɛ
t
B
ɛ
B
2
=
ɛ
t
B
ɛ
ɛ
B
3
=
t
A
ɛ
ɛ
t
A
C
1
=
ɛ
ɛ
t
C
ɛ
C
2
=
ɛ
t
A
ɛ
ɛ
C
3
=
t
B
ɛ
ɛ
t
B
C
4
=
ɛ
ɛ
ƞ
2
ɛ
C
5
=
ɛ
ɛ
t
C
ƞ
2
ƞ
1
ƞ
2As we can see, for a very simple reentrant manufacturing system (one product – one machine used twice), the Max-Plus algebra involved can be cumbersome to the common average industrial practitioner, where, on the other hand, the proposed four-step methodology allows the depiction of fairly more complex reentrant manufacturing systems (two products – two machines used twice). For this reason, we consider the proposed methodology as a value-adding original contribution of this paper.


## 3.3. Simulation Study

As we are interested in studying the impact of varying the processing time (at eachP
i
M
i,) has on the rest of the manufacturing system, we apply the recursive formulation (see (3)) or use the four-step algorithm proposed in the previous section to obtain evolution equations of the manufacturing system as follows:(17)
X
1
k
+
1
=
5
X
2
k
⊕
3
X
7
k
⊕
U
1
k
+
1
⊕
U
5
k
+
1
(18)
X
2
k
+
1
=
1
X
1
k
+
1
⊕
3
X
5
k
⊕
U
6
k
+
1
(19)
X
3
k
+
1
=
3
X
5
k
⊕
4
X
6
k
⊕
U
2
k
+
1
⊕
U
4
k
+
1
(20)
X
4
k
+
1
=
1
X
1
k
+
1
⊕
3
X
3
k
+
1
(21)
X
5
k
+
1
=
5
X
2
k
+
1
⊕
2
X
4
k
+
1
(22)
X
6
k
+
1
=
3
X
3
k
+
1
⊕
3
X
7
k
⊕
U
3
k
+
1
(23)
X
7
k
+
1
=
2
X
4
k
+
1
⊕
4
X
6
k
+
1Table 3 shows the starting time of each product P
i at each machine M
i, for a total processing of fifteen units per P
i (that is, fifteen units of P1, P2, and P3, respectively). Table 4 shows the seven different tested scenarios, consisting of varying the processing time of each P
i
M
i, one at the time, while maintaining the original processing times of rest of the P
i
M
i. For example, Scenario 1 consists of varying t1 from 1 to 5, while maintaining unchanged the original processing times of t2 through t7. Moreover, in order to test the validity of (17)–(23), a discrete-event simulation model of the case study was developed using the simulation software ARENA [65], validated, and used to collect sample data, following an approach similar to [66]. The running conditions of the seven different tested scenarios were the following: fifteen units per P
i were simulated; a new unit of P
i can enter the manufacturing system, only after the previous unit of P
i departed already; the running time was long enough to allow all the fifteen units per P
i to be processed and depart from the manufacturing system.Table 3
Starting time of each productP
i at each machine M
i.
1st unit
……..
15th unit
(K+1)
(K+1)
(K)
X1
X2
X3
X4
X5
X6
X7
(K)
X1
X2
X3
X4
X5
X6
X7
X2
ɛ
ɛ
ɛ
ɛ
ɛ
X2
130
131
131
136
133
X5
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
X5
133
133
136
138
136
140
X6
ɛ
ɛ
ɛ
ɛ
ɛ
X6
130
133
135
133
135
X7
ɛ
ɛ
ɛ
ɛ
ɛ
ɛ
X7
133
134
134
139
133
136
U1(K+1)
0
1
1
6
3
U1(K+1)
0
1
1
6
3
U2(K+1)
0
3
5
3
5
U2(K+1)
0
3
5
3
5
U3(K+1)
0
4
U3(K+1)
0
4
U4(K+1)
0
3
5
3
5
U4(K+1)
0
3
5
3
5
U5(K+1)
0
1
1
6
3
U5(K+1)
0
1
1
6
3
U6(K+1)
0
5
U6(K+1)
0
5
starting time
0
1
0
3
6
3
5
starting time
133
134
133
136
139
136
140Table 4
Processing time variations of the tested scenarios.
Processing time per P
i
M
i
Scenario
t1
t2
t3
t4
t5
t6
t7
1
1-5
5
3
2
3
4
3
2
1
1-5
3
2
3
4
3
3
1
5
1-5
2
3
4
3
4
1
5
3
1-5
3
4
3
5
1
5
3
2
1-5
4
3
6
1
5
3
2
3
1-5
3
7
1
5
3
2
3
4
1-5The impact of varying the processing time (at eachP
i
M
i) has on the rest of the manufacturing system was evaluated in terms of he following:(i)
Cycle time (C
t
i) of each P
i: Table 5 shows the calculation of the C
t
i (for the case of Scenario 1 with t1 = 1), based on Table 3. For example, Ct1 = starting time of P1M3 + processing time t6 – starting time of P1M2; Ct2 = starting time of P2M3 + processing time t7 – starting time of P2M1; Ct3 = starting time of P3M2 + processing time t5 – starting time of P3M1, and so on.
(ii)
Utilization (U
i) of each M
i: Table 5 shows the calculation of the U
i (for the case of Scenario 1 with t1 = 1), based on Table 3. For all three U
i, the utilization was calculated by adding the processing time of each P
i
M
i using the same machine, multiplied by the fifteen processed units, and divided by the highest departure time of the 15
t
h processed unit. For example, in Scenario 1 with t1 = 1, M1 consumed 6 minutes (t1 = 1 + t2 = 5), for a total of 90 minutes (15 ∗ 6), and was divided by 143 (starting time of X7 = 140 plus t7 = 3), for a utilization of 90/143 = 0.629.
(iii)
Waiting time (w.t.) before eachP
i
M
i: Table 6 shows the calculation of the w.t. (for the case of Scenario 1 with t1 = 1) before each P
i
M
i, according the following rules:
(a)
XI; IF t7 X7 (K+1) ≥ t2 X2(K+1) THEN w.t. = 0 ELSE w.t. = t2 X2(K+1) - t7 X7 (K+1).
(b)
X2; IF t5 X5 (K+1) ≥ t1 X1(K+2) THEN w.t. = 0 ELSE w.t. = t1 X1(K+2) - t5 X5 (K+1).
(c)
X3; IF t6 X6 (K+1) ≥ t5 X5(K+1) THEN w.t. = 0 ELSE w.t. = t5 X5(K+1) - t6 X6 (K+1).
(d)
X4; IF t1 X1 (K+1) ≥ t3 X3(K+1) THEN w.t. = 0 ELSE w.t. = t3 X3(K+1) - t1 X1 (K+1).
(e)
X5; IF t2 X2 (K+1) ≥ t4 X4(K+1) THEN w.t. = 0 ELSE w.t. = t4 X4(K+1) - t2 X2 (K+1).
(f)
X6; IF t3 X3 (K+2) ≥ t7 X7(K+1) THEN w.t. = 0 ELSE w.t. = t7 X7(K+1) - t3 X3 (K+2).
(g)
X7; IF t4 X4 (K+1) ≥ t6 X6(K+1) THEN w.t. = 0 ELSE w.t. = t6 X6(K+1) - t4 X4 (K+1).Table 5
Scenario 1 (with t1 = 1),C
t
i, and U
i calculations.
t1 = 1
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
X1
0
10
19
29
38
48
57
67
76
86
95
105
114
124
133
X2
1
11
20
30
39
49
58
68
77
87
96
106
115
125
134
X3
0
9
19
28
38
47
57
66
76
85
95
104
114
123
133
X4
3
12
22
31
41
50
60
69
79
88
98
107
117
126
136
X5
6
16
25
35
44
54
63
73
82
92
101
111
120
130
139
X6
3
12
22
31
41
50
60
69
79
88
98
107
117
126
136
X7
7
16
26
35
45
54
64
73
83
92
102
111
121
130
140
Ct
P1
7
7
7
7
7
7
7
7
7
7
7
7
7
7
7
7
7
P2
10
9
10
9
10
9
10
9
10
9
10
9
10
9
10
9.4
9.53
P3
8
8
8
8
8
8
8
8
8
8
8
8
8
8
8
8
8
U
M1
0.63
0.62
M2
0.84
0.83
M3
0.73
0.73Table 6
Waiting times before each PiMi, for the case of Scenario 1 (with t1 = 1).
P2M1 (XI), t1 = 1
0
10
19
29
38
48
57
67
76
86
95
105
114
124
133
w.t.
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
P3M1 (X2), t2 = 5
1
11
20
30
39
49
58
68
77
87
96
105
115
125
134
w.t.
2
1
2
1
2
1
2
1
2
1
2
1
2
1
1.500
1.466
P1M2 (X3), t3 = 3
0
9
19
28
38
47
57
66
76
85
95
104
114
123
133
w.t.
2
3
2
3
2
3
2
3
2
3
2
3
2
3
2.500
2.333
P2M2 (X4), t4 = 2
3
12
22
31
41
50
60
69
79
88
98
107
117
126
136
w.t.
2
1
2
1
2
1
2
1
2
1
2
1
2
1
1.500
1.400
P3M2 (X5), t5 = 3
6
16
25
35
44
54
63
73
82
92
101
111
120
130
139
w.t.
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0.000
0.000
P1M3 (X6), t6 = 4
3
12
22
31
41
50
60
69
79
88
98
107
117
126
136
w.t.
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0.000
0.000
P2M3 (X7), t7 = 3
7
16
26
35
45
54
64
73
83
92
102
111
121
130
140
w.t.
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2.000
1.866Regarding Table5, the row with the C
t
i calculation presents two extra values, i.e., for Ct2, the 9.4 and 9.53 values. The first value refers to the average Ct2 value obtained from Table 5 and the second value from the simulation run output. For both cases, the results are practically the same. The same applies to the U
i calculation. Regarding Table 6, the row with the w.t. calculation presents two extra values, i.e., for X2 and t2 = 5, the 1.5 and 1.466 values. The first value refers to the average w.t. value obtained from Table 6 and the second value from the simulation run output. For both cases, the results are practically the same.


## 3.4. Experimental Results

Figure4 shows the impact of varying the processing times of P1M1, from t1 = 1 to t1 = 5, have on the waiting times right in front of P3M1. Figures 5(a)–5(g) show a similar analysis for the tested scenarios presented in Table 4.Figure 4
Impact of varying the processing time of P1M1, on the waiting time of P3M1.Figure 5
Scenarios 1 – 7.
(a)
Scenario 1
(b)
Scenario 2
(c)
Scenario 3
(d)
Scenario 4
(e)
Scenario 5
(f)
Scenario 6
(g)
Scenario 7Three things must be noted regarding Figures5(a)–5(g):(i)
EachP
i
M
i is located according to the layout presented in Figure 1.
(ii)
TheC
t
i is presented at the end of each column, and the U at the end of each row.
(iii)
All the m.a.t. relationships are always presented on the horizontal direction, as each row represents a certain machineM
i.
(iv)
All the p.a.t. relationships are always on the vertical direction, as each column represents a certain partP
i.Figure6 shows an analysis of the impact of varying the processing time of P2M1 has on the rest of the manufacturing system (Scenario 1 of Table 4, Figure 5(a)), where such impact is typified in terms of degree link types and INC/CON effects:(i)
1
s
t, 2
n
d, 3
r
d, and 4
t
h degree link type refers to the degree of closeness between the affecting P
i
M
i (that is, the source of processing time variation) and the affected P
i
M
i. In this sense, a 1
s
t degree link type means the affecting and the affected P
i
M
is are adjacent, where, in a 2
n
d link type, there is a P
i
M
i in between the affecting and the affected P
i
M
is. This idea of using 1
s
t, 2
n
d, 3
r
d, and 4
t
h connectivity degree links is taken from the supply network risk propagation model proposed by [43].
(ii)
INC (increase) and CON (constant) effect refers to the effect the affectingP
i
M
i has on the m.a.t/p.a.t. of the affected P
i
M
i and the t
i value at which these last start increasing, or the value at which these last remain constant.Figure 6
Waiting time values behind Scenario 1.The analysis of Figure6 shows the following:(i)
X1 (P2M1): its m.a.t increases due to a 2
n
d degree link from P3M1; its p.a.t. remains constant (and then starts increasing) due to a 3
r
d degree link from P2M3. This results in the pattern A of Figure 7.
(ii)
X2 (P3M1): both its m.a.t and its p.a.t. increase due to a 1
s
t degree link from P2M1 and a 3
r
d degree link from P3M2. This results in pattern B of Figure 7.
(iii)
X3 (P1M2): its m.a.t. increases due to a 3
r
d degree link from P3M2; its p.a.t. remains always constant due to a 4
t
h degree link from P1M3. This results in the pattern C of Figure 7.
(iv)
X4 (P2M2): its m.a.t. remains constant due to a 4
t
h degree link from P1M2; its p.a.t. increases due to a 1
s
t degree link from P2M1. This results in the pattern D of Figure 7.
(v)
X5 (P3M2): its m.a.t. remains constant (and then starts increasing) due to a 2
n
d degree link from P2M2; its p.a.t. increases due to a 2
n
d degree link from P3M1. This results in the pattern E of Figure 7.
(vi)
X6 (P1M3): its m.a.t. remains constant (and then starts increasing) due to a 3
r
d degree link from P2M3; its p.a.t. increases due to a 4
t
h degree link from P1M2. This results in the pattern F of Figure 7.
(vii)
X7 (P2M3): its m.a.t. remains constant due to a 4
t
h degree link from P1M3: its p.a.t. remains constant due to a 2
n
d degree link from P2M2. This results in the pattern G of Figure 7.Figure 7
Scenario 1 in the context of the m.a.t./p.a.t. analysis.


## 3.5. Results Analysis and Derived Conclusions

The previous analysis seems to suggest that the further a certain affectedP
i
M
i is from the affecting P2M1, the less impact its varying processing time has on it. As the ultimate objective of this research effort is to study the conditions under which the manufacturing disturbance propagation takes place, we performed a similar analysis on the seven tested scenarios of Table 4 (Table 7 presents a summary of the results obtained from this analysis, while the Appendix presents the full details of the results obtained). For this matter, we defined the type of impact the affecting P
i
M
i has on the m.a.t. and p.a.t. of the affected P
i
M
i: primary (TOT), the m.a.t./p.a.t. of the affected P
i
M
i increases; secondary (PAR), the m.a.t./p.a.t. of the affected P
i
M
i remains constant and then increases; or null (NON), the m.a.t./p.a.t. of the affected P
i
M
i remains constant. Also, we classify each affected P
i
M
i in terms of its location: level I, for X1 (P2M1) and X2 (P3M1); level II, for X3 (P1M2), X4 (P2M2), and X5 (P3M2); level III for X6 (P1M3) and X7 (P2M3). An analysis of the results presented in Table 7 reveals the following facts:(1)
Independently at which level the affectingP
i
M
i is located (i.e., I, II, and III), it has a primary impact on the m.a.t. of its immediate affected P
i
M
is (1
s
t degree link), 85.71% of the times. From here, its frequency starts decreasing: 21.42% for a 2
n
d degree link, 11.76% for a 3
r
d degree link, and 9.09% for a 4
t
h degree link. Moreover, 50% of the times there is a primary impact on an affected P
i
M
is, is due to a 1
s
t degree link affecting P
i
M
i, and only 25% of the times from 2
n
d degree link, 16.66% from a 3
r
d degree link, and 8.33% from a 4
t
h degree link. A similar situation presents the p.a.t. case. These results lead us to conclude that the primary impact from the affecting P
i
M
i is mostly contained by the immediate, 1
s
t degree link affected P
i
M
is.
(2)
Independently at which level the affectingP
i
M
i is located (i.e., I, II, and III), it has a secondary impact on the m.a.t. of its immediate affected P
i
M
is (1
s
t degree link), 14.28% of the times. From here, its frequency increases (42.85% for a 2
n
d degree link) and then decreases again (23.52% for a 3
r
d degree link). Moreover, 14.28% of the times there is a secondary impact on an affected P
i
M
is, is due to a 1
s
t degree link affecting P
i
M
i, 42.85% of the times from 2
n
d degree link, and 23.52% from a 3
r
d degree link. A similar situation presents the p.a.t. case. These results lead us to think that there must be some residues of the primary impact (absorbed by the 1
s
t degree link affected P
i
M
is) that keeps propagating further. This issue is reflected on the fact that the secondary impact affecting the 2
n
d degree link affected P
i
M
is is almost three times higher than the frequency of the secondary impact affecting the 1
s
t degree link affected P
i
M
is (42.85% versus 14.28%).
(3)
Independently at which level the affectingP
i
M
i is located (i.e., I, II, and III), it has a null impact on the m.a.t. of its immediate affected P
i
M
is (1
s
t degree link), 0% of the times. From here, its frequency starts increasing: 35.71% for a 2
n
d degree link, 64.70% for a 3
r
d degree link, and 90.9% for a 4
t
h degree link. Moreover, 0% of the times there is a null impact on an affected P
i
M
is, due to a 1
s
t degree link affecting P
i
M
i, 19.23% of the times from 2
n
d degree link, 42.30% from a 3
r
d degree link affecting, and 38.46% for a 4
t
h degree link. A similar situation presents the p.a.t. case. As these results go on the opposite direction of those regarding the frequency of primary impact, and the fact that a null impact is the complement of a primary impact, we arrive to the same conclusion already stated: the primary impact from the affecting P
i
M
i is mostly contained by the immediate 1
s
t degree link affected P
i
M
is.Table 7
m.a.t./p.a.t. analysis, degree link and impact types context.
m.a.t. analysis
p.at. analysis
TOT
PAR
NON
TOT
PAR
NON
TOT
PAR
NON
TOT
PAR
NON
1st
85.71
14.28
100
1st
50
9.09
1st
85.71
14.28
100
1st
46.15
10
2nd
21.42
42.85
35.71
100
2nd
25
54.54
19.23
2nd
28.57
42.85
28.57
100
2nd
30.76
60
15.38
3rd
11.76
23.52
64.70
100
3rd
16.66
36.36
42.30
3rd
12.5
12.5
75
100
3rd
15.38
30
46.15
4th
9.09
90.90
100
4th
8.33
38.46
4th
8.33
8.33
83.33
100
4th
7.69
38.46
100
100
100
100
100
100Figure8 shows the same results of Table 7, where Series1 refer to the TOT impact type, Series2 to NON impact type, an Series3 to the PAR impact type; 1
s
t, 2
n
d, 3
r
d, and 4
t
h refer to the connectivity link type.Figure 8
m.a.t. and p.a.t. values behavior (frequency), according to Table7.


## 3.6. Implications for Robustness Improvement

According to [1, 19, 20, 37], the improvement of robustness in a production system is provided by the presence of inventory (to make up production short-falls) and capacity buffers (slack available at resources). On the other hand, [67] mentions that, in order to improve the robustness of a production/supply chain network, its topology needs to be modified. Author in [68] expands this idea when he states that a dynamic network reconfiguration bring robustness, as new machines can be plugged to provide redundancy. Reference [37] claims that the robustness of the system depends on the resource with the minimum ability to recover. Reference [40] proposes a dynamic network of the manufacturing system, where this critical resource is the one with the highest correlation between performance and system disturbance. For [43], there is the node that provokes two or more downstream nodes to be disrupted. It is the authors belief that, within a context of a dynamic reconfiguration of the network topology, typical of an Industry 4.0 environment, where the production system reconfiguration takes place through the easy integration of plug-and-produce, fully automated, digitized, highly cost-efficient, smart new manufacturing units [69], following an approach like the proposed in this paper would allow(i)
identifying under which conditions a resource will become critical
(ii)
establishing the appropriate stock/inventory and time/capacity buffers countermeasures to minimize it


## 4. Conclusions and Future Research

One of the major challenges faced by the implementation of the Industry 4.0 concept is the robustness of the manufacturing/production system, that is, the ability to absorb manufacturing disruptions without failing or breaking. In this paper, we proposed the use of a Max-Plus algebra approach to study the impact of manufacturing disturbances (i.e., processing time variations) on the performance of a manufacturing system, as they propagate through it. For this purpose, an algorithm to derive the evolution equations (of a manufacturing system) was proposed and applied to a case study. Seven different scenarios were analyzed and the findings show that, independently where the processing time variation source is,(
1
) the primary disturbance effect is contained mostly by the immediate elements and (
2
) the secondary disturbance effect initially increases (as it propagates) and then vanishes. Future research will consider using the Max-Plus modelling approach proposed in this paper, to identify, within a context of a dynamic reconfiguration of the network topology, under which conditions a resource will become critical, so the appropriate robustness improvement strategies can be taken.



---

*Source: 1932361-2018-12-20.xml*