System reliability for a multistate flow network with multiple joint minimal paths under time constraint

Simulation Modelling Practice and Theory 29 (2012) 78–92

Contents lists available at SciVerse ScienceDirect

Simulation Modelling Practice and Theory journal homepage: www.elsevier.com/locate/simpat

System reliability for a multistate flow network with multiple joint minimal paths under time constraint Yi-Kuei Lin ⇑, Po-Hsun Kuo, Chin-Chia Chang Department of Industrial Management, National Taiwan University of Science & Technology, Taipei 106, Taiwan

a r t i c l e

i n f o

Article history: Received 27 March 2012 Received in revised form 13 July 2012 Accepted 21 July 2012 Available online 3 September 2012 Keywords: Joint minimal paths Monte Carlo Simulation Multistate flow network Time constraint Congestion

a b s t r a c t The uncertainty of transmission time for a multistate flow network is the main issue when the quickest path problem is applied to such a network. A variant of quickest path problem focusing on the probability of complete transmission instead of finding the quickest path is studied in this paper. In transmission process, a congestion phenomenon may happen when there are joint arcs in minimal paths (MPs). That is, the some arcs are shared by some of MPs. Based on Monte Carlo Simulation, an algorithm is developed to assign the flows through MPs for evaluating the system reliability, the probability that the data can be sent within a time constraint. Experiments are performed to verify the approximate system reliability value obtained from the proposed algorithm is close to the real one. Furthermore, the experiments about CPU time of the algorithm and the characteristics of the variants of quickest path problem are also displayed. Ó 2012 Elsevier B.V. All rights reserved.

1. Introduction In operations research, computer science, and other related fields, one of well-known and important problems is the shortest path problem. For such a problem, the data or goods are transmitted from one node to another node in a flow network. Based on this basic concept, there are many developed versions of the shortest path problem which adopt the shortest path, least-cost path, largest capacity path, shortest delay path, or some path that fulfills multiple criteria [1–4]. In addition to the above versions, another version which considers two attributes, capacity and lead time for each arc in a flow network is called the quickest path problem (QPP) that has been proposed by Chen and Chin [5]. In QPP, the purpose is to find a path (named the quickest path) with minimum transmission time for a given amount of data. Notice that two attributes for each arc are deterministic [5–7]. Both attributes are known before testing the reliable communications in a flow network. Several variants of the QPP have been proposed, including the constrained quickest path problem [8,9], k quickest paths problem [9–11], and all-pairs quickest path problem [12,13]; however, these variants of the QPP were constructed in a deterministic flow network. In general, a real flow network consists of a series of components (nodes or arcs) whose capacities may be uncertain due to failure, maintenance, etc. That is, the components are multistate; as a result, such a network is called a multistate flow network (MFN) [14–24]. The traditional QPP is solved for a deterministic flow network. Focused on MFN, this version of QPP is called the multistate quickest path problem (MQPP) herein. In MQPP, the quickest path is not fixed because of non-deterministic capacity for each component. Hence, the aim in MQPP is to evaluate the system reliability, the probability that the MFN can send a given amount of data from the source to the sink through multiple minimal paths (MPs) within a time constraint. An MP is a series of arcs and nodes that does not form a cycle. All MPs on a network can be efficiently derived from those algorithms discussed in [25,26]. The algorithm of Al-Ghanim [25] showed an approximate linear time response ⇑ Corresponding author. Tel.: +886 2 27303277; fax: +886 2 27376344. E-mail address: [email protected] (Y.-K. Lin). 1569-190X/$ - see front matter Ó 2012 Elsevier B.V. All rights reserved. http://dx.doi.org/10.1016/j.simpat.2012.07.008

Y.-K. Lin et al. / Simulation Modelling Practice and Theory 29 (2012) 78–92

79

Nomenclature QPP MFN MP MQPP MCS DDF TANET n; n0 q; q0 ai ; a0i A; A0 N; N0 Mi ; M0i M; M0 0 li ; li L; L0 G; G0 Pj ci fj C F Wj g Ns Nr d P0j dj si S t T dj(t)

quickest path problem multistate flow network minimal path multistate quickest path problem Monte Carlo Simulation distribution of data flow Taiwan Academic Network number of arcs in MFN; number of arcs in the revised MFN number of nodes in MFN; number of nodes in the revised MFN component i in MFN, i = 1, 2, . . . , n + q; component i in the revised MFN, i = 1, 2, . . . , n0 + q0   {aij1 6 i 6 n}: the set of arcs in MFN; a0i j1 6 i  6 n0 : the set of arcs  in the revised MFN 0 0 {aijn + 1 6 i 6 n + q}: the set of nodes in MFN; ai jn þ 1 6 i 6 n0 þ q0 : the set of nodes in the revised MFN maximal capacity of ai, i = 1, 2, . . . , n; maximal capacity of a0i ; i ¼ 1; 2; . . . ; n0   fM i j1 6 i 6 ng; M0i j1 6 i 6 n0 lead time of ai, i = 1, 2, . . . , n; lead time of a0i , i = 1, 2, . . . , n0 0  fli j1 6 i 6 ng; li j1 6 i 6 n0 (A, N, L, M): a MFN; (A0 , N0 , L0 , M0 ): a revised MFN jth MP in MFN, j = 1, 2, . . . , m current capacity of ai current flow on Pj (c1, c2, . . . , cn): a capacity vector (f1, f2, . . . , fm): a flow vector. maximum capacity of Pj number of assigned MP cumulative number of successes cumulative number of runs total amount of data jth MP in the revised MFN, j = 1, 2, . . . , g amount of data assigned to P0j ; j ¼ 1; 2; . . . ; g current capacity of a0i ; i ¼ 1; 2; . . . ; n0 (s1, s2, . . . , sn0 ): a system state point in time time constraint amount of data assigned to P0j in data buffer at t, j = 1, 2, . . . , g, t = 1, 2, . . . , T

mj

number of components on P0j

ajk

kth component on P0j ; j ¼ 1; 2; . . . ; g; k ¼ 1; 2; . . . ; mj

sjk j lk X jk ðtÞ ðtÞ Y j;c k

capacity of ajk ; j ¼ 1; 2; . . . ; g; k ¼ 1; 2; . . . ; mj lead time of ajk ; j ¼ 1; 2; . . . ; g; k ¼ 1; 2; . . . ; mj amount of data assigned to P0j on ajk at t, t = 1, 2, . . . , T, j = 1, 2, . . . , g, k = 1,3, . . . , mj (odd number)

