Analysis of an N-policy MX/M/1 Two-phase Queueing System with State-dependent Arrival Rates and Unreliable Server

In many queueing systems that are observed in telecommunication networks, production systems etc. the arrival rates are not constant, but depends on the state of the server. Various authors have studied queueing problems with state dependent arrival rates under different queue control policies.

One of the critical works done towards this path was performed by Yechiali and Naor [1]. Haris and Marchal [2] studied the M/G/1 queue with server vacations whose distributions can be state dependent. Yijun and Quanlin [3] examined the two-stage queueing system with state dependent vacation policy. They derived the stationary distributions of the queue length and the cycle time for the closed state-dependent vacation model.

Kumar and Chandan [4], analyzed the performance of two-phase M^X/E_K/1 queuing system with server startups and N-policy. They derived optimal threshold value of N by computing average queue size of the system and average cost considering three batch size distributions. Al Hanbali, & Boxma [5] studied Busy period analysis of the state dependent M/M/1/K queue. Vasanta Kumar et al. [8, 9] examined the performance of two-phase M^X/M/1 and M^X/E_k/1 queueing systems with N-policy and unreliable server, respectively.

Madhu Jain and Agarval [6] dealt with a state dependent batch arrival queueing system with modified Bernoulli vacation under N-policy. They proposed a method to find the optimal value of the threshold parameter to minimize the total expected cost. Singh et al. [10] considered a single server, state dependent Poisson arrival system. They used supplementary variable technique to obtain probability generating function of number of units in the system. In addition, some special cases are also provided. Charan Jeet Singh and Binay Kumar [11] investigated a batch arrival queueing system with unreliable server, state dependent arrival rates and two stages of heterogeneous service. They studied the transient and steady-state behaviour of the queue length distribution.

Banik [7] assessed state-dependent arrival in GI/BMSP/1/∞ queue. Rashmita [12] evaluated M^X/G/1 queueing model. They considered the state-dependent server vacation and derived the explicit expressions for the system size. Charan Jeet Singh et al. [15] investigated a single repairable server queueing system with bulk input and state dependent rates considering the general distributions for the repair, delay to repair and service processes.

Hanumantha Rao et al. [13] studied two-phase M/M/1 queueing system with server breakdowns, delayed repair, and impatient customers. They found expected loss due to balking and reneging. Numerical illustrations are presented to support the model. Recently, Hanumantha Rao et al. [14] examined the M/M/1 two -phase queueing system with state dependent arrival rates under N policy.

As observed from the review of literature two-phase M^X/M/1 queueing system with state dependent arrival rates, N-policy and unreliable server has not been studied so far. Thus the present study is aimed at the analysis of this queueing system. The remainder of this paper is organized as follows: Section 2 describes the model and its assumptions. Section 3 describes the analysis under steady state, Section 4 presents the performance analysis of the system, Section 5 describes the cost function and the optimal operating policy, Section 6 describes the sensitivity analysis and summary is presented in Section 7.

In the current paper, the following notations have been utilized as below:

Π_0,m,0=Steady state probability when m customers are in the batch queue and the server is on vacation state, m= 0,1, 2, ...(N-1).

Π_1,m,0=Steady state probability when m customers are in the batch queue while the server is in startup state, where m = N, N+1, N+2, …

Π_2,m,0=Steady state probability when m customers are in the batch which is in batch service state, m = 1, 2, 3…

Π_3,m,0= Steady state probability when m customers are in the batch which is in batch service, while the server is found to be broken down and waiting for repair state, m = 1, 2, 3, …

Π_4,m,0= Steady state probability when m customers are in the batch which is in batch service, while the server is undergoing repair, m = 1, 2, 3, …

Π_5,m,n= Steady state probability when m customers are in the batch queue service and n customers in the individual service while the server is in individual service state, m = 0, 1, 2…, and n = 1, 2, 3, …

Π_6,m,n= Steady state probability when m customers are in the batch service and n customers in the individual service state, while the server is in individual service but found to be broken down and waiting for repair, m = 0,1,2 …, and n = 1,2,3….

Π_7,m,n = Steady state probability when m customers are in the batch service and n customers in the individual service queue while the server is in individual service, but undergoing repair, m =0, 1, 2…, and n =1, 2, 3….

The steady state equations for the queue length distribution

$\lambda_{1} \Pi_{0,0,0}=\mu \Pi_{5,0,1}$    (1)

$\lambda_{1} \Pi_{0, \mathrm{m}, 0}=\lambda_{1} \sum_{l=1}^{\mathrm{m}} \mathrm{a}_{1} \Pi_{0, \mathrm{m}-1,0}, 1 \leq \mathrm{m} \leq(\mathrm{N}-1)$    (2)

