Introduction

This research is devoted to verification by simulation a few performance results obtained in previous work [17,16]. Moreover, we focus on verification of the obtained stability conditions of the system with coupled orbits. In general, sufficient and necessary stability conditions are different, and the study of the "gap" between these conditions is important in practice. Indeed, as we show, sufficient condition is redundant, and possible extension of stability region seems to be useful to increase the efficiency (throughput) of the system. In this work we pay the main attention to the so-called symmetric model, in which all classes of customers have the same parameters. This scenario allows to simplify analysis and detect some properties which turn out to be useful in a more general setting. In simulation, we focus on the symmetric system with three classes of customers which follow independent Poisson inputs. It is worth mentioning that this research complements and develops previous works [17,16].

To study classical retrial queues, we mention the books [12,1], and the survey papers [2,13]. Also the stability analysis of a multi-class retrial queue with constant retrial rates, which do not depend on the states of other orbits, has been developed in [4,5].

As to application of the system with coupled orbits, we mention the modelling of wireless multiple access systems, in particular, relay-assisted cognitive cooperative wireless systems [18]. Moreover, in the modern cognitive radio [15] there exists a possibility to dynamically adjusts retransmission rates to improve spectrum utilization [6,8,11]. Furthermore, as it has been mentioned in [16], this model is suitable to describe dynamics of cellular networks, in which the transmission rate in a particular cell decreases as the number of users in the neighboring cells increase [6]. A similar effect is observed in the processor sharing models [7,14].

This paper is organized as follows. In Section 2 we describe the basic model. In Section 3 we summarize the main theoretical results obtained in [17] which we verify by simulation in Section 4. In particular, we introduce the symmetric model with coupled orbits. In section 4 we verify theoretical results simulating 3class symmetric model with exponential and Pareto service times. In particular, we estimate the stationary probability that a fixed orbit is busy and server is idle, and demonstrate the correctness of the lower and upper bounds of this probability. Also the accuracy of stability conditions is studied, using the "gap" between the necessary and sufficient stability conditions. Actually simulation shows that the necessary stability condition is in fact stability criterion of the system.

Description of the model

We study a single-server with no buffer for the waiting customers, and with three classes of customers. Nevertheless, it is worth mentioning that the theoretical results, which we will verify, hold for an arbitrary number N of customer classes. By this reason we will formulate below theoretical results for this general setting.

The class-i customers form Poisson input with rate λ i , i = 1, . . . , N . Because all inputs are assumed to be independent, then the summary input is Poisson as well, with rate λ := i λ i , in which an arbitrary arrival belongs to class i with the probability p i =: λ i /λ, i = 1, . . . , N. Then interarrival times of the input are exponential with generic element τ and expectation Eτ = 1/λ ∈ (0, ∞). We also assume service times of class-i customers, {S (i) n , n ≥ 1}, to be independent identically distributed (iid) with service rate

γ i =: 1 ES (i) ∈ (0, ∞), i = 1, . . . , N.

A customer, meeting server busy joins a class-dependent virtual orbit, where joins the end of the orbit queue. At that the head customer waiting in orbit i makes retrial attempts until he finds server idle to occupy it. The distance between attempts are exponentially distributed with rate µ J(i) , where J(i) is the current configuration of the orbits: busy or empty. Thus each orbit acts as a FIFO queueing system with state-dependent "service" rate µ J(i) . This dependence is a key new property of the model.

To be more precise, for each i, we define (N − 1)-dimensional vectors

J(i) = {j 1 , . . . , j i−1 , j i+1 , . . . , j N }

with binary components j k ∈ {0, 1}, where the ith component is omitted. If the k-th orbit is currently busy, we put j k = 1, otherwise, j k = 0. Each vector J(i) is called configuration. For each i, we introduce the set G(i) = {J(i)} of possible configurations. It is assumed that there is given constant µ J(i) , retransmission rate from orbit i, if current configuration is J(i). We denote M i the set of rates for all configurations belonging to G(i). This construction, proposed in [17], considerably generalizes the setting studied in previous works [16,11,10,9]. In these works, it is assumed that orbit i has rate µ i if at least one (other) orbit is busy, otherwise, the rate is µ * i , i = 1, . . . , N . Thus, in setting in [16], each set

M i = {µ * i , µ i }.

In this work we continue to study th new general setting from [17] with focus on verification some bounds and stability conditions proved in [17], for the system with three classes of customers.

Before to give main stationary performance measures to be verified by simulation, we mention that the main stochastic processes describing the dynamics of the system, such like accumulated work (workload) orbit size, ect., are regenerative, with regeneration instants T n . A regeneration occurs when a new customer meets an idle system [3]. The distances T n+1 − T n are iid regeneration periods, and we denote T the generic period.