qjk

amount of data assigned to P0j at cth position inside ajk at t, t = 1, 2, . . . , T, j = 1, 2, . . . , g, k = 2,4, . . . , mj  1 (even j number), c ¼ 1; 2; . . . ; lk  1 maximum available capacity of ajk

Rz MP inflow outflow

approximate system reliability for repetition z, z = 1, 2, . . . an MP is a path whose proper subsets are no longer paths the amount of something entering into somewhere at point in time the amount of something leaving from somewhere at point in time

versus the number of network nodes. And for a network of 101 nodes, including 10 branch nodes, the computational time per path is 0.518 second (run on a PC 486 machine). Kobayashi and Yamamoto [26] showed that to generate all minimal paths for a random network with 30 nodes and 100 arcs takes no more than 1300 s. In the beginning of development of MQPP, only time factor is concentrated [27,29]. Subsequently, cost (budget) factor and spare routing are involved and discussed [28]. However, all of the studies on the MQPP [27–29] are under the strict assumption that the MPs are disjoint. In real world, the joint situation is frequently occurred due to management consideration and risk reduction. Hence, a solution to MQPP with either joint or disjoint MPs is worth developing. For joint MPs, the transmission time of the system is difficult to calculate, because the data that comes from different MPs may pass through an arc at the same time such that a congestion phenomenon may occur. No algorithm has been developed yet to evaluate the system

80

Y.-K. Lin et al. / Simulation Modelling Practice and Theory 29 (2012) 78–92

reliability for an MQPP with some joint MPs. Either analytical approaches or mathematical means have difficulty on dealing with this phenomenon. Hence, the purpose of this paper is to evaluate an approximate system reliability through simulation techniques. Monte Carlo Simulation (MCS) is a useful approach for modeling a system to explore a phenomenon and for evaluating an approximate value of a complicated problem. The main idea in MCS is to rely on repeated random sampling to find an approximate value. Therefore, different sampling plans are further introduced and compared to each other in [30,31]. Other issues about MCS are how to determine the suitable length of experiment and to improve the quality of the approximate value [32]. The studies related to MCS are copious and the various fields are covered, because it is a straightforward and general technique compared to other methodologies. Focusing on a binary-state network, Zenklusen et al. [33] put efforts to keep the Monte Carlo based algorithm accurate for evaluating system reliability. In particular, the system reliability is the probability to connect source with sink. Lin and Yang [34] also concentrated on a binary-state network to propose a sampling estimator to speed up the efficiency of MCS under highly reliable environment. With respect to MFN, Ramirez-Marquez and Coit [35] proposed an algorithm based on MCS to estimate the system reliability without considering time attribute. Hence, this paper addresses the evaluation of system reliability for a MFN with multiple joint MPs under time constraint. The remainder of this paper is organized as follows: Section 2 covers the properties of MFN, the simplification of MFN and the assumptions about the MQPP. In Section 3, a series of equations are created in order to simulate the process of data flow. In Section 4, the structure of the algorithm is displayed and a simple example is provided. Section 5 lists not only the credibility of the results obtained from the proposed algorithm but also the dependence of system reliability on the amount of data and the time constraint. Finally, Section 6 draws conclusions about the present contribution and future developments. 2. Model building A MFN which owns a source and a sink node is considered a classical MFN. Such a network is denoted as G = (A, N, L, M) in which A = {aij1 6 i 6 n} is the set of arcs, N = {aijn + 1 6 i 6 n + q} is the set of nodes, L = {lij1 6 i 6 n} with li being the lead time of ai, and M = {Mij1 6 i 6 n} with Mi being the maximal capacity of ai. Besides, an MP is denoted as Pj consisting of arcs and nodes in the G, and does not involve any cycle in it. Furthermore, such a G is assumed to further satisfy the following assumptions: (1) (2) (3) (4) (5)

The data can be transmitted through multiple joint MPs simultaneously. Each node is perfectly reliable. The capacity of each arc is an integer-valued random variable with a given probability distribution. The capacities of different arcs are statistically independent. Data flow must satisfy the flow-conservation law [36] in transmission process.

2.1. MFN In MFN, each arc has several states (or capacities) with a probability distribution. Suppose that P1, P2, . . . , Pm are all MPs from source to sink node. The MFN is descried in terms of the capacity vector C = (c1, c2, . . . , cn) and flow vector F = (f1, f2, . . . , fm) where ci denotes the current capacity of ai and fj denotes the current flow on Pj. Let Wj denote the maximum capacity of Pj, i.e., Wj = min{Mijai 2 Pj}. A feasible flow vector F satisfies the following inequalities:

fj 6 W j for j ¼ 1; 2; . . . ; m; m X ffj jai 2 Pj g 6 M i for i ¼ 1; 2; . . . ; n:

ð1Þ ð2Þ

j¼1

Inequality (1) expresses that the flow on each Pj cannot exceed the maximum capacity of Pj. Inequality (2) expresses that the total flow through ai under F cannot exceed the maximum capacity of ai. In this paper, the MQPP is constructed in an MFN. The purpose of MQPP is to evaluate the system reliability which is defined as the probability a given amount of data can be sent from the source to the sink through the multiple joint MPs within a time constraint. 2.2. The revised MFN Since the data is determined to transmit through g known MPs simultaneously, the components (arcs or nodes) not on these MPs can original network G0 = (A0 , N0 , L0 , M0 ),  be ignored.  For this reason, the S  G can be simplified to the revised network S where A0 ¼ a0i j1 6 i 6 n0 is the set of arcs on gj¼1 Pj and N 0 ¼ a0i jn0 þ 1 6 i 6 n0 þ q0 is the set of nodes on gj¼1 Pj . There are 0  two attributes for each a0i 2 A0 :one is the lead time, denoted by L0 ¼ li j1 6 i 6 n0 ; and the other is the maximal capacity, 0 0 0 0 0 denoted by M ¼ M i j1 6 i 6 n0 . Furthermore, P1 ; P2 ; . . . ; P g are used to replace P1, P2, . . . , Pg after G is simplified. They consist of a series of renamed components (arcs or nodes). Due to failure, maintenance, etc., each arc has multiple possible capacities. The current capacity of an arc a0i is denoted by si. To make above transformation clear, considering a G shown in Fig. 1 if

Y.-K. Lin et al. / Simulation Modelling Practice and Theory 29 (2012) 78–92

81

Fig. 1. A simple network.

Fig. 2. A revised network.

