Time Dependent System State Probabilities of Single Server Queuing System with Infinite Queue

The analyses of queuing systems, based on Markov processes, are in most cases restricted to steady state behavior i.e. to systems in equilibrium. The reason for that lies in the fact that for obtaining time dependent system state probabilities of these queuing systems, a system of linear first order differential equations has to be solved. Unfortunately, analytical solutions rarely exist and if they exist, their obtaining tends to be quite difficult and complicated. Several authors have obtained the results time dependent system state probabilities (properties) of some queuing systems in analytical form. These analytical expressions are usually obtained by use of generating functions and transforms such as Laplace transform, z-transform etc. The derived expressions are very complicated and require alternative computational techniques by the fact that they often refer to Bessel functions. The transient analysis of single server queuing systems with infinite queue capacity, based on Markov processes i.e. M/M/1, can be found, for example, in the works of: Morse [1,2], Greenberg and Greenberg [3], Heathcote and Winer [4], Gross and Harris [5], Kleinrock [6], Cooper [7], Takacs [8]. The main idea of this paper is to transform, using elementary mathematical operations, the expression for time dependent system state probabilities of single server queuing system with finite queue capacity M/M/1/K, in order to find its limit value when the number of places in the queue (m) tends to infinity. In another words, main idea is to obtain analytical expression for time dependent system state probabilities of single server queuing system with infinite queue capacity M/M/1 indirectly without solving adequate system of linear differential equations, in the case that the system is empty at the beginning. 2. M/M/1/K QUEUING SYSTEM


INTRODUCTION
The analyses of queuing systems, based on Markov processes, are in most cases restricted to steady state behavior i.e. to systems in equilibrium.The reason for that lies in the fact that for obtaining time dependent system state probabilities of these queuing systems, a system of linear first order differential equations has to be solved.Unfortunately, analytical solutions rarely exist and if they exist, their obtaining tends to be quite difficult and complicated.
Several authors have obtained the results time dependent system state probabilities (properties) of some queuing systems in analytical form.These analytical expressions are usually obtained by use of generating functions and transforms such as Laplace transform, z-transform etc.The derived expressions are very complicated and require alternative computational techniques by the fact that they often refer to Bessel functions.
The main idea of this paper is to transform, using elementary mathematical operations, the expression for time dependent system state probabilities of single server queuing system with finite queue capacity M/M/1/K, in order to find its limit value when the number of places in the queue (m) tends to infinity.In another words, main idea is to obtain analytical expression for time dependent system state probabilities of single server queuing system with infinite queue capacity M/M/1 indirectly without solving adequate system of linear differential equations, in the case that the system is empty at the beginning.

M/M/1/K QUEUING SYSTEM
The corresponding system of linear first order differential equations, defining the time dependent state probabilities of the M/M/1/K queuing system, is [6]: where: λ=const. is mean arrival rate, µ=const. is mean service rate, K=m+1 is maximal possible number of units (customers) in the system, p i (t) is time dependent probability that the system at the time t will be in the state i, and Analytical solution of differential equations system (1) may be found in works by several authors, such as: Morse [2], Takacs [8], Sharma and Gupta [9], Tarabia and El-Baz [10], Bugaric [11] etc.
The case when the system is empty at the beginning means that there is no customers being present in the system initially i.e. at t=0.In such case the initial values of the system state probabilities are: Integration constants, determined upon the initial values of the system state probabilities (3), are [11]: By substituting (4) and ( 5) into (2) the final expression for time dependent system state probabilities, whose limit has to be found when m→∞, has following the form:

