Computationally efficient and numerically stable reliability bounds for repairable fault-tolerant systems

Computationally Efficient and Numerically Stable Reliability Bounds for Repairable Fault-Tolerant Systems Juan A. Carrasco, Member,IEEE AbstractÐThe transient analysis of large continuous time Markov reliability models of repairable fault-tolerant systems is computationally expensive due to model stiffness. In this paper, we develop and analyze a method to compute bounds for a measure defined on a particular, but quite wide, class of continuous time Markov models, encompassing both exact and bounding continuous time Markov reliability models of fault-tolerant systems. The method is numerically stable and computes the bounds with wellcontrolled and specifiable-in-advance error. Computational effort can be traded off with bounds accuracy. For a class of continuous time Markov models, class C00, including typical failure/repair reliability models with exponential failure and repair time distributions and repair in every state with failed components, the method can yield reasonably tight bounds at a very small computational cost. The method builds upon a recently proposed numerical method for the transient analysis of continuous time Markov models called regenerative randomization. Index TermsÐFault-tolerant systems, repairable systems, reliability, continuous time Markov models, bounds, randomization. æ 1INTRODUCTION INCREASING demand for system dependability has created great interest in fault-tolerant systems. In many applications, e.g., critical applications, an appropriate measure to quantify a system's dependability is the reliability, defined as the probability that the system has not failed by time t, or, alternatively, the complementary unreliability measure, urt, defined as the probability that the system has failed by time t. Homogeneous continuous time Markov chain (CTMC) models are commonly used to predict the unreliability of fault-tolerant systems, particularly when the system is repairable. Computation of the unreliability then requires the transient analysis of the CTMC model. Available numerical methods to perform that transient analysis include ODE (ordinary differential equation) solvers and randomization (also called uniformization) [12], [13], [19]. The randomization method is attractive because it is numerically stable and the computation error is well controlled and can be specified in advance. However, the performance of randomization is seriously affected by model stiffness. For CTMC models, a practical measure of stiffness is t[19], where is the maximum output rate of the model. For large t, randomization requires a number of steps tand will be highly inefficient if the model is large. CTMC reliability models of repairable fault-tolerant systems tend to be very stiff when the mission time of interest is large. To illustrate the point, Fig. 1 shows a small CTMC reliability model XfXt;t0gof a repairable fault-tolerant system using the pair-and-spare technique [9] in which active modules have failure rate M, the spare module does not fail, the failure of an active module is ªsoftº with probability SMand ªhardº with probability 1SM, and, whether soft or hard, the failure of an active module is covered with probability CM. Modules in soft failure are independently recovered at rate Sand modules in hard failure are repaired by a single repairman at rate H. The unreliability of the system is urtPXtf. For the model, 2SM120 h1and, for a mission time t1 year 8;760 h,t1;051;200. Several variants of the (standard) randomization method have been proposed to improve its efficiency: selective randomization [14], [15], multistepping [19, Section 3.1.2], adaptive uniformization [16], adaptive/standard uniformization [17], uniformization with steady-state detection [12], [21], and regenerative randomization [5], [6]. For large CTMC reliability models of repairable fault-tolerant systems and long mission times, regenerative randomization seems to be the best of them. The method has the same good properties as the standard randomization method (numerical stability, well-controlled computation error, and ability to specify the computation error in advance) and can be much faster than standard randomization. The regenerative randomization method covers CTMC models XfXt;t0gwith state space S[ff1;f 2;...;f Ag;jSj2;A0; where fiare absorbing states and either 1) all states in S are transient or 2) Shas a single trapping component 1 and the chosen regenerative state r2Sbelongs to that component, and all states are reachable from some state with nonnull initial probability. It is also assumed that X 254 IEEE TRANSACTIONS ON COMPUTERS, VOL. 51, NO. 3, MARCH 2002 .The author is with the Department d'Enginyeria Electro Ánica, Universitat Polite Ácnica de Catalunya, Diagonal 647, plta. 9, 08028 Barcelona, Spain. E-mail: [email protected]. Manuscript received 1 Apr. 2000; revised 21 Feb. 2001; accepted 7 Apr. 2001. For information on obtaining reprints of this article, please send e-mail to: [email protected], and reference IEEECS Log Number 111155. 1. Two states i,jof a CTMC are strongly connected if there are paths in the state transition diagram of the CTMC from ito jand from jto i; a state is strongly connected with itself; a component is a maximal subset of strongly connected states; a component is trapping if no state of the component has transition rates to states outside the component. 0018-9340/02/$17.00 ß2002 IEEE has some transition rate from rto S0frg, although that condition can be easily circumvented in practice [5]. The generic measure considered in [5] is mtX A i1 rfiPXtfi; A1, where rfiare different reward rates 0(in [6], more general measures are considered and A0is allowed). Those models with A1and the generic measure mt cover both exact and bounding CTMC reliability models of fault-tolerant systems (bounding models are useful when an exact model would have an unmanageable size). In an exact reliability model, Awould be equal to 1, Swould include all operational states, entry in f1would represent the failure of the system, rf1would be equal to 1, the initial probability of f1would be equal to the probability of the system being initially failed, and mtwould be the unreliability urt(an example of such an exact reliability model is the model given in Fig. 1 with Sf1;2;3;4;5;6g and f1f). In a lower bounding reliability model, Awould be equal to 2, Swould be a proper subset of the set of operational states, entry in f1would represent the failure of the system from a state in S, entry in f2would represent entry in an operational state outside S,rf1would be equal to 1, rf2would be equal to 0, the initial probability of f1would be equal to the probability of the system being initially failed, the initial probability of f2would be the probability of the system being initially in an operational state outside S, and mtwould be a lower bound for urt. Finally, in an upper bounding reliability model, Awould be equal to 1, Swould be a proper subset of the set of operational states, entry in f1would represent exit from S, rf1would be equal to 1, the initial probability of f1would be the probability of the system being initially either failed or in an operational state outside S, and mtwould be an upper bound for urt. The regenerative randomization method requires the selection of a regenerative state r2S. The performance of the method depends on that selection. In this paper, we consider CTMC models with the same structure and properties as the models considered in regenerative randomization with A1and develop a method called bounding regenerative randomization to obtain bounds for the generic measure mt. The method yields a lower bound for mt, an upper bound for mt, or both. The lower bound is obtained by solving, by regenerative randomization, a lower bounding CTMC, Xlb. The upper bound is obtained by solving, by regenerative randomization, an upper bounding CTMC, Xub. Both Xlb and Xub are obtained from Xby scaling some of its transition rates. The method has the same good properties as standard randomization. Although not restricted to them, the bounding regenerative randomization method is intended to be used for a class of models C00. Let i;j denote the transition rate of X from state ito state j,i6 j, let iPj2figi;j denote the output rate from state i, and let i;B Pj2Bi;j, Bfig. Class C00 includes the models Xwith the properties assumed in the regenerative randomization method with A1for which there exists a partition S0[ S1[[SNCfor Ssatisfying the following three properties: P1. S0fog(i.e., jS0j1). P2. max0kNCmaxi2Ski;Skfig[Sk1[[SNCis significantly smaller than min 0<kNC min i2Sk i;S0[[Sk1[ff1;...;fAg>0: P3. omini2Sfogi. The class covers failure/repair reliability models with exponential failure and repair time distributions and repair in every state with failed components when failure rates are significantly smaller than repair rates (the typical case), such as the model given in Fig. 1. For those models, a partition for which properties P1, P2, and P3 are satisfied is Sk{states in Swith kfailed components}. The class also CARRASCO: COMPUTATIONALLY EFFICIENT AND NUMERICALLY STABLE RELIABILITY BOUNDS FOR REPAIRABLE FAULT-TOLERANT... 255 Fig. 1. CTMC reliability model of a repairable fault-tolerant system using the pair-and-spare technique. covers failure/repair reliability models with exponential failure time distributions, repair times with acyclic phasetype distributions [18] (which can be used to fit distributions of nonexponential positive random variables [3]), and repair in every state with failed components, provided that the transition rates of the transient CTMCs defining the phase-type distributions are sufficiently large compared with failure rates. For those models, the proposed method can be extremely efficient and, yet, provide quite tight bounds. Tighter bounds can be obtained at the cost of increased computational effort. An approach to deal with stiffness is the aggregation technique proposed in [2]. For class C00 models, that technique could be used to aggregate the states in Sfog, yielding an aggregated CTMC model with a transient state and Aabsorbing states with symbolic solution. The aggregation would be done by replacing each state in Sfogby a switch and is equivalent to scale the transition rates i;j,i2Sfogwith i!1, keeping the relative values of the transition rates from a given state. The aggregated model would give an upper bound for the measure mtlooser than the upper bounds that are computed by the bounding regenerative randomization method proposed in this paper. The rest of the paper is organized as follows: Section 2 presents a brief review of both the standard randomization and the regenerative randomization methods (the latter particularized to the computation of the measure mt), including algorithmic descriptions for both methods. Section 3 describes the proposed bounding regenerative randomization method, proves that it yields bounds for the measure mt, and gives theoretical results assessing the efficiency of the method for class C00 models. Section 4 analyzes the performance of the bounding regenerative randomization method using a large reliability model belonging to class C00 and compares the computational cost of the method with that of regenerative randomization and standard randomization. Finally, Section 5 concludes the paper. 2REVIEW OF STANDARD AND REGENERATIVE RANDOMIZATION The review of the standard randomization method will be made for arbitrary rewarded CTMC models XfXt;t 0gwith finite state space and for the expected transient reward rate measure ETRRtErXtX i2 riPXti; where ri0,i2is the reward rate associated with state i. The quantity rihas the meaning of ªrateº at which reward is earned while Xis in state i. The measure mtis a particular case of ETRRt. The standard randomization method is based on the following result (see, for instance, [10, Theorem 4.19]). Consider any maxi2iand define the homogeneous discrete time Markov chain (DTMC) ^ X f^ Xk;k0;1;2;...gwith same state space and initial probability distribution as Xand transition probabilities Pi;j i;j=,i6 j,Pi;i 1i=. The DTMC ^ Xis called the randomized DTMC of Xwith randomization rate . The CTMC Xis said to be the derandomized CTMC of ^ X with randomization rate . Let QfQt;t0gbe a Poisson process with arrival rate independent of ^ X (PQtkettk=k!). Then, XfXt;t0gis probabilistically identical to f^ XQt;t0g.Thatresult allows expressing ETRRtin terms of the transient regime of ^ Xas: ETRRtX i2 riX 1 k0 P^ XkiPQtk X 1 k0X i2 riP^ Xkiettk k! X 1 k0 dkettk k!;1 with dkPi2riP^ Xki. Let  PX0ii2be the initial probability row vector of Xand let qk P^ Xkii2be the probability row vector of ^ Xat step k. We have q0. From q0,qk,k>0can be obtained using qk1qkP, where PPi;ji;j2is the transition probability matrix of ^ X. An approximate value for ETRRt,ETRRa Nt, can be obtained by truncating series (1): ETRRa NtX N k0 dkettk k!:2 Using dkrmax maxi2ri, the truncation error can be upper bounded as ETRRtETRRa Ntrmax X 1 kN1 ettk k!:3 Then, "being the allowed error for the computation of ETRRt, in the standard randomization method Nis chosen as Nminnm0:rmax X 1 km1 ettk k!"o; and ETRRtis approximated with error "by the ETRRa Ntgiven by (2). The computational cost of standard randomization is essentially the cost of performing the Nvector-matrix multiplications qk1qkP, k0;1;...;N1.Qthas, for t!1, an asymptotic normal distribution with mean and variance t[20], and, for large tand "1, the required Nis t, making standard randomization computationally very expensive if both Xand tare large. Since the performance of standard randomization degrades as increases, is usually taken equal to maxi2i. An algorithmic description of the standard randomization method is given in Fig. 2. The algorithm has as inputs the CTMC X, the reward rates ri,i2, the initial probability row vector ,the allowed error ", the number of time points nat which ETRRthas to be computed, and the time points t1;t 2;...;t n. The algorithm has as outputs the computed 256 IEEE TRANSACTIONS ON COMPUTERS, VOL. 51, NO. 3, MARCH 2002 values of ETRRt,g ETRRt1;g ETRRt2;...;g ETRRtn. The truncation error bound given by (3) increases with t and, therefore, that error is controlled for tmax maxft1;t 2;...;t ng: We review next the regenerative randomization method for the CTMC models Xconsidered in this paper and the measure mt. Let S0Sfrgand let iPX0i, i2. We will use the notation BPi2Bi. In the method, the behavior of Xfrom S0up to state ror a state fiand from r until the next hit of ror a state fiis approximately characterized by a truncated transformed model from which an approximate value with bounded error for mt can be computed and that approximate value is computed solving the truncated transformed model by the standard randomization method. To build the truncated transformed model, two DTMCs, Zand Z0[5], obtained from the randomized DTMC ^ Xof Xwith rate and a version ^ X0of ^ Xin which the initial probability distribution is concentrated in state r, have to be stepped in general. The randomization rate is taken slightly larger than maxi2Si (i.e., 1maxi2Si,being a small value, say 104). This simplifies considerably the description and implementation of the method and has negligible impact on its performance. The transition probability matrix of ^ Xwill be denoted as before by PPi;ji;j2.TheDTMCZ fZk;k0;1;2;...gfollows ^ Xfrom rtill reentry in r.Z has state space S[ff1;f 2;...;f A;ag, where fiand aare absorbing states and all states in Sare transient, initial state r, and its (possibly nonnull) transition probabilities are: PZk1jjZkiPi;j;i2S;j 2S0[ff1;f 2;...;f Ag; PZk1ajZkiPi;r;i2S; PZk1fijZkfiPZk1ajZka1;1iA: The DTMC Z0fZ0 k;k0;1;2;...gfollows ^ Xuntil its first visit to state r.Z0has state space S0[ff1;f 2;...;f A;ag, where fiand aare absorbing states and all states in S0are transient. The initial probability distribution of Z0is PZ0 0ii,i2S0[ff1;f 2;...;f Ag,PZ0 0ar, and its (possibly nonnull) transition probabilities are: PZ0 k1jjZ0 kiPi;j;i2S0;j2S0[ff1;f 2;...;f Ag; PZ0 k1ajZ0 kiPi;r;i2S0; PZ0 k1fijZ0 kfiPZ0 k1ajZ0 ka1;1iA: Let ikPZki,0 ikPZ0 kiand consider the row vectors kiki2Sand 0k0 iki2S0. Let PZ be the transition probability matrix of Zrestricted to S and let PZ0be the transition probability matrix of Z0 restricted to S0. From 0,k,k>0can be obtained using k1kPZ. From 00,0k,k>0can be obtained using 0k10kPZ0. Let akX i2S ik; vj kX i2S ikPi;fj=ak; qkX i2S ikPi;r=ak; wkX i2S ikPi;S0=ak; and, if S0>0, let a0kX i2S0 0 ik; v0j kX i2S0 0 ikPi;fj=a0k; q0 kX i2S0 0 ikPi;r=a0k; w0 kX i2S0 ikPi;S0=a0k; where Pi;S0Pj2S0Pi;j (akand, if S0>0,a0kare guaranteed to be >0(see [5])). Then, for the case S0>0, CARRASCO: COMPUTATIONALLY EFFICIENT AND NUMERICALLY STABLE RELIABILITY BOUNDS FOR REPAIRABLE FAULT-TOLERANT... 257 Fig. 2. Algorithmic description of standard randomization. the truncated transformed model is the CTMC VK;L  fVK;Lt;t0gwith state space fsk;0kKg[fs0 k;0kLg[ff1;f 2;...fA;ag; initial probability distribution PVK;L0s0r; PVK;L0s0 0S0; PVK;L0fifi; PVK;L0i0;i62fs0;s 0 0;f 1;f 2;...;f Ag; and the state transition diagram illustrated in Fig. 3 for A1. For the case S00, the truncated transformed model is the CTMC VKfVKt;t0gwith initial probability distribution PVK0s0S,PVK0fifi, PVK0i0,i62fs0;f 1;f 2;...;f Agand a state transition diagram identical to the state transition diagram of VK;L, but without states s0 k. For the case S0>0, the approximate value for mt given by VK;L is ma K;LtX A i1 rfiPVK;Ltfi and we have mtma K;Ltrmaxa0LX 1 kL1 ettk k! rmaxSaKX 1 kK1 kKettk k!; 4 where rmax max1iArfi. For the case S00, the approximate value for mtgiven by VKis: ma KtX A i1 rfiPVKtfi and we have mtma KtrmaxSaKX 1 kK1 kKettk k!:5 The model truncation error bounds given by (4) and (5) decrease for increasing Kand Land can be made arbritrarily small by choosing large enough values for K and L. In regenerative randomization, "being the allowed error for the computation of mt, suitable truncation parameters K,Lare chosen so that the model truncation error bounds are smaller than "=2and, then, an approximate value for mtis obtained by computing ma K;Lt (ma Kt) by solving the truncated transformed model VK;L (VK) by standard randomization with error upper bounded by "=2. An algorithmic description of the regenerative randomization method is given in Fig. 4, where Icdenotes the indicator function returning the value 1 if condition cis satisfied and the value 0 otherwise. The algorithm has as inputs the CTMC X, the number Aof absorbing states fi, the reward rates rf1;r f2;...;r fA, an initial probability distribution vector  ii2with S>0, the regenerative state r, the allowed error ", the number of time points nat which mthas to be computed, and the time points t1;t 2;...;t n. The algorithm has as outputs the computed values of mt,e mt1;e mt2;...;e mtn. Since the model truncation error bounds increase with t, they are controlled for tmax maxft1;t 2;...;t ng. For the case S0>0, the "=2 allocated for the model truncation error bound is divided equally between its two contributions. The truncation error bound associated with the solution of the truncated transformed model by standard randomization also increases with tand that error is controlled for tmax. The method requires stepping the randomized DTMC ^ VK;L (^ VK) of VK;L (VK) with rate . The state transition diagram of ^ VK;L is illustrated in Fig. 5 for A1. The state transition diagram of ^ VKis identical, but without the states s0 k. The regenerative randomization method (as the standard randomization method) requires the computation of the Poisson probabilities ettk=k!. Stable and efficient computation of those Poisson probabilities, avoiding overflows and intermediate underflows, is a delicate issue and 258 IEEE TRANSACTIONS ON COMPUTERS, VOL. 51, NO. 3, MARCH 2002 Fig. 3. State transition diagram of the CTMC VK;L for A1. several alternatives have been proposed [4], [8], [11], [17]. The method described in [11, pp. 1028±1029] (see also [1]) has good numerical stability and is the one we follow in our implementations. The regenerative randomization method involves the computation of Sm X 1 km1 etmax tmaxk=k! and S0m X 1 km1 kmetmax tmaxk=k! for increasing values of m(the standard randomization method also requires the computation of Smfor increasing values of m). Our implementations use the algorithms described in [5], which are numerically stable and efficient. The computational cost of regenerative randomization has two components: cost associated with the construction of the truncated transformed model and cost associated with the solution of the truncated transformed model by standard randomization. The first is roughly proportional to the number of steps on the DTMCs Z,Z0,KLif S0>0 and Kif S00, with a cost per step which, for large X, will typically be slightly larger than the cost per step in standard randomization. The second component is roughly proportional to the truncation parameter N(approximately equal to the truncation parameter Nof standard randomization) and to the size of the truncated transformed model. It is shown in [5] that the required Kis Ologt=" and, if S0>0,therequiredLis Olog1=".Thatiscalled ªbenignº behavior and implies that, for large enough X and large enough t, regenerative randomization will be significantly faster than standard randomization. The performance of regenerative randomization depends, of course, on the selection of the regenerative state r. That selection should be made so that akand a0kdecrease as fast as possible and the required Kand Lare as small as possible. Since class C00 models move fast to either state oor an absorbing state fi, a natural selection for those models is ro. Let R0maxi2Sfogi=mini2Sfogi. Then, we can state the following result: Theorem 1. For class C00, models with selection ro,ak hkand a0kS0h0k, where, for k!1, CARRASCO: COMPUTATIONALLY EFFICIENT AND NUMERICALLY STABLE RELIABILITY BOUNDS FOR REPAIRABLE FAULT-TOLERANT... 259 Fig. 4. Algorithmic description of regenerative randomization. hkBk p1  k; and h0kB0k p01  0k; with B>0,B0>0,p,p0integers 1,11=R0, and 011=R0. 2 Proof. The model class C00 is a subset of the model class C considered in [5]. In [5], it is considered the parameter Rmaxi2Si=mini2S0i. For class C00 models, it follows from Property P3 that, with selection ro,RR0and the result follows from Theorem 4 of [5] and the discussion following it. tu Theorem 1 asserts that, for class C00 models, the performance of regenerative randomization with the natural selection roshould be mainly determined by the parameter R0: the larger R0, the more costly the method. In particular, for R01,0and 00and the method should be very efficient. Those observations motivate the bounding regenerative randomization method. 3THE BOUNDING REGENERATIVE RANDOMIZATION METHOD The bounding regenerative randomization method obtains a lower bound for mt, an upper bound for mt, or both. The bounds are computed with an error upper bounded by a"given by the user. Depending on the nature of the CTMC model X, one or the other bound or both bounds could be of interest. Thus, if Xis an exact reliability model, both bounds would be of interest to have an assessment of the error on urt. However, if an exact reliability model cannot be used because its size would be unmanageable, then we could use a lower bounding reliability model and an upper bounding reliability model and use bounding regenerative randomization to compute a lower bound for the lower bound for the unreliability given by the first model and an upper bound for the upper bound for the unreliability given by the second model: the exact unreliability would be bracketed by those values. The bounding regenerative randomization method requires the selection of a regenerative state r2Sand has an input parameter Dcontrolling the accuracy of the bounds. Let min mini2S0iand max maxi2S0i. The method assumes that the controlling parameter Dis restricted by 1D< max=min. 3 To obtain the lower bound for mt, the method modifies the CTMC Xto obtain a CTMC Xlb. The CTMC Xlb is obtained from Xby scaling the transition rates from states in S0so that, calling lb ithe output rates of Xlb,lb ii,i2S0and maxi2S0lb i=mini2S0lb iD. That scaling is defined by lb i;j i;jlb i=i,lb iminfi;D ming,i2S0,wherelb i;j are the transition rates in Xlb. The lower bound for mtis given by mlbtX A i1 rfiPXlbtfi: That lower bound is obtained by solving Xlb by regenerative randomization with regenerative state r. To obtain the upper bound for mt, the method modifies the CTMC Xto obtain a CTMC Xub. The CTMC Xub is obtained from Xby scaling the transition rates from states in S0so that, calling ub i,i2S0the output rates of Xub,ub ii,i2S0and maxi2S0ub i=mini2S0ub iD. That scaling is defined by ub i;j i;jub i=i,ub imaxfi; max=Dg,i2S0. The upper bound for mtis given by mubtX A i1 rfiPXubtfi: That upper bound is obtained by solving Xub by regenerative randomization with regenerative state r. The particular case in which both bounds are to be computed, D1and min rallows a more efficient implementation of the bounding regenerative randomization method than that described in the previous paragraph. To justify that particular implementation, we will use the following result: 260 IEEE TRANSACTIONS ON COMPUTERS, VOL. 51, NO. 3, MARCH 2002 Fig. 5. State transition diagram of the DTMC ^ VK;L for A1. 2. ckdkfor k!1denotes limk!1 ck=dk1. 3. For class C00 models with the selection ro, when max min and no selection for Dis possible, regenerative randomization should be very efficient because of Theorem 1 and the fact that R01, obviating the need for the bounding regenerative randomization method. Lemma 1. For x>0,K0, and R>1, 1 RX 1 kK1 kKeRx Rxk k!>X 1 kK1 kKexxk k!: Proof. See the Appendix. tu For that particular case, denoting by superscripts lb and ub the terms referred to, respectively, Xlb and Xub and the objects involved in their solution by regenerative randomization, and letting R00 max=min >1(because D1and D< max=min), we have lb r;j ub r;j and, noting that lb i min and ub imax,i2S0, we have lb i;j ub i;j =R00,i2S0.We also have lb ub=R00 (because lb rub rrmin, lb imin,i2S0,ub imax,i2S0, and max > min). It follows that the transition probabilities of the randomized DTMCs ^ Xlb and ^ Xub are related as Plb r;j R00Pub r;j ,j6 r, Plb r;r 1R001Pub r;r and, for i2S0,Plb i;j Pub i;j . Then, taking into account that the transition probabilities of the DTMCs Zlb and Zub from rto j2S0are, respectively, Plb r;j and Pub r;j , that the transition probabilities of Zlb and Zub within S0are, respectively, Plb i;j and Pub i;j and that the initial state of Zlb and Zub is r, we have, for k1,lb ikR00ub ikand albkR00aubk. Then, using alb0aub01, we have vjlb 0Plb r;fjR00Pub r;fjR00vjub 0; qlb 0Plb r;r 1R001Pub r;r 1R001qub 0; wlb 0Plb r;S0R00Pub r;S0R00wub 0; and, taking into account that lb rkub rk0,k1, for k1, we have vjlb kPi2S0lb ikPlb i;fj albkPi2S0R00ub ikPub i;fj R00aubk Pi2S0ub ikPub i;fj aubkvjub k and, similarly, qlb kqub kand wlb kwub k. On the other hand, the DTMCs Z0lb and Z0ub have identical initial probability distributions and transition probabilities and, then, a0lbka0ubk,and,usingPlb i;j Pub i;j ,i2S0, v0jlb kv0jub k,q0lb kq0ub k,w0lb kw0ub k. For m1, rmaxSaubmX 1 km1 kmeubtmax ubtmaxk=k! rmaxSalbm=R00X 1 km1 kmeR00lbtmax R00lbtmaxk=k! >r maxSalbmX 1 km1 kmelbtmax lbtmaxk=k!; by Lemma 1 with Km,xlbtmax and RR00, implying Kub Klb.IfS0>0, for m1,ub >lb implies 4 rmaxa0ubmX 1 km1 eubtmax ubtmaxk=k! rmaxa0lbmX 1 km1 eubtmax ubtmaxk=k! >r maxa0lbmX 1 km1 elbtmax lbtmaxk=k!; implying Lub Llb. Then, if we start by computing the upper bounds mubtusing Xub and save aubk,vjub k,qub k,wub k, and, if S0>0,a0ubk,v0jub k,q0ub k,w0ub k,wecanusethe relationships between those parameters and the corresponding parameters for Xlb to avoid stepping Zlb and Z0lb when computing the lower bounds mlbt. An algorithmic description of the bounding regenerative randomization method, including the previously discussed particular implementation, is given in Fig. 6. The algorithm has as inputs the CTMC X, the number Aof absorbing states fi, the reward rates rf1;r f2;...;r fA,aninitial probability distribution vector  ii2with S>0, parameters lb and ub indicating, respectively, whether the lower and upper bounds for mtare desired or not, the regenerative state r, the controlling parameter D,the allowed error ", the number of time points nat which mlbt,mubthave to be computed, and the time points t1;t 2;...;t n. The algorithm has as outputs the computed values of mlbt,e mlbt1;e mlbt2;...;e mlbtnand of mubt, e mubt1;e mubt2;...;e mubtn. The algorithmic description makes reference to DTMCs ^ VK;L (S0>0) and ^ VK(S00). Those DTMCs are the randomized DTMCs with randomization rate lb of the truncated transformed models of Xlb used in the solution of Xlb by regenerative randomization. For the case S0>0,^ VK;L has state space fsk;0kKg[fs0 k;0kLg[ff1;f 2;...;f A;ag; initial probability distribution P ^ VK;L0s0r; P ^ VK;L0s0 0S0; P ^ VK;L0fifi; P ^ VK;L0i0;i62fs0;s 0 0;f 1;f 2;...;f Ag and the state transition diagram illustrated in Fig. 5 for the case A1, where K,L,vj k,qk,wk,v0j k,q0 k, and w0 khave the values computed in the algorithm. For the case S00,^ VK is the DTMC with initial probability distribution P ^ VK0s0S; P ^ VK0fifi; P ^ VK0i0;i62fs0;f 1;f 2;...;f Ag CARRASCO: COMPUTATIONALLY EFFICIENT AND NUMERICALLY STABLE RELIABILITY BOUNDS FOR REPAIRABLE FAULT-TOLERANT... 261 4. The inequality comes from the fact that P1 km1ettk=k!is the probability that the number of arrivals in the interval 0;tin a Poisson process with arrival rate is m1, which is increasing with t, and, therefore, with . and state transition diagram identical to the state transition diagram of ^ VK;L, but without states s0 k. In the following, we prove the correctness of the method, i.e., mlbtmtmubt. To that end, we consider the embedded DTMC of X,fk;k0;1;2;...g:has same state space and initial probability distribution as X and transition probabilities i;j i;j=i,i2S,j2fig, i;i 0,i2S, fi;fi1,1iA, fi;j 0,j6 fi. The behavior of Xcan be described in terms of the DTMC by saying that the sequence of states visited by Xis given by  with sojourn times in each state iof Xexponentially distributed with parameter i, independently on the path followed by . That interpretation is usually referred to as the ªstructureº of X[7, Section 8.3]. Since Xlb and Xub have 262 IEEE TRANSACTIONS ON COMPUTERS, VOL. 51, NO. 3, MARCH 2002 Fig. 6. Algorithmic description of bounding regenerative randomization.