P1 = {a7, a1, a8, a2, a10} and P2 = {a7, a5, a9, a4, a8, a2, a10} are assigned, G will be revised to G0 shown in Fig. 2. Obviously, two components are deleted are renumbered. So P1 is replaced by   in the process, and the  remaining components  P 01 ¼ a05 ; a01 ; a06 ; a04 ; a08 and P2 is replaced by P 02 ¼ a05 ; a02 ; a07 ; a03 ; a06 ; a04 ; a08 . For most situations, we can delete several components from original network, and then the revised one is simpler like the mentioned example. As a result, the process of the algorithm can be speeded up. 2.3. System state S The capacity of each a0i 2 A0 is influenced by its own probability distribution. Once the capacity of a0i is determined, the current capacity is denoted as si. The system state is a vector which involves all a0i 2 A0 based on G0 . It is further denoted by a vector S ¼ ðs1 ; s2 ; . . . ; sn0 Þ and directly affects transmission ability of each MP. In order to estimate system reliability through MCS, the status of success or failure for the system state S should be clarified. Therefore, we define success as the generated S triggers all data arriving the destination within time constraint T; failure as the generated S triggers the part of data arriving the destination beyond time constraint T. 3. Simulation for distribution of data flow The distribution of data flow (DDF) represents the amount of data located at different positions in the revised network G0 . The DDF changes over time due to inflow and outflow at positions in the network. Inflow represents the amount of data that is expected to enter somewhere (a position) at a certain point in time. Conversely, outflow represents the amount of data that is expected to leave somewhere (a position) at a certain point in time. During the progress of transmission, data can be located at positions of different types. We define three types of position: data buffer, node, and arc. A simple example is shown in Fig. 3 to illustrate these three types of position and the DDF at different points in time. Initially (t = 0), all data is in the data buffer, and then data start moving from the data buffer to another position in the network. The outflow is limited by arc capacity, so only three units of data are allowed to reach the source (node) at t = 1. Since the lead time is three units, the arc is divided into three segments and there are two positions within the arc. At t = 3, all data is located within the arc, but we can distinguish three units of data from the remaining two units of data based on the different positions. 3.1. Data assignment on minimal paths A given amount of data d should be transmitted through g MPs simultaneously. Before transmitting, data assignments that represent how much data is carried through each MP are determined. All possible data assignments are non-negative P integer solutions which are denoted by (d1, d2, . . . , dg) and satisfy the condition gj¼1 dj ¼ d, where dj denotes the amount of 0 data assigned to P j . Under a generated system state S, each (d1, d2, . . . , dg) will be used to build the initial DDF itself. The termination condition is that t > T or the amount of data at the sink node has already increased to d. If the latter condition is satisfied, the generated S is regarded as a success; otherwise, another (d1, d2, . . . , dg) is used to restart the DDF simulation. If the DDF built from every (d1, d2, . . . , dg) fails to result in the amount d reaching the sink node within the time constraint T, the generated S is regarded as a failure.

82

Y.-K. Lin et al. / Simulation Modelling Practice and Theory 29 (2012) 78–92

(a) A simple example.

(b) The evolution of DDF. Fig. 3. Illustration of the evolution of DDF.

3.2. Data at different positions Let dj(t) denote the amount of data assigned to P0j in the data buffer at time t = 0,1, . . . , T. The kth component on P 0j is denoted by ajk where k = 1, 2, . . . , mj. The components on P 0j can be a node or an arc. When k is an odd number, the component ajk is a node; otherwise, an arc. No doubt, the component aj1 and ajmj are the source and the sink, respectively. Notice that   there exists several ajk indicating the same component a0i . For example, considering P 01 ¼ a05 ; a01 ; a06 ; a04 ; a08 and   P02 ¼ a05 ; a02 ; a07 ; a03 ; a06 ; a04 ; a08 in Fig. 2, all of components, a05 ; a06 ; a04 and a08 can be represented by a1k and a2k . Only the compoj nent a04 is a joint arc, and the others are nodes. Let sjk and lk denote the capacity and the lead time of ajk , respectively. When k = 2,4, . . . , mj  1 (even number), the type of position is an arc. Let X jk ðtÞ denote the amount of data assigned to P0j on ajk at time t = 0,1, . . . , T where k = 1,3, . . . , mj (odd number) that means the type of position is a node. Furthermore, Y j;c k ðtÞ denotes j the amount of data assigned to P 0j at cth position inside ajk at time t = 0,1, . . . , T, where 1 6 c < lk and k = 2,4, . . . , mj  1 (even number) that means the type of position is an arc. 3.3. Equations for implementing simulation For implementing the simulation for DDF, a series of equations will be introduced here. Firstly, there are four terms which form the structure of each equation as follows: term term term term

(0): (1): (2): (3):

How much data will be stored at t + 1; the amount of data is already stored at t; inflow at point in time t; outflow at point in time t.

The structure of equation can be represented as the form: term (0) = term (1) + term (2)  term (3). Without loss of generality for P0j , let P 01 own the highest priority, then P 02 and so on. According to differences of position types, the equations in simulation for DDF are divided into three parts to be discussed: in data buffer, at nodes, and on arcs. 3.3.1. Simulation in data buffer The following equation is utilized to obtain dj(t) in data buffer for t = 0,1, . . . , T, with j = 1, 2, . . . , g:

n o dj ðt þ 1Þ ¼ dj ðtÞ  min dj ðtÞ; qj2 ;

ð3Þ

Y.-K. Lin et al. / Simulation Modelling Practice and Theory 29 (2012) 78–92

n o 8 j > s2 ; if j ¼ 1 or hjaj2 2 P0h ; h < j ¼ £; > > > > < 8 9 j > > where q2 ¼ < = X > j j > > d  min ðtÞ; s s h > 2 2 ; otherwise: > > > : : ; j h2fhja 2P 0 ;h
83

ð4Þ

h

P The data buffer is the original position which stores all data, so dj(t = 0) = dj and gj¼1 dj ¼ d. Furthermore, the term (2) representing the inflow has no possibility to happen in data buffer, so it is removed in Eq. (3). Owing to arc’s capacity limit, only n o the amount min dj ðtÞ; qj2 can be transmitted where qj2 denotes the maximum available (or residual) capacity of aj2 . For each aj2 , the qj2 is equal to sj2 if there is not any higher priority P0h which shares the capacity of aj2 . Otherwise, the qj2 is replaced by nP o sj2  min d ðtÞ; sj2 . j h2fhja 2P 0 ;h
h