LIMIT OF THE SYSTEM STATE PROBABILITIES WHEN m→ → → →∞ ∞ ∞ ∞
Applying known trigonometric formulas, expression (6) can be transformed as: ; 0,1, 2,..., 1 On the other side, part of expression ( 7) Previous definite integral, using the following formulas: ( ) can be transformed into the following form: where ( ) is modified Bessel function of the first kind of order n (n -integer number).[14] Solution of definite integral: (10) is Gamma function.Previous expression depends only from n (n = 0,1,...) and in further text it will be denoted as R(n).Analytical values of R(n) in dependence of λ and µ will be calculated later.
Finally, the solution of definite integrals given by expression (9) i.e. ( 8) is: with expression (10), expression (7) obtains the following form: For further analysis, previous expression should be written in expanded form i.e. decomposed into individual sums per k.After multiplying and dividing expression (12) by π, elementary mathematical transformations and application of known trigonometric formulas, a new convenient form for finding the limit of expression (11) is following: According to definition of definite integral, integral sum becomes definite integral if the following limit exists [15]: Each of the sums in expression (13) according to (14) can be transformed into definite integral: when m→∞ in the following way: Previous definite integral is, in fact, a modified Bessel function of the first kind of order i of argument i.e. [14,16]: ) The limit of all sums per k in expression ( 13) when m→∞ and 0 ) ( → ∆ m k θ can be determined in the same way. The limit of the first addend in expression (13) for ρ < 1 when m→∞ is: Replacing all sums per k in ( 13) with its limit and first addend of (13) with expression (15), the limit of time dependent system state probabilities when m→∞, is obtained as: Previous expression, using known relation between modified Bessel functions of first kind , can be reduced to [13]: Analytical values of R(n) can be determined indirectly from expression (17).First, it is necessary to replace initial conditions (3) into (17) i.e. for t=0, 0 (0) 1 p = , (0) 0 i p = , i=1,2, ..., ∞.In order to determine the values of expression (17) for t=0, for each system state probability the L'Hospital-Bernoulli rule has to be applied to solve undefined forms of type 0/0 which will appear.Also, the formula for differentiation of modified Bessel functions of the first kind:  as well as the fact that 0 (0) 1 I = and (0) 0 n I = , n=1, 2, ..., have to be applied.Applying previous to expression (17) for various values of parameter i (different system state probabilities), the system of recurrent formulas suitable for determining analytical values of coefficients R(n), n=0,1,2, ... , is obtained as: In order to solve the system (18) it is necessary to know first two values of coefficients R(n) i.e.R(0) and R (1).Those values can be obtained by solving definite integral (10) for n=0 and n=1.
General analytical solution of the system (18) is: Finally, by replacing (19) in (17), analytical expression for the time dependent system state probabilities of single server queuing system with infinite queue capacity M/M/1, obtained as the limit value of the system state probabilities of single server queuing system with finite queue capacity M/M/1/K when the number of places in the queue (m) tends to infinity, has the following form:

SIMULATION MODEL OF M/M/1 QUEUING SYSTEM
Validation of the expression for the time dependent system state probabilities of single server queuing system with infinite queue capacity M/M/1 (20) will be done by using adequate simulation model (discrete event simulation).The reason for that is the fact that the system of infinite differential equations can not be exactly solved with known numerical methods.The simulation model used for validation of expression (20) as output has only change of system state probabilities in time in dependence of λ -mean arrival rate and µ -mean service rate i.e. utilisation factor ρ. The algorithm of the developed simulation model is shown in figure 1.Initial conditions for simulation model (experiment) are: − number of units in the system as well as in the queue in equal to zero, − servicing channel is free (in the state "waiting"), and − first unit comes to the system at t=0.
Inter arrival times of units to the system are generated according to exponential distribution with parameter λ, while unit servicing times are generated according to exponential distribution with parameter µ.
For every time unit during simulation time, in each simulation experiment, time and system state are written in separate database (depending on system state).
When given number of simulation experiments are finished, from created databases, absolute frequencies of system states in given time intervals are calculated (histogram of absolute frequencies).After that relative frequencies histogram i.e. system state probabilities, in given time intervals, are calculated.
Overview of labels used in the simulation model: − State -state of servicing channel: "Wo" -work, "Wa" -waiting, − No sim -number of simulation experiments, − N w -current number of units in the queue, − N ws -current number of units in the system, − t -current system time, − t ar -moment of arrival of new unit to the system, − t sim -duration of simulation experiment, − t sc -moment of service completition i.e. change of servicing channel state, − t ser -servicing duration, − d int -interval length for system state absolute and relative histograms, − X sist -time dependent system state vector for each simulation experiment, − X j -time dependent absolute frequencies of j-th system state, − p j -time dependent relative frequencies (probabilities) of j-th system state, − RN -random number generated according to uniform distribution in interval from 0 to 1. Duration of each simulation experiment is 6000 time units, the number of executed simulation experiments is 4000 for one pair of values for λ and µ i.e. utilisation factor ρ.

COMPARISON OF ANALYTICAL AND SIMULA-TION RESULTS
Diagrams in pictures 2, 3, 4 and 5 show change in time of the first four system state probabilities (p 0 , p 1 , p 2 and p 3 ) of single server queuing system with infinite queue capacity M/M/1, for different values of utilisation factor ρ, such as 0.35, 0.5, 0.65 and 0.8 respectively.Time dependent system state probabilities obtained as a result of simulation are marked with symbols (*, +, , ), while time dependent system state probabilities obtained using expression (20) are not marked.Analysis of results presented in figures 2÷5, shows that the values of the system state probabilities obtained as a result of simulation match values of the system state probabilities calculated using expression (20).This leads to the conclusion that expression (20) for time dependent system state probabilities of single server queuing system with infinite queue capacity M/M/1 is correct.

CONCLUSION
In this paper, analytical expression for time dependent system state probabilities of single server queuing system with infinite queue capacity M/M/1 is derived as a limit value of expression for the time dependent sys-tem state probabilities of single server queuing system with finite queue capacity M/M/1/K.The limit value is found in the case when the number of places in the queue tends to infinity and the case that the system is empty at the beginning.
Validation of derived expression for the time dependent system state probabilities of single server queuing system with infinite queue capacity M/M/1 is done by developed simulation model.
The application of derived expression can be found in the analysis of non-stationary working regimes of transportation devices in industry.For example, the source of transportation units (pallets) is an output from production (final goods), while service consists of storing pallets into warehouse using one AS/RS device.All pallets have to be stored in the warehouse, so the limitation of input zone (queue) of the warehouse, theoretically, does not exists.