Model of the leaky bucket ATM generic flow control mechanism: a case study on solving large cyclic models

Model of the leaky bucket ATM generic flow control mechanism: a case study on solving large cyclic models J.A.Carrasco, V.Sutie, S.Mahevas and G.Rubino Abstract: The authors describe and solve a Markov model of the leaky bucket ATM generic flow control mechanism. The model has a space cardinality which grows quickly with its parameters and is challenging to solve. Exploiting the cyclic nature of the model, the authors develop a methodology which allows them to efficiently solve instances of the model with 3905 134 states and 53869532 transitions using 29.8Mbyte of memory and 222Mbyte of disc storage. The CPU utilisation is high (between 70% and 90%). The methodology is new and can be easily extended to any kind of finite cyclic Markov models. 1 Introduction Communication networks based upon the asynchronous transfer mode (ATM) provide high performance user/ network interfaces [l]. The high bit rates made available by ATM allow several ATM connections to share a common medium. Examples of such media are a bus, a dual bus, or a ring. Generic flow control (GFC) mechanisms are needed to arbitrate the access to the medium. Several GFC mechanisms/protocols have heen proposed [24. A simple and popular one is the leaky bucket method [5]. In [6] the performances of several GFC protocols are evaluated by simulation. In this paper we evaluate the performances of the leaky bucket protocol by numerically solving a discretetime Markov chain model. For typical model parameters, the resulting Markov chains are large. In addition the characterisation of the irreducible closed set in which the steady-state regime is established is digcult. That characterisation is necessary to solve the model efficiently. In this paper we first characterise the irreducible closed set of the Markov chain model. Then, we develop an ‘on-the-fly’ model generation methodology and an associated model solution methodology, which exploit the cyclic structure of the model to reduce memory requirements to a minimum. This allows us to solve instances of the model with 3 905 134 states and 73 869 532 transitions using 29.8 Mbyte of memory and 222Mbyte of disc storage with a CPU utilisation between 70”/0 and 90% depending on model parameters. Our model generation and model solution methodologies can be used to solve any class of finite cyclic Markov models. 0 IEE, 2001 ZEE Proceeditgs online no. 20010285 DO1 10.1049/ipcom:20010285 Paper first received 4th November 1999 and in revised form 13th December 2000 J.A. Carram and V. Sulk are with the Lkpartament d‘Enginyeria Elactrbnica, Universitat Politknica de Catalunya, Diagonal 647, plta. 9,08028 Barcelona, Spain S. Mahevas is with FREMER, MAERHA Laboratory, Rue de l’ile &Yen, BP 21105,44311 Nantes Cedex 03, France G. Rubino is with IRISA, Campus de Beaulieu, 35042 Rennes Cedex, France 2 Model description The purpose of the leaky bucket GFC protocol is to regulate traffic by limiting throughput and cell clumping. The leaky bucket protocol has a parameter K equal to the maximum number of credits given to local cells. At the beginning, the number of credits is equal to K. The number of credits is decremented each time a local cell is put in the medium. The number of credits is incremented up to Kat a rate l/A (once every A slots). A local cell can only be put in a free time slot (a slot not containing a cell from some other station upstream) if there are credits. If there is a waiting local cell when a free slot is received and there are not credits available, that local cell is either lost, or contained and it will wait until a credit is available, or it is accepted but marked, meaning that it may be handled differently from unmarked cells in the network. We will focus on the first case. An appropriate selection of K then allows a balance between local cell loss probability which decreases with increasing K, and cell clamping, which increases with increasing K. To reduce clamping, one should use the smallest K giving a cell loss probability smaller than or equal to the required value (for instance, Fig. 1 shows a model of the leaky bucket protocol (proposed to the authors by Guillemin and Dupuis [7]). A first queue with infinite capacity and slotted deterministic service time 1 models contention for the medium. The queue has two incoming streams of cells: a global stream G modelling the traffic generated by the upstream stations and a local stream L modelling the traffic generated by the station under consideration. The global stream is ‘discrete Poisson’ with load parameter p (the probability that the number of arriving cells during a time slot x is P(x) = @”/ x!)e-p). The local stream is assumed to be deterministic with interarrival time equal to D slots. This corresponds to assuming that the station under consideration has negotiated a traffic with peak rate 1/D and putting ourselves in the worst case scenario in which the station under consideration generates a cell to be put in the medium every D slots. Cells from the global stream (called global) have priority over cells from the local stream (called local). This models the fact that local cells can be put in the medium only in free time slots. A cell from the global stream 10-6). IEE Proc -Commun., Vol 148. No 3, June 2001 188 departing from the first queue models the fact that the current slot is filled with a cell coming from the upstream stations. A local cell departing from the upstream models the fact that the current slot is empty and an existing local cell in the buffer of the station under consideration can be put in the medium. However, whether the local cell will be put or not in the medium is determined by the leaky bucket protocol. The leaky bucket protocol is modelled by the second queue. The number of credits available at a given time in the leaky bucket protocol is K minus the number of cells in that queue. The queue has a server with deterministic service time A, modelling the fact that the number of credits is incremented up to K once every A slots. Local cells departing from the first queue are diverted to the second queue. If the second queue is not full, this means that there are credits available and the local cell will be put in the medium. The resulting decrease by one in the number of credits is modelled by putting the local cell in the second queue. If the second queue is full there are not credits available and the local cell will be thrown away. Thus, cell losses in the local stream diverted to the second queue due to that queue being full are cell losses in the local stream caused by the leaky bucket protocol, and local cells entering the second queue correspond to accepted local cells. The parameter A establishes an upper bound l/A for the average local traffic put into the medium and, therefore, l/A should be as close as possible to the negotiated peak rate 1/D. In addition, to have a small cell loss probability in the local stream, A should be strictly less than D. This leads to a typical selection A = D - 1. K . G - L Fig. 1 Queue model ojlenky bucket bused GFCprotocol The behaviour of the model can be described by a discrete time Markov chain (DTMC), whose generic state is the state of the queueing network defined at the beginning Table 1: Some generation rules of II, of time slots. That state can be described by the vector (ATG, NL, d, k, s), where NG and NL are the number of cells of, respectively, the global stream and the local stream in the first queue, d is the counter of the local stream source (1 I d 5 D and a cell arrives when d = D), k is the number of cells in the second queue, and S is the counter of the service time of the second queue (6 = 0 when k = 0, 1 I 6 5 A when k > 0, and a cell leaves the queue when 6 = A). Since the first queue has infinite capacity, the model described previously has infinitely many states. This makes an exact solution of the model extremely difficult, if not impossible (solution techniques for models with infinite state spaces require a certain regularity in the state space and the model lacks that structure). Therefore, we decided to give to the model an approximate solution obtained by solving for increasing values of C, models with the first queue of finite capacity C, monitoring the convergence of the solution, and stopping when some appropriate convergence condition is satisfied. We decided to start with C = 4 and increase C by increments of two. We found this an appropriate tradeoff between reducing the number of models which are solved and minimising the size of the largest model solved. Let IIc be the DTMC modelling the queueing network with a first queue of finite capacity C. To generate TIc we use a function which for a given state description (NG, NL, d, k, s) of a state s gives the state descriptions of the successors of s and the associated transition probabilities. This function is eficiently implemented using 28 generation rules. Each generation rule has a precondition on the state variables and specifies the successors and transition probabilities which apply when the precondition is satisfied. As an illustration, we give several generation rules in Table 1. Transitions with probabilities smaller than are discarded. 3 Model generation and solution The measure we want to compute is the cell loss probability. We took three criteria for convergence with respect to C. Let E be a small relative tolerance parameter. The first such criterion is to have a relative difference between the cell loss probabilities computed in successive iterations smaller than or equal to E; the second criterion is to have a relative error for the probability that the server of the first queue is serving a local cell (which is equal to 1/D) smaller than or equal to E; the last criterion is to have a relative error for the average number of global cells in the first queue at the beginning of time slots, which is equal to p(2 - p)/(2( 1 - p)) as it easily follows using well-known results for Rule Precondition Transition probability Successors 1 NG=OANL=OA~#D 7 NG=OA NL=OA~#D Ak=O A k= KA 6#A 13 NG>OAd#DAk=l AS=A 18 NG>OAd=D/\k>l A~=A 20 NG=OA NL>OA~= D A k>OA k<KA6#A (x, 0, d, + 1,0,0) (C, 0, d+ 1,0,0) (x, NL1, d+ 1, K, S+ 1) (CNL, NL1, d+ 1, K, 6+ 1) (NG+x-~, N~,d+1,0,0) (CNL1, NL, d+ 1,0,0) (N~+~-l,N~+l,l,k-l,l) (cNL1, NL, 1, k1, 1) (x, NL, 1, k+ 1,6+ 1) (CNL, NL1,1, k+ 1,6+ 1) -- fix) = (.d/k!)e-p is the probability that xglobal cells arrive in a time slot Q(y) = X&, fiz) is the probability that yor more global cells arrive in a time slot The integer variable x in f(x) extends from 0 to y1 (for instance, for rule 1 we have y= C) IEE Proc.-Comwzun.. Vol. 148. No. 3, June 2001 189 the M/G/l queue [8], smaller than or equal to E. The last two criteria were added to increase confidence in the procedure, since we are approximating an infinite model with a finite one. In the following we will describe the methodology used to generate and solve the models l7,. The space of feasible states of I& is SZ, = {(NG, NL, d, k, 0 if k = 0, 1 5 6 1. A if k > O}. Detailed analysis of nc reveals that some of the states in QC are transient. Intuitively, it is clear that Hc has a single irreducible closed set (we use the theory of the classification of states of finite DTMCs as presented in [9]), since otherwise the system would not show steady-state behaviour independent of its initial state. However, we should prove that formally and, furthermore, characterise the irreducible closed set. Given the structure of H, an explicit characterisation of its irreducible closed set is extremely difficult. However, the following theorem characterises it implicitly: Theorem 1: Assume D > 1 and A < D. Let Sc be the subset of SZ, including the state U = (0, 1, 1, 0, 0) and the states reachable in llc from U. Then, Sc is the only irreducible closed set of n,. Prooj See the Appendix. In the following we will consider llc restricted to Se Then, theorem 1 allows the restricted DTMC nc to be obtained easily by simply generating the model from state U. The DTMC IIc is periodic with period D. This follows easily by considering that all successors of states (NG, NL, d, k, 4, 1 I d < D are of the form (VG, VL, d + 1, k’, 8) and that all successors of states (NG, NL, D, k, 8) are of the form (Nk, Ni, 1, k’, 8). Let Sd denote the subset of Sc including the states with d = i + 1. Note that the successors of S& belong to Sd” mod D. For values of interest for C, K, D and A, the number of states of l7, can be extremely large (of the order of tens of millions). To reduce memory and disc storage requirements to a minimum, we decided to use ‘on-the-fly’ model generation and model solution techniques [lo]. Using ‘on-the-fly’ techniques, the successors of a given state and the corresponding transition probabilities are obtained dynamically as required by the model generation and model solution algorithms, thus avoiding the storage in memory or disc of the transition probability matrix of the model. A standard way of reducing memory requirements when dealing with large models is to use keys instead of state descriptions. To support the use of keys, two functions are required an encoding function computing the key from a state description and a decoding function obtaining the state description from a key. Since 0 I NG 2 C, 0 I NL 5 C, 1 I d 5 D, 0 I k I K and 0 I 6 I A, a suitable encoding function is: b( NG, NL , d, IC, 6) = NG(C+ 1) (D+ 1) (K+ 1) (A + 1) 4, NG 2 0, NI 2 0, NG + NL I C, 1 I d I D, 0 1.k I K, 6= + NL(D + 1)(K + l)(A + 1) + d(K+l)(A+l)+k(A+l)+S The decoding function can be implemented by: S = b mod(A + 1), bl = Lb/(A + l)], k = bl mod(K + 1)) bz = [bl/(K + l)], d = b2 mod(D + 1)) b3 = Lb2/(D + 1)J) NL = b3 mod(C + l), NG = lbS/(C + 1)J. We have represented keys using 32-bit-long unsigned integer variables. This sets a maximum value for a key of 232 - 1 = 4294967 295, which has been enough in all the instances of l7, we have tried. 190 Model generation requires a procedure of test of existence and insertion of a state with a given key into a subset of already generated states. If the key exists, the procedure should return the index of the state. If the key does not exist, the key should be inserted with a state index one unit greater than the maximum index in the subset. Holding the keys and the associated state indices in a 2-3 tree [Ill indexed by keys, the operation can be done in O(1og n) time, where n is the number of states held in the subset. We are now ready to present the model generation algorithm. The purpose of the algorithm is to identify Sc. According to theorem 1, Sc can be obtained by generating all its states starting from the state U = (0, 1, 1, 0, 0). The memory requirements of ‘on-the-fly’ techniques are typically vectors or data structures (e.g. search trees) with size equal to the cardinality of the generated state space. However, the cyclic nature of the model allows us to translate those requirements to disc and reduce the memory requirements to the components associated with two subsets s&. This is achieved by generating the state space following a breadth-first approach. We start by state U which belongs to Sco and generate all its successors, which will belong to Sc‘. Then, we expand the successors of U, obtaining states in Sc2, and so on. It is clear that the generation process can be organised in steps. At each step we expand the generated and unexpanded states of a subset SCt and obtain states of S,l+l mod which, if new, are added. The generation process can be finished when a complete cyclic sequence of subsets S,l has been visited without expanding any state. If data structures are partitioned according to the subsets S,l, only the data structures associated with the subset S& from which states are expanded and the subset S,l+‘ mod have to be kept in memory, holding in disc the data structures associated with the other subsets. We give in Fig. 2 a description of the generation algorithm. We use sequential files file-i, 0 2 i 2 D - 1, where file-i holds the data structures associated with a subset Sd. These data structures include the number of generated states, the last expanded state, an array of keys giving for each state index the corresponding key, and a 2-3 tree holding the keys of the generated states and, for each key, the corresponding state index. The function look-and-insert(key-tree, b, &n) looks for the key b in the search tree key-tree, returning YES and, in n, the index of the state, if the key exists, and returning NO and inserting the key b into the tree with index one unit greater than the maximum index in the tree, if the key does not exist. The variable previous-done is set to YES when no state has been expanded in a sequence of cyclically consecutive subsets S,l ending in the previously visited subset. When previous-done is equal to YES, the variablefirst-done is equal to the index i of the subset Scz starting the cyclic sequence of ‘done’ subsets. We next describe the solution procedure. Since llc is periodic, it does not have a steady-state probability distribution. However (see, for instance [9]) it has an invariant measure Y = (vJzESc, where vz can be interpreted as the long term average frequency of visits to state i. The invariant measure v is the only normalised vector (vTl = 1, where 1 is a vector with all its components equal to 1) satisfying the linear system: where Q is the transition probability matrix of nC. Defining A = Q - Z, eqn. 1 can be written as where 0 is a vector of the appropriate dimension with all its elements equal to 0. The cell loss probability can be vTQ = uT (1) vTA = OT (2) IEE Proc -Commtin ~ Vol 148, No 3, June 2001 n-states = 1; last-exp = 0; key[l] = encode(D,C',K,n,(O,l,1,0,0)); make empty tree key-tree; answer = look-andinsert(key-tree, 6, Scn); write nutates, /ast-exp, keyu and key-tree in file file-0; for (i =: 1 ; d 5 D - 1; i++) { n-states = 0; last-exp = 0; ma,ke empty array key[]; make empty tree key-tree; write n-states, last-exp, key1 and key-tree in file .tile-i; 1 e9id= NO; previous-done = NO; f = 0; while (!end){ read n-statesf and last-exp-f from file filef; if (East-exp-f == n-states-j) { if (!preoiozls-done) { prmvious-done = YES; ,fir&done = 1; 1 else if (firstdone == f + 1 mod D) end = YES; I else previous-done = NO; if (!end){ reid key-jj and key-tree-f from file filef; t= f-i-1rnodL); read n-states-t, last-exp-t, key-tu and key-tree-t from file J1e-l; for (i = last-exp-f + 1; i 5 n-states-J is+){ decode(D,C,ZC,A,ke~-~~], &(NG, NL, d, h, 6)); obtain the set S of descriptions of the SuccessorN of (NG, NL, d, IC, 6); for (each (Nc:, NL, d, k, 6) E S){ 6 = encode(D,C,K,A,(NG, NL,d, k, 6)); if (!look_andinsert(key-tree-t, b, &n)){ n-stutes-t++; /* new state */ key-t[ 12-states-t] = 6; 1 1 1 East-eq-f = ra-states-f; wrik n-states-f, East-exp-f, kcy-A] and key-tree-f in file file-f; write a-stakx.t, last-exp-t, key$ and key-tree-t in file file& f = t; 1 1 Fig. 2 Model generation algorithm obtained as lp = lrl(llD), where lr is the cell loss rate. The quantity lr can be computed by adding the vis of the states in which a cell is lost. These states are those with NG = 0, NL > 0 and k = K. The other quantities to be computed are P, the probability that the server of the first queue is serving a local cell, and B, the average number of global cells in the first queue at the beginning of a time slot. They can also be easily obtained from v and the corresponding state descriptions: P can be obtained by adding the vis of the states with NG = 0 and NL > 0; B can be computed by adding the products of the vis by NG. Sorting the states according to the succession Sco, Sc' , . . ., ScD-', we obtain for A the block structure: -10 Qo,~ 0 e.. 0 -11 Qi,2 a.. A = [ 0 0 -12 ... 0 ... Qo-l.0 0 0 ... ID-^ where 4 denote identity matrices of appropriate dimensions and Qij is the restriction of Q to pairs in Sd x 5'2. Such an ordering guarantees the convergence of the Gauss-Seidel method ([12], Section 7.4.2). SOR is an iterative method which is usually used to speed up Gauss-Seidel iteration. For matrices with the block structure of A, under the condition that all eigenvalues of JD, where J is the Jacobi iteration matrix, are real and non-negative, a theory is available ([ 121, Section 7.7.1) providing the optimum relaxation parameter value. Computation of the optimum value for the relaxation parameter requires the knowledge of the convergence factor for Jacobi p(J), which, using the known relationship ([12], Section 7.7.1) between the eigenvalues of the Jacobi and Gauss-Seidel iteration matrices, can be computed from the convergence factor of the Gauss-Seidel, p(G), which can be estimated using: IEE Proc.-Commun., Vol. 148, No. 3, June 2001 191 We performed numerical experiments and found that: (i) stabilkation of the estimate for p(C) takes, for large models, a very substantial portion of the iterations needed for the convergence under the Gauss-Seidel method (about 20%), (ii) the computed os are far from the optimum ones (due to the fact that JD has negative or complex eigenvalues). Thus, the conditions behind the SOR optimisation theory do not seem to hold for the matrices of our models and, given the large number of iterations required for the stabilkation of the estimates of the convergence factors, we decided that trying to develop a sophisticated adaptive SOR optimisation procedure such as that proposed in [13] was not worthwhile. Thus, we decided to use the GaussSeidel method. Denoting the elements of Q by qd, the Gauss-Seidel iterative step applied to the solution of eqn. 2 can be described as: (k+l) = c4 ja pfl) 3 +Cgjavik), 'a 3 <a 9>% 2=1,2,..,,JScJ (3) The trivial implementation of eqn. 3 requires access to the elements of Q in a column-by-column manner. This is inconvenient, since it requires functions running backwards in the model and these functions are difficult to obtain. To solve the problem, an algorithm has been proposed in [lo] that, at the cost of requiring two vectors instead of one, implements the Gauss-Seidel iteration with access to the elements of the matrix Q in a row by row manner. However, the cyclic nature of matrix Q allows a natural implementation of the Gauss-Seidel iteration in whch the elements of the matrix are also accessed in a row by row manner. To see that, it is enough to note that all predecessors of S& are in Sd-' mod and that this subset has only successors in Scl. Then, to perform, for instance, the iterative step on the states of S$, it is enough to generate successors of statesj E SCS1 and, for each successor i E SCo, to cumulate in v, the contribution q119. This observation leads to the algorithm given in Fig. 3. First, the files file-i left by the algorithm given in Fig. 2 are read, the iteration vectors initialised with all components equal to 1 and the files written with the information required by the soluold-lr = 0; old-P = 0; old-B = 0; for (i = 0; i 5 II - 1; if+){ read n-states-f, last..exp-,f, key-fn aid key-tree-f from file Jlc-i; for (j = 1; j 5 n-stutes; j++){ U&'] = 1: decode(D,C',I~,A,ke~~~~], &(NC;, NL, d, k, 6)); upda,t,P-contributioris((Nr;, NL, d, k, d), v-fi], &old-lr, &old-P, &old-B); } write n-staks-f, key-A], key-trecf and U-,@ in file ,fil~i; 1 read n-slates-& kcy-tu, key-tree-t and u-tO from file filclD - 1; end= NO; while (!end){ Ir' = 0; P = 0; B = 0; sum-I/ = 0; for (t = 0; t 5 D - 1; t++){ ra..,sidcsf = nstates-t; key-a = key-.tU; key-tree-f = key-tree-l; v-fl = u-tu; read n-states-t, key-tu, key-tree-t ancl v-tn from file file-t ; for (i = I; i 5 n-stutes-t; i++) v.44 z 0; for (i = 1; i i: n-states-f; i++){ decode(D,C,K,A,key-~Ai], &(Nr;, NL,~, k, 6)); obtain the set. S of descriptions of the successors of (NG, NL, d, IC, 6) and the corresponding transition probabilities; for (each (XG,NZ,,~,~,S) E S){ let q he the transition probability associated with (NG, NL, d, k, 6); b = enc.ode(L),C,K,A,( NG, NL, d, k:, 6)); answer = look-andinsert( key-tree-t, b, gLj) : v-s[j] += y x u-14; } SZ197tLU += v-t[i]; 1 for (i = 1; i 5 n-states-t; i++){ decode(U,C',l~,A,kear-t[i3, &(NG, NI,, d, Fc, 8)); update-conkcibutions((NG, NL, d, k, 6), u-t[i], &k, &P, &B); 1 write n-stu.tes4, key-to, keg-tree-t and v-20 in file filed; 1 if (IEr - old..lrl/k. <: E' && (P - old-PI/P 5 E' k& IB - old-UI/R 5 E') end = YES; else end = NO; old-lr = lr; old-P = P; old-B = B; 1 1p = Er x D/surn.-u; P = P/sum-u; B = Blsum-u; Fig. 3 Model solution algorithm using Gauss-Seidel method with access to elements of Q in row-by-row manner 192 IEE Proc-Commun., Vol. 148, No. 3, June 2001 tion algorithm: number of states, array of keys, key tree and iteration vector. The function update-contributions0 updates the unnormalised contributions to Ir, P and B of a state with given description. The iteration step involves looping through the subsets Scf from t = 0 to D - 1. In each iteration of the loop, n-states f, keyx], key-tree f and v_f[] hold the data structures of the predecessor subset sct-l mod D and n-states-t, key-t[ 1, key-tree-t and v-t[] hold the data structures of the current subset Set. The iteration vector is kept unnormalised during the Gauss-Seidel iterations. Convergence is considered achieved when the relative tolerance in all the unnormalised quantities Ir, P and B is smaller than or equal to a given small quantity E', significantly smaller than E. Once convergence is achieved, the cell loss probability Ip and the normalised P and B are computed using the sum of the components of the iteration vector which is cumulated in sum-v. Note that, at a given time, only the data structures associated with subsets Sci and S2-I mod are held in memory. 4 Results Table 2 illustrates the more relevant parameters related to the computational effort of the solution of the model. For several values of p, D, and K (A = D - 1 in all cases) and E = lt5, E' = 10-lo, we give the value of C for which the solution converged. For that value of C, we also give the number of states, number of transitions, memory and disc storage requirements in Mbytes, and number of required Gauss-Seidel iterations. The amount of disc storage is approximately the amount of memory which would be required were our techniques to reduce it based on the cyclic nature of the model not implemented. The required value of C increases with K and is very sensitive to the utilisation factor of the first queue (p + UD). The required number of Gauss-Seidel steps is also sensitive to the utilisation factor of the first queue and increases slightly with K. We can note that memory requirements are very small relative to the number of states and transitions of the models. Our code reads arrays from disc and writes arrays into disc using, respectively, the ANSI standard C functions fread0 and fwriteo. CPU utilisation was high (between 70% and 90%, depending basically on the value of C). That high CPU utilisation is due to the fact that the average number of transitions per state is high (the contribution to the wall clock time due to disc accesses is approximately proportional to the number of states while the CPU processing time is approximately proportional to the number of transitions). Thus, the price paid by our techniques to reduce the memory requirements exploiting the cyclic nature of the model is small. The flop rate our code achieved in a 167 MHz UltraSPARCl workstation was about 174Kflops, which is reasonable considering the use of the 2-3 trees to find the indices of the states from their keys. 4 12 20 28 36 44 52 Fi Relative tolerance in loss rate (criterion I) and relative errorsfor pr%*ility that server of irst queue is servin a local cell (criterion 2) mi expected number of globahh in first queue (%iterion 3) for incremmg C, for case p= 0.7, D = IO, A= 9, K= 4 -0criterion 1 - -+- - criterion 2 -0.- criterion 3 C 4 Fig. 4 gives the relative tolerance for the loss rate and the relative errors for the probability that the server of the first queue is serving a local cell and the expected number of global cells in the first queue for increasing C, for the case p = 0.7, D = 10, A = 9, K = 4 and E' = We can note first that both the relative tolerance and the relative errors decrease at a high rate when C increases. Then, we should expect the relative error in the loss rate to be smaller than the relative tolerance. Also, of the two other criteria, the strongest is the third one. For this particular set of parameters, the required C would be determined by the first criterion. Figs. 5-7 display the loss probability as a function of K for, respectively, p = 0.3, p = 0.5 and p = 0.7, for D = 5, D = 10 and D = 15, and A = D - 1 in all cases. All results were obtained using E = and E' = 10-lo. The K required to achieve a given loss probability increases with the load (i.e. with p, and with l/D). Except for the cases in which the load of the system p + 1/D is very high, the minimum K required to achieve a loss probability smaller than or equal to a typical value is moderate. Thus, for instance, Table 3 gives the minimum required K to achieve a loss Table 2: Crequired for convergence in results and number of states, number of transitions, memory and disc storage requirements in Mbytes, and number of required iterations for that value of Cfor several sets of model parameter values ~ ~ ~~ ~~ ~ Disc p D K C States Transitions Memory storage Iterations 0.5 10 2 18 29 513 223 968 0.341 1.64 25 0.5 10 8 32 365224 4557871 4.24 21.1 42 0.7 5 2 46 42 167 708 045 1.01 2.50 513 0.7 5 8 56 249129 4896546 5.66 14.0 631 0.7 5 14 66 598935 13159009 14.0 34.7 731 0.7 15 2 28 166578 1819 108 1.30 9.62 49 0.7 15 8 42 1479612 23397586 11.3 ' 84.3 70 0.7 15 14 52 3905134 73869532 29.8 222 86 A = D1 in all cases IEE Proc.-Commun.. Vol. 148, No. 3, June 2001 193 probability smaller than or equal to 10-6 for all pairs of p, D values considered in Figs. 5-7. Since the first queue has infinite capacity, the traffic entering the second queue has average rate 1/D. The average service rate of that queue is l/A = 1/(D - 1). Thus, the utilisation of the second queue is K Fig. 5 D-I -0D=5 .... 0 .... D= 15 Loss probability as ajinction of Kfor p = 0.3, D = 5, 10, 15 and A = --+-- D=lO 'F 2 4 6 8 10 12 14 K Lossprobability as afiozctwn of Kfor p = 0.5, D = 5, 10, 15 and A = Fig. 6 D-I -0D=5 --+-- D = 10 .... 0 .... D= 15 1081 I I I I I I I I1 I I, 2 4 6 8 10 12 14 K Loss probability as a function of Kfor p = 0.7, D'= 5, 10, 1.5 and A = Fig. 7 -0D= 5 .... 0 .... D= 15 D-1 --+-- D= 10 (D - l)/D, which increases with D. However, the results given in Figs. 5-7 show that the loss probability decreases with D, i.e. with the utilisation. The counter-intuitive behaviour can be explained by the fact that the burstiness of the traffic diverted to the second queue decreases with D, compensating the larger utilisation for larger D. Table 3: Minimum required K to achieve a loss probability smaller than or equal to IO4 for several values of p and D P D K 0.3 15 3 0.3 10 3 0.3 5 5 0.5 15 5 0.5 10 6 0.5 5 11 0.7 15 13 0.7 10 >I4 0.7 5 >I4 5 Conclusions We have developed a methodology for the solution of finite large cyclic DTMC models and have applied it to the performance analysis of an ATM leaky bucket GFC protocol. The methodology includes 'on-the-fly' model generation techniques and exploits the cyclic nature of the model to reduce to a minimum the memory requirements. In addition, we have shown that the class of ATM leaky bucket models considered in the paper have a single irreducible closed set and have characterised this set implicitly by showing that there exists a path from any feasible state of the model to a given particular state U included in the irreducible closed set. The characterisation allows the irreducible closed set to be generated very easily (it is enough to generate states from U). The approach we have followed to show that the models have a single irreducible closed set and characterise it could be used for other classes of models. Using the methodology developed in the paper, we have been able to solve very large models (3 905 134 states and 73 869 532 transitions) using moderate amounts of memory (29.8 Mbyte) and disc storage (222Mbyte). The results obtained can be used to select an appropriate value for the parameter K as a function of the negotiated portion (l/D) of the bandwidth of the medium and the monitored upstream load @) depending on the required loss probability (typically between and lod). CPU time requirements are, however, large. In the future we are planning to develop and test approximations requiring the solution of smaller models so that more efficient computations can be performed. In that sense, the methodology proposed in the paper will allow the goodness of such approximations to be tested. 6 Acknowledgments The work of the two first authors has been supported by the Comisi6n Interministerial de Ciencia y Tecnologia (CICYT) under the research grant TIC95-0707-C02-02. The authors are grateful to the anonymous reviewers, whose comments have allowed the improvement of a former version of the paper. IEE Proc.-Commun.. Vol. 148, No. 3, June 2001 194 7 References 1 RITTLER, M.: ‘Multirate models for dimensioning and performance evaluation of ATM networks’. COST-242 report, June 1994 2 ADAMS, J.: ‘A multiservice GFC protocol’, Znt. J. Digit. Analog Commun. Syst., 1991,4, pp. 123-130 3 BUDRIKIS, Z.I., MERCANKOU, G., BLASIKIEWIECZ, M., ZUKERMAN, M., YAO, L., and POTTER, P.: ‘A generic flow control for B-ISDN. Proceedings of IEEE Infocom’92, 1992, pp. 895-904 4 FALCONER, R., and ADAMS, J.: ‘Orwell: a protocol for an integrated services local network‘, British Telecommun. Tech. 1, 1985, 3, 5 NIESTEGGE, G.: ‘The leaky bucket policing method in ATM networks’, Int. J. Digital Analog Commun. Syst., 1990,3, pp. 187-197 6 GUILLEMIN, F., and MONIN, W.: ‘Analysis of cell clumping caused by ATM network GFC protocols’, Comput. Commun., 1994, 17, (9), pp, 637646 7 GUILLEMIN, F., and DUPUIS, A.: private communication, 1996 8 KLEINROCK, L.: ‘Queueing systems. Volume 1: Theory’ (John Wiley, 1975) 9 CINLAR, E.: ‘Introduction to stochastic processes’ (Prentice-Hall, 1975) 10 DEAVOURS, D.D., and SANDERS, W.H.: ‘On-the-fly solution techniaues for stochastic aetri nets and extensions’. Proceedinns of the (4), pp. 27-35 7th IEkE international iorkshop on Petri nets and perform&e models (PNPM97), 1997, pp. 132-141 11 AHO, A.V.. HOPCROFT. J.E., and ULLMAN, J.D.: ‘Data structures and algorithms’ (Addison-Wesley, 1985) 12 STEWART, W.J.: ‘Introduction to numerical solution of Markov chains’ (Princeton University Press, 1994) 13 CIARDO, G., BLAKEMORE, A., CHIMENTO, P.F., MUPPALA, J.K., and TRIVEDI, K.S.: ‘Automated generation and analysis of Markov reward models using stochastic reward nets’ in MEYER, C., and PLEMMONS, R. (Eds.): ‘Linear algebra, Markov chains and queumg models’ (IMA, 1992) 8 We will use the theory of classification of the states of finite DTMCs as presented in [9]. It is enough to prove that there exists a path in nC from any state s E Rc - {U} to U. To see that, note first that U cannot be transient, since, if it were, some state in Qc - (U) would be recurrent and the existence of a path from such state to U would contradict that U is transient (there cannot be paths from recurrent states to transient states). Not being transient, U has to belong to some irreducible closed set. Let Sc be the irreducible closed set of Ilc including U. States which are not reachable from U cannot belong to Se Also, every state s reachable from U belongs to Sc, since there exist paths between s and U in both directions. Thus, Sc includes precisely U and the states reachable from U. To show that Sc is the only irreducible closet set of nC, assume the existence of another irreducible closed set Sl, and pick any state SE S‘,. The existence of a path from s to U e Sc is in contradiction with SC being a closed set. We split the proof of the existence of a path from any state s E Qc - {U) to U in two parts. In the first part, we prove the existence of a path from any state (NG, NL, d, k, S, to a state (0, 0, d‘, k‘, 6’). In the second part, we prove the existence of a path from any state (0, 0, d, k, S, to U. We use the notation Appendix: Proof of theorem 1 s -3 s’ to indicate that state s‘ is the successor of s reached when x arrivals from the global stream occur. First part: Consider an arbitrary state (NG, NL, d, k, s) and the path made up of transitions associated with 0 arrivals of the global stream. Such a path reaches a state (0, 0, d‘, k’, 6’). To prove that, consider first the case NG > 0. Each transition will reduce NG by 1 and after NG transitions the path will reach a state (0, NL, d‘, k’, 6’). Then, it is enough to prove the result for states (0, NL, d, k, S, with NL > 0. Consider a path of transitions associated with 0 arrivals from the global stream. Each transition will decrement NL if d < D and will leave unchanged (or decrement if NL = C) IEE Proc -Commun , Vol. 148, No 3, June 2001 NL if d = D. Since D > 1, a d cycle will decrement NL and a state (0, 0, d‘, k’, 6’) will be eventually reached. Second part: We will analyse nine cases for the values of d, k, 6, D, A and K. For each case we will find a path to either U or a state of the form (0, 0, d’, k’, 6’) in such a way that the existence of a path to U will be guaranteed. Case 1: k = 0 Note that k = 0 implies 6 = 0. Therefore, the state will be of the form (0, 0, d, 0, 0). A path to U is: (O,O,d,O,O) -%(O,O,d+l,O,O) A... o\ (O,O, D, 0,O) -% (0,1, I, 0,O) Case2:A-6>D-d+l,O<k<K,K>l Note that D - d + 6 + 2 I A and that, with d I D and (sincek>0)62 1,A-D+d-6+ lIA-6+ 1 IA<D. These two inequalities and k < K, D 2 2 guarantee the existence of the path: (O,O, d, IC, 6) 3 . I. o\ (O,O, D, IC, D - d + 6) o\ (0,1,1, k,D - d+ 6 + 1) -% (0,0,2,k + l,D - d + 6 + 2) 5.. , o\ (O,O,A - D + d - 6, k + 1,A) -% (O,O,A - D + d6 + l,k, 1) Letd’=A-D+d-6+ land#= l.WehaveA-$-(D -d’)=d-6-2(D-A)<d-6-(0-A)=A-&(D-d). Then, by concatenating an enough number of such paths we will reach a state corresponding to case 5, 8 or 9. Case 3: A6> D-d+ 1, k = K > 1 Note that D - d + 6 + 2 5 A and, as shown in the discussion of case 2, A - D +d - 6 + 1 < D. These two inequalities with D 2 2 guarantee the existence of the path: (0, 0, d, K, 6) o\ . * * -% (O,O, D, K, D - d + 6) 3 (0,1,1,K,D-d+6+1) 3 (0,0,2,K1D-d+6+2) 3 (O,O,A - D +d6+ l,K - 1,l) 0 -+ . . . -% (0,O,A - D + d - 6,K,A) and we have reached a state corresponding to case 2, 5, 8 or 9. Case 4: A - 6> D - d + 1, k = K= 1 Note that D - d + 6 + 2 I A and, as shown in the discussion of case 2, A - D + d - 6 + 1 < D. Since D 2 2, the following path exists: (0,0,d,1,6) A... 5(0,0,D,11D-d+6) -% (0,1,1,1,D - d+ 6 + 1) 3- (0,0,2,1,D - d+6+2) 0 --+ . . . 5 (O,O, A - D + d - 6, 1,A) 5 (O,O, A - D + d - 6 + l,o,o) and we have reached a state corresponding to case 1. CaseS:A-6=D-d+l,O<k<K,K>l Since D 2 2 we have the path: (O,O, d, k,6) 3 * . . 5 (O,O, D, k, a - 1) -% (0,1,1, IC, A) -% (O,O, 2, IC, 1) = s 195 Calling the 'd' and '6' parameters of s, respectively, d' and 8, we have A - 6' = A - 1 5 D - 2 = D - d'. Notice also that k > 0. Then, we have reached a state corresponding to case 8 or case 9. Case 6: A - 6 = D - d + 1, k = K > 1 Since D 2 2 we have the path: (O,O,d,K,S) A... -%(O,O,D,K,A-1) 3 (O,l,l,K,A) 5 (0,0,2,K1,l) = s Calling the 'd' and '6' parameters of s, respectively, d' and 8, we have A - 6' = A - 1 I D - 2 = D - d'. Notice also that K - 1 > 0. Then, we have reached a state corresponding to case 8 or case 9. Case7: A-6 = D-d+ 1, k = K= 1 Since D 2 2 we have the path: 0 0 (0,0)d,1,6) -+... --+ (0,0,D,l,A1) -% (O,l,l,l,A) -% (0,0,2,0,0) and we have reached a state corresponding to case 1. Case 8: A6= D -d, 0 < k I K For k = 1 we have the path to U: (O,O,d,1,6)3.. .&(O,O)D, l,A)&(O, l,l,O,O) For k > 1 and A = 1, since D 2 2, we have the path (O,O,d,k,S) A... 3(0,0,D,k,A) 3 (0,1,1, IC - 1,l) -% (O,O, 2, IC - 1,l) (0,O,d,Ic,6) ~...~( O,o,D,Ic,n) 3 (O,l,l,k1,l) 3 (0,0,2,k,2) -+... 0 &(O,O,A,IC,A) ~(O,O,A+l,IC-l,l) For k > 1 and A > 1, since D 2 2, we have the path: In both last cases we reach a state with 'k' one less. Case 9: A6 < D - d, 0 < k I K Note that A - 6 + d + 1 I D. For k = 1 we have the path: (O,O,d,l,S) A... &(O,O,A-S+d,l,A) -% (O,O, A - 6 + d + 1,0,0) and we reach a state corresponding to case 1. For k > 1 the previous path leads to (0, 0, A - 6+d + 1, k - 1, l), a state with 'k' one less. To summarise, cases 4 and 7 are reduced to case 1; cases 2, 3, 5 and 6 are reduced to case 8 or case 9. In case 1 there exists a path to U. In cases 8 and 9 either there exists a path to U or a path to a state with 'k' one less, guaranteeing that the case 1 (k = 0) will be eventually reached. Thus from any state (0, 0, d, k, s) there exists a path to U and the proof of the theorem is finished. 196 IEE Proc-Commun.. Vol. 148. No. 3. June 2001