3.3.2. Simulation at nodes Here the position type what we focus on is a node, so all ajk ; j ¼ 1; 2; . . . ; g and k = 1,3, . . . , mj  2 (odd number) are discussed and all X jk ðt ¼ 0Þ ¼ 0. According to network structure, these nodes are further divided into a source, intermediated nodes, and a sink. If the node is the source (k = 1), the equations for obtaining X jk ðtÞ for t = 0,1, . . . , T, with j = 1, 2, . . . , g are as follows:

n o n o X j1 ðt þ 1Þ ¼ X j1 ðtÞ þ min dj ðtÞ; qj2  min X j1 ðtÞ; qj2 ; 8 n o > > sj2 ; if j ¼ 1 or hjaj2 2 P0h ; h < j ¼ £; > > > 8 9 < > > < = where qj2 ¼ X j > > dh ðtÞ; sj2 ; otherwise: > s2  min > > > > :h2fhjaj 2P0 ;h
ð5Þ

ð6Þ

h

In Eq. (5), the term (2) is equal to the term (3) in Eq. (3). The source is directly connected to the data buffer, so the outflow in data buffer is the inflow on source. If the node are intermediated nodes, the equations for obtaining X jk ðtÞ for t = 0,1, . . . , T, with j = 1, 2, . . . , g and k = 3, . . . , mj  2 (odd number) are as follows:

n o 8 j j j;l1 j j > < X k ðtÞ þ Y k1 ðtÞ  min X k ðtÞ; qkþ1 ; if lk1 > 1; j X k ðt þ 1Þ ¼ n o n o > j j : X j ðtÞ þ min X j ðtÞ; qj k k2 k1  min X k ðtÞ; qkþ1 ; otherwise;

where qjkþ1

n o 8 j > skþ1 ; if j ¼ 1 or hjajkþ1 2 P 0h ; h < j ¼ £; > > > > < 8 9 > > ¼ < = X > j h j > >  min X ðtÞ; s s ; otherwise; j h a > kþ1 kþ1 w:a > w > > k : : ; j h2fhja 2P 0 ;h
and qjk1

k1

Notice that in Eq. (7)

ð8Þ

h

n o 8 j j 0 > s ; if j ¼ 1 or hja 2 P ; h < j ¼ £; > h k1 k1 > > > < 8 9 > > ¼ < = X > j h j > >  min X ðtÞ; s s ; otherwise: h aj > k1 k1 w:a > w > > k2 : :h2fhjaj 2P0 ;h
ð7Þ

ð9Þ

h

j

¼ 0Þ ¼ 0, where l denotes lk1 for simplification and (l  1)th position is the last position inside

j

ajk1 . If lk1 ¼ 1, it means that there is no position inside ajk1 , so the previous position data located is ajk2 which is a node. The n o ðtÞ to min X jk2 ðtÞ; qjk1 . term (2) changes from Y j;l1 k1 If the node is the sink (k = mj), the equations for obtaining X jk ðtÞ for t = 0,1, . . . , T, with j = 1, 2, . . . , g are as follows:

8 j j j;l1 > < X mj ðtÞ þ Y mj 1 ðtÞ; if lmj 1 > 1; j X mj ðt þ 1Þ ¼ n o > j : X j ðtÞ þ min X j mj mj 2 ðtÞ; qmj 1 ; otherwise;

ð10Þ

84

Y.-K. Lin et al. / Simulation Modelling Practice and Theory 29 (2012) 78–92

where qjmj 1

n o 8 j > smj 1 ; if j ¼ 1 or hjajmj 1 2 P0h ; h < j ¼ £; > > > 8 9 > < > > > > < = ¼ X j j h > >  min X ðtÞ; s s > 2 m mj 1 ; otherwise: > mj 1 h > > > > > : :h2fhjajm 1 2P0h ;h
ð11Þ

j

The sink is the final position, so the term (3) indicating outflow is not involved in Eq. (10). It also has two types, because of lead time of the specific arcs. 3.3.3. Simulation on arcs Here the position type what we focus on is an arc, and all the following equations are utilized to obtain Y j;c ðtÞ for k j t = 0,1, . . . , T, with j = 1, 2, . . . , g, k = 2,4, . . . , mj  1 (even number) and c ¼ 1; 2; . . . ; lk  1:

Y j;c k ðt þ 1Þ ¼

8 n < min X j

j k1 ðtÞ; qk

:

o ; if c ¼ 1;

ðtÞ; otherwise; Y j;c1 k

8 n o > > sjk ; if j ¼ 1 or hjajk 2 P0h ; h < j ¼ £; > > > 8 9 < > > < = where qjk ¼ X j > > X hw:ah aj ðtÞ; sjk ; otherwise: > sk  min > w > > k1 > :h2fhjaj 2P0 ;h
ð12Þ

ð13Þ

h

Y j;c k ðt

In Eq. (12), all ¼ 0Þ ¼ 0. Corresponding to the structure of equation, only the term (2) is retained in Eq. (12). The amount of data located at positions inside an arc is already determined when it prepares to enter an arc, so the outflow do not need to consider arc’s capacity limit again; therefore, the term (3) is always the same to the term (1). Besides, there arentwo types in o Eq. (12) according to the positions inside an arc. If c = 1, the previous position is a node, so the inflow is min X jk1 ðtÞ; qjk . Otherwise, the previous position is the other positions inside the same arc, so the inflow is Y j;c1 ðtÞ. k 4. The proposed algorithm In this section, we summarize all concepts mentioned before. First of all, the original network is simplified into a revised network. A system state S is then generated randomly for each run and the simulation for DDF is performed accordingly. In each run, the system state is recorded as being a success or failure. Until the final run is finished, the approximate system reliability is obtained as Ns/Nr, where Ns and Nr denote the cumulative number of successes and the cumulative number of runs, respectively. 4.1. Algorithm to approximate the system reliability Input: G = (A, N, L, M); the number of runs; the amount of data d; time constraint T; g MPs. Step 0. [Initialization] Ns = 0, Nr = 0. Step 1. Build G0 = (A0 , N0 , L0 , M0 ). P Step 2. Find all non-negative integer solutions (d1, d2, . . . , dg) such that gj¼1 dj ¼ d. Step 3. For each run, generate a system state S ¼ ðs1 ; s2 ; . . . ; sn0 Þ randomly. Step 4. Perform simulation for DDF for the system state S. (4.1) For each non-negative integer solution (d1, d2, . . . , dg), start from t = 0. (4.2) Use Eq. (3)–(13) to compute dj ðt þ 1Þ; X jk ðt þ 1Þ for k = 1,3, . . . , mj (odd number), and Y j;c k ðt þ 1Þ for j k = 2,4, . . . , mj  1 (even number), j ¼ 1; 2; . . . ; g; c ¼ 1; 2; . . . ; lk  1. P (4.3) If gj¼1 X jmj ðt þ 1Þ ¼ d, then Ns = Ns + 1 and go to step 5. (4.4) t = t + 1. (4.5) If t 6 T, then go to step 4.2. (4.6) Next (d1, d2, . . . , dg). Go to step 4.1. Step 5. Nr = Nr + 1. Step 6. If Nr is equal to the number of runs, then go to Step 8. Step 7. Next run. Go to step 3. Step 8. The approximate system reliability = Ns/Nr. Output: The approximate system reliability. In step 4.2, the amount of data at different positions in revised MFN is calculated per unit time. Besides, the step 4.3 checks whether all amount of data stored at the latest position for each MP are the same to the given amount d of data or not.