A queuing process is called positive recurrent if the mean generic period is finite, that is ET < ∞ [3]. Under positive recurrence, there exists the stationary regime of the system [3].

Preliminary results

Now we define the main stationary performance metrics and give necessary and sufficient stability conditions proved in [17] by regenerative method. Let I(t) be the summary idle time of the server in interval [0, t], then busy time is defined as B(t) = t − I(t).

If the system is positive recurrent, then there exist the limits, with probability (w.p.) 1,

lim t→∞ I(t) t = P 0 = 1 − P b ,

where P 0 is the stationary idle probability of the server, and

P b = lim t→∞ B(t) t

is the stationary busy probability of the server. Denote B i (t) the summary time, in interval [0, t], when the server is occupied by class-i customers. It is proved in [17] that the stationary probability that the server is occupied by class-i customer is defined as

lim t→∞ B i (t) t = P (i) b = ρ i , i = 1, . . . , N.

Denote the traffic intensity for each class,

ρ i = λ i /γ i , i = 1, . . . , N,

and summary traffic intensity ρ = ρ i .

Because

P b = i P (i)

b then it follows that P b = ρ. Now we introduce the maximal possible rate from orbit i:

μi = max J(i)∈G(i) µ J(i) .

The following statement, which has been proved in [17], contains the necessary stability (positive recurrence) condition of our system. Theorem 1. If the N -class retrial system with coupled orbits is positive recurrent, then

P b = N i=1 ρ i = ρ ≤ min 1≤i≤N μi λ i + μi < 1.(1)

To formulate the next statement and sufficient stability conditions, we denote, for each class i, the minimal retrial rate

µ 0 i = min J(i)∈G(i) µ J(i) .

Also let P (i) 0 be the stationary probability that server is idle and orbit i is busy. The following statement is proved in [17]. Theorem 2. The following inequalities hold

λ i μi ρ ≤ P (i) 0 ≤ λ i µ 0 i ρ, i = 1, . . . , N.(2)

Now we formulate the sufficient stability condition [17]. Theorem 3. The sufficient stability condition of N -class retrial system with coupled orbits is

N i=1 ρ i + max 1≤i≤N λ µ 0 i + λ < 1.

This condition implies a negative drift of the workload process, and as a result, positive recurrence of the system, and can be written as

ρ = i ρ i < min i µ 0 i λ + µ 0 i .(3)

To compare (3) with the necessary stability condition (1), we define the gap between the necessary and sufficient conditions,

∆ =: min i μi λ i + μi − min i µ 0 i λ + µ 0 i > 0. (4)

An important special class constitute the symmetric systems. To explain symmetry for the system with coupled orbits, we first note that this system becomes standard retrial system with constant rate provided µ J(i) ≡ µ i , implying equality μi = µ 0 i = µ i . The latter standard retrial system with constant rate, becomes symmetrical, if the corresponding parameters are equal, that is

λ i ≡ λ, γ i ≡ γ, µ i ≡ µ.

However, the notion symmetrical system with coupled orbits is more flexible. To explain it in more detail, we define, for each i, the set of vectors describing retrial rates for each configuration J(i) = {j 1 , . . . , j N }. Because in our case N = 3, then the capacity |J(i)| = 4 for each i = 1, 2, 3. To compose J(i), we use lexicographical order, that is, for each orbit i, and two remaining orbits j < k, with k, j = i, the following four configurations J(i) are possible:

M i =: {(i j = 0, i k = 0), (i j = 1, i k = 0), (i j = 0, i k = 1), (i j = 1, i k = 1)}. (5)

(Recall that i k = 1 means that orbit k is busy, while i k = 0 means that it is empty.) We denote µ i 00 , µ i 10 , µ i 01 , µ i 11 , the retrial rates corresponding to each configuration in (5). Then we obtain that, in the symmetrical system, all sets

G(i) = {µ i 00 , µ i 10 , µ i 01 , µ i 11 }, i = 1, 2, 3,

are identical, although rates within G(i) may differ. One can expect that this structure leads to similar behavior of the orbits, and it is confirmed below by simulation.

Note that for the symmetric coupled orbits, the difference (4) becomes

∆ := μ λ/N + μ − µ 0 λ + µ 0 ,

where μ = μi , µ 0 = µ 0 i . Finally, for symmetric classical (non-coupled) orbits, μ = µ 0 and (4) becomes

∆ = µ λ/N + µ − µ λ + µ .

In the next section, containing simulation results, we analyze the symmetric model and leave studying a more general model for a further work.

Remark. For non-coupled orbits, µ J(i) = µ i for all configurations J(i), and each busy orbit i has a fixed retrial rate µ i . Then relation (2) becomes

λ i P b = µ i P (i) 0 ,

and we obtain the following explicit formula for the stationary probability that orbit queue i is busy and server is idle, see [17]:

P (i) 0 = λ i µ i N k=1