$\left(\lambda_{1}+\theta\right) \Pi_{1, \mathrm{N}, 0}=\lambda_{1} \sum_{l=1}^{\mathrm{N}} \mathrm{a}_{1} \Pi_{0, \mathrm{N}-1,0}$    (3)

$\left(\lambda_{1}+\theta\right) \Pi_{1, \mathrm{m}, 0}=\lambda_{1} \sum_{\mathrm{l}=1}^{\mathrm{m}-\mathrm{N}} \mathrm{a}_{1} \Pi_{1, \mathrm{m}-1,0}+\lambda_{1} \sum_{l=\mathrm{m}-(\mathrm{N}-1)}^{\mathrm{m}} \mathrm{a}_{1} \Pi_{0, \mathrm{m}-1,0}, \mathrm{m} \geq \mathrm{N}+1$    (4)

$\left(\lambda_{2}+\beta+\alpha_{1}\right) \Pi_{2, \mathrm{m}, 0}=\lambda_{2} \sum_{l=1}^{\mathrm{m}} \mathrm{a}_{1} \Pi_{2, \mathrm{m}-1,0}+\mu \Pi_{5, \mathrm{m}, 1}+\gamma \Pi_{4, \mathrm{m}, 0}, 1 \leq \mathrm{m} \leq(\mathrm{N}-1)$    (5)

$\left(\lambda_{2}+\beta+\alpha_{1}\right) \Pi_{2, \mathrm{m}, 0}=\lambda_{2} \sum_{l=1}^{\mathrm{m}} \mathrm{a}_{1} \Pi_{2, \mathrm{m}-1,0}+\mu \Pi_{5, \mathrm{m}, 1}+\gamma \Pi_{4, \mathrm{m}, 0}+\theta \Pi_{1, \mathrm{m}, 0}, \mathrm{m} \geq \mathrm{N}$    (6)

$\left(\lambda_{3}+\delta\right) \Pi_{3,1,0}=\alpha_{1} \Pi_{2,1,0}$    (7)

$\left(\lambda_{3}+\delta\right) \Pi_{3, \mathrm{m}, 0}=\alpha_{1} \Pi_{2, \mathrm{m}, 0}+\lambda_{3} \sum_{\mathrm{l}=1}^{\mathrm{m}} \mathrm{a}_{1} \Pi_{3, \mathrm{m}-1,0}, \quad \mathrm{m}>1$    (8)

$\left(\lambda_{3}+\gamma\right) \Pi_{4,1,0}=\delta \Pi_{3,1,0}$    (9)

$\left(\lambda_{3}+\gamma\right) \Pi_{4, \mathrm{m}, 0}=\delta \Pi_{3, \mathrm{m}, 0}+\lambda_{3} \sum_{l=1}^{\mathrm{m}} \mathrm{a}_{1} \Pi_{3, \mathrm{m}-1,0}, \mathrm{m}>1$    (10)

$\left(\lambda_{2}+\alpha_{2}+\mu\right) \Pi_{5,0, n}=\mu \Pi_{5,0, n+1}+\beta \Pi_{2, n, 0}+\gamma \Pi_{7,0, n}, n \geq 1$    (11)

                 $\left(\lambda_{2}+\alpha_{2}+\mu\right) \Pi_{5, \mathrm{m}, \mathrm{n}}=\lambda_{2} \sum_{\mathrm{l}=1}^{\mathrm{m}} \mathrm{a}_{1} \Pi_{\mathrm{5}, \mathrm{m}-1, \mathrm{n}}+\mu \Pi_{\mathrm{5}, \mathrm{m}, \mathrm{n}+1}+\gamma \Pi_{7, \mathrm{m}, \mathrm{n}}, \mathrm{m} \geq 1, \mathrm{\&n} \geq 1$    (12)

$\left(\lambda_{3}+\delta\right) \Pi_{6,0, \mathrm{n}}=\alpha_{2} \Pi_{5,0, \mathrm{n}}, \mathrm{n} \geq 1$    (13)

$\left(\lambda_{3}+\delta\right) \Pi_{6, \mathrm{m}, \mathrm{n}}=\alpha_{2} \Pi_{5, \mathrm{m}, \mathrm{n}}+\lambda_{3} \sum_{l=1}^{\mathrm{m}} \mathrm{a}_{1} \Pi_{6, \mathrm{m}-1, \mathrm{n}}, \mathrm{m} \geq 1, \mathrm{n} \geq 1$    (14)