85

Y.-K. Lin et al. / Simulation Modelling Practice and Theory 29 (2012) 78–92

4.2. An illustrative example We use the MFN shown in Fig. 1 to illustrate the proposed algorithm. Table 1 lists the probability distribution and lead time for each arc. Two MPs are P1 = {a7, a1, a8, a2, a10} and P2 = {a7, a5, a9, a4, a8, a2, a10}. For d = 12, time constraint T = 14, and 6000 runs, the system reliability can be calculated by the following steps: Step Step Step Step Step

0. Ns = 0, Nr = 0.     1. The revised network is shown in Fig. 2. Two MPs are P01 ¼ a05 ; a01 ; a06 ; a04 ; a08 and P 02 ¼ a05 ; a02 ; a07 ; a03 ; a06 ; a04 ; a08 . 2. All non-negative integer solutions are (0, 12), (1, 11), . . . , (12, 0). 3. For a randomly generated S = (s1, s2, s3, s4) = (2,3,1,2), do the following steps. 4. Perform simulation for DDF. (4.1) For (d1, d2) = (5, 7) and start from t = 1.   (4.2) d1 ð1Þ ¼ d1 ð0Þ  min d1 ð0Þ; s12 ¼ 5  minf5; 2g ¼ 3. n o   X 11 ð1Þ ¼ X 11 ð0Þ þ min d1 ð0Þ; s12  min X 11 ð0Þ; s12 ¼ 0 þ minf5; 2g  minf0; 2g ¼ 2. n o 1 1 Y 1;1 ¼ 0. 2 ð1Þ ¼ min X 1 ð0Þ; s2 1;1 Y 1;2 2 ð1Þ ¼ Y 2 ð0Þ ¼ 0. 1;2 Y 1;3 2 ð1Þ ¼ Y 2 ð0Þ ¼ 0.

n o ð0Þ  min X 13 ð0Þ; s14 ¼ 0. X 13 ð1Þ ¼ X 13 ð0Þ þ Y 1;3 2 n o 1 1 ¼ 0. Y 1;1 4 ð1Þ ¼ min X 3 ð0Þ; s4 1;1 Y 1;2 4 ð1Þ ¼ Y 4 ð0Þ ¼ 0.

X 15 ð1Þ ¼ X 15 ð0Þ þ Y 1;2 4 ð0Þ ¼ 0.   d2 ð1Þ ¼ d2 ð0Þ  min d2 ð0Þ; s22 ¼ 7  minf7; 3g ¼ 4.   X 21 ð1Þ ¼ X 21 ð0Þ þ min d2 ð0Þ; s22  minfX 21 ð0Þ; s22 g ¼ 0 þ minf7; 3g  0 ¼ 3. n o 2 2 ¼ 0. Y 2;1 2 ð1Þ ¼ min X 1 ð0Þ; s2 2;1 Y 2;2 2 ð1Þ ¼ Y 2 ð0Þ ¼ 0.

n o 2 2 X 23 ð1Þ ¼ X 23 ð0Þ þ Y 2;3 ¼ 0. 2 ð0Þ  min X 3 ð0Þ; s4 n o 2;1 2 2 Y 4 ð1Þ ¼ min X 3 ð0Þ; s4 ¼ 0. n n oo X 25 ð1Þ ¼ X 25 ð0Þ þ Y 2;1 ¼0 ð0Þ  min X 25 ð0Þ; s26  min X 13 ð0Þ; s26 n 4 o 2 1 2 2 Y 2;1 ð1Þ ¼ min X ð0Þ; s  minfX ð0Þ; s g ¼ 0. 6 6 5 3 6 2;1 Y 2;2 6 ð1Þ ¼ Y 6 ð0Þ ¼ 0.

X 27 ð1Þ ¼ X 27 ð0Þ þ Y 2;2 6 ð0Þ ¼ 0. Table 1 The arc data of Fig. 1.

a

Arc

Revised arc

Capacity

Probability

Lead time

a1

a01

2a 1 0

0.97 0.02 0.01

4

a2

a04

3 2 1 0

0.95 0.02 0.02 0.01

3

a3

N/A

2 1 0

0.95 0.04 0.01

1

a4

a03

2 1 0

0.97 0.02 0.01

2

a5

a02

3 2 1 0

0.95 0.02 0.02 0.01

3

a6

N/A

1 0

0.95 0.05

2

Pr {the capacity of a1 is 2} = 0.97.

86

Y.-K. Lin et al. / Simulation Modelling Practice and Theory 29 (2012) 78–92

Fig. 4. A benchmark network.

Table 2 The arc data of Fig. 4. Arc

Capacity

Probability

Lead time

Arc

Capacity

Probability

Lead time

a1

3 2 1 0 3 2 1 0 3 2 1 0 1 0 2 1 0

0.80 0.10 0.05 0.05 0.80 0.10 0.05 0.05 0.80 0.10 0.05 0.05 0.90 0.10 0.80 0.10 0.10

2

a6

4 3 2 1 0 5 4 3 2 1 0 4 3 2 1 0

0.60 0.20 0.10 0.05 0.05 0.55 0.10 0.10 0.10 0.10 0.05 0.70 0.10 0.10 0.05 0.05

2

a3

a4 a5

1 a7

3

a8 3 1

Approximate system reliability

a2

Cumulative number of runs Fig. 5. The trend of variation of approximate system reliability.

2

1

87

Y.-K. Lin et al. / Simulation Modelling Practice and Theory 29 (2012) 78–92 Table 3 The comparisons in different path combinations.