ρ k , i = 1, . . . , N.

Simulation results

In this section we verify by simulation some obtained above theoretical results for three classes of customers considering symmetric model. To this end, we define the following variables,

Γ 1 := min i μi λ i + μi − ρ, Γ 2 := min i µ 0 i λ + µ 0 i − ρ.

which delimit the boundary of stability region. More exactly, if Γ i > 0, i = 1, 2, then both stability conditions (1) and ( 3) are satisfied. If Γ i < 0, i = 1, 2, then instability of the orbits is expected. However, if the intermediate case Γ 1 > 0, Γ 2 < 0 holds, then, as simulation shows, the orbits are stable in all experiments. It indicates that the necessary stability condition is in fact stability criterion, while sufficient stability condition is redundant. (However we can not prove it strictly for the system with general service times.) These observations are illustrated by simulation below. We emphasize again that in this work we pay attention simulation symmetric model with coupled orbits. Now we present numerical results obtained by simulation, to verify theoretical analysis of stationary regime and necessary and sufficient conditions (1), (3). Everywhere we use the black, grey and dotted curve to demonstrate the dynamics of the 1st, 2nd and 3rd orbit, respectively. (The axis t counts the number of discrete events: arrivals, departures, attempts, in the applied discrete-event simulation algorithm.)

First we perform some experiments for the completely symmetric model and demonstrate stability/instability of all orbits. This analysis also shows the redundancy of sufficient condition (3) for stability because simulation stays all orbits stable even for Γ 1 < 0.

In the 1st experiment, we use the following input and service rates,

λ 1 = λ 2 = λ 3 = 3, γ 1 = γ 2 = γ 3 = 15,

and the following retrial rates:

M 1 = {µ 1 00 = 20, µ 1 10 = 30, µ 1 01 = 15, µ 1 11 = 25}, M 2 = {µ

Thus, with those parameters we receive ρ = 0.6 and Γ 1 = 0.3, Γ 2 = 0. Therefore, both stability conditions are satisfied and all orbits are stable, as expected, and we can see that at 3) hold: Γ1 > 0, Γ2 > 0; all orbits are stable.

In the following experiments we show the results for Pareto service times (Pareto model), and not equal service rates, ceteris paribus. That is we still keep a symmetry in the input rates.

Fig. 4 shows the dynamics of the orbits in Pareto model with service time distribution

F i (x) = 1 − ( x i 0 x ) α , x ≥ x i 0 (F i (x) = 0, x ≤ x i 0 ),

and expectation

ES (i) = α x i 0 α − 1 , α > 1, x i 0 > 0, i = 1, 2, 3.

We select α = 2 and the following values of the shape parameter x i 0 for orbit i = 1, 2, 3, respectively:

x i 0 = 1 24 , i = 1, 2, 3.

This choice gives the following service rates

γ 1 = γ 2 = γ 3 = 12,

ceteris paribus. Here ρ = 0.76, implying Γ 1 = 0.14 and Γ 2 = −0.16. Thus condition (1) holds while condition (3) is violated. As we see, all three orbits remain stable, however the stability is reached at a higher level.

In the 3rd experiment, shown on the Fig. 3, we further increase service rate

is violated: Γ1 > 0, Γ2 < 0; all orbits are stable at a higher level.

(ceteris paribus)

γ 1 = γ 2 = γ 3 = 10.

Here we obtain ρ = 0.9 and it gives Γ 1 = 0, Γ 2 = −0.3. Thus both conditions (1) and (3) are violated. (Γ 1 = 0 is called boundary case.) As we see on Fig. 3, all orbits become now unstable.

Thus in the these experiments, a gradual decreasing service rates (which implies reduction Γ 1 and Γ 2 ) makes orbits unstable only if the necessary condition (1) is violated. So we suggest that condition (3) is redundant, and moreover that the necessary condition (1) is indeed stability criterion.

In the following experiment we demonstrate the estimation the stationary probability P

(1) 0 for the symmetric model with exponential service time. In that experiment we use the following input and service rates,

λ 1 = λ 2 = λ 3 = 3γ 1 = γ 2 = γ 3 = 15,

while the retrial rates remain (6). Thus in this case µ 0 i and μi are different, and exact value of the target probability is unknown. However, by the positive recurrence, the sample mean estimate still converges to a limit. Fig. 4 shows that the sample mean estimator of P 3) are violated: Γ1 < 0, Γ2 < 0; all orbits are unstable.

Conclusion

In this work, we simulate a 3 -class symmetric retrial system with independent Poisson inputs and the coupled orbits to verify some theoretical results found earlier. In this system, a new customer meeting server busy joins the corresponding infinite capacity orbit. The retrial rate from orbit i depends on the current configuration of other orbits: busy or idle. We verify by simulation some stationary performance measures and the accuracy of the found earlier stability conditions of this model.