$\left(\lambda_{3}+\gamma\right) \Pi_{7,0, \mathrm{n}}=\delta \Pi_{6,0, \mathrm{n}}, \mathrm{n} \geq 1$    (15)

$\left(\lambda_{3}+\gamma\right) \Pi_{7, \mathrm{m}, \mathrm{n}}=\delta \Pi_{6, \mathrm{m}, \mathrm{n}}+\lambda_{3} \sum_{\mathrm{l}=1}^{\mathrm{m}} \mathrm{a}_{1} \Pi_{7, \mathrm{m}-1, \mathrm{n}}, \mathrm{m} \geq \& \mathrm{n} \geq 1$    (16)

In the next step PGF’s (Probability generating functions) of queue size at an arbitrary time epoch are derived for different states of the system. Prior to that PGF’s are defined below:

$\mathrm{F}_{0}(\mathrm{s})=\sum_{\mathrm{m}=0}^{\infty} \Pi_{0, \mathrm{m}, 0} \mathrm{s}^{\mathrm{m}}$, $\mathrm{F}_{1}(\mathrm{s})=\sum_{\mathrm{m}=0}^{\infty} \Pi_{1, \mathrm{m}, 0} \mathrm{s}^{\mathrm{m}}$, $\mathrm{F}_{2}(\mathrm{s})=\sum_{\mathrm{m}=0}^{\infty} \Pi_{2, \mathrm{m}, 0} \mathrm{s}^{\mathrm{m}}$, $\mathrm{F}_{3}(\mathrm{s})=\sum_{\mathrm{m}=0}^{\infty} \Pi_{3, \mathrm{m}, \mathrm{0}} \mathrm{s}^{\mathrm{m}}$, $\mathrm{F}_{4}(\mathrm{s})=\sum_{\mathrm{m}=0}^{\infty} \Pi_{4, \mathrm{m}, 0} \mathrm{s}^{\mathrm{m}}$, $\mathrm{F}_{5}(\mathrm{s}, \mathrm{y})=\sum_{\mathrm{m}=0}^{\infty} \sum_{\mathrm{n}=1}^{\infty} \Pi_{\mathrm{5}, \mathrm{m}, \mathrm{n}} \mathrm{s}^{\mathrm{m}} \mathrm{y}^{\mathrm{n}}$, $\mathrm{F}_{6}(\mathrm{s}, \mathrm{y})=\sum_{\mathrm{m}=0}^{\infty} \sum_{\mathrm{n}=1}^{\infty} \Pi_{6, \mathrm{m}, \mathrm{n}} \mathrm{s}^{\mathrm{m}} \mathrm{y}^{\mathrm{n}}$, $\mathrm{F}_{7}(\mathrm{s}, \mathrm{y})=\sum_{\mathrm{m}=0}^{\infty} \sum_{\mathrm{n}=1}^{\infty} \Pi_{7, \mathrm{m}, \mathrm{n}} \mathrm{s}^{\mathrm{m}} \mathrm{y}^{\mathrm{n}}$, $S_{n}(s)=\sum_{m=0}^{\infty} \Pi_{6, m, n} s^{m}, T_{n}(s)=\sum_{m=0}^{\infty} \Pi_{7, m, n} s^{m}$ and $R_{n}(s)=\sum_{m=0}^{\infty} \Pi_{5, m, n} s^{m},|s| \leq 1,|y| \leq 1$

Let $B(s)=\sum_{m=1}^{\infty} a_{m} s^{m}$ be the probability generating function of the arrival batch size random variable X and $\mathrm{B}^{\prime}(\mathrm{s}), B^{\prime \prime}(\mathrm{s})$ represents the first and second order derivatives of  B(s) respectively.

From equation (1) to (16), using the PGFs, we will get $\mathrm{F}_{0}(\mathrm{s})=\Pi_{0,0,0} \mathrm{Y}_{\mathrm{N}}(\mathrm{s})$,

where 

$Y_{N}(s)=\sum_{m=0}^{N-1} y_{m} s^{m}, Y_{N}(1)=\sum_{m=0}^{N-1} y_{m} \& Y_{N}^{\prime}(1)=\sum_{m=1}^{N-1} m y_{m}$    (17)