a

Combination

Mean of approximate system reliability

St.dev.

Actual system reliability

Differencea

P1 and P5 P1 and P7 P1 and P8 P1 and P9 P1 and P10 P2 and P5 P2 and P7 P2 and P8 P2 and P9 P2 and P10 P3 and P6 P3 and P7 P3 and P8 P3 and P9 P3 and P11 P4 and P5 P4 and P6 P4 and P10 P4 and P11 P5 and P8 P5 and P9 P5 and P11 P6 and P8 P6 and P9 P6 and P10 P7 and P10 P7 and P11 P1, P5 and P8 P1, P5 and P9 P1, P7 and P10 P2, P5 and P8 P2, P5 and P9 P2, P7 and P10 P3, P6 and P8 P3, P6 and P9 P3, P7 and P11 P4, P5 and P11 P4, P6 and P10

0.87821 0.51216 0.50015 0.48864 0.78433 0.85303 0 0.43348 0.38919 0.60945 0.61125 0.644 0.62647 0.60355 0.57868 0.86564 0 0.63928 0.38919 0.8916 0.87849 0.86293 0.58632 0.52609 0.6784 0.71905 0.52609 0.94728 0.941573 0.8782 0.942293 0.934293 0.74551 0.85288 0.83021 0.83503 0.932467 0.73459

0.00662 0.01115 0.01081 0.01065 0.00722 0.00756 0 0.01072 0.00913 0.00819 0.01044 0.01005 0.0098 0.01014 0.00937 0.00642 0 0.00857 0.00913 0.005462 0.00574 0.00701 0.00897 0.00971 0.00993 0.00764 0.00971 0.004188 0.004063 0.00608 0.004254 0.004592 0.00933 0.00715 0.00791 0.00831 0.004651 0.01052

0.8784 0.5136 0.5015 0.4871 0.7855 0.8532 0 0.4322 0.3889 0.6097 0.6125 0.6441 0.626 0.6042 0.5779 0.8669 0 0.6395 0.3889 0.8923 0.8796 0.8643 0.5852 0.5267 0.6781 0.7185 0.5267 0.9476 0.9407 0.8783 0.9428 0.9338 0.7464 0.8539 0.8298 0.8351 0.9321 0.7358

0.00019 0.00144 0.00135 0.00154 0.00117 0.00017 0 0.00128 0.00029 0.00025 0.00125 1E04 0.00047 0.00065 0.00078 0.00126 0 0.00022 0.00029 0.0007 0.00111 0.00137 0.00112 0.00061 0.0003 0.00055 0.00061 0.00032 0.000873 1E04 0.000507 0.000493 0.00089 0.00102 0.00041 7E05 0.000367 0.00121

Difference = (mean of approximate system reliability)  (actual system reliability).

(4.3) X 15 ð1Þ þ X 27 ð1Þ ¼ 0–12. (4.4) t = 1 + 1 = 2. (4.5) t 6 14 and go to step 4.2. .. . (the loop to repeat steps 4.2–4.5) (4.3) X 15 ð14Þ þ X 27 ð14Þ ¼ 5 þ 6 ¼ 11 – 12. (4.4) t = 14 + 1 = 15. (4.5) t > 14. (4.6) Next (d1, d2) = (6, 6), go to step 4.1. .. . (the loop to repeat steps 4.1–4.6) (4.3) X 15 ð14Þ þ X 27 ð14Þ ¼ 6 þ 6 ¼ 12, then Ns = 0 + 1 = 1 and go to step 5. Step 5. Nr = 0 + 1 = 1. Step 6. Nr = 1 – 6000. . Step 7. Next run. Go to step 3...(the loop to repeat steps 3–7) Step 6. Nr = 6000, then go to Step 8. Step 8. The approximate system reliability is Ns/Nr = 5884/6000 = 0.98067. 5. Experiment procedures and results In this section, two scenarios are involved to validate the proposed algorithm and to observe the CPU time. In the first scenario, we determine whether the approximate system reliability obtained by the proposed algorithm is close to the actual system reliability. In the second scenario, a large-scale MFN is utilized to show how approximate system reliability and CPU time change with different amount of data and time constraints.

88

Y.-K. Lin et al. / Simulation Modelling Practice and Theory 29 (2012) 78–92

Fig. 6. The Taiwan Academic network [37].

5.1. Approximate system reliability versus actual system reliability For the purpose of demonstrating the credibility of the approximate system reliability obtained by the proposed algorithm, we use the network shown in Fig. 4 and all disjoint combinations of MP to be taken in account. The arc data is listed in Table 2. Considering MQPP with these disjoint MPs, the actual system reliability can be calculated by the approaches from [28,29] which find candidates of lower boundary vector to each data assignment, and then compare them in order to obtain lower boundary vectors. Once such vectors are determined, the probability (system reliability) that the system state is greater than them can be calculated. The comparisons can be performed once the actual system reliability is determined. Suppose that the given amount of data d = 12 and time constraint T = 14. In Fig. 4, the network has eleven MPs: P1 = {a13, a1, a14, a2, a15, a3, a19}, P2 = {a13, a1, a14, a4, a16, a5, a15, a3, a19}, P3 = {a13, a1, a14, a4, a16, a7, a19}, P4 = {a13, a1, a14, a4, a16, a9, a18, a12, a19}, P5 = {a13, a6, a16, a7, a19}, P6 = {a13, a6, a16, a5, a15, a3, a19}, P7 = {a13, a6, a16, a9, a18, a12, a19}, P8 = {a13, a10, a17, a11, a18, a12, a19}, P9 = {a13, a10, a17, a8, a16, a9, a18, a12, a19}, P10 = {a13, a10, a17, a8, a16, a7, a19}, and P11 = {a13, a10, a17, a8, a16, a5, a15, a3, a19}. The total of different combinations of disjoint MP is up to 38. We observe the trend of the variation of approximate system reliability with the increase of number of runs. Fig. 5 displays that the variation of approximate system reliability becomes stable as number of runs increases. Besides, when the cumulative number of runs is over 2000, the curve is quite smooth. Hence, according to Fig. 5 we prefer to set the number greater than 2000 as the number of runs. Here, we set 2500 as the number of runs for each combination. The algorithm adopting 2500 runs is performed for 30 repetitions and approximate system reliabilities R1, R2, . . . , R30 are obtained. For each path combination, the experimental results are listed in Table 3. The differences (Mean of approximate system reliability  Actual system reliability) show that the error is less than 0.02. In addition, the standard deviation is also less than 0.02. It implies that the proposed algorithm is robust.

