A two-queue JSQ-SLQ service system with absolute preemptive priority and non-zero switch-over times

Full Text

Introduction In queueing theory, an important role is played by models that incorporate control mechanisms for customer arrivals and the service process. Such models are widely used in applications where high throughput, efficient load balancing, and low delays are required, in particular in telecommunication networks, data centers, and transportation systems [1-5]. As the load on such systems grows and their structure becomes more complex, dynamic policies that take the current state of the system into account become especially important. These policies make it possible to distribute the customer flow more rationally and to use service resources more efficiently. © 2026 E. A. Khityaev, I. S. Zaryadov, T. A. Milovanova This work is licensed under a Creative Commons “Attribution-NonCommercial 4.0 International” license. One of the classical approaches is the Join-the-Shortest-Queue (JSQ) routing discipline, under which an arriving customer is routed to the subsystem with the shortest queue. This model was first studied in 1958 for a system with two parallel queues and was described by a two-dimensional Markov process [6]. The classical setting assumes Poisson arrivals and exponentially distributed service times, which allows the system to be described in Markovian terms and makes Markov-process methods applicable to its analysis. Since then, the JSQ discipline has been studied in a wide range of queueing models [7- 14], including more recent studies [15-19]. The main advantage of the JSQ discipline is that it balances the load by distributing customers among the queues, thereby reducing the risk of overflow in systems with finite buffers and decreasing the mean waiting time [3, 5, 7]. Moreover, it has been shown that in large systems with parallel queues and servers, the mean waiting time under JSQ may converge to zero as the number of queues and servers increases, which illustrates the high efficiency of this discipline [5]. At the same time, JSQ requires information about the current lengths of all queues, which may introduce additional overhead in real systems [20]. Another important approach for improving system performance is the Serve-the-Longest-Queue (SLQ) service discipline, which was considered by Cohen in 1987 [21]. Under this discipline, the server selects the longest queue for service at each decision epoch, which helps equalize queue lengths and prevents individual queues from growing disproportionately. The SLQ discipline was later studied in several related queueing models [22-25]. As a result, SLQ reduces the probability that individual queues become overloaded, which helps to decrease the mean waiting time and, in systems with finite buffers, the loss probability. This is especially important under heavy traffic and in systems with multiple servers, where an uneven distribution of customers may overload some queues even when other servers are idle. In this context, the SLQ discipline has also been studied in the Halfin-Whitt regime [26]. Thus, SLQ can effectively complement routing approaches that balance customer arrivals. At the same time, it requires information about the current queue lengths to select a queue for service. Combining the JSQ and SLQ disciplines is of particular interest because it allows the routing of arriving customers to be coordinated with the dynamic control of the service process. Models of this type were studied by Perel et al., including a preemptive two-queue system with zero switch-over times [27], a non-preemptive two-queue system with zero and non-zero switch-over times [28], and a preemptive three-queue system with zero switch-over times [29]. In these works, the rule by which the server selects a queue for service is deterministic and based on the comparison of queue lengths: the server selects the longest queue for service. In the preemptive models, the server switches as soon as a non-attended queue becomes longer than the attended one, whereas in the non-preemptive model switching decisions are made only upon service completions. Ties between queues are treated separately. In the two-queue models, when the queue lengths become equal, the server remains at the current queue. In the three-queue model, the same deterministic longest-queue rule is used when there is a unique longest queue; however, if two non-attended queues share the longest length while the attended queue is shorter, the choice between them is made according to a probabilistic tie-breaking rule. Thus, systems in which the server may switch queues even when the queue lengths are equal, according to a probabilistic tie-breaking rule for selecting the queue to be served, are of interest and have not yet been considered. In systems with non-zero switch-over times, preemptive priority may be relative, meaning that an ongoing switch-over cannot be interrupted, or absolute, meaning that the server may interrupt the current switch-over and change its direction of motion. The present paper considers non-zero switch-over times between the queues and, in contrast to previously studied models, introduces a probabilistic rule for selecting the queue to be served when the queue lengths are equal. This generalization broadens the scope of the model and makes it better suited for systems in which server switching is not purely deterministic. The proposed model can be used to analyze and optimize systems of this type, including data centers, server clusters, telecommunication networks, and distributed computing systems, where both load distribution and server switch-over times play an important role. The remainder of the paper is organized as follows. Section 2 describes the model, including the routing and service rules, the switch-over mechanism, and tie-breaking probabilities. In Section 3, a two-dimensional Markov process is introduced and the corresponding stationary balance equations are derived. Section 4 introduces the conditional probability generating functions and discusses the stability condition for the existence of a stationary regime. Section 5 presents expressions for the main stationary performance measures. Section 6 contains numerical results showing how the main system parameters affect the stationary characteristics and phase probabilities. Section 7 provides a discussion and interpretation of the obtained numerical results. Finally, Section 5 summarizes the results and discusses possible extensions of the model. 2. System description The system consists of two infinite queues Qi, i = 1, 2, a single server, and non-zero server switch-over times between the queues. It operates under an absolute preemptive priority rule, under which the server may interrupt an ongoing switch-over and change its direction. Customers arrive according to a Poisson process with rate . Under the JSQ discipline, each arriving customer joins the currently shortest queue. If the queue lengths are equal, a probabilistic tie-breaking rule is applied: the customer joins queue Qi with probability pi, i = 1, 2, where p1 + p2 = 1. Throughout the paper, we use notation partly based on Perel et al. [27-29], adapted to the present model. Let Li(t) denote the number of customers in queue Qi, i = 1, 2, at time t > 0. According to the SLQ discipline, the server selects the longest queue for service. The service time of a customer from queue Qi is exponentially distributed with rate i, i = 1, 2. Customers in each queue are served according to the First-Come, First-Served (FCFS) discipline. When the server switches over to another queue, the customers remaining in the previous queue stay in the system and wait until the server returns. The dynamic priority service rule depends on the current queue lengths and is defined as follows. o If, while the server is serving a customer from queue Qi, a new customer arrives at queue Qj, j ≠ i, and after this arrival Lj(t) > Li(t), then the server interrupts service in queue Qi and starts a switch-over to queue Qj. If after this arrival Li(t) = Lj(t), then the server continues the current service in queue Qi with probability qi or interrupts it and starts a switch-over to queue Qj with probability qj = 1 - qi. o If, after completing the service of a customer from queue Qi, the condition Li(t) < Lj(t), j ≠ i, is satisfied, then the server starts a switch-over to queue Qj. If Li(t) = Lj(t), then the server either remains at the current queue Qi with probability qi or switches over to queue Qj with probability qj = 1 - qi. o Since the priority is absolute, if a new customer arrives at queue Qi while the server is switching over from queue Qi to queue Qj, then, if Li(t) > Lj(t), the server interrupts the current switch-over and starts switching back from Qj to Qi. If Li(t) = Lj(t), the server either continues switching over to Qj with probability qj = 1 - qi or returns to queue Qi with probability qi. In all cases, whenever the server starts a switch-over to queue Qj, the corresponding switch-over time is exponentially distributed with rate j , j = 1, 2. 3. Two-dimensional Markov process and balance equations Let D(t) = L1(t) - L2(t) denote the difference between the numbers of customers in the first and second queues. To describe the position and switch-over state of the server, we introduce the phase variable J(t), which encodes both the value of D(t) and the position or switch-over state of the server. Then the state of the system at time t can be described by the Markov process (t) = {L1(t), J(t)}. The probability distribution of this process is defined as follows: Pn,d(t) = ℙ{L1(t) = n, J(t) = d}, n ≥ 0, d ∈ , where only admissible pairs (n, d) are considered, and = {1, 1s1 , 01, 0s1 , 0s2 , 02, -1s2 , -1}. In each phase, the numerical part of the notation indicates the value of D(t), while the subscript, when present, specifies the position or switch-over state of the server. The phases d are defined as follows: o 1 means that D(t) = 1, the server is located at Q1, and a new customer joins Q2. o 1s1 means that D(t) = 1, the server is switching over to Q1, and a new customer joins Q2. o 01 means that D(t) = 0, the server is located at Q1, and a new customer joins Qi with probability pi, i = 1, 2, where p1 + p2 = 1. o 0s1 means that D(t) = 0 and the server is switching over to Q1. o 0s2 means that D(t) = 0 and the server is switching over to Q2. o 02 means that D(t) = 0, the server is located at Q2, and a new customer joins Qi with probability pi, i = 1, 2, where p1 + p2 = 1. o -1s2 means that D(t) = -1, the server is switching over to Q2, and a new customer joins Q1. o -1 means that D(t) = -1, the server is located at Q2, and a new customer joins Q1. Under the stability condition, which ensures the existence of a stationary regime, the stationary probabilities are defined by t Pn,d = lim Pn,d(t), n ≥ 0, d ∈ . (1) →∞ Then the stationary probabilities satisfy the following balance equations. For the state (n, 1): ( + 1)Pn,1 = 1Pn,1s1 + p1Pn-1,01 , n ≥ 1. For the state (n, 1s1): ( + 1)Pn,1s1 = p1Pn-1,0s1 + p1Pn-1,0s2 + p1Pn-1,02 + 2Pn,02 , n ≥ 1. For the state (n, 01): For the state (n, 0s1): P0,01 = 1q1P1,1 + 1P0,0s1 , n = 0, ( + 1)Pn,01 = q1Pn,1 + 1q1Pn+1,1 + 1Pn,0s1 , n ≥ 1. ( + 1)P0,0s1 = 2q1P0,-1, n = 0, ( + 1)Pn,0s1 = q1Pn,1s1 + q1Pn-1,-1s2 + q1Pn-1,-1 + 2q1Pn,-1, n ≥ 1. For the state (n, 0s2): ( + 2)P0,0s2 = 1q2P1,1, n = 0, ( + 2)Pn,0s2 = q2Pn,1 + 1q2Pn+1,1 + q2Pn,1s1 + q2Pn-1,-1s2 , n ≥ 1. For the state (n, 02): P0,02 = 2P0,0s2 + 2q2P0,-1, n = 0, ( + 2)Pn,02 = 2Pn,0s2 + q2Pn-1,-1 + 2q2Pn,-1, n ≥ 1. For the state (n, -1s2 ): ( + 2)Pn,-1s2 = p2Pn,01 + 1Pn+1,01 + p2Pn,0s1 + p2Pn,0s2 , n ≥ 0. For the state (n, -1): ( + 2)Pn,-1 = p2Pn,02 + 2Pn,-1s2 , n ≥ 0. As a result, we obtain the following system of balance equations: ⎧ ( + 1)Pn,1 = 1Pn,1s1 + p1Pn-1,01 , n ≥ 1, ⎪ ⎪( + 1)Pn,1s1 = p1Pn-1,0s1 + p1Pn-1,0s2 + p1Pn-1,02 + 2Pn,02 , n ≥ 1, ⎪ ⎪ P0,01 = 1q1P1,1 + 1P0,0s1 , n = 0, ⎪ ⎪( + 1)Pn,01 = q1Pn,1 + 1q1Pn+1,1 + 1Pn,0s1 , n ≥ 1, ⎪ ⎪( + 1)P0,0s1 = 2q1P0,-1, n = 0, ⎪( + 1)Pn,0s1 = q1Pn,1s1 + q1Pn-1,-1s2 + q1Pn-1,-1 + 2q1Pn,-1, n ≥ 1, ⎪ ⎨( + 2)P0,0s2 = 1q2P1,1, n = 0, ⎪ ⎪( + 2)Pn,0s2 = q2Pn,1 + 1q2Pn+1,1 + q2Pn,1s1 + q2Pn-1,-1s2 , n ≥ 1, ⎪ ⎪ P0,02 = 2P0,0s2 + 2q2P0,-1, n = 0, ⎪ ⎪( + 2)Pn,02 = 2Pn,0s2 + q2Pn-1,-1 + 2q2Pn,-1, n ≥ 1, ⎪ ⎪( + 2)Pn,-1s2 = p2Pn,01 + 1Pn+1,01 + p2Pn,0s1 + p2Pn,0s2 , n ≥ 0, ⎩( + 2)Pn,-1 = p2Pn,02 + 2Pn,-1s2 , n ≥ 0. 4. Conditional probability generating functions For each phase (2) d ∈ = {1, 1s1 , 01, 0s1 , 0s2 , 02, -1s2 , -1}, we introduce the corresponding conditional probability generating function (PGF) Gd(z): Gd(z) = ∑ Pn,d zn, |z| ≤ 1. (3) n Taking into account the admissible values of n for each phase, definition (3) gives G1(z) = ∑ Pn,1zn, ⎪ ⎧ n≥1 ⎪G1s (z) = ∑ Pn,1s zn, ⎪ ⎪ 1 n≥1 1 ⎪G01 (z) = ∑ Pn,01 zn, ⎪ ⎪ n≥0 ⎪G0s (z) = ∑ Pn,0s zn, ⎪ 1 n≥0 1 (4) ⎨G0s (z) = ∑ Pn,0s zn, ⎪ ⎪ 2 n≥0 2 ⎪G02 (z) = ∑ Pn,02 zn, ⎪ ⎪ n≥0 ⎪G-1s (z) = ∑ Pn,-1s zn, ⎪ ⎪ 2 n≥0 2 ⎪G-1(z) = ∑ Pn,-1zn. ⎩ n≥0 ( + 1)G1(z) - 1G1s1 (z) - p1z G01 (z) ⎪( + 1)G1s1 (z) - p1z G0s1 (z) - p1z G0s2 (z) - ( p1z + 2)G02 (z) ⎪-( z + 1)q1G1(z) + ( + 1)z G01 (z) - 1z G0s1 (z) ⎪- q1G1s1 (z) + ( + 1)G0s1 (z) - q1z G-1s2 (z) - ( z + 2)q1G-1(z) = 0, = - 2P0,02 , = 1z P0,01 , = 0, ⎨-( z + 1)q2G1(z) - q2z G1s (z) + ( + 2)z G0s (z) - q2z2G-1s (z) ⎪ 1 2 2 = 0, ⎪- 2G0s2 (z) + ( + 2)G02 (z) - ( z + 2)q2G-1(z) ⎪-( p2z + 1)G01 (z) - p2z G0s1 (z) - p2z G0s2 (z) + ( + 2)z G-1s2 (z) ⎩- p2G02 (z) - 2G-1s2 (z) + ( + 2)G-1(z) = 2P0,02 , = - 1P0,01 , = 0, which can be written in matrix form as A(z) G⃗(z) = P⃗(z), (5) Multiplying the corresponding equations of system (2) by zn, summing over n, and using relations (4), we obtain the following system of equations for the generating functions: ⎧ ⎪ ⎪ ⎪ ⎪ where G⃗(z) = (G1(z), G1 (z), G0 (z), G0 (z), G0 (z), G0 (z), G-1 (z), G-1(z))T s1 1 s1 is the vector of generating functions, s2 2 s2 2 1 2 1 P⃗(z) = (0, - 2P0,0 , 1zP0,0 , 0, 0, 2P0,0 , - 1P0,0 , 0)T is the right-hand side vector containing the unknown boundary probabilities P0,01 and P0,02 , and the matrix A(z) is given by ⎛ + 1 - 1 - p1z 0 0 0 0 0 ⎞ ⎜ 0 + 1 0 - p1z - p1z -( p1z + 2) 0 0 ⎟ ⎜ ⎟ ⎜-( z + 1)q1 0 ( + 1)z - 1z 0 0 0 0 ⎟ ⎜ ⎟ ⎜ 0 - q1 0 + 1 0 0 - q1z -( z + 2)q1⎟ . ⎜-( z + 1)q2 - q2z 0 0 ( + 2)z 0 - q2z2 0 ⎟ 0 0 0 0 - 2 + 2 0 0 0 -( p2z + 1) - p2z - p2z 0 ( + 2)z 0 0 0 0 0 - p2 - 2 ⎜ -( z + 2)q2⎟ 0 ⎜ ⎟ ⎜ ⎟ ⎝ + 2 ⎠ By Cramer’s rule, the generating functions, which are the solutions of system (5), are given by |A(z)| Gd(z) = |Ad(z)| , d ∈ , where |A(z)| is the determinant of the matrix A(z), and |Ad(z)| is the determinant of the matrix obtained by replacing the corresponding column of A(z) with the right-hand side vector P⃗(z). We note that the generating functions Gd(z), d ∈ , depend on the boundary probabilities P0,01 and P0,02, whereas the remaining probabilities at n = 0 are not independent and can be expressed in terms of these boundary probabilities using the stationary balance equations (2). The determinant |A(z)| is a sixth-degree polynomial in z and can be written as |A(z)| = z(1 - z)h(z), where h(z) is a fourth-degree polynomial in z whose coefficients depend on the system parameters. The explicit expression for h(z) is omitted because of its considerable length and complex structure. Setting the determinant |A(z)| equal to zero, we obtain that its roots are z = 0, z = 1, and the roots of the polynomial h(z). The root z = 1 corresponds to the normalization condition for the stationary probabilities: ∑ Gd(1) = G1(1) + G01 (1) + G1s1 (1) + G0s1 (1) + G02 (1) + G0s2 (1) + G-1s2 (1)+ d∈ + G-1(1) = ∑ lim |Ad(z)| = 1, (6) d∈ z→1 |A(z)| whereas the factor z corresponding to the root z = 0 appears in both the numerator and the denominator of the expressions for the generating functions Gd(z). Since the multiplicity of this root is the same in both, this factor cancels out, and no singularity arises at z = 0. Therefore, the root z = 0 does not affect the analytic properties of these functions. Since the generating functions Gd(z) must be analytic and finite for |z| ≤ 1, all roots of the equation |A(z)| = 0 that lie inside the unit circle and are different from z = 1 must also be roots of |Ad(z)|. Thus, the boundary probabilities are determined from two equations: one follows from the normalization condition (6), while the other follows from the requirement that, after substituting a root of the denominator |A(z)| lying inside the unit circle, the numerators |Ad(z)| also vanish. This ensures that the generating functions Gd(z) remain finite and allows the boundary probabilities P0,01 and P0,02 to be determined uniquely. We further note that a stationary regime exists only if the polynomial h(z) has exactly one root inside the unit circle |z| < 1. The location of the remaining roots outside the unit circle does not affect the existence of the stationary regime. Since, in general, an analytical study of the number of roots of the polynomial h(z) inside the unit circle is difficult due to the complicated dependence of its coefficients on the system parameters, the subsequent analysis is carried out under the assumption that a stationary regime exists (see Section 5). To determine the stability boundary, we consider the value of polynomial h(z) at z = 1 and introduce the notation g( ) = h(1), where g( ) is an eighth-degree polynomial in the parameter : 8 7 7 g( ) = ∑ i i = 8 8 + ∑ i i + 0 = (p1q2 + q1p2 + p1q2q1p2) 8 + ∑ i i - g2g2 2 2. (7) i=0 i=1 i=1 1 2 1 2 Since the stability boundary corresponds to a change in the sign of the polynomial g( ), it is necessary to consider the positive roots of the equation g( ) = 0. (8) However, an analytical determination of these roots is also difficult in the general case due to the complicated dependence of the coefficients of g( ) on the system parameters. At the same time, it follows from (7) that g(0) < 0, lim →+∞ g( ) = +∞. Therefore, by the Bolzano-Cauchy theorem, equation (8) has at least one positive root. Since the arrival rate is a non-negative parameter, the stability boundary corresponds to the smallest positive root of this equation: ∗ = min{ > 0 ∶ g( ) = 0}. Thus, the interval 0 < < ∗ defines the range of arrival rates for which the existence of a stationary regime is considered. Moreover, for all in this interval, the condition g( ) < 0 (9) is satisfied. In the numerical experiments presented in Section 6, we considered parameter sets for which condition (9) is satisfied. For each such set, we additionally checked that the polynomial h(z) has exactly one root inside the unit circle |z| < 1, which ensures the existence of a stationary regime. 5. Stationary performance measures We use the stationary probabilities defined earlier in (1) and the generating functions introduced in (3). Then, at z = 1, we obtain Gd(1) = ∑ Pn,d, d ∈ . n≥0 We determine the values of the first and second derivatives of the generating functions Gd(z) at z = 1: G′ (z) = d (∑ Pn,dzn) = ∑ nPn,d z n-1 ⇒ G′ (1) = ∑ nPn,d, d dz n≥0 n≥0 d n≥0 G″(z) = d (∑ nPn,d zn-1) = ∑ n(n - 1)Pn,d z n-2 ⇒ G″(1) = ∑ n(n - 1)Pn,d. d dz n≥0 n≥0 d n≥0 The stationary marginal probability of phase d is given by P⋅,d = lim ℙ{J(t) = d} = ∑ Pn,d = Gd(1), d ∈ . t→∞ n≥0 The stationary marginal probability that queue Q1 contains n customers is given by Pn,⋅ = lim ℙ{L1(t) = n} = ∑ Pn,d, n ≥ 0. t→∞ d∈ We now compute the mean queue lengths (L1) and (L2), as well as the mean difference (D). Since the numerical part of each phase d ∈ = {1, 1s1 , 01, 0s1 , 0s2 , 02, -1s2 , -1} corresponds to the value of the queue-length difference, by the definition of expectation (D) = ∑ d Gd(1) = ∑ d P⋅,d. (10) d∈ d∈ Expanding (10) in terms of the generating functions, we obtain (D) = 1 ⋅ G1(1) + 1 ⋅ G1s1 (1) + 0 ⋅ G01 (1) + 0 ⋅ G0s1 (1) + 0 ⋅ G0s2 (1) + 0 ⋅ G02 (1)- - 1 ⋅ G-1s2 (1) - 1 ⋅ G-1(1) = G1(1) + G1s1 (1) - G-1s2 (1) - G-1(1). Equivalently, expanding (10) in terms of the marginal phase probabilities, we obtain (D) = 1 ⋅ P⋅,1 + 1 ⋅ P⋅,1s1 + 0 ⋅ P⋅,01 + 0 ⋅ P⋅,0s1 + 0 ⋅ P⋅,0s2 + 0 ⋅ P⋅,02 - - 1 ⋅ P⋅,-1s2 - 1 ⋅ P⋅,-1 = P⋅,1 + P⋅,1s1 - P⋅,-1s2 - P⋅,-1. The mean length of the first queue is given by d (L1) = ∑ n Pn,⋅ = ∑ n ∑ Pn,d = ∑ ∑ n Pn,d = ∑ G′ (1). n≥0 n≥0 d∈ d∈ n≥0 d∈ The mean length of the second queue can be obtained from the relation between the queue lengths. Since D = L1 - L2, we have L2 = L1 - D. Hence, by the linearity of expectation, (L2) = (L1) - (D) = ∑ G′ (1) - ∑ d Gd(1) = ∑ (G′ (1) - d Gd(1)). d d∈ d∈ d d∈ The covariance of the queue lengths Cov(L1, L2) is defined as Cov(L1, L2) = (L1L2) - (L1) (L2). Since L2 = L1 - D, we have (L1L2) = (L1(L1 - D)) = (L2 - L1D) = (L2) - (L1D). 1 1 Therefore, where Cov(L1, L2) = (L2) - (L1D) - (L1) (L2), 1 1 (L2) = ∑ n2Pn,⋅ = ∑ ∑ n2Pn,d = ∑ ∑ n(n - 1)Pn,d + ∑ ∑ nPn,d = n≥0 d∈ n≥0 d∈ n≥0 d∈ n≥0 = ∑ G″(1) + ∑ nPn,⋅ = ∑ G″(1) + (L1) = ∑ G″(1) + ∑ G′ (1), d d∈ n≥0 d d∈ d d∈ d d∈ (L1D) = ∑ ∑ nd Pn,d = ∑ nPn,1 + ∑ nPn,1s1 - ∑ nPn,-1s2 - ∑ nPn,-1 = d∈ n≥0 n≥0 n≥0 n≥0 n≥0 = G′ (1) + G′ (1) - G′ (1) - G′ (1). 1 1s1 -1s2 -1 The variances of the numbers of customers in the first and second queues, Var(L1) and Var(L2), are computed as follows: 2 Var(L1) = (L2) - ( (L1))2 = ∑ G″(1) + ∑ G′ (1) - ( ∑ G′ (1)) = 1 d d∈ d d∈ d d∈ = ∑ G″(1) + (1 - ∑ G′ (1)) ∑ G′ (1), d d∈ d d∈ d d∈ Var(L2) = (L2) - ( (L2))2 = ((L1 - D)2) - ( (L1 - D))2 = (L2 - 2L1D + D2)- 2 Moreover, since 1 1 - ( (L1) - (D))2 = (L2) - 2 (L1D) + (D2) - ( (L1) - (D))2. D2 = 1, d ∈ {1, 1s1 , -1s2 , -1}, {0, d ∈ {01, 0s , 0s , 02}, 1 2 we have (D2) = G1(1) + G1s (1) + G-1s (1) + G-1(1). 1 2 The correlation coefficient Cor(L1, L2) is computed as Cor(L1, L2) = Cov(L1, L2) . √Var(L1) Var(L2) i Next, we determine the effective arrival rates eff, i = 1, 2, of customers joining queue Qi, i = 1, 2. For queue Q1: If L1 < L2, then an arriving customer joins queue Q1 with probability 1, and the corresponding contribution is (P⋅,-1s2 + P⋅,-1). If L1 = L2, then the customer joins queue Q1 with probability p1, and the corresponding contribution is p1(P⋅,01 + P⋅,0s1 + P⋅,0s2 + P⋅,02 ). Consequently, the effective arrival rate to queue Q1 is eff = (P⋅,-1 + P⋅,-1 + p1(P⋅,0 + P⋅,0 + P⋅,0 + P⋅,0 )). 1 s2 1 s1 s2 2 For queue Q2: If L1 > L2, then an arriving customer joins queue Q2 with probability 1, and the corresponding contribution is (P⋅,1 + P⋅,1s1 ). If L1 = L2, then the customer joins queue Q2 with probability p2, and the corresponding contribution is p2(P⋅,01 + P⋅,0s1 + P⋅,0s2 + P⋅,02 ). Consequently, the effective arrival rate to queue Q2 is eff = (P⋅,1 + P⋅,1 + p2(P⋅,0 + P⋅,0 + P⋅,0 + P⋅,0 )). 2 s1 1 s1 s2 2 Since the sum of the effective arrival rates to queues Q1 and Q2 equals the total arrival rate , we have eff + eff = . 1 2 The effective load factor of queue Qi, interpreted as the fraction of time during which the server serves this queue, is given by eff eff = i , i = 1, 2. i i The mean fraction of time during which the server is in a switch-over state is given by Psw = P⋅,1s1 + P⋅,0s1 + P⋅,0s2 + P⋅,-1s2 . The probability that the server is idle, that is, neither serving nor switching over, is defined as follows: Pidle = P0,0 + P0,0 = 1 - eff - eff - Psw. Since 1 2 1 2 ∞ Gd(z) = ∑ Pn,d zn = P0,d + P1,dz + P2,dz2 + ⋯ , |z| ≤ 1, n=0 substituting z = 0 into the expression for Gd(z) gives Gd(0) = P0,d + P1,d ⋅ 0 + P2,d ⋅ 02 + ⋯ = P0,d. It follows that the mean fraction of time during which the server is in a switch-over state while the system is empty is given by Psw0 = P0,0s1 + P0,0s2 = G0s1 (0) + G0s2 (0). Then the probability that the system is empty, regardless of the server state, is given by P0 = P0,01 + P0,0s1 + P0,0s2 + P0,02 = Pidle + Psw0 . 6. Numerical results This section presents the numerical results in two groups of tables. Tables 1-4 show how changes in the customer arrival rate, service rates, server switch-over rates, and tie-breaking probabilities affect the main stationary characteristics of the system. These characteristics include the mean queue lengths, the variances of the queue lengths, the correlation coefficient, the effective customer arrival rates, and the fractions of time during which the server is serving, switching over, or idle for different parameter values. Tables 5-7 present the stationary probabilities of the process phases, each of which corresponds to a particular value of the queue-length difference D and a particular server state, and show how these probabilities change as the customer arrival rate and the tie-breaking probabilities vary. We introduce the load parameter ∗ = , ∈ (0, 1), and set = ∗ for = 0.3, 0.5, 0.7, 0.85, 0.95 to evaluate the effect of the arrival rate on the stationary characteristics of the system (see Table 1). In Tables 2-7, the results are presented for a moderate load, = 0.7, that is, for = 0.7 ∗, where ∗ is the corresponding stability boundary. This choice makes it possible to illustrate the effect of the model parameters on the system characteristics more clearly. We recall that the parameter = ∗ i characterizes the load level relative to the stability boundary, whereas the quantities eff represent the fractions of time during which the server serves the corresponding queues Qi, i = 1, 2. Dependence of the main performance measures on the load parameter for ∗ ≈ 1.491719, 1 = 2 = 3, 1 = 2 = 3, p1 = p2 = q1 = q2 = 0.5. Table 1 (L1) (L2) Var(L1) Var(L2) Cor(L1, L2) Pidle 0.30 0.448 0.172 0.172 0.169 0.169 0.404 0.698 0.50 0.746 0.416 0.416 0.470 0.470 0.667 0.496 0.70 1.044 0.996 0.996 1.630 1.630 0.875 0.295 0.85 1.268 2.459 2.459 7.510 7.510 0.969 0.146 0.95 1.417 8.326 8.326 74.062 74.062 0.997 0.049 Table 2 Influence of the service rates on the main performance measures at = 0.7, 1 = 2 = 2, p1 = p2 = q1 = q2 = 0.5. 1 2 ∗ (L1) (L2) Cor(L1, L2) eff eff Psw Pidle 0.1 0.2 0.12815 1.189 1.088 0.897 0.403 0.247 0.053 0.297 1 3 0.87537 1.073 0.969 0.881 0.274 0.113 0.320 0.293 5 5 1.51539 1.037 1.037 0.887 0.106 0.106 0.505 0.283 0.5 10 0.76393 1.256 1.017 0.907 0.407 0.033 0.286 0.274 1 2 Table 2 shows the effect of the service rates 1 and 2 on the stationary characteristics of the system. The influence of the server switch-over rates 1 and 2 on the system characteristics is shown in Table 3. Table 4 presents the system characteristics for different combinations of the tie-breaking probabilities p1 and q1. Next, we consider the stationary probabilities of the process phases for different values of the load parameter (see Table 5). Due to rounding, the sums of probabilities in the rows may differ slightly from one. The effect of the tie-breaking probability p1 on the stationary distribution over the process phases is shown in Table 6. Finally, Table 7 shows how the tie-breaking probability q1 affects the stationary distribution over the process phases. The numerical results presented in this section are further discussed and interpreted in Section 7. 7. Discussion The results presented in Table 1 show that, as the load parameter increases, the mean queue lengths (L1) and (L2) grow significantly. As approaches one, the mean queue lengths increase sharply and nonlinearly. This is because the arrival rate becomes close to the stability boundary ∗, so that the server spends less time idle and more time serving customers or switching over between the queues. Under such conditions, even a small increase in the load may lead to customer accumulation in the queues and, consequently, to a sharp increase in the queue lengths. At the same time, the variances of the queue lengths grow faster than the corresponding mean values, especially for high values of . This indicates that the queue lengths become more variable as the load increases. Together with the growth of the mean queue lengths, this indicates a decrease in the efficiency of the system. Table 3 Influence of the switch-over rates on the main performance measures at = 0.7, 1 = 2 = 2, p1 = p2 = q1 = q2 = 0.5. 1 2 ∗ (L1) (L2) Cor(L1, L2) eff eff Psw Pidle 0.1 0.2 0.15689 1.269 1.174 0.919 0.025 0.030 0.688 0.257 1 3 0.89238 1.063 0.978 0.882 0.143 0.170 0.401 0.287 5 5 1.38724 1.020 1.020 0.877 0.243 0.243 0.213 0.301 0.5 10 0.78691 1.228 1.027 0.904 0.110 0.165 0.463 0.262 = 2 = 2, 1 = 2 = Table 4 2. p1 q1 ∗ (L1) (L2) Cor(L1, L2) eff eff Psw Pidle 0.2 0.2 1.02190 0.931 0.932 0.868 0.114 0.244 0.323 0.320 0.5 0.99603 0.931 1.057 0.880 0.134 0.214 0.355 0.296 0.8 0.93088 0.877 1.127 0.895 0.146 0.180 0.393 0.281 0.5 0.2 0.97267 1.037 0.909 0.878 0.148 0.192 0.362 0.298 0.5 0.99447 0.996 0.996 0.875 0.174 0.174 0.357 0.295 0.8 0.97267 0.909 1.037 0.878 0.192 0.148 0.362 0.298 0.8 0.2 0.93088 1.127 0.877 0.895 0.180 0.146 0.393 0.281 0.5 0.99603 1.057 0.931 0.880 0.214 0.134 0.355 0.296 0.8 1.02190 0.932 0.931 0.868 0.244 0.114 0.323 0.320 1 2 Influence of tie-breaking probabilities on the main performance measures at = 0.7, 1 1 2 Table 5 Stationary phase probabilities for different values of at ∗ ≈ 1.491719, 1 = 2 = 3, 1 = 2 = 3, p1 = p2 = q1 = q2 = 0.5. P⋅,1 P⋅,1s1 P⋅,01 P⋅,0s1 P⋅,0s2 P⋅,02 P⋅,-1s2 P⋅,-1 0.30 0.0596 0.0413 0.3640 0.0351 0.0351 0.3640 0.0413 0.0596 0.50 0.0854 0.0710 0.2867 0.0569 0.0569 0.2867 0.0710 0.0854 0.70 0.1031 0.1009 0.2185 0.0776 0.0776 0.2185 0.1009 0.1031 0.85 0.1120 0.1229 0.1725 0.0925 0.0925 0.1725 0.1229 0.1120 0.95 0.1164 0.1373 0.1441 0.1022 0.1022 0.1441 0.1373 0.1164 Table 6 Stationary phase probabilities for different values of p1 at = 0.7, 1 = 2 = 3, 1 = 2 = 3, q1 = 0.5. p1 ∗ P⋅,1 P⋅,1s1 P⋅,01 P⋅,0s1 P⋅,0s2 P⋅,02 P⋅,-1s2 P⋅,-1 0.2 1.459 0.1200 0.1455 0.0901 0.0268 0.1364 0.3444 0.0533 0.0835 0.5 1.491 0.1031 0.1009 0.2185 0.0776 0.0776 0.2185 0.1009 0.1031 0.8 1.459 0.0835 0.0533 0.3444 0.1364 0.0268 0.0901 0.1455 0.1200 Stationary phase probabilities for different values of q1 at = 0.7, 1 = 2 = 3, 1 = 2 = 3, p1 = 0.5. Table 7 q1 ∗ P⋅,1 P⋅,1s1 P⋅,01 P⋅,0s1 P⋅,0s2 P⋅,02 P⋅,-1s2 P⋅,-1 0.2 1.494 0.0660 0.0745 0.2081 0.0962 0.0589 0.2298 0.1257 0.1407 0.5 1.491 0.1031 0.1009 0.2185 0.0776 0.0776 0.2185 0.1009 0.1031 0.8 1.494 0.1407 0.1257 0.2298 0.0589 0.0962 0.2081 0.0745 0.0660 The correlation coefficient between the queue lengths also increases with and becomes close to one as the stability boundary is approached. This suggests that the queue lengths vary in a more coordinated way, indicating effective balancing due to the combined use of the JSQ and SLQ disciplines. The server idle probability, on the contrary, decreases as increases, since the server is less likely to be idle when the arrival rate becomes higher. Table 2 shows that, in the symmetric case 1 = 2 = 5, the mean queue lengths are equal, (L1) = (L2) = 1.037, as are the effective load factors, eff = eff = 0.106. This is consistent with 1 2 the symmetry of the model parameters and reflects the symmetric behavior of the two queues. 1 2 When the service rates are asymmetric, the load is redistributed between the queues. If the first queue is served more slowly than the second one, its mean length increases. This effect is clearly seen for 1 = 0.5, 2 = 10, where (L1) = 1.256, whereas (L2) = 1.017. Along with the increase in the mean queue length, the effective load factor of the first queue also increases to eff = 0.407, while the effective load factor of the second queue is much smaller, eff = 0.033. Thus, the queue with the lower service rate becomes more heavily loaded. The values of the correlation coefficient Cor(L1, L2) remain rather high in all the cases considered, which indicates that a positive dependence between the queue lengths is preserved regardless of whether the service rates are symmetric or asymmetric. This effect may be explained by the combined use of the JSQ and SLQ disciplines, which tend to balance the queue lengths by routing new customers to the shortest queue and giving service priority to the longest queue. It is also important to note how the fraction of time during which the server is in a switch-over state, Psw, changes. In the symmetric case ( 1 = 2 = 5), it attains the largest value among all the considered cases, Psw = 0.505. This is because, for equal service rates, the queue lengths become equal more often, and the server therefore switches between the queues more frequently. In the asymmetric case, for example, when 1 = 0.5 and 2 = 10, the fraction of switch-over time decreases to Psw = 0.286, since the server mostly remains at the slower and, consequently, more heavily loaded queue, and switches over less frequently than in the symmetric case. The idle probability Pidle also depends on the combination of service rates, although in the considered cases it remains approximately at the same level, in the range 0.274-0.297. The smallest value, Pidle = 0.274, is observed under strong asymmetry, 1 = 0.5 and 2 = 10. In this case, the first queue is served much more slowly than the second one and therefore becomes more heavily loaded, as reflected by eff = 0.407, while eff = 0.033. Thus, in the considered examples, changes in the 1 2 service rates affect the idle probability only moderately. Table 3 shows that, for small values of 1 and 2, the server spends most of its time in a switch-over state. For example, when 1 = 0.1 and 2 = 0.2, we have Psw = 0.688. This is due to the long switch- over process. Moreover, because the system operates under absolute preemptive priority, longer switch-overs increase the probability of arrivals during a switch-over, which further increases Psw, the fraction of time spent in switch-over. As a result, in this case the effective load factors of the queues are small, eff = 0.025 and eff = 0.030, while the mean queue lengths are relatively large. 1 2 As the switch-over rates increase, the fraction of time spent in switch-over generally decreases. For example, it decreases from Psw = 0.688 for 1 = 0.1, 2 = 0.2 to Psw = 0.213 for 1 = 2 = 5. This indicates that faster switch-overs reduce the time spent in switch-over and allow the server to serve customers more efficiently. Consequently, the mean queue lengths also decrease: from (L1) = 1.269 and (L2) = 1.174 for 1 = 0.1, 2 = 0.2 to (L1) = (L2) = 1.020 for 1 = 2 = 5. When the switch-over rates are asymmetric, the fraction of time spent in switch-over may remain rather high, especially when one of the switch-over rates is much smaller than the other. For example, when 1 = 0.5 and 2 = 10, we have Psw = 0.463, which is considerably higher than in the symmetric case 1 = 2 = 5, where Psw = 0.213. In this case, since 1 < 2, switch-overs to Q1 are slower, and queue Q1 becomes more heavily loaded in terms of its mean length: (L1) = 1.228, whereas (L2) = 1.027. The effective load factors are also unequal, eff = 0.110 and eff = 0.165. 1 2 The idle probability Pidle changes only moderately as the switch-over rates vary. In the considered cases, it remains in the range 0.257-0.301. The largest value, Pidle = 0.301, is observed in the symmetric case 1 = 2 = 5, where the fraction of time spent in switch-over is the smallest. Thus, in these examples, changes in the switch-over rates mainly affect the fraction of time spent in switch-over, whereas the idle probability varies only slightly. The correlation coefficient Cor(L1, L2), as in Table 2, remains high, which indicates that the queue lengths vary in a coordinated way. This can similarly be explained by the tendency of the JSQ and SLQ disciplines to balance the queue lengths even in the presence of non-zero switch-over times. When the queue lengths are equal, the parameters p1 and q1 define two different tie-breaking mechanisms. The parameter p1 controls the routing of arriving customers, whereas q1 controls the server’s choice of the queue to be served. As Table 4 shows, when these parameters differ from 0.5, the load tends to shift toward one of the queues. For example, for (p1, q1) = (0.2, 0.8), the second queue is more heavily loaded: (L2) = 1.127 and eff = 0.180, whereas (L1) = 0.877 and eff = 0.146. Since 2 1 the remaining parameters are symmetric, the case (p1, q1) = (0.8, 0.2) exhibits mirror symmetry: the same numerical values are interchanged between the first and second queues, namely (L1) = 1.127, eff = 0.180, (L2) = 0.877, and eff = 0.146. 1 2 The correlation coefficient Cor(L1, L2) remains high and positive in all considered cases. This indicates a coordinated behavior of the queue lengths and can be explained in the same way as in the previous tables. The mean fraction of time during which the server is in a switch-over state, Psw, also depends on the choice of p1 and q1, but varies within a moderate range. In the considered cases, it ranges from 0.323 to 0.393. The idle probability Pidle is less sensitive to changes in these parameters and remains between 0.281 and 0.320. Table 5 shows that, for a low load value ( = 0.3), the largest probabilities correspond to the phases 01 and 02. This means that the system spends most of the time in a balanced state, where the queue lengths are equal. As the load increases, the probabilities of the phases 01 and 02, where the queue lengths are equal and the server is at one of the queues, decrease. At the same time, the probabilities of the phases 1 and -1, as well as those of the switch-over phases 1s1 , 0s1 , 0s2 , and -1s2 , increase. This is because, as the arrival rate increases, the system has less time to restore equality between the queue lengths. It is important to note that, even as the load increases, a substantial part of the probability mass remains concentrated in the phases with D = 0. This reflects the effect of the JSQ and SLQ disciplines, which tend to equalize the queue lengths: arriving customers are routed to the shortest queue, while the server gives priority to the longest queue. As a result, the system regularly returns to phases with equal queue lengths. As the load parameter approaches one, the probability distribution becomes more evenly spread across the phases, while the probability that the server is in a switch-over phase increases. This is due to more frequent changes in the queue-length difference, which cause the server to switch between the queues more often. It should also be noted that, due to the symmetry of the model parameters, the probabilities of the corresponding phases for the first and second queues coincide. Table 6 shows that changing the parameter p1 leads to a redistribution of probabilities among the process phases. For p1 = 0.5, the distribution is symmetric: the probabilities of the corresponding phases for the first and second queues coincide. When p1 deviates from the symmetric value, the probability distribution shifts toward one of the queues. For example, for p1 = 0.2, the total probability of the phases associated with the second queue increases, while the first queue becomes less loaded. For p1 = 0.8, the opposite pattern is observed: the distribution shifts toward the first queue, while the second queue becomes less loaded. It can also be seen that changing p1 affects the probabilities of the switch-over phases. As p1 increases, the probabilities of the phases 0s1 and -1s2 become larger, whereas the probabilities of the phases 1s1 and 0s2 decrease. When p1 decreases, the opposite pattern is observed: the probabilities of the phases 0s1 and -1s2 become smaller, whereas the probabilities of the phases 1s1 and 0s2 increase. Thus, the parameter p1 has a significant effect on the distribution of phase probabilities and determines the direction of the load shift in the system. Table 7 shows that the parameter q1 affects the distribution of probabilities over the process phases, although this effect is weaker than that of the parameter p1 in Table 6. For q1 = 0.5, as in Table 6, the distribution is symmetric: the probabilities of the corresponding phases for the first and second queues coincide. When q1 differs from 0.5, the probability mass shifts toward one of the queues. For q1 = 0.2, the system is more often in phases associated with the second queue, whereas for q1 = 0.8, the shift is toward the first queue. These changes are most noticeable in the phases with unequal queue lengths. Thus, compared with p1, the parameter q1 plays a secondary role in the load distribution between the queues. 8. Conclusion This paper studies a two-queue queueing system with the JSQ-SLQ disciplines, absolute preemptive priority, and non-zero server switch-over times. A Markov model is constructed in terms of the length of the first queue and a phase variable that encodes both the queue-length difference and the server state. This makes it possible to apply the probability generating function method and to obtain analytical expressions for the main stationary characteristics, as well as a condition for the existence of a stationary regime. The obtained analytical results were used to numerically evaluate how the arrival rate, service rates, server switch-over rates, and tie-breaking parameters affect the stationary characteristics and phase probabilities. The numerical results show that the behavior of the system is shaped by its load, service rates, tie- breaking parameters, and server switch-over times. As the load increases, the mean queue lengths and their variances grow rapidly, which is related to the arrival rate approaching the stability boundary. At the same time, the positive correlation between the queue lengths increases, reflecting the balancing effect of the jointly applied JSQ and SLQ disciplines. It is shown that the service rates noticeably affect the load distribution between the queues. The queue with the lower service rate becomes more heavily loaded, which is reflected in an increase in its mean length and effective load factor. However, in the considered examples, the idle probability changes only moderately as the service rates vary. The server switch-over times also have a significant effect on the system characteristics. When switching is slow, the server spends a substantial fraction of time in switch-over, which leads to larger mean queue lengths. Faster switching, in contrast, reduces the time spent in switch-over and improves service efficiency. The tie-breaking parameters also affect the system behavior. The parameter p1 plays an important role in shifting the load between the queues and noticeably affects the distribution of phase probabilities. The parameter q1 also affects the phase probabilities, but its effect is weaker than that of p1, and it plays a secondary role in the load distribution between the queues. Overall, the combined use of the JSQ and SLQ disciplines provides stable equalization of queue lengths even in the presence of parameter asymmetry. However, the system performance is sensitive to the combination of the routing, service, and switch-over parameters. Further research may be directed toward generalizing the proposed model to systems with a larger number of queues and servers. Such a generalization would make it possible to analyze multi- server systems with a complex interaction structure, which is typical of distributed computing and telecommunication systems. In addition, it is of particular interest to study other rules for selecting a queue when the queue lengths are equal, as well as systems with different priority mechanisms, including preemptive and non-preemptive variants. Another relevant direction is the study of systems with more general distributions of service and server switch-over times, as well as non-Poisson arrival processes. This would make it possible to go beyond the Markovian assumptions and extend the applicability of the model. The model considered in this paper can be used to analyze and optimize queueing systems in which the routing of arriving customers and the selection of a queue for service depend on the current state of the system, for example, in telecommunication networks, data centers, and distributed computing systems.

About the authors

Evgeny A. Khityaev

RUDN University

Email: evgenykhityaev@yandex.ru
ORCID iD: 0000-0002-2375-1873
ResearcherId: PYW-1479-2026

Postgraduate Student of the Department of Probability Theory and Cyber Security of RUDN University

6 Miklukho-Maklaya St, Moscow, 117198, Russian Federation

Ivan S. Zaryadov

RUDN University; Institute of Informatics Problems, Federal Research Center “Computer Science and Control” of the Russian Academy of Sciences

Email: zaryadov-is@rudn.ru
ORCID iD: 0000-0002-7909-6396
Scopus Author ID: 35294470000
ResearcherId: B-8154-2018

Candidate of Physical and Mathematical Sciences, Assistant Professor of the Department of Probability Theory and Cyber Security of RUDN University; Senior Researcher of the Institute of Informatics Problems of the Federal Research Center “Computer Science and Control” of the Russian Academy of Sciences

6 Miklukho-Maklaya St, Moscow, 117198, Russian Federation; 44-2 Vavilova St, Moscow, 119333, Russian Federation

Tatiana A. Milovanova

RUDN University

Author for correspondence.
Email: milovanova-ta@rudn.ru
ORCID iD: 0000-0002-9388-9499
ResearcherId: X-6718-2018

Candidate of Physical and Mathematical Sciences, Assistant Professor of the Department of Probability Theory and Cyber Security of RUDN University

6 Miklukho-Maklaya St, Moscow, 117198, Russian Federation

References

Supplementary files