$\left[\lambda_{1}(1-\mathrm{B}(\mathrm{s})+\theta] \mathrm{F}_{1}(\mathrm{s})=\lambda_{1} \Pi_{0,0,0}+\lambda_{1}(\mathrm{B}(\mathrm{s})-1) \mathrm{F}_{0}(\mathrm{s})\right.$    (18)

$\left[\lambda_{2}\left(1-\mathrm{B}(\mathrm{s})+\beta+\alpha_{1}\right] \mathrm{F}_{2}(\mathrm{s})=\mu \mathrm{R}_{1}(\mathrm{s})+\theta \mathrm{F}_{1}(\mathrm{s})-\lambda_{1} \Pi_{0,0,0}+\gamma \mathrm{F}_{4}(\mathrm{s})\right.$    (19)

$\left[\lambda_{3}(1-\mathrm{B}(\mathrm{s})+\delta] \mathrm{F}_{3}(\mathrm{s})=\alpha_{1} \mathrm{F}_{2}(\mathrm{s})\right.$    (20)

$\left[\lambda_{3}(1-\mathrm{B}(s)+\gamma] \mathrm{F}_{4}(\mathrm{s})=\delta \mathrm{F}_{3}(\mathrm{s})\right.$    (21)

$\left[\lambda_{2}\left(1-\mathrm{B}(s)+\alpha_{2}+\mu\right] \mathrm{R}_{\mathrm{n}}(\mathrm{s})=\mu \mathrm{R}_{\mathrm{n}+1}(\mathrm{s})+\gamma \mathrm{T}_{\mathrm{n}}(\mathrm{s})+\beta \Pi_{2, \mathrm{n}, 0}\right.$    (22)

$\left[\lambda_{2} y\left(1-B(s)+\alpha_{2} y+\mu(y-1)\right] F_{5}(s, y)=-\mu y R_{1}(s)+\gamma y F_{7}(s, y)+\beta y F_{2}(y)\right.$    (23)

$\left[\lambda_{3}(1-\mathrm{B}(\mathrm{s})+\delta] \mathrm{S}_{\mathrm{n}}(\mathrm{s})=\alpha_{2} \mathrm{R}_{\mathrm{n}}(\mathrm{s})\right.$    (24)

$\left[\lambda_{3}(1-\mathrm{B}(\mathrm{s})+\delta] \mathrm{F}_{6}(\mathrm{s}, \mathrm{y})=\alpha_{2} \mathrm{F}_{5}(\mathrm{s}, \mathrm{y})\right.$    (25)

$\left[\lambda_{3}(1-\mathrm{B}(\mathrm{s})+\gamma] \mathrm{T}_{\mathrm{n}}(\mathrm{s})=\delta \mathrm{S}_{\mathrm{n}}(\mathrm{s})\right.$    (26)

$\left[\lambda_{3}(1-\mathrm{B}(s)+\gamma] \mathrm{F}_{7}(\mathrm{s}, \mathrm{y})=\delta \mathrm{F}_{6}(\mathrm{s}, \mathrm{y})\right.$    (27)

Put the value of $\mathrm{F}_{6}(\mathrm{s}, \mathrm{y})$ in equation (27) we obtain

$\left[\lambda_{3}(1-\mathrm{B}(\mathrm{s})+\gamma] \mathrm{F}_{7}(\mathrm{s}, \mathrm{y})=\frac{\delta \alpha_{2} \mathrm{F}_{5}(\mathrm{s}, \mathrm{y})}{\left[\lambda_{3}(1-\mathrm{A}(s)+\delta]\right.}\right.$    (28)

Put y = s in Eq. (23), we get

$\left[\lambda_{2} s(1-B(s))+\alpha_{2} s+\mu(s-1)\right] F_{5}(s, s)=-\mu s R_{1}(s)+\gamma s F_{7}(s, s)+\beta s F_{2}(s)$

Put the values of $\mathrm{F}_{2}(\mathrm{s})$ and $\mathrm{F}_{7}(\mathrm{s}, \mathrm{s})$ obtained through equations (19) and (28) respectively, we get

$\left[\left(\lambda_{2} s\left(1-B(s)+\alpha_{2} s+\mu(s-1)\left(\lambda_{3}(1-B(s)+\delta)\left(\lambda_{3}(1-B(s)+\gamma)-\gamma \delta \alpha_{2} s\right] F_{5}(s, s)\right.\right.\right.\right.=\left(\beta s F_{2}(s)-\mu s R_{1}(s)\right)\left[\lambda_{3}(1-B(s)+\right.\gamma]\left[\lambda_{3}(1-\mathrm{B}(\mathrm{s})+\delta]\right.$    (29)

Put s=1 and y=1 in Eqns. (17), (18), (19), (20), (21), (25), (28), and (29), we get

$\mathrm{F}_{0}(1)=\mathrm{Y}_{\mathrm{N}}(1) \Pi_{0,0,0}$    (30)

$\mathrm{F}_{1}(1)=\frac{\lambda_{1} \Pi_{0,0,0}}{\theta}$    (31)

$\mathrm{F}_{2}(1)=\frac{\mu}{\beta} \mathrm{R}_{1}(1)$    (32)

$\mathrm{F}_{3}(1)=\frac{\alpha_{1}}{\delta} \mathrm{F}_{2}(1)$    (33)

$\mathrm{F}_{4}(1)=\frac{\alpha_{1}}{\gamma} \mathrm{F}_{2}(1)$    (34)

$\mathrm{F}_{5}(1,1)=\frac{\frac{\lambda_{2} A^{\prime}(1)}{\mu}\left(1+\frac{\lambda_{3}}{\lambda_{2}}\left(\frac{\alpha_{1}}{\delta}+\frac{\alpha_{1}}{\gamma}\right)\right) \mathrm{F}_{2}(1)+\frac{\theta}{\mu} F_{1}^{\prime}(1)}{1-\frac{\lambda_{2} A^{\prime}(1)}{\mu}\left(1+\frac{\lambda_{3}}{\lambda_{2}}\left(\frac{\alpha_{2}}{\delta}+\frac{\alpha_{2}}{\gamma}\right)\right)}$    (35)

$\mathrm{F}_{6}(1,1)=\frac{\alpha_{2}}{\delta} \mathrm{F}_{5}(1,1)$    (36)

$\mathrm{F}_{7}(1,1)=\frac{\alpha_{2}}{\gamma} \mathrm{F}_{5}(1,1)$    (37)

Thequeue length distribution is given by

$\mathrm{F}(\mathrm{s}, \mathrm{s})=\mathrm{F}_{0}(\mathrm{s})+\mathrm{F}_{1}(\mathrm{s})+\mathrm{F}_{2}(\mathrm{s})+\mathrm{F}_{3}(\mathrm{s})+\mathrm{F}_{4}(\mathrm{s})+\mathrm{F}_{5}(\mathrm{s}, \mathrm{s})+\mathrm{F}_{6}(\mathrm{s}, \mathrm{s})+\mathrm{F}_{7}(\mathrm{s}, \mathrm{s})$

Probability that the server is not doing any service is given by

$\mathrm{F}_{0}(1)+\mathrm{F}_{1}(1)=1-\frac{\lambda_{2} B^{\prime}(1)}{\beta}-\frac{\lambda_{3} B^{\prime}(1)}{\beta}\left(\frac{\alpha_{1}}{\gamma}+\frac{\alpha_{1}}{\delta}\right)-\frac{\lambda_{2} \mathrm{B}^{\prime}(1)}{\mu}-\frac{\lambda_{3} B^{\prime}(1)}{\mu}\left(\frac{\alpha_{1}}{\gamma}+\frac{\alpha_{1}}{\delta}\right)$

This gives

$\left(\frac{\lambda_{1}+\theta Y_{N}(1)}{\theta}\right) \Pi_{0,0,0}=1-\rho_{1}-\rho_{2}$    (38)

where  $\rho_{1}=\frac{\lambda_{2} B^{\prime}(1)}{\beta}\left(1+\frac{\lambda_{3}}{\lambda_{2}}\left(\frac{\alpha_{1}}{\gamma}+\frac{\alpha_{1}}{\delta}\right)\right)$, $\rho_{2}=\frac{\lambda_{2} B^{\prime}(1)}{\mu}\left(1+\frac{\lambda_{3}}{\lambda_{2}}\left(\frac{\alpha_{2}}{\gamma}+\frac{\alpha_{2}}{\delta}\right)\right)$

and $\rho=\rho_{1}+\rho_{2}$ is the utilizing factor of the system.

Hence, $\Pi_{0,0,0}=(1-\rho) \frac{\theta}{\left(\lambda_{1}+\theta Y_{N}(1)\right)}$.

The Normalizing condition is

$\mathrm{F}(1,1)=\mathrm{F}_{0}(1)+\mathrm{F}_{1}(1)+\mathrm{F}_{2}(1)+\mathrm{F}_{3}(1)+\mathrm{F}_{4}(1)+\mathrm{F}_{5}(1,1)+\mathrm{F}_{6}(1,1)+\mathrm{F}_{7}(1,1)=1$    (39)

Using the normalizing condition, we get

$R_{1}(1)=\frac{\beta\left[\left(\rho_{1}+\rho_{2}\right)\left(1-\rho_{2}\right) \mu-\left(1+\frac{\alpha_{2}}{\gamma}+\frac{\alpha_{2}}{\delta}\right) \theta \mathrm{F}_{1}^{\prime}(1)\right]}{\mu\left[\left(1+\frac{\alpha_{1}}{\gamma}+\frac{\alpha_{1}}{\delta}\right) \mu+\left(\alpha_{2}-\alpha_{1}\right)\left(\frac{1}{\gamma}+\frac{1}{\delta}\right)\left(\lambda_{2}-\lambda_{3}\right) B^{\prime}(1)\right]}$    (40)

Let F(1,1)=1 in the forms $\lim _{y \rightarrow 1} F(1, y)=1$ and $\lim _{s \rightarrow 1} \mathrm{F}(s, 1)=1$ we find

$R_{1}^{\prime}(1)=\frac{\rho_{2}}{\left(1-\rho_{2}\right)}\left[\left(\frac{\lambda_{2} B^{\prime}(1)}{\mu}+\frac{\lambda_{3} B^{\prime}(1)}{\beta}\left(\frac{\alpha_{1}}{\delta}+\frac{\alpha_{1}}{\gamma}\right)\right) \mathrm{F}_{2}(1)+\frac{\lambda_{1} B^{\prime}(1)}{\mu}\left(1-\rho_{1}-\rho_{2}\right)\right]$    (41)

System managers are keen on limiting the total expenditure. For this reason, we establish the cost function in terms of adequate performance measures and interrelated cost components to ascertain the optimal threshold parameters.

Let C_A (N) be the average cost per unit of time. Then

$C_{A}(\mathrm{N})=\mathrm{C}_{\mathrm{h}} \mathrm{L}(\mathrm{N})+\mathrm{C}_{\circ}\left(\frac{\mathrm{W}_{\mathrm{b}}}{\mathrm{W}_{\mathrm{c}}}+\frac{\mathrm{W}_{\mathrm{i}}}{\mathrm{W}_{\mathrm{c}}}\right)+\mathrm{C}_{\mathrm{m}}\left(\frac{\mathrm{W}_{\mathrm{s}}}{\mathrm{W}_{\mathrm{c}}}\right)+\mathrm{C}_{\mathrm{b}}\left(\frac{\mathrm{W}_{\mathrm{bb}}+\mathrm{W}_{\mathrm{db}}+\mathrm{W}_{\mathrm{bi}}+\mathrm{W}_{\mathrm{di}}}{\mathrm{W}_{\mathrm{c}}}\right)+\mathrm{C}_{\mathrm{s}}\left(\frac{1}{\mathrm{W}_{\mathrm{c}}}\right)-\mathrm{C}_{\mathrm{r}}\left(\frac{\mathrm{W}_{\mathrm{v}}}{\mathrm{W}_{\mathrm{c}}}\right)$    (60)

Neglecting the terms independent of N in C_A(N), from (61) we get the new cost function

$\mathrm{T}_{A}(\mathrm{N})=\left(\frac{\lambda_{1} B^{\prime}(1)}{\mu\left(1-\rho_{2}\right)}+1\right) \mathrm{Y}_{\mathrm{N}}^{\prime}(1) \Pi_{0,0,0} \frac{\mathrm{C}_{\mathrm{h}}}{k}+\mathrm{C}_{\mathrm{m}} \frac{\lambda_{1}}{\theta} \Pi_{0,0,0}+\lambda_{1} \Pi_{0,0,0} \mathrm{C}_{\mathrm{s}}-\mathrm{C}_{\mathrm{r}} \mathrm{Y}_{\mathrm{N}}(1) \Pi_{0,0,0}$

$=\frac{1}{\left(1-\rho_{2}\right)}\left[\left(\lambda_{1} B^{\prime}(1)+\mu\left(1-\rho_{2}\right)\right) \mathrm{Y}_{\mathrm{N}}^{\prime}(1) \frac{\mathrm{c}_{\mathrm{h}}}{k}+\mu\left(1-\rho_{2}\right)\left(C_{m} \frac{\lambda_{1}}{\theta}+\lambda_{1} C_{s}-C_{r} \mathrm{Y}_{\mathrm{N}}(1)\right)\right] \Pi_{0,0,0}$

$=\left[\left(\lambda_{1} B^{\prime}(1)+\mu\left(1-\rho_{2}\right)\right) \mathrm{Y}_{\mathrm{N}}^{\prime}(1) \frac{\mathrm{C}_{\mathrm{h}}}{k}+\mu\left(1-\rho_{2}\right)\left(\mathrm{C}_{\mathrm{m}} \frac{\lambda_{1}}{\theta}+\lambda_{1} \mathrm{C}_{\mathrm{s}}\right.\right.$

$\left.\left.-Y_{N}(1) \mathrm{C}_{\mathrm{r}}\right)\right]\left(\frac{\theta}{\lambda_{1}+\theta \mathrm{Y}_{\mathrm{N}}(1)}\right)\left(1-\rho_{1}-\rho_{2}\right)$    (61)

Subsequently for determination of the superior operating N-policy, minimizing C_A(N) in (60) is equivale to minimizing $\mathrm{T}_{A}(\mathrm{N})$ in (61).

It is tight to prove that TA(N) is convex but now we presented a technique to determine the optimal threshold N*.

Result

Utilizing the long run expected average cost criterion, the optimal threshold N* is given by

$\mathrm{N}^{*}=\min \left\{\mathrm{k} \geq 1 /\left(\sum_{\mathrm{n}=0}^{\mathrm{k}-1}(\mathrm{k}-\mathrm{n}) \mathrm{y}_{\mathrm{n}}+\frac{\mathrm{k} \lambda_{1}}{\theta}\right)>\frac{\lambda_{1}}{\mathrm{c}_{\mathrm{h}}}\left(\frac{\mathrm{C}_{\mathrm{m}}+\mathrm{C}_{\mathrm{r}}}{\theta}+\mathrm{C}_{\mathrm{s}}\right)\left(1-\rho_{2}\right)\right\}$    (62)

Proof: Let $\mathrm{J}(\mathrm{k})=\Sigma_{\mathrm{m}=1}^{\mathrm{k}} \mathrm{y}_{\mathrm{m}}+\frac{\lambda_{1}}{\theta}$ and $\mathrm{I}(\mathrm{k})=\sum_{\mathrm{m}=1}^{\mathrm{k}} \mathrm{my}_{\mathrm{m}}$.

Consider the following difference

$\Delta \mathrm{T}_{\mathrm{A}}(\mathrm{k})=\mathrm{T}_{\mathrm{A}}(\mathrm{k}+1)-\mathrm{T}_{\mathrm{A}}(\mathrm{k})=\frac{\mathrm{MY}_{\mathrm{k}}}{\left(1-\rho_{2}\right)} \frac{\mathrm{H}(\mathrm{k})}{\mathrm{J}(\mathrm{k}) \mathrm{J}(\mathrm{k}-1)}$

where $\Delta$ is the difference operator, and

$\mathrm{H}(\mathrm{k})=\mathrm{C}_{\mathrm{h}}(\mathrm{k} \mathrm{J}(\mathrm{k})-\mathrm{I}(\mathrm{k}))-\lambda_{1}\left(\frac{\mathrm{C}_{\mathrm{m}}+\mathrm{C}_{\mathrm{r}}}{\theta}+\mathrm{C}_{\mathrm{s}}\right)\left(1-\rho_{2}\right)$

$=\mathrm{C}_{\mathrm{h}}\left[\sum_{\mathrm{n}=0}^{\mathrm{k}-1}(\mathrm{k}-\mathrm{n}) \mathrm{y}_{\mathrm{k}}+\frac{\mathrm{k} \lambda_{1}}{\theta}\right]-\lambda_{1}\left(\frac{\mathrm{C}_{\mathrm{m}}+\mathrm{C}_{\mathrm{r}}}{\theta}+\mathrm{C}_{\mathrm{s}}\right)\left(1-\rho_{2}\right)$    (63)

By definition, $\mathrm{C}_{\mathrm{h}}\left[\sum_{\mathrm{n}=0}^{\mathrm{k}-1}(\mathrm{k}-\mathrm{n}) \mathrm{y}_{\mathrm{k}}+\frac{\mathrm{k} \lambda_{1}}{\theta}\right]>0$ and $\frac{\mathrm{MY}_{\mathrm{k}}}{\left(1-\rho_{2}\right) \mathrm{J}(\mathrm{k}) \mathrm{J}(\mathrm{k}-1)}>0$

Then it follows that $\Delta \mathrm{T}_{\mathrm{A}}(\mathrm{k})>0$.

Thus, the sign of H(k) determines whether T_A(N) increases or decreases,

Let m be the first k such that $H(k)>0$, then we have

$H(m+1)=\mathrm{C}_{\mathrm{h}}[(\mathrm{m}+1) \mathrm{J}(\mathrm{m}+1)-\mathrm{P}(\mathrm{m}+1)]-\lambda_{1}\left(\frac{\mathrm{C}_{\mathrm{m}}+\mathrm{C}_{\mathrm{r}}}{\theta}+\mathrm{C}_{\mathrm{s}}\right)\left(1-\rho_{2}\right)=\mathrm{H}(\mathrm{m})+\mathrm{C}_{\mathrm{h}} P(m)$  

It follows that H(m+1)> H(m).

Hence for some n>m, we have $\mathrm{T}_{\mathrm{A}}(\mathrm{n})>\mathrm{T}_{\mathrm{A}}(\mathrm{m})$.

Let N^*be the optimal value of N, which minimizes $\mathrm{T}_{\mathrm{A}}(\mathrm{N})$, Then from equation (63) we have

$\mathrm{N}^{*}=\min \left\{\mathrm{k} \geq 1 /\left(\sum_{\mathrm{n}=0}^{\mathrm{k}-1}(\mathrm{k}-\mathrm{n}) \mathrm{y}_{\mathrm{n}}+\frac{\mathrm{k} \lambda_{1}}{\theta}\right)>\frac{\lambda_{1}}{\mathrm{c}_{\mathrm{h}}}\left(\frac{\mathrm{C}_{\mathrm{m}}+\mathrm{C}_{\mathrm{r}}}{\theta}+\mathrm{C}_{\mathrm{s}}\right)\left(1-\rho_{2}\right)\right\}$    (64)

therefore, the premier threshold of N may be evaluated from equation (64), through selecting the pleasant value of k, that is one of the integers surrounding ‘N’.

Also, note that if  $\frac{\mathrm{C}_{\mathrm{h}}}{\left(\frac{\mathrm{C} \mathrm{m}+\mathrm{Cr}}{\theta}+\mathrm{C}_{\mathrm{s}}\right)}>\lambda_{1}\left(1-\rho_{2}\right)$, the optimal threshold N^*should be 1.

To carry out the sensitivity analysis, we assume that the arrival batch size follows geometric distribution. Then $\mathrm{b}_{\mathrm{k}}=\mathrm{P}(\mathrm{X}=\mathrm{k})=\mathrm{p}(1-\mathrm{p})^{\mathrm{k}-1}, 0<p<1, k=1,2, \ldots$

with probability generating function $\mathrm{B}(\mathrm{s})=\frac{1-\mathrm{p}}{(1-\mathrm{ps})}$. $E(X)=B^{\prime}(1)=\frac{1}{p}$ and $\mathrm{E}(\mathrm{X}(\mathrm{X}-1))=B^{\prime \prime}(1)=\frac{2(1-\mathrm{p})}{\mathrm{p}^{2}}$. Then

$L(\mathrm{N})=\mathrm{Y}_{\mathrm{N}}^{\prime}(1) \Pi_{0,0,0}+\frac{\lambda_{1}}{\mathrm{p}} \frac{\left(\lambda_{1}+\theta \mathrm{Y}_{\mathrm{N}}(1)\right)}{\theta^{2}} \Pi_{0,0,0}+\left(1+\frac{\alpha_{1}}{\gamma}+\frac{\alpha_{1}}{\delta}\right) \mathrm{F}_{2}^{\prime}(1)+\frac{\lambda_{3}}{\mathrm{p}}\left(\frac{\alpha_{1}}{\delta^{2}}+\frac{\alpha_{1}}{\delta \gamma}+\frac{\alpha_{1}}{\gamma^{2}}\right) \mathrm{F}_{2}(1)+\frac{\lambda_{3}}{\mathrm{p}}\left(\frac{\alpha_{2}}{\delta^{2}}+\frac{\alpha_{2}}{\delta \gamma}+\frac{\alpha_{2}}{\gamma^{2}}\right) \mathrm{F}_{5}(1,1)+\left(1+\frac{\alpha_{2}}{\delta}+\frac{\alpha_{2}}{\gamma}\right) \mathrm{F}_{5}^{\prime}(1,1)$

where,  $\Pi_{0,0,0}=\left(1-\rho_{1}-\rho_{2}\right)\left(\frac{\theta}{\lambda_{1}+\theta Y_{N}(1)}\right)$, $\rho_{1}=\frac{\lambda_{2}}{p \beta}\left(1+\frac{\lambda_{3}}{\lambda_{2}}\left(\frac{\alpha_{1}}{\gamma}+\frac{\alpha_{1}}{\delta}\right)\right)$, and $\rho_{2}=\frac{\lambda_{2}}{p \mu}\left(1+\frac{\lambda_{3}}{\lambda_{2}}\left(\frac{\alpha_{2}}{\delta}+\frac{\alpha_{2}}{\gamma}\right)\right)$.