89

Y.-K. Lin et al. / Simulation Modelling Practice and Theory 29 (2012) 78–92 Table 4 The arc data of Fig. 6. Arc

a1 a2 a3 a4 a5 a6 a7 a8 a9 a10 a11 a12 a13 a14 a15 a16 a17 a18 a19 a20 a21 a22 a23 a24 a25 a26 a27 a28 a29 a30 a31

Lead time (s)

1 1 2 1 1 1 2 1 1 1 2 1 1 3 2 1 1 2 2 1 3 1 1 1 2 2 1 3 2 3 1

Capacity (Gbps) 0

1

2

3

4

5

0.000064 0.000343 0.000100 0.000041 0.000216 0.000125 0.001000 0.000008 0.000064 0.000027 0.000064 0.000343 0.000125 0.000010 0.000125 0.000024 0.000008 0.000024 0.000003 0.000003 0.000003 0.008100 0.002500 0.001600 0.003600 0.004900 0.004900 0.003600 0.000006 0.000006 0.000041

0.004608 0.013671 0.003600 0.001884 0.010152 0.007125 0.027000 0.001176 0.004608 0.002619 0.004608 0.013671 0.007125 0.000450 0.007125 0.001276 0.001176 0.001276 0.000188 0.000188 0.000188 0.163800 0.095000 0.076800 0.112800 0.130200 0.130200 0.112800 0.000475 0.000475 0.001884

0.110592 0.181629 0.048600 0.032502 0.159048 0.135375 0.243000 0.057624 0.110592 0.084681 0.110592 0.181629 0.135375 0.008100 0.135375 0.025428 0.057624 0.025428 0.004334 0.004334 0.004334 0.828100 0.902500 0.921600 0.883600 0.864900 0.864900 0.883600 0.013538 0.013538 0.032502

0.884736 0.804357 0.291600 0.249180 0.830584 0.857375 0.729000 0.941192 0.884736 0.912673 0.884736 0.804357 0.857375 0.072900 0.857375 0.225220 0.941192 0.225220 0.049836 0.049836 0.049836 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.171475 0.171475 0.249180

0.000000 0.000000 0.656100 0.716393 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.328050 0.000000 0.748052 0.000000 0.748052 0.286557 0.286557 0.286557 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.814506 0.814506 0.716393

0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.590490 0.000000 0.000000 0.000000 0.000000 0.659082 0.659082 0.659082 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000

5.2. Taiwan Academic Network In the second scenario, the large-scale MFN that we employ is the Taiwan Academic Network (TANET) with 31 arcs and 25 nodes [37], as shown in Fig. 6. The arc information is listed in Table 4. As three MPs we select P1 = {a32, a1, a33, a2, a34, a3, a35, a4, a36, a5, a37, a6, a38, a7, a39, a8, a40, a18, a48, a19, a49, a20, a58}, P2 = {a32, a14, a45, a15, a46, a16, a47, a17, a40, a18, a48, a19, a49, a20, a58} and P3 = {a32, a22, a50, a23, a51, a24, a52, a25, a53, a26, a54, a27, a55, a28, a58}, where a18, a19, and a20 are joint arcs shared by P1 and P2. According to these MPs, a revised network with 22 arcs and 21 nodes is built in Fig. 7. Suppose that the number of runs 10,000 is set for the proposed algorithm, a series of procedures in it are further performed. All system reliability and CPU time corresponding to different arrangements of amount of data d and time constraint T are gathered in Table 5. With respect to system reliability, the decrease is insensitive with the increase of d considering a loose T. However, it is quite sensitive with the increase of d considering a strict T. This phenomenon is clear in Fig. 8. With respect to CPU time, the increase is insensitive with the increase of d considering a loose T. However, it is quite sensitive with the increase of d considering a strict T. This phenomenon is shown in Fig. 9.

6. Conclusion and further research Considering a QPP for an MFN, we call the kind of problem an MQPP. For such an MQPP there are some discussions in the literature [27–29]; however, these assume that the MPs are disjoint from each other. In the real world, we may use relatively few arcs in order to easily manage them. As a result, it is worth discussing that data transmitted through multiple joint MPs. However, a feasible algorithm for an MQPP with multiple joint MPs had never been addressed. In order to deal with the congestion phenomenon, an algorithm based on MCS has been developed. The main contribution of this paper is the development of an algorithm to evaluate the system reliability for an MQPP with multiple joint MPs. Although the system reliability is an approximate value, we do the comparison between approximate system reliability and actual system reliability in Section 5.1 which shows that the error is less than 0.02. Moreover, the results for system reliability and CPU time as influenced by two parameters, d and T, are shown in Table 5, Figs. 8 and 9 in

90

Y.-K. Lin et al. / Simulation Modelling Practice and Theory 29 (2012) 78–92

Fig. 7. A revised network of Fig. 6.

Table 5 The approximate system reliability and CPU time of TANET. Amount of data

Time constraint

9 10 11 12 13 14 15 16

0.9858 0.9628 0.8307 0.8222 0.8264 0.3969 0.2994 0.2953

15

a

16 (6.359)a (8.938) (21.969) (26.484) (29.953) (89.593) (116.047) (130.969)

0.9988 0.9954 0.9897 0.9878 0.9804 0.9463 0.9047 0.8333

17 (4.813) (5.703) (7.5) (9) (11.531) (19.563) (29.531) (43.031)

0.9999 0.9995 0.9986 0.9982 0.9946 0.9951 0.9904 0.9822

(4.047) (4.812) (5.781) (6.828) (8.797) (10.484) (12.547) (15.391)

18

19

0.9999 (3.375) 1 (4.203) 0.9997 (5.094) 0.9998 (6.016) 0.9984 (7.141) 0.9983 (8.344) 0.9984 (10.063) 0.9947 (12.219)

0.9999 0.9997 0.9998 0.9999 0.9998 0.9999 0.9986 0.9978

20 (2.687) (3.61) (4.547) (5.406) (6.406) (7.25) (8.64) (9.859)

1 (1.969) 1 (3.156) 0.9999 (3.875) 1 (4.828) 0.9998 (5.797) 1 (6.625) 0.9996 (7.656) 0.9999 (8.61)

CPU time which is measured in seconds.

Section 5.2. In the future, we may consider more complicated situations such as the loss of data due to a delay, unknown g MPs, the policy of spare routing, multi-commodity allocation, or budget constraint in MQPP.

91

Approximate system reliability

