1 Introduction
The Erlang queue is a queueing system in which the arrival of customers is modeled according to a Poisson process with rate λ and its service system has Erlang distribution with shape parameter k and mean $1/\mu $. Its inter-arrival times are exponentially distributed with parameter λ and its service system consists of k phases, where the service time of each phase is exponentially distributed with parameter $k\mu $. In the field of telecommunication, the Erlang queues are quite useful. In call centers, the Erlang queue is used to predict number of agents required to handle incoming calls (see [29]).
In [21, 22] a single channel queue with Poisson arrivals and a general class of service-time distributions is studied. The transient solution of Erlang queue is obtained in [16]. Cahoy et al. [7] proposed the first fractional generalization of $M/M/1$ queue where they obtained its state probabilities and an algorithm to simulate fractional $M/M/1$ queue. Di Crescenzo et al. [9] studied $M/M/1$ queue with catastrophes and its fractional variant is discussed in [3]. Giorno et al. [14] studied a single-server queueing system with Poisson arrivals and state-dependent service mechanism characterized by logarithmic steady-state distribution. The Erlang queue time-changed with inverse stable subordinator is studied by Ascione et al. [4], where its state probabilities, mean queue length along with the distributions its of busy period, inter-arrival, inter-phase and service times are derived.
Let $\mathbb{N}$ denote the set of positive integers. The Erlang queue ${\{\mathcal{Q}(t)\}_{t\ge 0}}$ with state space ${\mathcal{H}_{0}}=\mathcal{H}\cup \{(0,0)\}$, where $\mathcal{H}=\{(n,s)\in \mathbb{N}\times \mathbb{N}:s\le k\}$ can be described as follows: $\mathcal{Q}(t)=(\mathcal{N}(t),\mathcal{S}(t))$. Here, $\mathcal{N}(t)$ denotes the number of customers in the system at time $t\ge 0$ and $\mathcal{S}(t)$ is the phase of customer being served at time t provided $\mathcal{N}(t)\gt 0$. For $\mathcal{N}(t)=0$, we take $\mathcal{S}(t)=0$.
For $t\ge 0$, let us denote the transient state probabilities of Erlang queue as follows: ${p_{0}}(t)=\mathrm{Pr}(\mathcal{Q}(t)=(0,0)|\mathcal{Q}(0)=(0,0))$ and ${p_{n,s}}(t)=\mathrm{Pr}(\mathcal{Q}(t)=(n,s)|\mathcal{Q}(0)=(0,0))$, $(n,s)\in \mathcal{H}$. That is, ${p_{0}}(t)$ is the probability that there are no customers in the system at time t and ${p_{n,s}}(t)$ is the probability that there are n customers in the system at time t and the customer being served is at phase s.
Equivalently, the Erlang queue can be represented as a queue length process in terms of phase count (see [4], p. 3252). Consider a bijective map ${g_{k}}:{\mathcal{H}_{0}}\to {\mathbb{N}_{0}}$ defined as follows:
where ${\mathbb{N}_{0}}=\mathbb{N}\cup \{0\}$. Its inverse map $({a_{k}}(m),{b_{k}}(m))$ is such that
and
Then, $\mathcal{L}(t)={g_{k}}(\mathcal{Q}(t))$ is the length of Erlang queue in terms of number of phases at time $t\ge 0$. Its state probabilities are denoted by ${p_{n}}(t)=\mathrm{Pr}(\mathcal{L}(t)=n|\mathcal{L}(0)=0)$, $n\ge 0$. Note that ${p_{n}}(t)={p_{{a_{k}}(n),{b_{k}}(n)}}(t)$.
(1)
\[ {g_{k}}(n,s):=\left\{\begin{array}{l@{\hskip10.0pt}l}k(n-1)+s,\hspace{0.1667em}\hspace{0.1667em}\hspace{0.1667em}\hspace{0.1667em}& (n,s)\in \mathcal{H},\\ {} 0,& (n,s)=(0,0),\end{array}\right.\](2)
\[ {b_{k}}(m)=\left\{\begin{array}{l@{\hskip10.0pt}l}\min \{s\gt 0:s\equiv m\hspace{0.3em}(\mathrm{mod} \hspace{0.3em}k)\},\hspace{0.1667em}& m\gt 0,\\ {} 0,\hspace{0.1667em}& m=0,\end{array}\right.\](3)
\[ {a_{k}}(m)=\left\{\begin{array}{l@{\hskip10.0pt}l}\frac{m-{b_{k}}(m)}{k}+1,\hspace{0.1667em}\hspace{0.1667em}& m\gt 0,\\ {} 0,& m=0.\end{array}\right.\]It is well known that the occurrence of catastrophic events such as earthquakes, tsunami, etc. exhibits long memory. It is empirically proven that these real life situations can be modeled by time-changed processes as they exhibits long range behaviour. Recently, the time-changed processes are extensively studied by many researchers, for example, the time fractional Poisson process (see [6, 23]), space fractional Poisson process (see [27]), time-changed birth-death processes (see [19, 26]), etc. Beghin [5] studied random-time processes governed by differential equations of fractional distributed order, Kataria and Khandakar [18] studied mixed fractional risk process which are based on the inverse mixed stable subordinator. Recently, Pote and Kataria [28] studied a time-changed variant of the Erlang queue with multiple arrivals where the time-changing component used is the inverse stable subordinator. For some more recent literature in this direction, we refer the reader to the works on multiparameter generalized counting process and its time-changed variants (see [10]), tempered space-time fractional negative binomial process (see [11]), bivariate gamma subordination for a Poisson shock model (see [25]), non-homogeneous generalized fractional Skellam process (see [30]), etc.
In a queue time-changed with inverse stable subordinator, the waiting times depends on a single fractional parameter, and hence has a fixed memory scale (see [4]). However, whenever a process is time-changed via an inverse mixed stable subordinator it mixes different waiting times. This makes it potentially applicable for modeling systems with distributed-order fractional dynamics, which cannot be captured by an inverse stable subordinator (see [5]).
In this paper, we introduce and study a time-changed variant of the Erlang queue where the time changing component used is the first hitting time of a mixed stable subordinator. We call it the mixed time-changed Erlang queue. The paper is organized as follows:
First, we give some preliminary results and known definitions of special functions, fractional derivative, subordinators and their inverse processes, etc. Then, we define the mixed time-changed Erlang queue by time-changing it with an inverse mixed stable subordinator which is the first hitting time of a mixed stable subordinator. The system of fractional differential equations that governs its state probabilities is derived. The Laplace transform of its state probabilities is obtained whose inversion yields its distribution. Also, the mixed time-changed Erlang queue is represented as a queue length process which characterizes it in terms of number of phases. The fractional differential equation that governs its mean queue length is derived, and an explicit expression for the mean queue length is obtained. Moreover, we derive the distribution of inter-arrival times, inter-phase times, service times and busy period of mixed time-changed Erlang queue, and discuss its conditional waiting time. Some plots of sample paths simulation are given.
2 Preliminaries
Here, we collect some known definitions and results that will be used.
2.1 Mittag-Leffler function
The three-parameter Mittag-Leffler function also known as the Prabhakar function is defined as (see [20])
where $\mathbb{R}$ denotes the set of real numbers.
(4)
\[ {E_{\alpha ,\beta }^{\gamma }}(t):={\sum \limits_{r=0}^{\infty }}\frac{\Gamma (\gamma +r){t^{r}}}{r!\Gamma (\gamma )\Gamma (\alpha r+\beta )},\hspace{0.1667em}t\in \mathbb{R},\hspace{0.1667em}\alpha ,\hspace{0.1667em}\beta ,\hspace{0.1667em}\gamma \gt 0,\]The following result is useful in the theory of stable subordinators (see [17]):
The following result will be used (see [17], Eq. (17.6)):
where $\alpha \gt \beta \gt 0$, $\alpha -\rho +1\gt 0$ and $|a{z^{\beta }}/({z^{\alpha }}+b)|\lt 1$.
(6)
\[ {\mathbb{L}^{-1}}\Big(\frac{{z^{\rho -1}}}{{z^{\alpha }}+a{z^{\beta }}+b};t\Big)={t^{\alpha -\rho }}{\sum \limits_{h=0}^{\infty }}{(-a)^{h}}{t^{(\alpha -\beta )h}}{E_{\alpha ,\alpha +(\alpha -\beta )h-\rho +1}^{h+1}}(-b{t^{\alpha }}),\]2.2 Fractional derivative
The Caputo fractional derivative of an absolutely continuous function $g(\cdot )$ is defined as follows (see [20]):
Its Laplace transform is given by
where $\tilde{g}(z)$ denotes the Laplace transform of the function $g(t)$.
(7)
\[ \frac{{\mathrm{d}^{\alpha }}}{\mathrm{d}{t^{\alpha }}}g(t)=\left\{\begin{array}{l}\frac{1}{\Gamma (1-\alpha )}{\textstyle\textstyle\int _{0}^{t}}{(t-y)^{-\alpha }}{g^{\prime }}(y)\mathrm{d}y,\hspace{0.1667em}0\lt \alpha \lt 1,\hspace{1em}\\ {} {g^{\prime }}(t),\hspace{0.1667em}\alpha =1.\hspace{1em}\end{array}\right.\](8)
\[ \mathbb{L}\Big(\frac{{\mathrm{d}^{\alpha }}}{\mathrm{d}{t^{\alpha }}}g(t)\Big)(z)={z^{\alpha }}\tilde{g}(z)-{z^{\alpha -1}}g(0),\]2.3 Stable subordinators and their inverse
An α-stable subordinator is a one dimensional Lévy process ${\{{D_{\alpha }}(t)\}_{t\ge 0}}$, $0\lt \alpha \lt 1$ whose Laplace transform is given by $\mathbb{E}({e^{-z{D_{\alpha }}(t)}})={e^{-t{z^{\alpha }}}}$, $z\gt 0$ (see [2]). Its first hitting time process ${\{{Y_{\alpha }}(t)\}_{t\ge 0}}$, where ${Y_{\alpha }}(t)=\inf \{s\ge 0:{D_{\alpha }}(s)\gt t\}$ is known as the inverse α-stable subordinator. The Laplace transform of its probability density function is given by (see [24])
The mixed stable subordinator ${\{{D_{{\alpha _{1}},{\alpha _{2}}}}(t)\}_{t\ge 0}}$ is a Lévy process whose Laplace transform is given by
where ${c_{1}}\ge 0$, ${c_{2}}\ge 0$, ${c_{1}}+{c_{2}}=1$, $0\lt {\alpha _{2}}\lt {\alpha _{1}}\lt 1$. Its first hitting time process ${\{{Y_{{\alpha _{1}},{\alpha _{2}}}}(t)\}_{t\ge 0}}$ is known as the inverse mixed stable subordinator, where
(9)
\[ \mathbb{E}({e^{-z{D_{{\alpha _{1}},{\alpha _{2}}}}(t)}})={e^{-t({c_{1}}{z^{{\alpha _{1}}}}+{c_{2}}{z^{{\alpha _{2}}}})}},\hspace{0.1667em}z\gt 0,\]
\[ {Y_{{\alpha _{1}},{\alpha _{2}}}}(t)=\text{inf}\{s\ge 0:{D_{{\alpha _{1}},{\alpha _{2}}}}(s)\gt t\},\hspace{0.1667em}t\gt 0\]
with ${Y_{{\alpha _{1}},{\alpha _{2}}}}(0)=0$. For ${c_{1}},{c_{2}}\in \{0,1\}$, it reduces to inverse stable subordinator. The probability density function ${f_{{\alpha _{1}},{\alpha _{2}}}}(t,u)=\frac{\mathrm{d}}{\mathrm{d}u}\mathrm{Pr}({Y_{{\alpha _{1}},{\alpha _{2}}}}(t)\le u)$ of inverse mixed stable subordinator has the following Laplace transform (see [1]):
3 Mixed time-changed Erlang queue
In this section, we introduce a time-changed Erlang queue where the time changing component is the first hitting time of mixed stable subordinator. We call it the mixed time-changed Erlang queue.
Let ${\{\mathcal{Q}(t)\}_{t\ge 0}}$ be an Erlang queue, that is, $\mathcal{Q}(t)=(\mathcal{N}(t),\mathcal{S}(t))$, $t\ge 0$, where $\mathcal{N}(t)$ denotes the number of customers in the system at time t and $\mathcal{S}(t)$ denotes the phase of customer being served at time t. We define the mixed time-changed Erlang queue ${\{{\mathcal{Q}^{{\alpha _{1}},{\alpha _{2}}}}(t)\}_{t\ge 0}}$, $0\lt {\alpha _{2}}\lt {\alpha _{1}}\lt 1$ as follows:
\[ {\mathcal{Q}^{{\alpha _{1}},{\alpha _{2}}}}(t):=\mathcal{Q}({Y_{{\alpha _{1}},{\alpha _{2}}}}(t)),\]
where the Erlang queue ${\{\mathcal{Q}(t)\}_{t\ge 0}}$ is independent of the inverse mixed stable subordinator ${\{{Y_{{\alpha _{1}},{\alpha _{2}}}}(t)\}_{t\ge 0}}$. For $t\ge 0$, let us denote its state probabilities by
\[\begin{aligned}{}{p_{0}^{{\alpha _{1}},{\alpha _{2}}}}(t)& =\mathrm{Pr}({\mathcal{Q}^{{\alpha _{1}},{\alpha _{2}}}}(t)=(0,0)|{\mathcal{Q}^{{\alpha _{1}},{\alpha _{2}}}}(0)=(0,0)),\\ {} {p_{n,s}^{{\alpha _{1}},{\alpha _{2}}}}(t)& =\mathrm{Pr}({\mathcal{Q}^{{\alpha _{1}},{\alpha _{2}}}}(t)=(n,s)|{\mathcal{Q}^{{\alpha _{1}},{\alpha _{2}}}}(0)=(0,0)),\hspace{0.1667em}(n,s)\in \mathcal{H},\end{aligned}\]
where $\mathcal{H}=\{(n,s)\in \mathbb{N}\times \mathbb{N}:s\le k\}$.In the next result, we obtain the system of fractional differential equations whose solution gives the state probabilities of mixed time-changed Erlang queue.
Theorem 1.
Let ${c_{1}}\ge 0$, ${c_{2}}\ge 0$ such that ${c_{1}}+{c_{2}}=1$ and $\frac{{\mathrm{d}^{\alpha }}}{\mathrm{d}{t^{\alpha }}}$ be the Caputo fractional derivative of order α as defined in (7). Then, the state probabilities of mixed time-changed Erlang queue solve the following system of differential equations:
\[\begin{aligned}{}\Big({c_{1}}\frac{{\mathrm{d}^{{\alpha _{1}}}}}{\mathrm{d}{t^{{\alpha _{1}}}}}+{c_{2}}\frac{{\mathrm{d}^{{\alpha _{2}}}}}{\mathrm{d}{t^{{\alpha _{2}}}}}\Big){p_{0}^{{\alpha _{1}},{\alpha _{2}}}}(t)& =-\lambda {p_{0}^{{\alpha _{1}},{\alpha _{2}}}}(t)+k\mu {p_{1,1}^{{\alpha _{1}},{\alpha _{2}}}}(t),\\ {} \Big({c_{1}}\frac{{\mathrm{d}^{{\alpha _{1}}}}}{\mathrm{d}{t^{{\alpha _{1}}}}}+{c_{2}}\frac{{\mathrm{d}^{{\alpha _{2}}}}}{\mathrm{d}{t^{{\alpha _{2}}}}}\Big){p_{1,s}^{{\alpha _{1}},{\alpha _{2}}}}(t)& =-(\lambda +k\mu ){p_{1,s}^{{\alpha _{1}},{\alpha _{2}}}}(t)+k\mu {p_{1,s+1}^{{\alpha _{1}},{\alpha _{2}}}}(t),\hspace{0.1667em}1\le s\le k-1,\\ {} \Big({c_{1}}\frac{{\mathrm{d}^{{\alpha _{1}}}}}{\mathrm{d}{t^{{\alpha _{1}}}}}+{c_{2}}\frac{{\mathrm{d}^{{\alpha _{2}}}}}{\mathrm{d}{t^{{\alpha _{2}}}}}\Big){p_{1,k}^{{\alpha _{1}},{\alpha _{2}}}}(t)& =-(\lambda +k\mu ){p_{1,k}^{{\alpha _{1}},{\alpha _{2}}}}(t)+k\mu {p_{2,1}^{{\alpha _{1}},{\alpha _{2}}}}(t)+\lambda {p_{0}^{{\alpha _{1}},{\alpha _{2}}}}(t),\\ {} \Big({c_{1}}\frac{{\mathrm{d}^{{\alpha _{1}}}}}{\mathrm{d}{t^{{\alpha _{1}}}}}+{c_{2}}\frac{{\mathrm{d}^{{\alpha _{2}}}}}{\mathrm{d}{t^{{\alpha _{2}}}}}\Big){p_{n,s}^{{\alpha _{1}},{\alpha _{2}}}}(t)& =-(\lambda +k\mu ){p_{n,s}^{{\alpha _{1}},{\alpha _{2}}}}(t)+k\mu {p_{n,s+1}^{{\alpha _{1}},{\alpha _{2}}}}(t)+\lambda {p_{n-1,s}^{{\alpha _{1}},{\alpha _{2}}}}(t),\\ {} & \hspace{130.88284pt}1\le s\le k-1,\hspace{0.1667em}n\ge 2,\\ {} \Big({c_{1}}\frac{{\mathrm{d}^{{\alpha _{1}}}}}{\mathrm{d}{t^{{\alpha _{1}}}}}+{c_{2}}\frac{{\mathrm{d}^{{\alpha _{2}}}}}{\mathrm{d}{t^{{\alpha _{2}}}}}\Big){p_{n,k}^{{\alpha _{1}},{\alpha _{2}}}}(t)& =-(\lambda +k\mu ){p_{n,k}^{{\alpha _{1}},{\alpha _{2}}}}(t)+k\mu {p_{n+1,1}^{{\alpha _{1}},{\alpha _{2}}}}(t)\\ {} & \hspace{8.5359pt}+\lambda {p_{n-1,k}^{{\alpha _{1}},{\alpha _{2}}}}(t),\hspace{0.1667em}n\ge 2,\end{aligned}\]
with ${p_{0}^{{\alpha _{1}},{\alpha _{2}}}}(0)=1$ and ${p_{n,s}^{{\alpha _{1}},{\alpha _{2}}}}(0)=0$ for $1\le s\le k$, $n\ge 1$.
Proof.
For $(n,s)\in \mathcal{H}$, we have
and
where ${p_{n,s}}(y)$ and ${p_{0}}(y)$ are the state probability and zero state probability of classical Erlang queue, respectively.
(11)
\[\begin{aligned}{}{p_{n,s}^{{\alpha _{1}},{\alpha _{2}}}}(t)& =\mathrm{Pr}(\mathcal{Q}({Y_{{\alpha _{1}},{\alpha _{2}}}}(t))=(n,s)|\mathcal{Q}(0)=(0,0))\\ {} & ={\int _{0}^{\infty }}\mathrm{Pr}(\mathcal{Q}(y)=(n,s)|\mathcal{Q}(0)=(0,0))\mathrm{Pr}({Y_{{\alpha _{1}},{\alpha _{2}}}}(t)\in \mathrm{d}y)\\ {} & ={\int _{0}^{\infty }}{p_{n,s}}(y)\mathrm{Pr}({Y_{{\alpha _{1}},{\alpha _{2}}}}(t)\in \mathrm{d}y),\hspace{0.1667em}t\ge 0\end{aligned}\](12)
\[ {p_{0}^{{\alpha _{1}},{\alpha _{2}}}}(t)={\int _{0}^{\infty }}{p_{0}}(y)\mathrm{Pr}({Y_{{\alpha _{1}},{\alpha _{2}}}}(t)\in \mathrm{d}y),\hspace{0.1667em}t\ge 0\]Let the generating function of ${\{{\mathcal{Q}^{{\alpha _{1}},{\alpha _{2}}}}(t)\}_{t\ge 0}}$ be given by
On substituting (11) and (12) in (13), we get
where $G(x,t)$ is the generating function of ${\{\mathcal{Q}(t)\}_{t\ge 0}}$. On taking the Laplace transform of (12) and (14), and using (10), we get
and
(13)
\[ {G^{{\alpha _{1}},{\alpha _{2}}}}(x,t)={p_{0}^{{\alpha _{1}},{\alpha _{2}}}}(t)+{\sum \limits_{n=1}^{\infty }}{\sum \limits_{s=1}^{k}}{p_{n,s}^{{\alpha _{1}},{\alpha _{2}}}}(t){x^{(n-1)k+s}},\hspace{0.1667em}t\ge 0,\hspace{0.1667em}|x|\le 1.\](14)
\[ {G^{{\alpha _{1}},{\alpha _{2}}}}(x,t)={\int _{0}^{\infty }}G(x,y)\mathrm{Pr}({Y_{{\alpha _{1}},{\alpha _{2}}}}(t)\in \mathrm{d}y),\hspace{0.1667em}t\ge 0,\]On multiplying the first equation by x, the second by ${x^{s+1}}$, the third by ${x^{k+1}}$, the fourth by ${x^{k(n-1)+s+1}}$ and the fifth by ${x^{nk+1}}$ in the system of differential equations given in Theorem 1, and then summing these equations over the range of n and s, it can be established that the state probabilities of ${\{{\mathcal{Q}^{{\alpha _{1}},{\alpha _{2}}}}(t)\}_{t\ge 0}}$ are solution of the system of differential equations given in Theorem 1 if and only if its generating function is the solution of following equation:
with initial condition ${G^{{\alpha _{1}},{\alpha _{2}}}}(x,0)=1$. On taking the Laplace transform on both sides of (17) and using (8), we obtain
By using (15) and (16) in (18), we have
Now, on substituting Eq. (6) of [16] in (19), we get the following identity:
(17)
\[\begin{aligned}{}& x\Big({c_{1}}\frac{{\mathrm{d}^{{\alpha _{1}}}}}{\mathrm{d}{t^{{\alpha _{1}}}}}+{c_{2}}\frac{{\mathrm{d}^{{\alpha _{2}}}}}{\mathrm{d}{t^{{\alpha _{2}}}}}\Big){G^{{\alpha _{1}},{\alpha _{2}}}}(x,t)\\ {} & \hspace{1em}=(1-x)\big((k\mu -\lambda (x+\cdots +{x^{k}})){G^{{\alpha _{1}},{\alpha _{2}}}}(x,t)-k\mu {p_{0}^{{\alpha _{1}},{\alpha _{2}}}}(t)\big),\end{aligned}\](18)
\[\begin{aligned}{}x({c_{1}}{z^{{\alpha _{1}}-1}}+& {c_{2}}{z^{{\alpha _{2}}-1}})(z{\tilde{G}^{{\alpha _{1}},{\alpha _{2}}}}(x,z)-1)\\ {} & =(1-x)\big((k\mu -\lambda (x+\cdots +{x^{k}})){\tilde{G}^{{\alpha _{1}},{\alpha _{2}}}}(x,z)-k\mu {\tilde{{p_{0}}}^{{\alpha _{1}},{\alpha _{2}}}}(z)\big),\hspace{0.1667em}z\gt 0.\end{aligned}\](19)
\[\begin{aligned}{}x\Big(& ({c_{1}}{z^{{\alpha _{1}}}}+{c_{2}}{z^{{\alpha _{2}}}}){\int _{0}^{\infty }}G(x,y){e^{-y({c_{1}}{z^{{\alpha _{1}}}}+{c_{2}}{z^{{\alpha _{2}}}})}}\mathrm{d}y-1\Big)\\ {} =& {\int _{0}^{\infty }}\hspace{-0.1667em}(1-x)\big((k\mu -\lambda (x+\cdots +{x^{k}}))G(x,y)-k\mu {p_{0}}(y)\big){e^{-y({c_{1}}{z^{{\alpha _{1}}}}+{c_{2}}{z^{{\alpha _{2}}}})}}\mathrm{d}y,\hspace{0.1667em}z\gt 0.\end{aligned}\]
\[\begin{aligned}{}x\Big(({c_{1}}{z^{{\alpha _{1}}}}+{c_{2}}{z^{{\alpha _{2}}}}){\int _{0}^{\infty }}G(x,y)& {e^{-y({c_{1}}{z^{{\alpha _{1}}}}+{c_{2}}{z^{{\alpha _{2}}}})}}\mathrm{d}y-1\Big)\\ {} & =x{\int _{0}^{\infty }}\frac{\partial }{\partial y}G(x,y){e^{-y({c_{1}}{z^{{\alpha _{1}}}}+{c_{2}}{z^{{\alpha _{2}}}})}}\mathrm{d}y,\hspace{0.1667em}z\gt 0.\end{aligned}\]
This proves the required result. □In Theorem 1, the Caputo fractional derivatives can be written as follows:
\[ \Big({c_{1}}\frac{{\mathrm{d}^{{\alpha _{1}}}}}{\mathrm{d}{t^{{\alpha _{1}}}}}+{c_{2}}\frac{{\mathrm{d}^{{\alpha _{2}}}}}{\mathrm{d}{t^{{\alpha _{2}}}}}\Big)f(t)={\int _{0}^{1}}\frac{{\mathrm{d}^{\alpha }}}{\mathrm{d}{t^{\alpha }}}f(t)c(\alpha )\mathrm{d}\alpha ,\]
where $c(\alpha )={c_{1}}\delta (\alpha -{\alpha _{1}})+{c_{2}}\delta (\alpha -{\alpha _{2}})$, ${c_{1}}\ge 0$, ${c_{2}}\ge 0$, ${c_{1}}+{c_{2}}=1$ and $0\lt {\alpha _{2}}\lt {\alpha _{1}}\lt 1$. So, the fractional derivative depends on the entire time span $[0,t]$ via the weight function involved in the definition of Caputo derivative given in (7). Also, the memory depends on the two fractional parameters ${\alpha _{1}}$ and ${\alpha _{2}}$. For instance, with the increase or decrease in the values of ${\alpha _{1}}$ and ${\alpha _{2}}$, the memory effects can increase or decrease.In the following result, we obtain the Laplace transform of zero state probability of mixed time-changed Erlang queue which we use to derive its explicit expression:
Proposition 1.
The Laplace transform of zero state probability of ${\{{\mathcal{Q}^{{\alpha _{1}},{\alpha _{2}}}}(t)\}_{t\ge 0}}$ is given by
where
and
(20)
\[ {\tilde{p}_{0}^{{\alpha _{1}},{\alpha _{2}}}}(z)={\sum \limits_{m=1}^{\infty }}{\sum \limits_{r=0}^{\infty }}{A_{m,r}^{0}}\frac{({c_{1}}{z^{{\alpha _{1}}-1}}+{c_{2}}{z^{{\alpha _{2}}-1}})}{{(\lambda +k\mu +{c_{1}}{z^{{\alpha _{1}}}}+{c_{2}}{z^{{\alpha _{2}}}})^{{a_{m,r}^{0}}}}},\hspace{0.1667em}z\gt 0\]Proof.
By using Eq. (5) of [4] in (15), we obtain
On substituting (21) and (22) in (23), we get the required result. □
(23)
\[\begin{aligned}{}{\tilde{{p_{0}}}^{{\alpha _{1}},{\alpha _{2}}}}(z)& ={\sum \limits_{m=1}^{\infty }}{\sum \limits_{r=0}^{\infty }}\frac{m{\lambda ^{r}}{(k\mu )^{m+rk-1}}}{r!(m+rk)!}({c_{1}}{z^{{\alpha _{1}}-1}}+{c_{2}}{z^{{\alpha _{2}}-1}}){\int _{0}^{\infty }}{y^{m+r(k+1)-1}}\\ {} & \hspace{142.26378pt}\cdot {e^{-(\lambda +k\mu +{c_{1}}{z^{{\alpha _{1}}}}+{c_{2}}{z^{{\alpha _{2}}}})y}}\mathrm{d}y\\ {} & ={\sum \limits_{m=1}^{\infty }}{\sum \limits_{r=0}^{\infty }}\frac{m{\lambda ^{r}}{(k\mu )^{m+rk-1}}}{r!(m\hspace{-0.1667em}+\hspace{-0.1667em}rk)!}\frac{({c_{1}}{z^{{\alpha _{1}}\hspace{-0.1667em}-\hspace{-0.1667em}1}}\hspace{-0.1667em}+\hspace{-0.1667em}{c_{2}}{z^{{\alpha _{2}}-1}})(m\hspace{-0.1667em}+\hspace{-0.1667em}r(k\hspace{-0.1667em}+\hspace{-0.1667em}1)\hspace{-0.1667em}-\hspace{-0.1667em}1)!}{{(\lambda \hspace{-0.1667em}+\hspace{-0.1667em}k\mu \hspace{-0.1667em}+\hspace{-0.1667em}{c_{1}}{z^{{\alpha _{1}}}}\hspace{-0.1667em}+\hspace{-0.1667em}{c_{2}}{z^{{\alpha _{2}}}})^{m+r(k+1)}}},\hspace{0.1667em}z\gt 0.\end{aligned}\]For a positive integer N, let us define
and
(24)
\[ {f_{N}}(t)={t^{{\alpha _{1}}-\big(\frac{{\alpha _{1}}-1}{N}\big)-1}}{\sum \limits_{h=0}^{\infty }}{\Big(-\frac{{c_{2}}}{{c_{1}}}\Big)^{h}}{t^{({\alpha _{1}}-{\alpha _{2}})h}}{E_{{\alpha _{1}},{\alpha _{1}}+({\alpha _{1}}-{\alpha _{2}})h+\frac{1-{\alpha _{1}}}{N}}^{h+1}}\Big(-\Big(\frac{\lambda +k\mu }{{c_{1}}}\Big){t^{{\alpha _{1}}}}\Big)\](25)
\[ {g_{N}}(t)={t^{{\alpha _{1}}-\big(\frac{{\alpha _{2}}-1}{N}\big)-1}}{\sum \limits_{h=0}^{\infty }}{\Big(-\frac{{c_{2}}}{{c_{1}}}\Big)^{h}}{t^{({\alpha _{1}}-{\alpha _{2}})h}}{E_{{\alpha _{1}},{\alpha _{1}}+({\alpha _{1}}-{\alpha _{2}})h+\frac{1-{\alpha _{2}}}{N}}^{h+1}}\Big(-\Big(\frac{\lambda +k\mu }{{c_{1}}}\Big){t^{{\alpha _{1}}}}\Big).\]The next result gives the probability of the system being empty at time t.
Theorem 2.
Let ${f^{\ast (n)}}$ denote the n-fold convolution of any function $f(\cdot )$ and ${E_{\alpha ,\beta }^{\gamma }}(\cdot )$ be the three parameter Mittag-Leffler function defined in (4). Then, the zero state probability of ${\{{\mathcal{Q}^{{\alpha _{1}},{\alpha _{2}}}}(t)\}_{t\ge 0}}$ is given by
where ${a_{m,r}^{0}}$, ${A_{m,r}^{0}}$, ${f_{{a_{m,r}^{0}}}}(\cdot )$ and ${g_{{a_{m,r}^{0}}}}(\cdot )$ are given in (21), (22), (24), and (25) respectively.
(26)
\[ {p_{0}^{{\alpha _{1}},{\alpha _{2}}}}(t)={\sum \limits_{m=1}^{\infty }}{\sum \limits_{r=0}^{\infty }}\frac{{A_{m,r}^{0}}}{{c_{1}^{{a_{m,r}^{0}}}}}\Big({c_{1}}{f_{{a_{m,r}^{0}}}^{\ast ({a_{m,r}^{0}})}}(t)+{c_{2}}{g_{{a_{m,r}^{0}}}^{\ast ({a_{m,r}^{0}})}}(t)\Big),\hspace{0.1667em}t\ge 0,\]Proof.
On taking the inverse Laplace transform of (20), we obtain
On using (6) in (27), we get the required result. □
(27)
\[\begin{aligned}{}{p_{0}^{{\alpha _{1}},{\alpha _{2}}}}(t)& ={\sum \limits_{m=1}^{\infty }}{\sum \limits_{r=0}^{\infty }}{A_{m,r}^{0}}\Bigg(\frac{1}{{c_{1}^{{a_{m,r}^{0}}-1}}}{\mathbb{L}^{-1}}\Bigg({\bigg(\frac{{z^{\frac{{\alpha _{1}}-1}{{a_{m,r}^{0}}}}}}{{z^{{\alpha _{1}}}}+\frac{{c_{2}}}{{c_{1}}}{z^{{\alpha _{2}}}}+\frac{\lambda +k\mu }{{c_{1}}}}\bigg)^{{a_{m,r}^{0}}}};t\Bigg)\\ {} & \hspace{56.9055pt}+\frac{{c_{2}}}{{c_{1}^{{a_{m,r}^{0}}}}}{\mathbb{L}^{-1}}\Bigg({\bigg(\frac{{z^{\frac{{\alpha _{2}}-1}{{a_{m,r}^{0}}}}}}{{z^{{\alpha _{1}}}}+\frac{{c_{2}}}{{c_{1}}}{z^{{\alpha _{2}}}}+\frac{\lambda +k\mu }{{c_{1}}}}\bigg)^{{a_{m,r}^{0}}}};t\Bigg)\Bigg).\end{aligned}\]Remark 1.
In (26), if ${c_{1}}=1$ then ${c_{2}}=0$. So,
and
On taking the Laplace transform of (29) and using (5), we have
Thus, the Laplace transform of ${a_{0}^{m,r}}$-fold convolution of (30) is given by
In (28), on taking the inverse Laplace transform of (31) and by using (5), we get the zero state probability of fractional Erlang queue (see [4], Eq. (36)). Similar result holds if ${c_{2}}=1$.
(28)
\[ {p_{0}^{{\alpha _{1}},{\alpha _{2}}}}(t)={\sum \limits_{m=1}^{\infty }}{\sum \limits_{r=0}^{\infty }}{A_{m,r}^{0}}{f_{{a_{m,r}^{0}}}^{\ast ({a_{m,r}^{0}})}}(t),\hspace{0.1667em}t\ge 0\](29)
\[ {f_{{a_{m,r}^{0}}}}(t)={t^{{\alpha _{1}}-\big(\frac{{\alpha _{1}}-1}{{a_{m,r}^{0}}}\big)-1}}{E_{{\alpha _{1}},{\alpha _{1}}+\frac{1-{\alpha _{1}}}{{a_{m,r}^{0}}}}^{1}}(-(\lambda +k\mu ){t^{{\alpha _{1}}}}),\hspace{0.1667em}t\ge 0.\](30)
\[ {\tilde{f}_{{a_{m,r}^{0}}}}(z)=\frac{{z^{\frac{{\alpha _{1}}-1}{{a_{m,r}^{0}}}}}}{{z^{{\alpha _{1}}}}+\lambda +k\mu },\hspace{0.1667em}z\gt 0.\](31)
\[ {\tilde{f}_{{a_{m,r}^{0}}}^{\ast ({a_{m,r}^{0}})}}(z)=\frac{{z^{{\alpha _{1}}-1}}}{{({z^{{\alpha _{1}}}}+\lambda +k\mu )^{{a_{m,r}^{0}}}}},\hspace{0.1667em}z\gt 0.\]In the subsequent results, we will be using the following notations:
where ${a_{m,r}^{0}}$ and ${A_{m,r}^{0}}$ are as given in (21) and (22), respectively.
(32)
\[ \left.\begin{aligned}{}{B_{i}^{n,s}}& =\frac{{\lambda ^{n+i}}{(k\mu )^{k(i+1)-s}}}{(k(i+1)-s)!(n+i)!}(n+k-s+i(k+1))!,\\ {} {\delta _{i}^{n,s}}& =n-s+(i+1)(k+1),\\ {} {C_{i,m,r}^{n,s}}& =k\mu {B_{i}^{n,s}}{A_{m,r}^{0}},\\ {} {\nu _{i,m,r}^{n,s}}& ={\delta _{i}^{n,s}}+{a_{m,r}^{0}},\\ {} {D_{i,m,r}^{n,s}}& =\left\{\begin{array}{l}k\mu {B_{i}^{n,s+1}}{A_{m,r}^{0}},\hspace{0.1667em}s=1,2,\dots ,k-1,\hspace{1em}\\ {} k\mu {B_{i}^{n+1,1}}{A_{m,r}^{0}},\hspace{0.1667em}s=k,\hspace{1em}\end{array}\right.\\ {} {\pi _{i,m,r}^{n,s}}& =\left\{\begin{array}{l}{\delta _{i}^{n,s+1}}+{a_{m,r}^{0}},\hspace{0.1667em}s=1,2,\dots ,k-1,\hspace{1em}\\ {} {\delta _{i}^{n+1,1}}+{a_{m,r}^{0}},\hspace{0.1667em}s=k,\hspace{1em}\end{array}\right.\end{aligned}\right\}\]In the following result, we obtain the Laplace transform of state probabilities of mixed time-changed Erlang queue which we use to derive its explicit expressions.
Proposition 2.
The Laplace transform of the state probabilities of ${\{{\mathcal{Q}^{{\alpha _{1}},{\alpha _{2}}}}(t)\}_{t\ge 0}}$ is given by
(33)
\[\begin{aligned}{}{\tilde{p}_{n,s}^{{\alpha _{1}},{\alpha _{2}}}}(z)& =({c_{1}}{z^{{\alpha _{1}}-1}}+{c_{2}}{z^{{\alpha _{2}}-1}})\Big({\sum \limits_{i=0}^{\infty }}\frac{{B_{i}^{n,s}}}{{(\lambda +k\mu +{c_{1}}{z^{{\alpha _{1}}}}+{c_{2}}{z^{{\alpha _{2}}}})^{{\delta _{i}^{n,s}}}}}\\ {} & \hspace{8.5359pt}+{\sum \limits_{i=0}^{\infty }}{\sum \limits_{m=1}^{\infty }}{\sum \limits_{r=0}^{\infty }}\frac{{C_{i,m,r}^{n,s}}}{{(\lambda +k\mu +{c_{1}}{z^{{\alpha _{1}}}}+{c_{2}}{z^{{\alpha _{2}}}})^{{\nu _{i,m,r}^{n,s}}}}}\\ {} & \hspace{8.5359pt}-{\sum \limits_{i=0}^{\infty }}{\sum \limits_{m=1}^{\infty }}{\sum \limits_{r=0}^{\infty }}\frac{{D_{i,m,r}^{n,s}}}{{(\lambda +k\mu +{c_{1}}{z^{{\alpha _{1}}}}+{c_{2}}{z^{{\alpha _{2}}}})^{{\pi _{i,m,r}^{n,s}}}}}\Big),\hspace{0.1667em}n\ge 1,\hspace{0.1667em}1\le s\le k,\hspace{0.1667em}z\gt 0.\end{aligned}\]Proof.
On taking the Laplace transform of (11) and using (10), we get
For $n\ge 1$ and $1\le s\le k-1$, by using Eq. (6) of [4] in (34), we obtain
Also, by using Eq. (7) of [4] in (34), we have
By using Lemma 5.4 of [4] and (32) in (35) and (36), we get the required result. □
(34)
\[ {\tilde{p}_{n,s}^{{\alpha _{1}},{\alpha _{2}}}}(z)=({c_{1}}{z^{{\alpha _{1}}-1}}+{c_{2}}{z^{{\alpha _{2}}-1}}){\int _{0}^{\infty }}{p_{n,s}}(y){e^{-y({c_{1}}{z^{{\alpha _{1}}}}+{c_{2}}{z^{{\alpha _{2}}}})}}\mathrm{d}y,\hspace{0.1667em}z\gt 0.\](35)
\[\begin{aligned}{}{\tilde{p}_{n,s}^{{\alpha _{1}},{\alpha _{2}}}}(z)& =({c_{1}}{z^{{\alpha _{1}}-1}}+{c_{2}}{z^{{\alpha _{2}}-1}})\Big({\sum \limits_{i=0}^{\infty }}\frac{{\lambda ^{n+i}}{(k\mu )^{k(i+1)-s}}}{(k(i+1)-s)!(n+i)!}{\int _{0}^{\infty }}{y^{n+k-s+i(k+1)}}\\ {} & \hspace{3.33333pt}\hspace{3.33333pt}\cdot {e^{-(\lambda +k\mu +{c_{1}}{z^{{\alpha _{1}}}}+{c_{2}}{z^{{\alpha _{2}}}})y}}\mathrm{d}y+{\sum \limits_{i=0}^{\infty }}\frac{{\lambda ^{n+i}}{(k\mu )^{k(i+1)-s+1}}}{(k(i+1)-s)!(n+i)!}\\ {} & \hspace{3.33333pt}\hspace{3.33333pt}\cdot {\int _{0}^{\infty }}{\int _{0}^{y}}{p_{0}}(u){(y-u)^{n+k-s+i(k+1)}}{e^{-(\lambda +k\mu )(y-u)}}{e^{-({c_{1}}{z^{{\alpha _{1}}}}+{c_{2}}{z^{{\alpha _{2}}}})y}}\mathrm{d}u\mathrm{d}y\\ {} & \hspace{3.33333pt}\hspace{3.33333pt}-{\sum \limits_{i=0}^{\infty }}\frac{{\lambda ^{n+i}}{(k\mu )^{k(i+1)-s}}}{(k(i+1)-s-1)!(n+i)!}{\int _{0}^{\infty }}{\int _{0}^{y}}{p_{0}}(u){(y-u)^{n+k-s-1+i(k+1)}}\\ {} & \hspace{3.33333pt}\hspace{3.33333pt}\cdot {e^{-(\lambda +k\mu )(y-u)}}{e^{-({c_{1}}{z^{{\alpha _{1}}}}+{c_{2}}{z^{{\alpha _{2}}}})y}}\mathrm{d}u\mathrm{d}y\Big),\hspace{0.1667em}z\gt 0.\end{aligned}\](36)
\[\begin{aligned}{}& {\tilde{p}_{n,k}^{{\alpha _{1}},{\alpha _{2}}}}(z)\\ {} & \hspace{1em}=({c_{1}}{z^{{\alpha _{1}}-1}}+{c_{2}}{z^{{\alpha _{2}}-1}})\Big({\sum \limits_{i=0}^{\infty }}\frac{{\lambda ^{n+i}}{(k\mu )^{ki}}}{(ki)!(n+i)!}{\int _{0}^{\infty }}{y^{n+i(k+1)}}{e^{-(\lambda +k\mu +{c_{1}}{z^{{\alpha _{1}}}}+{c_{2}}{z^{{\alpha _{2}}}})y}}\mathrm{d}y\\ {} & \hspace{2em}+{\sum \limits_{i=0}^{\infty }}\frac{{\lambda ^{n+i}}{(k\mu )^{ki+1}}}{(ki)!(n+i)!}{\int _{0}^{\infty }}{\int _{0}^{y}}{p_{0}}(u){(y-u)^{n+i(k+1)}}{e^{-(\lambda +k\mu )(y-u)}}\\ {} & \hspace{2em}\cdot {e^{-({c_{1}}{z^{{\alpha _{1}}}}+{c_{2}}{z^{{\alpha _{2}}}})y}}\mathrm{d}u\mathrm{d}y-{\sum \limits_{i=0}^{\infty }}\frac{{\lambda ^{n+i+1}}{(k\mu )^{k(i+1)}}}{(k(i+1)-1)!(n+i+1)!}\\ {} & \hspace{2em}\cdot {\int _{0}^{\infty }}{\int _{0}^{y}}{p_{0}}(u){(y-u)^{n+k+i(k+1)}}{e^{-(\lambda +k\mu )(y-u)}}{e^{-({c_{1}}{z^{{\alpha _{1}}}}+{c_{2}}{z^{{\alpha _{2}}}})y}}\mathrm{d}u\mathrm{d}y\Big),\hspace{0.1667em}z\gt 0.\end{aligned}\]Remark 2.
For ${c_{1}}=1$ or ${c_{2}}=1$ in (33), the Laplace transform of state probabilities of mixed time-changed Erlang queue reduces to that of fractional Erlang queue (see [4], Theorem 5.5). Observe that the state probabilities of fractional Erlang queue have dependence only on one fractional parameter, whereas that of mixed time-changed Erlang queue have dependence on two distinct fractional parameters. That is, past memory effect is due to two distinct parameters in mixed time-changed Erlang queue.
Next, we obtain the probability that there are $n\ge 1$ customers in the system and the customer being served is in phase $s\in \{1,2,\dots ,k\}$ at time t.
Theorem 3.
For $1\le s\le k$, the state probabilities of ${\{{\mathcal{Q}^{{\alpha _{1}},{\alpha _{2}}}}(t)\}_{t\ge 0}}$ are given by
\[\begin{aligned}{}{p_{n,s}^{{\alpha _{1}},{\alpha _{2}}}}(t)& ={\sum \limits_{i=0}^{\infty }}\frac{{B_{i}^{n,s}}}{{c_{1}^{{\delta _{i}^{n,s}}}}}\Big({c_{1}}{f_{{\delta _{i}^{n,s}}}^{\ast ({\delta _{i}^{n,s}})}}(t)+{c_{2}}{g_{{\delta _{i}^{n,s}}}^{\ast ({\delta _{i}^{n,s}})}}(t)\Big)\\ {} & \hspace{3.57777pt}\hspace{3.57777pt}+{\sum \limits_{i=0}^{\infty }}{\sum \limits_{m=1}^{\infty }}{\sum \limits_{r=0}^{\infty }}\frac{{C_{i,m,r}^{n,s}}}{{c_{1}^{{\nu _{i,m,r}^{n,s}}}}}\Big({c_{1}}{f_{{\nu _{i,m,r}^{n,s}}}^{\ast ({\nu _{i,m,r}^{n,s}})}}(t)+{c_{2}}{g_{{\nu _{i,m,r}^{n,s}}}^{\ast ({\nu _{i,m,r}^{n,s}})}}(t)\Big)\\ {} & \hspace{3.57777pt}\hspace{3.57777pt}-{\sum \limits_{i=0}^{\infty }}{\sum \limits_{m=1}^{\infty }}{\sum \limits_{r=0}^{\infty }}\frac{{D_{i,m,r}^{n,s}}}{{c_{1}^{{\pi _{i,m,r}^{n,s}}}}}\Big({c_{1}}{f_{{\pi _{i,m,r}^{n,s}}}^{\ast ({\pi _{i,m,r}^{n,s}})}}(t)+{c_{2}}{g_{{\pi _{i,m,r}^{n,s}}}^{\ast ({\pi _{i,m,r}^{n,s}})}}(t)\Big),\hspace{0.1667em}n\ge 1,\hspace{0.1667em}t\ge 0,\end{aligned}\]
where ${f_{N}}(\cdot )$ and ${g_{N}}(\cdot )$ are given in (24) and (25), respectively.
Proof.
On taking the inverse Laplace transform of (33), we obtain
On using (6) in (37), we get the required result. □
(37)
\[\begin{aligned}{}{p_{n,s}^{{\alpha _{1}},{\alpha _{2}}}}(t)& ={\sum \limits_{i=0}^{\infty }}{B_{i}^{n,s}}\Bigg(\frac{1}{{c_{1}^{{\delta _{i}^{n,s}}-1}}}{\mathbb{L}^{-1}}\Bigg({\bigg(\frac{{z^{\frac{{\alpha _{1}}-1}{{\delta _{i}^{n,s}}}}}}{{z^{{\alpha _{1}}}}+\frac{{c_{2}}}{{c_{1}}}{z^{{\alpha _{2}}}}+\frac{\lambda +k\mu }{{c_{1}}}}\bigg)^{{\delta _{i}^{n,s}}}};t\Bigg)\\ {} & \hspace{3.33333pt}\hspace{3.33333pt}+\frac{{c_{2}}}{{c_{1}^{{\delta _{i}^{n,s}}}}}{\mathbb{L}^{-1}}\Bigg({\bigg(\frac{{z^{\frac{{\alpha _{2}}-1}{{\delta _{i}^{n,s}}}}}}{{z^{{\alpha _{1}}}}+\frac{{c_{2}}}{{c_{1}}}{z^{{\alpha _{2}}}}+\frac{\lambda +k\mu }{{c_{1}}}}\bigg)^{{\delta _{i}^{n,s}}}};t\Bigg)\Bigg)\\ {} & \hspace{3.33333pt}\hspace{3.33333pt}+{\sum \limits_{i=0}^{\infty }}{\sum \limits_{m=1}^{\infty }}{\sum \limits_{r=0}^{\infty }}{C_{i,m,r}^{n,s}}\Bigg(\frac{1}{{c_{1}^{{\nu _{i,m,r}^{n,s}}-1}}}{\mathbb{L}^{-1}}\Bigg({\bigg(\frac{{z^{\frac{{\alpha _{1}}-1}{{\nu _{i,m,r}^{n,s}}}}}}{{z^{{\alpha _{1}}}}+\frac{{c_{2}}}{{c_{1}}}{z^{{\alpha _{2}}}}+\frac{\lambda +k\mu }{{c_{1}}}}\bigg)^{{\nu _{i,m,r}^{n,s}}}};t\Bigg)\\ {} & \hspace{3.33333pt}\hspace{3.33333pt}+\frac{{c_{2}}}{{c_{1}^{{\nu _{i,m,r}^{n,s}}}}}{\mathbb{L}^{-1}}\Bigg({\bigg(\frac{{z^{\frac{{\alpha _{2}}-1}{{\nu _{i,m,r}^{n,s}}}}}}{{z^{{\alpha _{1}}}}+\frac{{c_{2}}}{{c_{1}}}{z^{{\alpha _{2}}}}+\frac{\lambda +k\mu }{{c_{1}}}}\bigg)^{{\nu _{i,m,r}^{n,s}}}};t\Bigg)\Bigg)\\ {} & \hspace{3.33333pt}\hspace{3.33333pt}-{\sum \limits_{i=0}^{\infty }}{\sum \limits_{m=1}^{\infty }}{\sum \limits_{r=0}^{\infty }}{D_{i,m,r}^{n,s}}\Bigg(\frac{1}{{c_{1}^{{\pi _{i,m,r}^{n,s}}-1}}}{\mathbb{L}^{-1}}\Bigg({\bigg(\frac{{z^{\frac{{\alpha _{1}}-1}{{\pi _{i,m,r}^{n,s}}}}}}{{z^{{\alpha _{1}}}}+\frac{{c_{2}}}{{c_{1}}}{z^{{\alpha _{2}}}}+\frac{\lambda +k\mu }{{c_{1}}}}\bigg)^{{\pi _{i,m,r}^{n,s}}}};t\Bigg)\\ {} & \hspace{3.33333pt}\hspace{3.33333pt}+\frac{{c_{2}}}{{c_{1}^{{\pi _{i,m,r}^{n,s}}}}}{\mathbb{L}^{-1}}\Bigg({\bigg(\frac{{z^{\frac{{\alpha _{2}}-1}{{\pi _{i,m,r}^{n,s}}}}}}{{z^{{\alpha _{1}}}}+\frac{{c_{2}}}{{c_{1}}}{z^{{\alpha _{2}}}}+\frac{\lambda +k\mu }{{c_{1}}}}\bigg)^{{\pi _{i,m,r}^{n,s}}}};t\Bigg)\Bigg),\hspace{0.1667em}t\ge 0.\end{aligned}\]3.1 Queue length process
Here, we define the queue length process of mixed time-changed Erlang queue as follows:
\[ {\mathcal{L}^{{\alpha _{1}},{\alpha _{2}}}}(t):={g_{k}}({\mathcal{Q}^{{\alpha _{1}},{\alpha _{2}}}}(t)),\hspace{0.1667em}t\ge 0,\]
where ${g_{k}}(\cdot )$ is a bijective map as defined in (1). Let ${p_{n}^{{\alpha _{1}},{\alpha _{2}}}}(t)=\mathrm{Pr}({\mathcal{L}^{{\alpha _{1}},{\alpha _{2}}}}(t)=n|{\mathcal{L}^{{\alpha _{1}},{\alpha _{2}}}}(0)=0)$, $n\ge 0$ be its state probabilities such that ${p_{-n}^{{\alpha _{1}},{\alpha _{2}}}}(t)=0$ for all $n\ge 1$.The mean queue length of mixed time-changed Erlang queue is defined as
As ${g_{k}}(\cdot )$ is a bijective map, the mixed time-changed Erlang queue $\{{\mathcal{Q}^{{\alpha _{1}},{\alpha _{2}}}}$ ${(t)\}_{t\ge 0}}$ is equivalently represented by the queue length process ${\{{\mathcal{L}^{{\alpha _{1}},{\alpha _{2}}}}(t)\}_{t\ge 0}}$. Thus,
and
where ${a_{k}}(\cdot )$ and ${b_{k}}(\cdot )$ are as given in (3) and (2), respectively.
(39)
\[ {p_{n,s}^{{\alpha _{1}},{\alpha _{2}}}}(t)={p_{{g_{k}}(n,s)}^{{\alpha _{1}},{\alpha _{2}}}}(t),\hspace{0.1667em}t\ge 0\](40)
\[ {p_{n}^{{\alpha _{1}},{\alpha _{2}}}}(t)={p_{{a_{k}}(n),{b_{k}}(n)}^{{\alpha _{1}},{\alpha _{2}}}}(t),\hspace{0.1667em}t\ge 0,\]Now, by using (39) in Theorem 1, we get the following system of differential equations that governs the state probabilities of ${\{{\mathcal{L}^{{\alpha _{1}},{\alpha _{2}}}}(t)\}_{t\ge 0}}$:
with initial condition ${p_{0}^{{\alpha _{1}},{\alpha _{2}}}}(0)=1$. Here, $\frac{{\mathrm{d}^{{\alpha _{i}}}}}{\mathrm{d}{t^{{\alpha _{i}}}}}$ is a Caputo fractional derivative of order ${\alpha _{i}}$ defined in (7).
(41)
\[ \left.\begin{aligned}{}\Big({c_{1}}\frac{{\mathrm{d}^{{\alpha _{1}}}}}{\mathrm{d}{t^{{\alpha _{1}}}}}+{c_{2}}\frac{{\mathrm{d}^{{\alpha _{2}}}}}{\mathrm{d}{t^{{\alpha _{2}}}}}\Big){p_{0}^{{\alpha _{1}},{\alpha _{2}}}}(t)& =-\lambda {p_{0}^{{\alpha _{1}},{\alpha _{2}}}}(t)+k\mu {p_{1}^{{\alpha _{1}},{\alpha _{2}}}}(t),\\ {} \Big({c_{1}}\frac{{\mathrm{d}^{{\alpha _{1}}}}}{\mathrm{d}{t^{{\alpha _{1}}}}}+{c_{2}}\frac{{\mathrm{d}^{{\alpha _{2}}}}}{\mathrm{d}{t^{{\alpha _{2}}}}}\Big){p_{n}^{{\alpha _{1}},{\alpha _{2}}}}(t)& =-(\lambda +k\mu ){p_{n}^{{\alpha _{1}},{\alpha _{2}}}}(t)+k\mu {p_{n+1}^{{\alpha _{1}},{\alpha _{2}}}}(t),\hspace{0.1667em}1\le n\le k-1,\\ {} \Big({c_{1}}\frac{{\mathrm{d}^{{\alpha _{1}}}}}{\mathrm{d}{t^{{\alpha _{1}}}}}+{c_{2}}\frac{{\mathrm{d}^{{\alpha _{2}}}}}{\mathrm{d}{t^{{\alpha _{2}}}}}\Big){p_{n}^{{\alpha _{1}},{\alpha _{2}}}}(t)& =-(\lambda +k\mu ){p_{n}^{{\alpha _{1}},{\alpha _{2}}}}(t)+k\mu {p_{n+1}^{{\alpha _{1}},{\alpha _{2}}}}(t)+\lambda {p_{n-k}^{{\alpha _{1}},{\alpha _{2}}}}(t),\hspace{0.1667em}n\ge k\end{aligned}\right\}\]Next, we obtain the expected number of remaining phases to be served.
Theorem 4.
The mean queue length of mixed time-changed Erlang queue solves the following fractional differential equation:
with initial condition ${\mathcal{M}^{{\alpha _{1}},{\alpha _{2}}}}(0)=0$.
(42)
\[ \Big({c_{1}}\frac{{\mathrm{d}^{{\alpha _{1}}}}}{\mathrm{d}{t^{{\alpha _{1}}}}}+{c_{2}}\frac{{\mathrm{d}^{{\alpha _{2}}}}}{\mathrm{d}{t^{{\alpha _{2}}}}}\Big){\mathcal{M}^{{\alpha _{1}},{\alpha _{2}}}}(t)=k(\lambda -\mu )+k\mu {p_{0}^{{\alpha _{1}},{\alpha _{2}}}}(t),\hspace{0.1667em}t\ge 0\]Proof.
From (38), we have
Note that ${\mathcal{M}^{{\alpha _{1}},{\alpha _{2}}}}(0)=0$ as ${p_{0}^{{\alpha _{1}},{\alpha _{2}}}}(0)=1$. By using (41) and (43), we get
This proves the result. □
(43)
\[ {\mathcal{M}^{{\alpha _{1}},{\alpha _{2}}}}(t)={\sum \limits_{n=0}^{\infty }}n{p_{n}^{{\alpha _{1}},{\alpha _{2}}}}(t),\hspace{0.1667em}t\ge 0.\](44)
\[\begin{aligned}{}\Big({c_{1}}\frac{{\mathrm{d}^{{\alpha _{1}}}}}{\mathrm{d}{t^{{\alpha _{1}}}}}+{c_{2}}\frac{{\mathrm{d}^{{\alpha _{2}}}}}{\mathrm{d}{t^{{\alpha _{2}}}}}\Big){\mathcal{M}^{{\alpha _{1}},{\alpha _{2}}}}(t)& =-(\lambda +k\mu ){\mathcal{M}^{{\alpha _{1}},{\alpha _{2}}}}(t)+k\mu {\sum \limits_{n=1}^{\infty }}n{p_{n+1}^{{\alpha _{1}},{\alpha _{2}}}}(t)\\ {} & \hspace{3.33333pt}\hspace{3.33333pt}+\lambda {\sum \limits_{n=k}^{\infty }}n{p_{n-k}^{{\alpha _{1}},{\alpha _{2}}}}(t)\\ {} & =-(\lambda +k\mu ){\mathcal{M}^{{\alpha _{1}},{\alpha _{2}}}}(t)\\ {} & \hspace{3.33333pt}\hspace{3.33333pt}+k\mu ({\mathcal{M}^{{\alpha _{1}},{\alpha _{2}}}}(t)+{p_{0}^{{\alpha _{1}},{\alpha _{2}}}}(t)-1)\\ {} & \hspace{3.33333pt}\hspace{3.33333pt}+\lambda ({\mathcal{M}^{{\alpha _{1}},{\alpha _{2}}}}(t)+k)\\ {} & =k(\lambda -\mu )+k\mu {p_{0}^{{\alpha _{1}},{\alpha _{2}}}}(t).\end{aligned}\]For any positive integer N, let us define the following function:
(45)
\[ {h_{N}}(t)={t^{{\alpha _{1}}-\big(\frac{N-1}{N}\big)}}{\sum \limits_{j=0}^{\infty }}{\Big(-\frac{{c_{2}}}{{c_{1}}}\Big)^{j}}{t^{({\alpha _{1}}-{\alpha _{2}})j}}{E_{{\alpha _{1}},{\alpha _{1}}+({\alpha _{1}}-{\alpha _{2}})j+\frac{1}{N}}^{j+1}}\Big(-\Big(\frac{\lambda +k\mu }{{c_{1}}}\Big){t^{{\alpha _{1}}}}\Big).\]Theorem 5.
Let ${E_{\alpha ,\beta }^{\gamma }}(\cdot )$ be the three parameter Mittag-Leffler function as defined in (4) and ${h_{{a_{m,r}^{0}}}}(\cdot )$ be as defined in (45). Then, the mean queue length of mixed time-changed Erlang queue is given by
where ${a_{m,r}^{0}}$ and ${A_{m,r}^{0}}$ are given in (21) and (22), respectively.
(46)
\[\begin{aligned}{}{\mathcal{M}^{{\alpha _{1}},{\alpha _{2}}}}(t)& =\frac{k(\lambda -\mu )}{{c_{1}}}{t^{{\alpha _{1}}}}{E_{{\alpha _{1}}-{\alpha _{2}},{\alpha _{1}}+1}^{1}}\Big(-\frac{{c_{2}}}{{c_{1}}}{t^{{\alpha _{1}}-{\alpha _{2}}}}\Big)\\ {} & \hspace{3.57777pt}\hspace{3.57777pt}+k\mu {\sum \limits_{m=1}^{\infty }}{\sum \limits_{r=0}^{\infty }}\frac{{A_{m,r}^{0}}}{{c_{1}^{{a_{m,r}^{0}}}}}{h_{{a_{m,r}^{0}}}^{\ast ({a_{m,r}^{0}})}}(t),\hspace{0.1667em}z\gt 0\end{aligned}\]Proof.
By using (11) and (40) in (43), we get
where the second equality follows from Eq. (21) of [22]. Here, ${p_{0}}(u)$ is the zero state probability of classical Erlang queue. On taking the Laplace transform on both sides of (47) and using (10), we get
where we have used (20) in the last step.
(47)
\[\begin{aligned}{}{\mathcal{M}^{{\alpha _{1}},{\alpha _{2}}}}(t)& ={\int _{0}^{\infty }}\mathcal{M}(y)\mathrm{Pr}({Y_{{\alpha _{1}},{\alpha _{2}}}}(t)\in \mathrm{d}y)\\ {} & =k(\lambda -\mu ){\int _{0}^{\infty }}y\mathrm{Pr}({Y_{{\alpha _{1}},{\alpha _{2}}}}(t)\in \mathrm{d}y)\\ {} & \hspace{3.33333pt}\hspace{3.33333pt}+k\mu {\int _{0}^{\infty }}{\int _{0}^{y}}{p_{0}}(u)\mathrm{d}u\mathrm{Pr}({Y_{{\alpha _{1}},{\alpha _{2}}}}(t)\in \mathrm{d}y),\hspace{0.1667em}t\ge 0,\end{aligned}\](48)
\[\begin{aligned}{}{\tilde{\mathcal{M}}^{{\alpha _{1}},{\alpha _{2}}}}(z)& =({c_{1}}{z^{{\alpha _{1}}-1}}+{c_{2}}{z^{{\alpha _{2}}-1}})\Big(k(\lambda -\mu ){\int _{0}^{\infty }}y{e^{-y({c_{1}}{z^{{\alpha _{1}}}}+{c_{2}}{z^{{\alpha _{2}}}})}}\mathrm{d}y\\ {} & \hspace{96.73918pt}+k\mu {\int _{0}^{\infty }}{\int _{0}^{y}}{p_{0}}(u){e^{-y({c_{1}}{z^{{\alpha _{1}}}}+{c_{2}}{z^{{\alpha _{2}}}})}}\mathrm{d}u\mathrm{d}y\Big)\\ {} & =({c_{1}}{z^{{\alpha _{1}}-1}}+{c_{2}}{z^{{\alpha _{2}}-1}})\Big(\frac{k(\lambda -\mu )}{{({c_{1}}{z^{{\alpha _{1}}}}+{c_{2}}{z^{{\alpha _{2}}}})^{2}}}+k\mu {\int _{0}^{\infty }}{p_{0}}(u)\\ {} & \hspace{99.58464pt}\cdot {\int _{u}^{\infty }}{e^{-y({c_{1}}{z^{{\alpha _{1}}}}+{c_{2}}{z^{{\alpha _{2}}}})}}\mathrm{d}y\mathrm{d}u\Big)\\ {} & =({c_{1}}{z^{{\alpha _{1}}-1}}+{c_{2}}{z^{{\alpha _{2}}-1}})\Big(\frac{k(\lambda -\mu )}{{({c_{1}}{z^{{\alpha _{1}}}}+{c_{2}}{z^{{\alpha _{2}}}})^{2}}}+\frac{k\mu }{({c_{1}}{z^{{\alpha _{1}}}}+{c_{2}}{z^{{\alpha _{2}}}})}\\ {} & \hspace{99.58464pt}\cdot {\int _{0}^{\infty }}{p_{0}}(u){e^{-u({c_{1}}{z^{{\alpha _{1}}}}+{c_{2}}{z^{{\alpha _{2}}}})}}\mathrm{d}u\Big)\\ {} & =\frac{k(\lambda -\mu )}{{c_{1}}{z^{{\alpha _{1}}+1}}+{c_{2}}{z^{{\alpha _{2}}+1}}}+\frac{k\mu }{z}{\int _{0}^{\infty }}{p_{0}}(u){e^{-u({c_{1}}{z^{{\alpha _{1}}}}+{c_{2}}{z^{{\alpha _{2}}}})}}\mathrm{d}u\\ {} & =\frac{k(\lambda -\mu )}{{c_{1}}{z^{{\alpha _{1}}+1}}+{c_{2}}{z^{{\alpha _{2}}+1}}}+\frac{k\mu }{z}\frac{{\tilde{{p_{0}}}^{{\alpha _{1}},{\alpha _{2}}}}(z)}{({c_{1}}{z^{{\alpha _{1}}-1}}+{c_{2}}{z^{{\alpha _{2}}-1}})}\hspace{0.1667em}\hspace{0.1667em}(\text{by using (15)})\\ {} & =\frac{k(\lambda -\mu )}{{c_{1}}{z^{{\alpha _{1}}+1}}+{c_{2}}{z^{{\alpha _{2}}+1}}}+\frac{k\mu }{z}{\sum \limits_{m=1}^{\infty }}{\sum \limits_{r=0}^{\infty }}\frac{{A_{m,r}^{0}}}{{(\lambda +k\mu +{c_{1}}{z^{{\alpha _{1}}}}+{c_{2}}{z^{{\alpha _{2}}}})^{{a_{m,r}^{0}}}}},\hspace{0.1667em}z\gt 0,\end{aligned}\]On taking the inverse Laplace transform on both sides of (48), we get
\[\begin{aligned}{}{\mathcal{M}^{{\alpha _{1}},{\alpha _{2}}}}(t)& =\frac{k(\lambda -\mu )}{{c_{1}}}{\mathbb{L}^{-1}}\Big(\frac{{z^{-{\alpha _{2}}-1}}}{{z^{{\alpha _{1}}-{\alpha _{2}}}}+\frac{{c_{2}}}{{c_{1}}}};t\Big)\\ {} & \hspace{28.45274pt}+k\mu {\sum \limits_{m=1}^{\infty }}{\sum \limits_{r=0}^{\infty }}\frac{{A_{m,r}^{0}}}{{c_{1}^{{a_{m,r}^{0}}}}}{\mathbb{L}^{-1}}\bigg({\bigg(\frac{{z^{\frac{-1}{{a_{m,r}^{0}}}}}}{{z^{{\alpha _{1}}}}+\frac{{c_{2}}}{{c_{1}}}{z^{{\alpha _{2}}}}+\frac{\lambda +k\mu }{{c_{1}}}}\bigg)^{{a_{m,r}^{0}}}};t\bigg).\end{aligned}\]
From (5) and (6), the required result follows. □3.2 Inter-arrival and inter-phase times of mixed time changed Erlang queue
Here, we view the mixed time-changed Erlang queue in terms of its inter-arrival and service times.
Let us denote the following:
\[\begin{aligned}{}{q_{n}^{{\alpha _{1}},{\alpha _{2}}}}(t)& =\mathrm{Pr}({\mathcal{N}^{{\alpha _{1}},{\alpha _{2}}}}(t)=n|{\mathcal{Q}^{{\alpha _{1}},{\alpha _{2}}}}(0)=(0,0)),\hspace{0.1667em}n\ge 1,\hspace{0.1667em}t\ge 0\\ {} {r_{s}^{{\alpha _{1}},{\alpha _{2}}}}(t)& =\mathrm{Pr}({\mathcal{S}^{{\alpha _{1}},{\alpha _{2}}}}(t)=s|{\mathcal{Q}^{{\alpha _{1}},{\alpha _{2}}}}(0)=(0,0)),\hspace{0.1667em}1\le s\le k,\hspace{0.1667em}t\ge 0.\end{aligned}\]
Then,
\[ {q_{n}^{{\alpha _{1}},{\alpha _{2}}}}(t)={\sum \limits_{s=1}^{k}}{p_{n,s}^{{\alpha _{1}},{\alpha _{2}}}}(t)\hspace{2.5pt}\text{and}\hspace{2.5pt}{r_{s}^{{\alpha _{1}},{\alpha _{2}}}}(t)={\sum \limits_{n=1}^{\infty }}{p_{n,s}^{{\alpha _{1}},{\alpha _{2}}}}(t),\hspace{0.1667em}t\ge 0.\]
For $1\le s\le k$, we add the differential equations given in Theorem 1 to obtain
\[\begin{aligned}{}\Big({c_{1}}\frac{{\mathrm{d}^{{\alpha _{1}}}}}{\mathrm{d}{t^{{\alpha _{1}}}}}+{c_{2}}\frac{{\mathrm{d}^{{\alpha _{2}}}}}{\mathrm{d}{t^{{\alpha _{2}}}}}\Big){p_{0}^{{\alpha _{1}},{\alpha _{2}}}}(t)& =-\lambda {p_{0}^{{\alpha _{1}},{\alpha _{2}}}}(t)+k\mu {p_{1,1}^{{\alpha _{1}},{\alpha _{2}}}}(t),\\ {} \Big({c_{1}}\frac{{\mathrm{d}^{{\alpha _{1}}}}}{\mathrm{d}{t^{{\alpha _{1}}}}}+{c_{2}}\frac{{\mathrm{d}^{{\alpha _{2}}}}}{\mathrm{d}{t^{{\alpha _{2}}}}}\Big){q_{1}^{{\alpha _{1}},{\alpha _{2}}}}(t)& =-\lambda {q_{1}^{{\alpha _{1}},{\alpha _{2}}}}(t)+k\mu ({p_{2,1}^{{\alpha _{1}},{\alpha _{2}}}}(t)-{p_{1,1}^{{\alpha _{1}},{\alpha _{2}}}}(t))\\ {} & \hspace{3.33333pt}\hspace{3.33333pt}+\lambda {p_{0}^{{\alpha _{1}},{\alpha _{2}}}}(t),\\ {} \Big({c_{1}}\frac{{\mathrm{d}^{{\alpha _{1}}}}}{\mathrm{d}{t^{{\alpha _{1}}}}}+{c_{2}}\frac{{\mathrm{d}^{{\alpha _{2}}}}}{\mathrm{d}{t^{{\alpha _{2}}}}}\Big){q_{n}^{{\alpha _{1}},{\alpha _{2}}}}(t)& =-\lambda {q_{n}^{{\alpha _{1}},{\alpha _{2}}}}(t)+k\mu ({p_{n+1,1}^{{\alpha _{1}},{\alpha _{2}}}}(t)-{p_{n,1}^{{\alpha _{1}},{\alpha _{2}}}}(t))\\ {} & \hspace{3.33333pt}\hspace{3.33333pt}+\lambda {q_{n-1}^{{\alpha _{1}},{\alpha _{2}}}}(t),\hspace{0.1667em}n\ge 2,\end{aligned}\]
with ${p_{0}^{{\alpha _{1}},{\alpha _{2}}}}(0)=1$ and ${q_{n}^{{\alpha _{1}},{\alpha _{2}}}}(0)=0$, $n\ge 1$.Similarly, for the classical Erlang queue, we have (see [4], p. 3257)
\[\begin{aligned}{}\frac{\mathrm{d}}{\mathrm{d}t}{p_{0}}(t)& =-\lambda {p_{0}}(t)+k\mu {p_{1,1}}(t),\\ {} \frac{\mathrm{d}}{\mathrm{d}t}{q_{1}}(t)& =-\lambda {q_{1}}(t)+k\mu ({p_{2,1}}(t)-{p_{1,1}}(t))+\lambda {p_{0}}(t),\\ {} \frac{\mathrm{d}}{\mathrm{d}t}{q_{n}}(t)& =-\lambda {q_{n}}(t)+k\mu ({p_{n+1,1}}(t)-{p_{n,1}}(t))+\lambda {q_{n-1}}(t),\hspace{0.1667em}n\ge 2,\end{aligned}\]
with ${p_{0}}(0)=1$ and ${q_{n}}(0)=0$, $n\ge 1$. Here,
\[\begin{aligned}{}{q_{n}}(t)& =\mathrm{Pr}(\mathcal{N}(t)=n|\mathcal{Q}(0)=(0,0)),\hspace{0.1667em}n\ge 1,\hspace{0.1667em}t\ge 0\\ {} {r_{s}}(t)& =\mathrm{Pr}(\mathcal{S}(t)=s|\mathcal{Q}(0)=(0,0)),\hspace{0.1667em}1\le s\le k,\hspace{0.1667em}t\ge 0.\end{aligned}\]
Note that ${q_{n}}(t)={\textstyle\sum _{s=1}^{k}}{p_{n,s}}(t)\hspace{2.5pt}\text{and}\hspace{2.5pt}{r_{s}}(t)={\textstyle\sum _{n=1}^{\infty }}{p_{n,s}}(t)$.Recall that the first arrival time is the random time at which the first arrival occurs, and the inter-arrival time is defined as the time between two successive arrivals.
In the following results, we obtain the distributions of inter-arrival, inter-phase and sojourn times of mixed time-changed Erlang queue. These results are useful in understanding the customers arrival, phase service and overall system pattern.
Theorem 6.
The inter-arrival times ${T^{{\alpha _{1}},{\alpha _{2}}}}$ of the mixed time-changed Erlang queue ${\{{\mathcal{Q}^{{\alpha _{1}},{\alpha _{2}}}}(t)\}_{t\ge 0}}$ are independent and identically distributed (iid) with the following distribution:
where ${E_{\alpha ,\beta }^{\gamma }}(\cdot )$ is the three parameter Mittag-Leffler function as defined in (4).
(49)
\[\begin{aligned}{}\mathrm{Pr}({T^{{\alpha _{1}},{\alpha _{2}}}}\gt t)& ={\sum \limits_{r=0}^{\infty }}{\Big(-\frac{{c_{2}}}{{c_{1}}}\Big)^{r}}{t^{({\alpha _{1}}-{\alpha _{2}})r}}{E_{{\alpha _{1}},({\alpha _{1}}-{\alpha _{2}})r+1}^{r+1}}\Big(-\frac{\lambda }{{c_{1}}}{t^{{\alpha _{1}}}}\Big)\\ {} & \hspace{8.5359pt}-{\sum \limits_{r=0}^{\infty }}{\Big(-\frac{{c_{2}}}{{c_{1}}}\Big)^{r+1}}{t^{({\alpha _{1}}-{\alpha _{2}})(r+1)}}{E_{{\alpha _{1}},({\alpha _{1}}-{\alpha _{2}})(r+1)+1}^{r+1}}\Big(-\frac{\lambda }{{c_{1}}}{t^{{\alpha _{1}}}}\Big),\hspace{0.1667em}t\ge 0\end{aligned}\]Proof.
Note that the inter-arrival times of mixed time-changed Erlang queue are independent. To obtain the distribution of inter-arrival times, let us define a process ${\{{\mathcal{N}_{\ast }}(t)\}_{t\ge 0}}$ whose state probabilities ${\textbf{q}_{n}}(t)=\mathrm{Pr}({\mathcal{N}_{\ast }}(t)=n)$, $n\ge 0$ solve the following Cauchy problem (see [4], p. 3258):
with initial conditions ${\textbf{q}_{m}}(0)=1$ and ${\textbf{q}_{n}}(0)=0$ for $n\ge 0$, $n\ne m$. Here, ${\{{\mathcal{N}_{\ast }}(t)\}_{t\ge 0}}$ is a counting process with the same arrival rate λ as that of ${\{\mathcal{N}(t)\}_{t\ge 0}}$. It starts at m and $m+1$ is its absorbing state. Observe that ${\textbf{q}_{m}}(t)+{\textbf{q}_{m+1}}(t)=1$, $t\ge 0$.
(50)
\[ \left.\begin{aligned}{}\frac{\mathrm{d}}{\mathrm{d}t}{\textbf{q}_{m}}(t)& =-\lambda {\textbf{q}_{m}}(t),\\ {} \frac{\mathrm{d}}{\mathrm{d}t}{\textbf{q}_{m+1}}(t)& =\lambda {\textbf{q}_{m}}(t),\\ {} \frac{\mathrm{d}}{\mathrm{d}t}{\textbf{q}_{n}}(t)& =0,\hspace{0.1667em}n\ge 0,\hspace{0.1667em}n\ne m,\hspace{0.1667em}m+1,\end{aligned}\right\}\]Let us denote T as the arrival time of first new customer in the classical Erlang queue ${\{\mathcal{Q}(t)\}_{t\ge 0}}$ starting from ${\mathcal{N}_{\ast }}(0)=m$. Then, its distribution function ${F_{T}}(\cdot )$ is given by
Let ${\mathcal{N}_{\ast }^{{\alpha _{1}},{\alpha _{2}}}}(t)={\mathcal{N}_{\ast }}({Y_{{\alpha _{1}},{\alpha _{2}}}}(t))$, $t\ge 0$ be the time-changed process whose state probabilities are denoted by ${\textbf{q}_{m}^{{\alpha _{1}},{\alpha _{2}}}}(t)=\mathrm{Pr}({\mathcal{N}_{\ast }^{{\alpha _{1}},{\alpha _{2}}}}(t)=m)\hspace{2.5pt}\text{and}\hspace{2.5pt}{\textbf{q}_{m+1}^{{\alpha _{1}},{\alpha _{2}}}}(t)$ $=\mathrm{Pr}({\mathcal{N}_{\ast }^{{\alpha _{1}},{\alpha _{2}}}}(t)=m+1)$. Here, ${\{{\mathcal{N}_{\ast }}(t)\}_{t\ge 0}}$ and ${\{{Y_{{\alpha _{1}},{\alpha _{2}}}}(t)\}_{t\ge 0}}$ are independent. So, ${\{{\mathcal{N}_{\ast }^{{\alpha _{1}},{\alpha _{2}}}}(t)\}_{t\ge 0}}$ is a process which starts at m with $m+1$ being its absorbing state. Let ${T^{{\alpha _{1}},{\alpha _{2}}}}$ be the arrival time of the first new customer in the mixed time-changed Erlang queue ${\{{\mathcal{Q}^{{\alpha _{1}},{\alpha _{2}}}}(t)\}_{t\ge 0}}$ starting from ${\mathcal{N}_{\ast }^{{\alpha _{1}},{\alpha _{2}}}}(0)=m$. Thus, ${F_{{T^{{\alpha _{1}},{\alpha _{2}}}}}}(t)=\mathrm{Pr}({T^{{\alpha _{1}},{\alpha _{2}}}}\le t)={\textbf{q}_{m+1}^{{\alpha _{1}},{\alpha _{2}}}}(t)$.
Let ${\{\mathcal{K}(t)\}_{t\ge 0}}$ be the counting process which counts the number of arrivals in the Erlang queue ${\{\mathcal{Q}(t)\}_{t\ge 0}}$. This counting process is obtained by assuming service times to be degenerate to ∞, that is, by setting $\mu =0$ in the Erlang queue. As ${\{\mathcal{K}(t)\}_{t\ge 0}}$ is a Markov process, its time-changed variant ${\{{\mathcal{K}^{{\alpha _{1}},{\alpha _{2}}}}(t)\}_{t\ge 0}}$, where ${\mathcal{K}^{{\alpha _{1}},{\alpha _{2}}}}(t)=\mathcal{K}({Y_{{\alpha _{1}},{\alpha _{2}}}}(t))$, $t\ge 0$ is a semi-Markov process. Let ${T_{n}}$ be its nth jump time. Then, discrete time process ${\mathcal{K}^{{\alpha _{1}},{\alpha _{2}}}}({T_{n}})$ is a time-homogeneous Markov chain (see [13]). Thus, we can say that the inter-arrival times of mixed time-changed Erlang queue ${\{{\mathcal{Q}^{{\alpha _{1}},{\alpha _{2}}}}(t)\}_{t\ge 0}}$ are distributed according to ${T^{{\alpha _{1}},{\alpha _{2}}}}$.
Now, we show that the state probabilities ${\textbf{q}_{m}^{{\alpha _{1}},{\alpha _{2}}}}(t)$ and ${\textbf{q}_{m+1}^{{\alpha _{1}},{\alpha _{2}}}}(t)$ of $\{{\mathcal{N}_{\ast }^{{\alpha _{1}},{\alpha _{2}}}}$ ${(t)\}_{t\ge 0}}$ solve the following differential equations:
with ${\textbf{q}_{m}^{{\alpha _{1}},{\alpha _{2}}}}(0)=1$ and ${\textbf{q}_{m+1}^{{\alpha _{1}},{\alpha _{2}}}}(0)=0$.
(51)
\[\begin{aligned}{}\Big({c_{1}}\frac{{\mathrm{d}^{{\alpha _{1}}}}}{\mathrm{d}{t^{{\alpha _{1}}}}}+{c_{2}}\frac{{\mathrm{d}^{{\alpha _{2}}}}}{\mathrm{d}{t^{{\alpha _{2}}}}}\Big){\textbf{q}_{m}^{{\alpha _{1}},{\alpha _{2}}}}(t)& =-\lambda {\textbf{q}_{m}^{{\alpha _{1}},{\alpha _{2}}}}(t),\end{aligned}\](52)
\[\begin{aligned}{}\Big({c_{1}}\frac{{\mathrm{d}^{{\alpha _{1}}}}}{\mathrm{d}{t^{{\alpha _{1}}}}}+{c_{2}}\frac{{\mathrm{d}^{{\alpha _{2}}}}}{\mathrm{d}{t^{{\alpha _{2}}}}}\Big){\textbf{q}_{m+1}^{{\alpha _{1}},{\alpha _{2}}}}(t)& =\lambda {\textbf{q}_{m}^{{\alpha _{1}},{\alpha _{2}}}}(t),\end{aligned}\]The Eq. (51) holds true if and only if the Laplace transform ${\tilde{\textbf{q}}_{m}^{{\alpha _{1}},{\alpha _{2}}}}(z)$ of ${\textbf{q}_{m}^{{\alpha _{1}},{\alpha _{2}}}}(t)$ solves:
where we have used (8).
(53)
\[ ({c_{1}}{z^{{\alpha _{1}}-1}}+{c_{2}}{z^{{\alpha _{2}}-1}})(z{\tilde{\textbf{q}}_{m}^{{\alpha _{1}},{\alpha _{2}}}}(z)-1)=-\lambda {\tilde{\textbf{q}}_{m}^{{\alpha _{1}},{\alpha _{2}}}}(z),\hspace{0.1667em}z\gt 0,\]On taking the Laplace transform on both sides of
On substituting (54) in (53), we obtain
By using (50) in (55) and on solving it, we get an identity. This establishes that (51) holds true. Similarly, it can be shown that (52) holds true.
\[ {\textbf{q}_{m}^{{\alpha _{1}},{\alpha _{2}}}}(t)=\mathrm{Pr}({\mathcal{N}_{\ast }}({Y_{{\alpha _{1}},{\alpha _{2}}}}(t))=m)={\int _{0}^{\infty }}{\textbf{q}_{m}}(y)\mathrm{Pr}({Y_{{\alpha _{1}},{\alpha _{2}}}}(t)\in \mathrm{d}y),\hspace{0.1667em}t\ge 0\]
and by using (10), we get
(54)
\[ {\tilde{\textbf{q}}_{m}^{{\alpha _{1}},{\alpha _{2}}}}(z)=({c_{1}}{z^{{\alpha _{1}}-1}}+{c_{2}}{z^{{\alpha _{2}}-1}}){\int _{0}^{\infty }}{\textbf{q}_{m}}(y){e^{-y({c_{1}}{z^{{\alpha _{1}}}}+{c_{2}}{z^{{\alpha _{2}}}})}}\mathrm{d}y,\hspace{0.1667em}z\gt 0.\](55)
\[ ({c_{1}}{z^{{\alpha _{1}}}}+{c_{2}}{z^{{\alpha _{2}}}}){\int _{0}^{\infty }}{\textbf{q}_{m}}(y){e^{-y({c_{1}}{z^{{\alpha _{1}}}}+{c_{2}}{z^{{\alpha _{2}}}})}}\mathrm{d}y-1=-\lambda {\int _{0}^{\infty }}{\textbf{q}_{m}}(y){e^{-y({c_{1}}{z^{{\alpha _{1}}}}+{c_{2}}{z^{{\alpha _{2}}}})}}\mathrm{d}y.\]Remark 3.
From (53), it follows that the Laplace transform of ${\textbf{\textit{q}}_{m}^{{\alpha _{1}},{\alpha _{2}}}}(t)=\mathrm{Pr}({T^{{\alpha _{1}},{\alpha _{2}}}}\gt t)$ is
The probability density function of inter-arrival times of mixed time-changed Erlang queue is given by (see [5], Eq. (2.38))
\[ {f_{{T^{{\alpha _{1}},{\alpha _{2}}}}}}(t)=\frac{\lambda }{{c_{1}}}{t^{{\alpha _{1}}-1}}{\sum \limits_{r=0}^{\infty }}{\Big(-\frac{{c_{2}}}{{c_{1}}}\Big)^{r}}{t^{({\alpha _{1}}-{\alpha _{2}})r}}{E_{{\alpha _{1}},{\alpha _{1}}+({\alpha _{1}}-{\alpha _{2}})r}^{r+1}}\Big(-\frac{\lambda {t^{{\alpha _{1}}}}}{{c_{1}}}\Big),\hspace{0.1667em}t\ge 0,\]
which has the following asymptotic relation (see [5], Eq. (2.47)): ${f_{{T^{{\alpha _{1}},{\alpha _{2}}}}}}(t)\sim \frac{\lambda {t^{{\alpha _{1}}-1}}}{{c_{1}}\Gamma ({\alpha _{1}})}$ as $t\to 0$. Here, $f(t)\sim g(t)$ as $t\to a$ means f and g are asymptotically equal as $t\to a$, that is, ${\lim \nolimits_{t\to a}}\frac{f(t)}{g(t)}=1$. So, the asymptotic behavior of inter-arrival times depends on fractional parameter ${\alpha _{1}}$ as $t\to 0$. Also, ${f_{{T^{{\alpha _{1}},{\alpha _{2}}}}}}(t)\sim \frac{{c_{2}}{\alpha _{2}}{t^{-{\alpha _{2}}-1}}}{\lambda \Gamma (1-{\alpha _{2}})}$ as $t\to \infty $ (see [5], (2.48)). Hence, the asymptotic behavior of inter-arrival times depends on the fractional parameter ${\alpha _{2}}$ as $t\to \infty $.The proof of the following result follows similar lines to that of Theorem 6. Thus, it is omitted.
Theorem 7.
The inter-phase times ${\mathcal{P}^{{\alpha _{1}},{\alpha _{2}}}}$ of mixed time-changed Erlang queue are iid with the following distribution:
where ${E_{\alpha ,\beta }^{\gamma }}(\cdot )$ is the three parameter Mittag-Leffler function as defined in (4).
(56)
\[\begin{aligned}{}& \mathrm{Pr}({\mathcal{P}^{{\alpha _{1}},{\alpha _{2}}}}\gt t)\\ {} & \hspace{1em}={\sum \limits_{r=0}^{\infty }}{\Big(-\frac{{c_{2}}}{{c_{1}}}\Big)^{r}}{t^{({\alpha _{1}}-{\alpha _{2}})r}}{E_{{\alpha _{1}},({\alpha _{1}}-{\alpha _{2}})r+1}^{r+1}}\Big(-\frac{k\mu }{{c_{1}}}{t^{{\alpha _{1}}}}\Big)\\ {} & \hspace{2em}\hspace{1em}-{\sum \limits_{r=0}^{\infty }}{\Big(-\frac{{c_{2}}}{{c_{1}}}\Big)^{r+1}}{t^{({\alpha _{1}}-{\alpha _{2}})(r+1)}}{E_{{\alpha _{1}},({\alpha _{1}}-{\alpha _{2}})(r+1)+1}^{r+1}}\Big(-\frac{k\mu }{{c_{1}}}{t^{{\alpha _{1}}}}\Big),\hspace{0.1667em}t\ge 0\end{aligned}\]The probability density function of inter-phase times of mixed time-changed Erlang queue is given by (see [5], Eq. (2.38))
\[ {f_{{\mathcal{P}^{{\alpha _{1}},{\alpha _{2}}}}}}(t)=\frac{k\mu }{{c_{1}}}{t^{{\alpha _{1}}-1}}{\sum \limits_{r=0}^{\infty }}{\Big(-\frac{{c_{2}}}{{c_{1}}}\Big)^{r}}{t^{({\alpha _{1}}-{\alpha _{2}})r}}{E_{{\alpha _{1}},{\alpha _{1}}+({\alpha _{1}}-{\alpha _{2}})r}^{r+1}}\Big(-\frac{k\mu {t^{{\alpha _{1}}}}}{{c_{1}}}\Big),\hspace{0.1667em}t\ge 0,\]
which has the following asymptotic relation (see [5], Eq. (2.47)): ${f_{{\mathcal{P}^{{\alpha _{1}},{\alpha _{2}}}}}}(t)\sim \frac{k\mu {t^{{\alpha _{1}}-1}}}{{c_{1}}\Gamma ({\alpha _{1}})}$ as $t\to 0$. So, the asymptotic behavior of inter-phase times depends on fractional parameter ${\alpha _{1}}$ as $t\to 0$. Also, ${f_{{\mathcal{P}^{{\alpha _{1}},{\alpha _{2}}}}}}(t)\sim \frac{{c_{2}}{\alpha _{2}}{t^{-{\alpha _{2}}-1}}}{k\mu \Gamma (1-{\alpha _{2}})}$ as $t\to \infty $ (see [5], (2.48)). Hence, the asymptotic behavior of inter-phase times depends on the fractional parameter ${\alpha _{2}}$ as $t\to \infty $.Remark 4.
Let ${Y_{1}},{Y_{2}},\dots ,{Y_{k}}$ be the inter-phase times of mixed time changed Erlang queue and ${f_{{Y_{1}}}}(\cdot )$ be the probability density function of ${Y_{1}}$, that is, ${f_{{Y_{1}}}}(t)=\frac{\mathrm{d}}{\mathrm{d}t}\mathrm{Pr}({Y_{1}}\le t)$. It is given by
Let $X={Y_{1}}+{Y_{2}}+\cdots +{Y_{k}}$. As ${Y_{i}}$’s are iid, we have
So, its Laplace transform is (see [5], Eq. (2.37))
Hence, the service times of mixed time-changed Erlang queue are iid with the probability density function as given in (57).
(57)
\[ {f_{X}}(t)=\frac{\mathrm{d}}{\mathrm{d}t}\mathrm{Pr}(X\le t)={f_{{Y_{1}}}^{\ast (k)}}(t),\hspace{0.1667em}t\ge 0.\](58)
\[ {\tilde{f}_{X}}(z)={({\tilde{f}_{{Y_{1}}}}(z))^{k}}={\Big(\frac{k\mu }{k\mu +{c_{1}}{z^{{\alpha _{1}}}}+{c_{2}}{z^{{\alpha _{2}}}}}\Big)^{k}},\hspace{0.1667em}z\gt 0.\]Recall that the sojourn time of a process is the time spent in a particular state before its transition to another state.
The proof of the following result follows along the similar lines to that of Theorem 6. Thus, it is omitted. Here, we need to work with the state probabilities of queue length process ${\{{\mathcal{L}^{{\alpha _{1}},{\alpha _{2}}}}(t)\}_{t\ge 0}}$.
Theorem 8.
The sojourn times ${S^{{\alpha _{1}},{\alpha _{2}}}}$ of queue length process ${\{{\mathcal{L}^{{\alpha _{1}},{\alpha _{2}}}}(t)\}_{t\ge 0}}$ in a state $n\gt 0$ are iid with the following distribution function:
where ${E_{\alpha ,\beta }^{\gamma }}(\cdot )$ is the three parameter Mittag-Leffler function defined in (4).
(59)
\[\begin{aligned}{}& \mathrm{Pr}({S^{{\alpha _{1}},{\alpha _{2}}}}\gt t)\\ {} & \hspace{1em}={\sum \limits_{r=0}^{\infty }}{\Big(-\frac{{c_{2}}}{{c_{1}}}\Big)^{r}}{t^{({\alpha _{1}}-{\alpha _{2}})r}}{E_{{\alpha _{1}},({\alpha _{1}}-{\alpha _{2}})r+1}^{r+1}}\Big(-\frac{(\lambda +k\mu )}{{c_{1}}}{t^{{\alpha _{1}}}}\Big)\\ {} & \hspace{2em}\hspace{1em}-{\sum \limits_{r=0}^{\infty }}{\Big(-\frac{{c_{2}}}{{c_{1}}}\Big)^{r+1}}{t^{({\alpha _{1}}-{\alpha _{2}})(r+1)}}{E_{{\alpha _{1}},({\alpha _{1}}-{\alpha _{2}})(r+1)+1}^{r+1}}\Big(-\frac{(\lambda +k\mu )}{{c_{1}}}{t^{{\alpha _{1}}}}\Big),\hspace{0.1667em}t\ge 0\end{aligned}\]The probability density function of sojourn times times for the queue length process of mixed time-changed Erlang queue is given by (see [5], Eq. (2.38))
\[\begin{array}{r}\displaystyle {f_{{S^{{\alpha _{1}},{\alpha _{2}}}}}}(t)=\frac{(\lambda +k\mu )}{{c_{1}}}{t^{{\alpha _{1}}-1}}{\sum \limits_{r=0}^{\infty }}{\Big(-\frac{{c_{2}}}{{c_{1}}}\Big)^{r}}{t^{({\alpha _{1}}-{\alpha _{2}})r}}{E_{{\alpha _{1}},{\alpha _{1}}+({\alpha _{1}}-{\alpha _{2}})r}^{r+1}}\Big(\hspace{-0.1667em}-\frac{(\lambda +k\mu ){t^{{\alpha _{1}}}}}{{c_{1}}}\Big),\\ {} \displaystyle t\ge 0,\end{array}\]
which has the following asymptotic relation (see [5], Eq. (2.47)): ${f_{{S^{{\alpha _{1}},{\alpha _{2}}}}}}(t)\sim \frac{(\lambda +k\mu ){t^{{\alpha _{1}}-1}}}{{c_{1}}\Gamma ({\alpha _{1}})}$ as $t\to 0$. So, the asymptotic behavior of sojourn times depends on fractional parameter ${\alpha _{1}}$ as $t\to 0$. Also, ${f_{{S^{{\alpha _{1}},{\alpha _{2}}}}}}(t)\sim \frac{{c_{2}}{\alpha _{2}}{t^{-{\alpha _{2}}-1}}}{(\lambda +k\mu )\Gamma (1-{\alpha _{2}})}$ as $t\to \infty $ (see [5], (2.48)). Hence, the asymptotic behavior of sojourn times depends on the fractional parameter ${\alpha _{2}}$ as $t\to \infty $. Remark 5.
Recall that the mixed time-changed Erlang queue is a generalization of fractional Erlang queue, for ${c_{1}}=1$ or ${c_{2}}=1$ it reduces to the fractional Erlang queue. As the inter-arrival and inter-phase times of fractional Erlang queue are not independent (see [4], Corollary 4.6), the inter-arrival and inter-phase times of mixed time-changed Erlang queue are not necessarily independent.
3.3 Distribution of busy period
The busy period of a queue is defined as a duration of time since the first customer enters the empty system till the system becomes empty again.
In the next result, we obtain the distribution of a busy period. It is useful in understanding for how much time system stays busy when customers arrive in an empty system.
Theorem 9.
Let ${\mathcal{B}^{{\alpha _{1}},{\alpha _{2}}}}$ be the busy period of mixed time-changed Erlang queue, and let ${F_{{\mathcal{B}^{{\alpha _{1}},{\alpha _{2}}}}}}(t)=\mathrm{Pr}({\mathcal{B}^{{\alpha _{1}},{\alpha _{2}}}}\le t)$ be its distribution function. Then
\[ {F_{{\mathcal{B}^{{\alpha _{1}},{\alpha _{2}}}}}}(t)={\sum \limits_{h=0}^{\infty }}\frac{k\mu {A_{k,h}^{0}}}{{c_{1}^{{a_{k,h}^{0}}}}}{h_{{a_{k,h}^{0}}}^{\ast ({a_{k,h}^{0}})}}(t),\hspace{0.1667em}t\ge 0\]
where ${a_{k,h}^{0}}$, ${A_{k,h}^{0}}$ and ${h_{{a_{k,h}^{0}}}}(\cdot )$ are as given in (21), (22) and (45), respectively.
Proof.
Let us consider a process ${\{{\mathcal{L}_{\ast }}(t)\}_{t\ge 0}}$ whose state probabilities ${\textbf{p}_{n}}(t)=\mathrm{Pr}({\mathcal{L}_{\ast }}(t)=n|{\mathcal{L}_{\ast }}(0)=k)$, $n\ge 0$ solve the following system of differential equations (see [22], Eq. (36)-Eq. (38)):
\[\begin{aligned}{}\frac{\mathrm{d}}{\mathrm{d}t}{\textbf{p}_{0}}(t)& =k\mu {\textbf{p}_{1}}(t),\\ {} \frac{\mathrm{d}}{\mathrm{d}t}{\textbf{p}_{n}}(t)& =-(\lambda +k\mu ){\textbf{p}_{n}}(t)+k\mu {\textbf{p}_{n+1}}(t),\hspace{0.1667em}1\le n\le k,\\ {} \frac{\mathrm{d}}{\mathrm{d}t}{\textbf{p}_{n}}(t)& =-(\lambda +k\mu ){\textbf{p}_{n}}(t)+k\mu {\textbf{p}_{n+1}}(t)+\lambda {\textbf{p}_{n-k}}(t),\hspace{0.1667em}n\ge k+1,\end{aligned}\]
with the initial condition ${\textbf{p}_{k}}(0)=1$.Note that the process ${\{{\mathcal{L}_{\ast }}(t)\}_{t\ge 0}}$ behaves similar to ${\{\mathcal{L}(t)\}_{t\ge 0}}$ except that it starts from k instead of 0, which indicates that the first customer entered the system, and 0 is its absorbing state. Thus, by construction, the distribution of busy period of classical Erlang queue is ${F_{\mathcal{B}}}(t)={\textbf{p}_{0}}(t)$, where ${F_{\mathcal{B}}}(t)=\mathrm{Pr}(\mathcal{B}\le t)$ is given in [22].
Let ${\{{\mathcal{L}_{\ast }}(t)\}_{t\ge 0}}$ and ${\{{Y_{{\alpha _{1}},{\alpha _{2}}}}(t)\}_{t\ge 0}}$ be independent. Also, let us define ${\mathcal{L}_{\ast }^{{\alpha _{1}},{\alpha _{2}}}}(t):={\mathcal{L}_{\ast }}({Y_{{\alpha _{1}},{\alpha _{2}}}}(t))$, $t\ge 0$ and denote the nth jump time of ${\{{\mathcal{L}^{{\alpha _{1}},{\alpha _{2}}}}(t)\}_{t\ge 0}}$ by ${T_{n}}$. Then, ${\mathcal{L}^{{\alpha _{1}},{\alpha _{2}}}}({T_{n}})$ is time-homogeneous Markov chain. So, ${F_{{\mathcal{B}^{{\alpha _{1}},{\alpha _{2}}}}}}(t)={\textbf{p}_{0}^{{\alpha _{1}},{\alpha _{2}}}}(t)$. Thus,
On taking the Laplace transform of (60), we get
On taking the inverse Laplace transform of (61), we obtain
Finally, on using (6) in (62), we get the required result. □
(60)
\[ {F_{{\mathcal{B}^{{\alpha _{1}},{\alpha _{2}}}}}}(t)={\int _{0}^{\infty }}{F_{\mathcal{B}}}(y)\mathrm{Pr}({Y_{{\alpha _{1}},{\alpha _{2}}}}(t)\in \mathrm{d}y).\](61)
\[\begin{aligned}{}{\tilde{F}_{{\mathcal{B}^{{\alpha _{1}},{\alpha _{2}}}}}}(z)& =({c_{1}}{z^{{\alpha _{1}}-1}}+{c_{2}}{z^{{\alpha _{2}}-1}}){\int _{0}^{\infty }}{F_{\mathcal{B}}}(y){e^{-y({c_{1}}{z^{{\alpha _{1}}}}+{c_{2}}{z^{{\alpha _{2}}}})}}\mathrm{d}y\\ {} & =({c_{1}}{z^{{\alpha _{1}}-1}}+{c_{2}}{z^{{\alpha _{2}}-1}}){\sum \limits_{h=0}^{\infty }}\frac{k{\lambda ^{h}}{(k\mu )^{k(h+1)}}}{h!(hk+k)!}{\int _{0}^{\infty }}\hspace{-0.1667em}\hspace{-0.1667em}{\int _{0}^{y}}{u^{k+h(k+1)-1}}{e^{-(\lambda +k\mu )u}}\\ {} & \hspace{119.50148pt}\cdot {e^{-y({c_{1}}{z^{{\alpha _{1}}}}+{c_{2}}{z^{{\alpha _{2}}}})}}\mathrm{d}u\mathrm{d}y\hspace{0.1667em}\hspace{0.1667em}\hspace{0.1667em}(\text{see [4], Eq. (13)})\\ {} & =({c_{1}}{z^{{\alpha _{1}}-1}}+{c_{2}}{z^{{\alpha _{2}}-1}}){\sum \limits_{h=0}^{\infty }}\frac{k{\lambda ^{h}}{(k\mu )^{k(h+1)}}}{h!(hk+k)!}{\int _{0}^{\infty }}{u^{k+h(k+1)-1}}{e^{-(\lambda +k\mu )u}}\\ {} & \hspace{119.50148pt}\cdot {\int _{u}^{\infty }}{e^{-y({c_{1}}{z^{{\alpha _{1}}}}+{c_{2}}{z^{{\alpha _{2}}}})}}\mathrm{d}y\mathrm{d}u\\ {} & =\frac{1}{z}{\sum \limits_{h=0}^{\infty }}\frac{k{\lambda ^{h}}{(k\mu )^{k(h+1)}}}{h!(hk+k)!}{\int _{0}^{\infty }}{u^{k+h(k+1)-1}}{e^{-(\lambda +k\mu +{c_{1}}{z^{{\alpha _{1}}}}+{c_{2}}{z^{{\alpha _{2}}}})u}}\mathrm{d}u\\ {} & ={\sum \limits_{h=0}^{\infty }}\frac{k{\lambda ^{h}}{(k\mu )^{k(h+1)}}(k\hspace{-0.1667em}+\hspace{-0.1667em}h(k\hspace{-0.1667em}+\hspace{-0.1667em}1)\hspace{-0.1667em}-\hspace{-0.1667em}1)!}{h!(hk\hspace{-0.1667em}+\hspace{-0.1667em}k)!}\frac{{z^{-1}}}{{(\lambda \hspace{-0.1667em}+\hspace{-0.1667em}k\mu \hspace{-0.1667em}+\hspace{-0.1667em}{c_{1}}{z^{{\alpha _{1}}}}\hspace{-0.1667em}+\hspace{-0.1667em}{c_{2}}{z^{{\alpha _{2}}}})^{k+h(k+1)}}}.\end{aligned}\](62)
\[ {F_{{\mathcal{B}^{{\alpha _{1}},{\alpha _{2}}}}}}(t)={\sum \limits_{h=0}^{\infty }}\frac{k\mu {A_{k,h}^{0}}}{{c_{1}^{{a_{k,h}^{0}}}}}{\mathbb{L}^{-1}}\bigg({\bigg(\frac{{z^{\frac{-1}{{a_{k,h}^{0}}}}}}{{z^{{\alpha _{1}}}}+\frac{{c_{2}}}{{c_{1}}}{z^{{\alpha _{2}}}}+\frac{\lambda +k\mu }{{c_{1}}}}\bigg)^{{a_{k,h}^{0}}}};t\bigg).\]3.4 Conditional waiting time distribution
Methods used for deriving the distribution of waiting time in classical single-server queue depends on its Markovian structure (see [12]). However, these methods are not applicable in the case of mixed time-changed Erlang queue as it is a semi-Markov process. But a partial information about the distribution of its waiting time can be discussed by using some sort of conditioning arguments, as done in [4].
The following distribution function will be used: ${\psi ^{{\alpha _{1}},{\alpha _{2}}}}(t)=\mathrm{Pr}({\mathcal{P}^{{\alpha _{1}},{\alpha _{2}}}}\gt t)$. It is given in (56). Let Y be a random variable with the following distribution function:
with parameters ${t_{0}}\gt 0$ and $k\mu $.
(63)
\[ {F_{Y}}(t)=\mathrm{Pr}({\mathcal{P}^{{\alpha _{1}},{\alpha _{2}}}}\le {t_{0}}+t|{\mathcal{P}^{{\alpha _{1}},{\alpha _{2}}}}\ge {t_{0}})=1-\frac{{\psi ^{{\alpha _{1}},{\alpha _{2}}}}({t_{0}}+t)}{{\psi ^{{\alpha _{1}},{\alpha _{2}}}}({t_{0}})},\hspace{0.1667em}t\ge 0\]On taking the Laplace transform of ${f_{\psi }}(t)={({t_{0}}+t)^{r}}$, $r\gt 0$, we obtain
where $\Gamma (y,z)$ is the analytic extension of the upper incomplete gamma function on ${\mathbb{C}^{2}}$ (see [15]). Also, let ${f_{Y}}(t)$ be the probability density function of Y. Then,
On taking the Laplace transform of (65) and by using (64), we get
where we have used (4) and (56).
(64)
\[\begin{aligned}{}{\tilde{f}_{\psi }}(z)& ={\int _{0}^{\infty }}{e^{-zt}}{({t_{0}}+t)^{r}}\mathrm{d}t\\ {} & ={e^{z{t_{0}}}}{\int _{{t_{0}}}^{\infty }}{e^{-zy}}{y^{r}}\mathrm{d}y\\ {} & =\frac{{e^{z{t_{0}}}}}{{z^{r+1}}}\Gamma (r+1,z{t_{0}}),\hspace{0.1667em}z\gt 0,\end{aligned}\](65)
\[ {f_{Y}}(t)=-\frac{1}{{\psi ^{{\alpha _{1}},{\alpha _{2}}}}({t_{0}})}\frac{\mathrm{d}}{\mathrm{d}t}{\psi ^{{\alpha _{1}},{\alpha _{2}}}}({t_{0}}+t).\](66)
\[\begin{aligned}{}{\tilde{f}_{Y}}(z)& =1-\Bigg({\sum \limits_{r=0}^{\infty }}{\sum \limits_{m=0}^{\infty }}{\Big(-\frac{{c_{2}}}{{c_{1}}}\Big)^{r}}{\Big(-\frac{k\mu }{{c_{1}}}\Big)^{m}}\frac{(m+r)!\Gamma (({\alpha _{1}}-{\alpha _{2}})r+{\alpha _{1}}m+1,z{t_{0}})}{r!m!\Gamma ({\alpha _{1}}m+({\alpha _{1}}-{\alpha _{2}})r+1)}\\ {} & \hspace{11.38092pt}\cdot {z^{({\alpha _{2}}-{\alpha _{1}})r-{\alpha _{1}}m}}-{\sum \limits_{r=0}^{\infty }}{\sum \limits_{m=0}^{\infty }}{\Big(-\frac{{c_{2}}}{{c_{1}}}\Big)^{r+1}}{\Big(-\frac{k\mu }{{c_{1}}}\Big)^{m}}\frac{(m+r)!}{r!m!}\\ {} & \hspace{11.38092pt}\cdot \frac{\Gamma (({\alpha _{1}}-{\alpha _{2}})(r+1)+{\alpha _{1}}m+1,z{t_{0}})}{\Gamma ({\alpha _{1}}m+({\alpha _{1}}-{\alpha _{2}})(r+1)+1)}{z^{({\alpha _{2}}-{\alpha _{1}})(r+1)-{\alpha _{1}}m}}\Bigg)\frac{{e^{z{t_{0}}}}}{{\psi ^{{\alpha _{1}},{\alpha _{2}}}}({t_{0}})},\hspace{0.1667em}z\gt 0,\end{aligned}\]Let ${W_{t}^{{\alpha _{1}},{\alpha _{2}}}}$ be the total time spend by the customer who enters the system at time $t\gt 0$. Consider the following random set:
\[ {\mathcal{E}_{t}^{{\alpha _{1}},{\alpha _{2}}}}(\omega )=\{0\le s\le t:{\mathcal{L}^{{\alpha _{1}},{\alpha _{2}}}}({s^{-}})(\omega )={\mathcal{L}^{{\alpha _{1}},{\alpha _{2}}}}(s)(\omega )+1\},\]
which denotes the set of time instants before t at which a phase completion occurs. Now, we define the following random variable:
\[ {\mathcal{N}_{t}^{{\alpha _{1}},{\alpha _{2}}}}=|{\mathcal{E}_{t}^{{\alpha _{1}},{\alpha _{2}}}}|,\hspace{0.1667em}t\ge 0\]
which gives the count of phases completed before time t. Also, consider a random variable
\[ {\tau _{t}^{{\alpha _{1}},{\alpha _{2}}}}=\left\{\begin{array}{l@{\hskip10.0pt}l}\sup {\mathcal{E}_{t}^{{\alpha _{1}},{\alpha _{2}}}},\hspace{1em}& {\mathcal{N}_{t}^{{\alpha _{1}},{\alpha _{2}}}}\gt 0,\\ {} t,\hspace{1em}& {\mathcal{N}_{t}^{{\alpha _{1}},{\alpha _{2}}}}=0,\end{array}\right.\]
namely, the last time instant before t at which a phase is completed. For a fixed t, the random time ${\tau _{t}^{{\alpha _{1}},{\alpha _{2}}}}$ is a Markov time which takes values among the jump times of ${\{{\mathcal{S}^{{\alpha _{1}},{\alpha _{2}}}}(t)\}_{t\ge 0}}$, except for the case when ${\tau _{t}^{{\alpha _{1}},{\alpha _{2}}}}=t$. Here, ${\{{\mathcal{S}^{{\alpha _{1}},{\alpha _{2}}}}(t)\}_{t\ge 0}}$ is a semi-Markov process and its semi-regenerative set coincides with its discontinuity points (see [8]). Moreover, if ${\tau _{t}^{{\alpha _{1}},{\alpha _{2}}}}={t_{0}}$ for some ${t_{0}}\lt t$, then the past and future of ${\{{\mathcal{S}^{{\alpha _{1}},{\alpha _{2}}}}(t)\}_{t\ge 0}}$ are independent. Therefore, for the waiting time, we condition on the events $\{{\mathcal{L}^{{\alpha _{1}},{\alpha _{2}}}}(t)=n\}$ and $\{{\tau _{t}^{{\alpha _{1}},{\alpha _{2}}}}={t_{0}}\}$.Let
\[ {w^{{\alpha _{1}},{\alpha _{2}}}}(x;t,{t_{0}},n)\hspace{0.1667em}\mathrm{d}x=\mathrm{Pr}({W_{t}^{{\alpha _{1}},{\alpha _{2}}}}\in \mathrm{d}x|{\tau _{t}^{{\alpha _{1}},{\alpha _{2}}}}={t_{0}},{\mathcal{L}^{{\alpha _{1}},{\alpha _{2}}}}(t)=n).\]
A customer has to wait for the completion of n independent phases (as we have conditioned on $\{{\mathcal{L}^{{\alpha _{1}},{\alpha _{2}}}}(t)=n\}$), one of the n phases started at ${t_{0}}$ that has not been completed till time t. Let ${T_{j}}$, $1\le j\le n$ be independent random variables, such that for $1\le j\le n-1$, ${T_{j}}$ denotes the duration of jth phase whose service is not yet started and ${T_{n}}$ be the remaining time duration of the phase which is in service at time t. Thus, we have
\[ \mathbb{E}({W_{t}^{{\alpha _{1}},{\alpha _{2}}}}|{\tau _{t}^{{\alpha _{1}},{\alpha _{2}}}}={t_{0}},{\mathcal{L}^{{\alpha _{1}},{\alpha _{2}}}}(t)=n)={\sum \limits_{j=1}^{n}}{T_{j}}.\]
Note that ${T_{j}}$, $1\le j\le n-1$ are iid from Theorem 7. Thus, $W={\textstyle\sum _{j=1}^{n-1}}{T_{j}}$ has the distribution given by (57). Also, ${T_{n}}$ is distributed as (63) with parameters $t-{t_{0}}$ and $k\mu $. This is because we have conditioned on $\{{\tau _{t}^{{\alpha _{1}},{\alpha _{2}}}}={t_{0}}\}$ which indicates that the last phase was completed at ${t_{0}}$. Then, we obtain
On taking the Laplace transform of (67), and using (58) and (66), we get
(68)
\[\begin{aligned}{}& \tilde{w}(z;t,{t_{0}},n)\\ {} & \hspace{1em}={\Big(\frac{k\mu }{k\mu +{c_{1}}{z^{{\alpha _{1}}}}+{c_{2}}{z^{{\alpha _{2}}}}}\Big)^{n-1}}\Bigg(1-\frac{{e^{z(t-{t_{0}})}}}{{\psi ^{{\alpha _{1}},{\alpha _{2}}}}(t-{t_{0}})}\Bigg({\sum \limits_{r=0}^{\infty }}{\sum \limits_{m=0}^{\infty }}{\Big(-\frac{{c_{2}}}{{c_{1}}}\Big)^{r}}{\Big(-\frac{k\mu }{{c_{1}}}\Big)^{m}}\\ {} & \hspace{2em}\hspace{1em}\cdot \frac{(m+r)!\Gamma (({\alpha _{1}}-{\alpha _{2}})r+{\alpha _{1}}m+1,z(t-{t_{0}}))}{r!m!\Gamma ({\alpha _{1}}m+({\alpha _{1}}-{\alpha _{2}})r+1)}{z^{({\alpha _{2}}-{\alpha _{1}})r-{\alpha _{1}}m}}\\ {} & \hspace{2em}\hspace{1em}-{\sum \limits_{r=0}^{\infty }}{\sum \limits_{m=0}^{\infty }}{\Big(-\frac{{c_{2}}}{{c_{1}}}\Big)^{r+1}}{\Big(-\frac{k\mu }{{c_{1}}}\Big)^{m}}\frac{(m+r)!}{r!m!}\\ {} & \hspace{2em}\hspace{1em}\cdot \frac{\Gamma (({\alpha _{1}}-{\alpha _{2}})(r+1)+{\alpha _{1}}m+1,z(t-{t_{0}}))}{\Gamma ({\alpha _{1}}m+({\alpha _{1}}-{\alpha _{2}})(r+1)+1)}{z^{({\alpha _{2}}-{\alpha _{1}})(r+1)-{\alpha _{1}}m}}\Bigg)\Bigg),\hspace{0.1667em}z\gt 0.\end{aligned}\]3.5 Sample paths simulation
Here, we give plots of sample path simulations of the mixed time-changed Erlang queue. The following results give the distributional properties of exponential random variable and mixed stable subordinator:
Lemma 1.
Let X be an exponential random variable with rate λ and ${\{{D_{{\alpha _{1}},{\alpha _{2}}}}(t)\}_{t\ge 0}}$ be mixed stable subordinator which is independent of X. Then, ${D_{{\alpha _{1}},{\alpha _{2}}}}(X)\stackrel{d}{=}{T^{{\alpha _{1}},{\alpha _{2}}}}$, where ${T^{{\alpha _{1}},{\alpha _{2}}}}$ is the inter-arrival time of the mixed time-changed Erlang queue.
Proof.
Consider the following:
On taking the Laplace transform of (69) and by using (10), we obtain
whose inversion yields the required result. □
(69)
\[\begin{aligned}{}\mathrm{Pr}({D_{{\alpha _{1}},{\alpha _{2}}}}(X)\gt t)& =\mathrm{Pr}(X\gt {Y_{{\alpha _{1}},{\alpha _{2}}}}(t))\\ {} & ={\int _{0}^{\infty }}\mathrm{Pr}(X\gt y)\mathrm{Pr}({Y_{{\alpha _{1}},{\alpha _{2}}}}(t)\in \mathrm{d}y)\\ {} & ={\int _{0}^{\infty }}{e^{-\lambda y}}\mathrm{Pr}({Y_{{\alpha _{1}},{\alpha _{2}}}}(t)\in \mathrm{d}y).\end{aligned}\]Lemma 2.
Let X be a non-negative random variable independent of the mixed stable subordinator ${\{{D_{{\alpha _{1}},{\alpha _{2}}}}(t)\}_{t\ge 0}}$. Then,
\[ {D_{{\alpha _{1}},{\alpha _{2}}}}(X)\stackrel{d}{=}{c_{1}^{1/{\alpha _{1}}}}{X^{1/{\alpha _{1}}}}{D_{{\alpha _{1}}}}(1)+{c_{2}^{1/{\alpha _{2}}}}{X^{1/{\alpha _{2}}}}{D_{{\alpha _{2}}}}(1),\]
where ${D_{{\alpha _{1}}}}(1)$ and ${D_{{\alpha _{2}}}}(1)$ are independent random variables.
Proof.
Note that,
\[\begin{aligned}{}\mathrm{Pr}({D_{{\alpha _{1}},{\alpha _{2}}}}(X)\gt t)& ={\int _{0}^{\infty }}\mathrm{Pr}({D_{{\alpha _{1}},{\alpha _{2}}}}(y)\gt t)\mathrm{Pr}(X\in \mathrm{d}y)\\ {} & ={\int _{0}^{\infty }}\mathrm{Pr}({c_{1}^{1/{\alpha _{1}}}}{D_{{\alpha _{1}}}}(y)+{c_{2}^{1/{\alpha _{2}}}}{D_{{\alpha _{2}}}}(y)\gt t)\mathrm{Pr}(X\in \mathrm{d}y)\\ {} & ={\int _{0}^{\infty }}\hspace{-0.1667em}\hspace{-0.1667em}\mathrm{Pr}({c_{1}^{1/{\alpha _{1}}}}{y^{1/{\alpha _{1}}}}{D_{{\alpha _{1}}}}(1)+{c_{2}^{1/{\alpha _{2}}}}{y^{1/{\alpha _{2}}}}{D_{{\alpha _{2}}}}(1)\gt t)\mathrm{Pr}(X\in \mathrm{d}y)\\ {} & =\mathrm{Pr}({c_{1}^{1/{\alpha _{1}}}}{X^{1/{\alpha _{1}}}}{D_{{\alpha _{1}}}}(1)+{c_{2}^{1/{\alpha _{2}}}}{X^{1/{\alpha _{2}}}}{D_{{\alpha _{2}}}}(1)\gt t),\hspace{0.1667em}t\ge 0,\end{aligned}\]
where we have used ${D_{{\alpha _{1}},{\alpha _{2}}}}(t)\stackrel{d}{=}{c_{1}^{1/{\alpha _{1}}}}{D_{{\alpha _{1}}}}(t)+{c_{2}^{1/{\alpha _{2}}}}{D_{{\alpha _{2}}}}(t)$, $t\ge 0$ (see [1], p. 704) such that ${\{{D_{{\alpha _{1}}}}(t)\}_{t\ge 0}}$ and ${\{{D_{{\alpha _{2}}}}(t)\}_{t\ge 0}}$ are independent stable subordinators. This proves the required result. □Fig. 1.
Plot (a) represents the sample path simulation of ${\mathcal{N}^{{\alpha _{1}},{\alpha _{2}}}}(t)$ and Plot (b) represents the sample path simulation of ${\mathcal{N}^{{\alpha _{1}},{\alpha _{2}}}}({T_{n}})$ such that ${T_{n}}$ are the jump times of ${\{{\mathcal{N}^{{\alpha _{1}},{\alpha _{2}}}}(t)\}_{t\ge 0}}$ for parameters ${c_{1}}=0.4$, ${c_{2}}=0.6$, ${\alpha _{1}}=0.5$, ${\alpha _{2}}=0.3$, $k=4$, $\lambda =6$ and $\mu =5$
Proposition 3.
Let X be an exponential random variable with rate $\lambda \gt 0$, $\{{D_{{\alpha _{1}},{\alpha _{2}}}}$ ${(t)\}_{t\ge 0}}$ be mixed stable subordinator independent of X, and Y be a random variable such that $Y\stackrel{d}{=}{T^{{\alpha _{1}},{\alpha _{2}}}}$. Then,
\[ Y\stackrel{d}{=}{c_{1}^{1/{\alpha _{1}}}}{X^{1/{\alpha _{1}}}}{D_{{\alpha _{1}}}}(1)+{c_{2}^{1/{\alpha _{2}}}}{X^{1/{\alpha _{2}}}}{D_{{\alpha _{2}}}}(1),\]
where ${D_{{\alpha _{1}}}}(1)$ and ${D_{{\alpha _{2}}}}(1)$ are independent.
Fig. 2.
Plot (a) represents the sample path simulation of $\mathcal{N}(t)$ with parameters $k=4$, $\lambda =6$ and $\mu =5$, and Plot (b) represents the sample path simulation of ${\mathcal{N}^{{\alpha _{1}},{\alpha _{2}}}}(t)$ with parameters ${c_{1}}=0.4$, ${c_{2}}=0.6$, ${\alpha _{1}}=0.5$, ${\alpha _{2}}=0.3$, $k=4$, $\lambda =6$ and $\mu =5$
Fig. 3.
Plot (a) represents the sample path simulation of ${\mathcal{N}^{{\alpha _{1}}}}(t)$ with parameters ${\alpha _{1}}=0.5$, $k=4$, $\lambda =6$ and $\mu =5$, and Plot (b) represents the sample path simulation of ${\mathcal{N}^{{\alpha _{1}},{\alpha _{2}}}}(t)$ with parameters ${c_{1}}=0.4$, ${c_{2}}=0.6$, ${\alpha _{1}}=0.5$, ${\alpha _{2}}=0.3$, $k=4$, $\lambda =6$ and $\mu =5$
We use the modified version of Gillespie algorithm (see [7]) for simulating the mixed time-changed Erlang queue. It is based on the fact that ${\mathcal{Q}^{{\alpha _{1}},{\alpha _{2}}}}({T_{n}})$ is again a Markov chain with the same transition probabilities as that of Erlang queue ${\{\mathcal{Q}(t)\}_{t\ge 0}}$, where ${T_{n}}$ is the nth jump time of ${\{{\mathcal{Q}^{{\alpha _{1}},{\alpha _{2}}}}(t)\}_{t\ge 0}}$. This algorithm is used in [4] for simulating the fractional Erlang queue. We use a similar algorithm to simulate the mixed time-changed Erlang queue in which a different distribution is utilized for the random variable I that appears in the algorithm used in [4].
The simulation algorithm is described as follows.
Step 1. Initialize the process by setting
Step 2. Let ${T_{n}}$ be the nth jump time of the mixed time-changed Erlang queue ${\{{\mathcal{Q}^{{\alpha _{1}},{\alpha _{2}}}}(t)\}_{t\ge 0}}$. Assume that the queue has been simulated till the nth event. If $\mathcal{Q}({T_{n}})=0$, generate a random variable $I\stackrel{d}{=}{T^{{\alpha _{1}},{\alpha _{2}}}}$; otherwise, generate $I\stackrel{d}{=}{S^{{\alpha _{1}},{\alpha _{2}}}}$, where the distributions of ${T^{{\alpha _{1}},{\alpha _{2}}}}$ and ${S^{{\alpha _{1}},{\alpha _{2}}}}$ are given in (49) and (59), respectively.
Step 4. Whenever ${\mathcal{Q}^{{\alpha _{1}},{\alpha _{2}}}}({T_{n}})=(0,0)$, set
If ${\mathcal{Q}^{{\alpha _{1}},{\alpha _{2}}}}({T_{n}})\in \mathcal{H}$, generate a random variable U such that it is uniformly distributed on $(0,1)$. For $U\lt \frac{\lambda }{\lambda +k\mu }$, set
\[ \mathcal{N}({T_{n+1}})=\mathcal{N}({T_{n}})+1,\hspace{2em}\mathcal{S}({T_{n+1}})=\mathcal{S}({T_{n}}).\]
For $U\ge \frac{\lambda }{\lambda +k\mu }$, phase is served and the state is updated according to the following cases:
Step 5. The sample path of the whole queue can be obtained by setting
Remark 6.
On substituting ${c_{1}}=1$ or ${c_{2}}=1$ in (20), (42), Theorem 1, Lemma 2 and Proposition 3, the results for mixed time-changed Erlang queue reduce to the corresponding results for fractional Erlang queue. On substituting ${c_{1}}=1$ in (46), (68) and Theorem 3, Theorem 6 - Theorem 9, the results for mixed time-changed Erlang queue reduce to the corresponding results for fractional Erlang queue (see [4]).