Y.-K. Lin et al. / Simulation Modelling Practice and Theory 29 (2012) 78–92

Amount of data

CPU time

Fig. 8. System reliability influenced by amount of data and time constraint.

Amount of data Fig. 9. CPU time influenced by amount of data and time constraint.

Acknowledgment This work was supported in part by the National Science Council, Taiwan, Republic of China, under Grant No. NSC 992221-E-011-066-MY3. References [1] R.K. Ahuja, Minimum cost-reliability ratio path problem, Computers and Operations Research 15 (1998) 83–89. [2] L.D. Bodin, B.L. Golden, A.A. Assad, M.O. Ball, Routing and scheduling of vehicles and crews: the state of the art, Computers and Operations Research 10 (1983) 63–211. [3] M.L. Fredman, R.E. Tarjan, Fibonacci heaps and their uses in improved network optimization algorithms, Journal of the ACM 34 (1987) 596–615. [4] B.L. Golden, T.L. Magnanti, Deterministic network optimization: a bibliography, Networks 7 (1977) 149–183. [5] Y.L. Chen, Y.H. Chin, The quickest path problem, Computers & Operations Research 17 (1990) 153–161. [6] Y.C. Hung, G.H. Chen, Distributed algorithms for the quickest path problem, Parallel Computing 18 (1992) 823–834. [7] E.D.Q.V. Martins, J.L.E.D. Santos, An algorithm for the quickest path problem, Operations Research Letters 20 (1997) 195–198. [8] G.H. Chen, Y.C. Hung, Algorithms for the constrained quickest path problem and the enumeration of quickest paths, Computers & Operations Research 21 (1994) 113–118. [9] J.B. Rosen, S.Z. Sun, G.L. Xue, Algorithms for the quickest path problem and the quickest paths, Computers & Operations Research 18 (1991) 579–584.

92 [10] [11] [12] [13] [14] [15] [16] [17] [18] [19] [20] [21] [22] [23] [24] [25] [26] [27] [28] [29] [30] [31] [32] [33] [34] [35] [36] [37]

Y.-K. Lin et al. / Simulation Modelling Practice and Theory 29 (2012) 78–92 Y.L. Chen, An algorithm for finding the k quickest paths in a network, Computers & Operations Research 20 (1993) 59–65. Y.L. Chen, Finding the k quickest simple paths in a network, Information Processing Letters 50 (1994) 89–92. G.H. Chen, Y.C. Hung, On the quickest path problem, Information Processing Letters 46 (1993) 125–128. D.T. Lee, E. Papadopoulou, The all-pairs quickest path problem, Information Processing Letters 45 (1993) 261–267. C.C. Jane, J.S. Lin, J. Yuan, On reliability evaluation of a limited-flow network in terms of minimal cutsets, IEEE Transactions on Reliability 42 (1993) 354–361. J.S. Lin, C.C. Jane, J. Yuan, On reliability evaluation of a capacitated-flow network in terms of minimal pathsets, Networks 25 (1995) 131–138. Y.K. Lin, J. Yuan, A new algorithm to generate d-minimal paths in a multistate flow network with noninteger arc capacities, International Journal of Reliability, Quality, and Safety Engineering 5 (1998) 269–285. Y.K. Lin, A simple algorithm for reliability evaluation of a stochastic-flow network with node failure, Computers & Operations Research 28 (2001) 1277–1285. Y.K. Lin, On reliability evaluation of a stochastic-flow network in terms of minimal cuts, Journal of Chinese Institute of Industrial Engineers 18 (2001) 49–54. Y.K. Lin, J. Yuan, Study on the multicommodity transportation reliability by containers under budget constraint and probabilistic capacities, Journal of Chinese Institute of Industrial Engineers 16 (1999) 639–647. J. Xue, On multistate system analysis, IEEE Transactions on Reliability 34 (1985) 329–337. W.C. Yeh, A revised layered-network algorithm to search for all d-Minpaths of a limited-flow acyclic network, IEEE Transactions on Reliability 47 (1998) 436–442. W.C. Yeh, A simple algorithm to search for all d-MPs with unreliable nodes, Reliability Engineering & System Safety 73 (2001) 49–54. J.E. Ramirez-Marquez, D.W. Coit, Multi-state component criticality analysis for reliability improvement in multi-state systems, Reliability Engineering & System Safety 92 (2007) 1608–1619. G. Levitin, A. Lisnianski, A new approach to solving problems of multi-state system reliability optimization, Quality and Reliability Engineering International 17 (2001) 93–104. A.M. Al-Ghanim, A heuristic technique for generating minimal paths and cutsets of a general network, Computers and Industrial Engineering 36 (1999) 45–55. K. Kobayashi, H. Yamamoto, A new algorithm in enumerating all minimal paths in a sparse network, Reliability Engineering & System Safety 65 (1999) 11–15. Y.K. Lin, Extend the quickest path problem to the system reliability evaluation for a stochastic-flow network, Computers & Operations Research 30 (2003) 567–575. Y.K. Lin, Spare routing reliability for a stochastic flow network through two minimal paths under budget constraint, IEEE Transactions on Reliability 59 (2010) 2–10. Y.K. Lin, System reliability of a stochastic-flow network through two minimal paths under time threshold, International Journal of Production Economics 124 (2010) 382–387. G.S. Fishman, A comparison of four Monte Carlo methods for estimating the probability of s-t connectedness, IEEE Transactions on Reliability 35 (1986) 145–155. G.S. Fishman, A Monte Carlo sampling plan for estimating network reliability, Operation Research 34 (1986) 581–594. M.Y. Ata, A convergence criterion for the Monte Carlo estimates, Simulation Modelling Practice and Theory 15 (2007) 237–246. R. Zenklusen, M. Laumanns, High-confidence estimation of small s-t reliabilities in directed acyclic networks, Networks 57 (2011) 376–388. C.H. Lin, W.N. Yang, A simple and efficient importance sampling scheme for stochastic network unreliability estimation, Simulation Modelling Practice and Theory 19 (2011) 924–935. J.E. Ramirez-Marquez, D.W. Coit, A monte-carlo simulation approach for approximating multi-state two-terminal reliability, Reliability Engineering & System Safety 87 (2005) 253–264. L.R. Ford, D.R. Fulkerson, Flows in Networks, Princeton University Press, Princeton, New Jersey, 1962. Y.K. Lin, C.T. Yeh, Using minimal cuts to optimize network reliability for a stochastic computer subject to assignment budget, Computers and Operations Research 38 (2011) 1175–1187.