# Multi-supplier Systems with Seasonal Demand

Ekaterina Bulinskaya, Larisa Afanasyeva

► To cite this version:

Ekaterina Bulinskaya, Larisa Afanasyeva. Multi-supplier Systems with Seasonal Demand. International Conference on Advances in Production and Management Systems (APMS), Sep 2009, Paris, France. pp.267-274, 10.1007/978-3-642-16358-6_34. hal-01060459

HAL Id: hal-01060459

https://inria.hal.science/hal-01060459

Submitted on 4 Sep 2014

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Multi-supplier Systems with Seasonal Demand*

Ekaterina Bulinskaya and Larisa Afanasyeva

Department of Mathematics and Mechanics,
Moscow State University, Moscow 119991, Russia,
ebulinsk@mech.math.msu.su

**Abstract.** Multi-supplier inventory systems with seasonal demand are investigated via queuing theory methods in the framework of cost approach. Optimal distribution of replenishment orders between suppliers is obtained under assumption of make-to-order strategy. To this end asymptotic analysis of systems behavior is performed. We apply various mathematical tools such as weak convergence and use properties of regenerative, doubly stochastic Poisson and Markov processes.

**Keywords:** Multi-supplier systems, seasonal demand, queuing theory approach.

## 1 Introduction

The aim of this research is twofold. On one hand, we would like to further develop the investigation of inventory systems via queuing theory methods. One of the first interpretations of inventory models as queuing ones was given in the book [4]. Since then many researchers have used such analogies, see, e.g., [2], [7]. On the other hand, we are interested in the study of seasonal demand impact on the inventory policies in systems with many suppliers. A single supplier is assumed to be available in many inventory control models. However, there exist such situations in which more than one supplier is necessary to reduce the total system cost or to sustain a desired service standard. Multi-supplier strategies can create a suppliers competition and ensure their providing faster delivery, see, e.g., [5]. The models considered below can be applied to logistics in agriculture or fuel supply in countries with continental climate. Due to lack of space we omit almost all proofs. Moreover, we consider only make-to-order (MTO) strategy.

## 2 Model Description

The system consists of one dealer (or producer) and *N* suppliers. Let customers (each demanding a unit of product) arrive according to a Poisson flow with intensity $\lambda(t)$. To take into account the seasonality of demand we assume $\lambda(t)$ to be periodic function with the period *T*. A customer arriving at time *t* is satisfied immediately if inventory on hand is positive, waiting otherwise. It is

<footnote>* The research was supported by RFBR grant 10-01-00266a.</footnote>
also supposed that each time when a product unit is demanded an order is sent to the i-th supplier with probability $\alpha_i$, $\alpha_i \ge 0$, $\sum_{i=1}^N \alpha_i = 1$. Such a procedure is the result of application of continuous-review ($S-1, S$) policy for any $S$.

Each supplier has a random delivery time. Suppose that arriving orders are served according to FIFO rule (First In First Out). The order processing time of the i-th supplier has a distribution function $B_i(x)$ with finite mean $b^{(i)}$ and second moment $b_2^{(i)}$. It is well known that under above assumptions the input flow of orders to the i-th supplier is a Poisson one with periodic intensity $\alpha_i \lambda(t)$. Thus, we have an open queuing network consisting of $N$ independent queuing systems of type $M(t)|GI|1|\infty$ working in parallel.

There exists a limit cyclic distribution $P_j^{(i)}(t) = \lim_{n \to \infty} P(q_i(nT+t) = j)$ iff $\rho_i = \alpha_i \lambda b^{(i)} < 1$, $i = 1, \dots, N$, where $\lambda = T^{-1} \int_0^T \lambda(y) dy$ and $q_i(t)$ is the orders number in the i-th node of network at time $t$, see, e.g., [1].

## 2.1 Objective Function

We consider optimization of the system performance from the dealer view-point. Namely, the objective function we are going to minimize is the long-run average costs per unit time under MTO strategy ($S=0$). We take into account only unit stock-out cost $c_s$ and want to choose $\alpha_i$, $i = 1, \dots, N$, providing the minimum of the function

$$W(\alpha_1, \dots, \alpha_N) = c_s \sum_{j=1}^{\infty} j P_j$$

with $P_j = T^{-1} \int_0^T P_j(y) dy$, here $P_j(y) = P(\sum_{i=1}^N q_i(y) = j)$.

A model of this type was considered in [2] under assumption that Poisson input has a constant intensity $\lambda$ and each of the nodes is $M|M|1$ system.

## 2.2 Seasonal Demand and Periodic Orders Processing

Taking into account the seasonal demand we get, in the simplest case, the system of the type $M(t)|M(t)|1|\infty$ in each node. Thus, the input flow to the network is Poisson with intensity $\lambda(t)$. Service intensity $\nu_i(t)$ in the i-th node also depends on time. Moreover, functions $\lambda(t)$ and $\nu_i(t)$ are periodic with period $T$. The cyclic limit distribution exists if $\alpha_i \lambda < \nu_i$, $i = 1, \dots, N$, with $\lambda = T^{-1} \int_0^T \lambda(y) dy$ and $\nu_i = T^{-1} \int_0^T \nu_i(y) dy$. A system $M(t)|M(t)|1|\infty$ was investigated in [6] under assumption that initially the system has $k$ orders. Introducing a time transformation $\tau =: \tau(t) = \int_0^t \nu(y) dy$ one can obtain the following system of differential equations for the probabilities $Q_j(\tau) = P(q(\tau) = j)$

$$Q_0'(\tau) = -Q_0(\tau) + \gamma(\tau)Q_1(\tau), \\ Q_j' = -(1+\gamma(\tau))Q_j(\tau) + Q_{j-1}(\tau) + \gamma(\tau)Q_{j+1}(\tau), \quad j > 0, \quad (1)$$

where $\gamma(\tau) = \nu(\tau)\lambda^{-1}(\tau)$.
If $\nu_i(t) = \beta_i \lambda(t)$ then for the i-th node we have $\gamma_i(t) = \beta_i$ and (1) is a system of differential equations with constant coefficients. It is easily shown that, for $\alpha_i < \beta_i$, $i = 1, \dots, N$, there exist

$$ \lim_{\tau \to \infty} Q_j^{(i)}(\tau) = \lim_{n \to \infty} P(q_i(nT + t) = j) = \rho_i^j(1 - \rho_i), \quad j \ge 0. $$

Thus, the cyclic distribution does not depend on time and coincides with that considered in [2]. If $\gamma_i(t) \neq \beta_i$ it is possible to propose an algorithm for solving (1), hence, for calculation of objective function and its minimization.

## 2.3 Arbitrary Service Times

To analyze a system of type $M(t)|GI|1|\infty$ with an arbitrary service distribution function $B(x)$ having two finite moments $b$ and $b_2$ it is useful to study the Markov process $(q(t), \zeta_t)$. Here $q(t)$ is the customers number in the system and $\zeta_t$ is the time a customer has been already served (if $q(t) > 0$). Introducing $P_j(t, x) = P(q(t) = j, \zeta_t \le x)$ and densities $p_j(t, x) = \partial P_j(t, x)/\partial x$ we get, as usual, the differential equations for the functions $g_0(t, x) = p_0(t, x)$ and $g_j(t, x) = p_j(t, x)/[1 - B(x)]$, $j > 0$. Solving the system one obtains the cyclic distribution.

Considering MTO strategy we can use the following procedure to calculate the objective function. Denote by $V(t)$ virtual waiting time. As shown in [1], in cyclic regime $m(t) = EV(t)$ has the following form

$$ m(t) = \frac{\lambda b_2}{2(1 - \rho)} - \frac{T(1 - \rho)}{2} + \pi(t) - \frac{1}{T(1 - \rho)} \int_t^{t+T} \pi(y) P_0(y) dy $$

with $\pi(t) = b\Lambda(t) - t$, $\Lambda(t) = \int_0^t \lambda(y) dy$.

According to Little's formula, see, e.g., [9]

$$ Eq = T^{-1} \int_{0}^{T} Eq(t) dt = \lambda T^{-1} \int_{0}^{T} m(t) dt + \rho. $$

In some cases it is possible to take $1 - \rho$ as a first approximation for $P_0(y)$.

Since it is impossible, with rare exceptions, to find the explicit form of cyclic distribution an important role belongs to the asymptotic analysis.

# 3 Asymptotic Analysis of Systems with Periodic Input

Considering approximations of cyclic distributions there arise a lot of questions to answer, among them the following:

- how one obtains the estimates of approximation precision using instead of initial periodic intensity $\lambda(t)$ a simpler one, for example, a step-function,

- whether it is possible to estimate the fluctuation of system characteristics, in particular $Eq(t)$, on interval $T$ (period) because that enables us to decide when we can disregard the dependence of intensity on time,

- it is useful to get some limit theorems concerning the heavy ($\rho \approx 1, \rho < 1$) or light ($\rho \approx 0$) traffic, moreover, we can take into account such properties of intensity as its being slow or quick varying function and so on.
## 3.1 Sensitivity Analysis

In order to find a class of functions $\lambda(t)$ for which one can use a certain approximation calculating the cyclic distribution it is necessary to answer the above mentioned questions. Moreover, for different values of parameters $\alpha_i$, $i = 1, \dots, N$, some of the network nodes may be in conditions of heavy traffic, the others having light traffic. Sensitivity analysis of queuing systems with respect to their input intensities is based on the paper [8] where it is established that the function $H(t, x) = P(V(t) \le x)$ is a probability of non-crossing a certain boundary by a compound Poisson process. Below we formulate one of the results using this function.

Consider a family of periodic intensities $\{\lambda_{\varepsilon}(t), \varepsilon > 0\}$. All the characteristics for the system with intensity $\lambda_{\varepsilon}(t)$ will be labeled by $\varepsilon$. Assume also that $\rho < 1$, $\rho_{\varepsilon} < 1$. For any $t \in [0, T]$, $x \in [0, \infty)$ and $\varepsilon > 0$ there exists

$$h(t, x, \varepsilon) = \sup_{v \ge 0, |x-v| \le \lambda_*^{-1}\varepsilon, 0 \le t-u \le \lambda_*^{-1}\varepsilon} \frac{\partial H(u, v)}{\partial v},$$

where $\lambda_* = \inf_{t \in [0,T]} \lambda(t)$.

**Theorem 1.** Let $\lambda_{\varepsilon} = \lambda$ and $\max_{t \in [0,T]} |\Lambda_{\varepsilon}(t) - \Lambda(t)| < \varepsilon$ then, for any $t \in [0,T]$ and $x \ge 0$, the following inequality is valid

$$|H(t, x) - H_{\varepsilon}(t, x)| \le 2\lambda_{*}^{-1}[h(t, x, 2\varepsilon) + \lambda^{*}]\varepsilon, \quad \lambda^{*} = \sup_{t \in [0,T]} \lambda(t).$$

In particular, for $\lambda(t) \equiv \lambda$ one has $|H_{\varepsilon}(t, x) - \Phi_{\lambda}(x)| \le 2\varepsilon$, here $\Phi_{\lambda}(x)$ is a stationary distribution of virtual waiting time $V(t)$ in a system $M|GI|1|\infty$ with input intensity $\lambda$.

Theorem 1 gives sufficient conditions for using the stationary distribution $\Phi_{\lambda}(x)$ as approximation of cyclic limit distribution.

## 3.2 Limit Theorems

**Heavy traffic.** Consider a family $S_{\varepsilon}$ of single-server queuing systems with service distribution function $B(x)$ and Poisson input with intensity $\lambda_{\varepsilon}(t)$ depending on time $t$. Heavy traffic means that $\rho_{\varepsilon} \nearrow 1$, as $\varepsilon \to \varepsilon_0$. Denote by $q_{\varepsilon}(t)$ the customers number in $S_{\varepsilon}$ at time $t$ and $G_{\varepsilon}(t, x) = P(q_{\varepsilon}(t) \le x)$ its cyclic distribution.

**Theorem 2.** Let $\lambda_{\varepsilon}(t) = \varepsilon\lambda(t)$ and $\rho = bT^{-1} \int_{0}^{T} \lambda(y) dy < 1$. Then

$$G_{\varepsilon}(t, x(1 - \varepsilon\rho)^{-1}) \to 1 - \exp(-2b^2x/b_2), \quad \text{as } \varepsilon \to \rho^{-1}.$$

*Proof.* Asymptotic behavior of $V_{\varepsilon}(t)$ and $bq_{\varepsilon}(t)$ is the same, as $\varepsilon \to \rho^{-1}$. Therefore it is sufficient to investigate the function $H_{\varepsilon}(t, x) = P(V_{\varepsilon}(t) \le x)$. This function satisfies the Takacs equation, see, e.g., [9]. Hence, using the periodicity,
it is possible to obtain for the Laplace transform $H_{\varepsilon}^*(t, s) = \int_{0}^{\infty} e^{-sx} H_{\varepsilon}(t, dx)$ the following expression

$$ \frac{s \int_{0}^{T} H_{\varepsilon}(u + t, 0) \exp\{-su + \varepsilon(1 - b^*(s))(\Lambda(t+u) - \Lambda(t))\} du}{1 - \exp\{(-s + \varepsilon\lambda(1 - b^*(s))T\}} $$

with $b^*(s) = \int_{0}^{\infty} e^{-sx} dB(x)$.

This enables us to establish that $H_{\varepsilon}^*(t, s(1 - \varepsilon\rho)) \to (1 + sb_2/2b)^{-1}$, whence easily follows the desired statement. □

The result of Theorem 2 can be used to obtain the representation of Eq$_{\varepsilon}(t)$, as $\varepsilon \nearrow \rho^{-1}$, in the form $\lambda b_2[(2b^2(1 - \varepsilon\rho) + u(t) + O(1 - \varepsilon\rho)]^{-1}$. Here $u(t)$ is given by

$$ c_1 \lambda^{-2} T^{-1} \int_{t}^{t+T} \lambda(y) [\Lambda(y) - \Lambda(t)] dy + (\rho T)^{-1} \int_{t}^{t+T} (y-t) \lambda(y) dy + c_2 $$

where $c_1 = (T - b + be^{-T/b})[Tb(1 - e^{-T/b})]^{-1}$ and $c_2$ is a constant that can be calculated. It is worth noting that points of extrema for $u(t)$ are obtained by solving the equation $\lambda(t) = \lambda$. This lets us estimate the fluctuation of the function Eq$_{\varepsilon}(t)$ on interval of length $T$ (period) and decide whether it is possible to approximate it by a constant when the traffic intensity is close to 1.

**Light traffic.** Let $\lambda_{\varepsilon}(t) = \varepsilon\lambda(t)$ and $\varepsilon \to 0$. It is not difficult to establish that

$$ P_{0}^{\varepsilon}(t) = P(q_{\varepsilon}(t) = 0) = H_{\varepsilon}(t, 0) = 1 - \varepsilon \int_{0}^{\infty} \lambda(t-x)[1-B(x)] dx + o(\varepsilon). $$

This is useful for obtaining the term of order $\varepsilon$ in expression of the virtual waiting time mean as well as the mean customers number. For getting the terms of higher order in $\varepsilon$ it is necessary to use the following theorem proved in [3].

**Theorem 3.** Let $\int_{0}^{\infty} x^{n^2+1+\delta} dB(x) < \infty$ for some $\delta > 0$ and positive integer $n$. Then there exist functions $F_i(t, y)$, $i = 1, \dots, n$, such that

$$ H_{\varepsilon}(t, y) = 1 - \sum_{i=1}^{n} \varepsilon^{i} F_{i}(t, y) + o(\varepsilon^{n}). $$

An algorithm for obtaining functions $F_i(t, y)$ is also provided.

**Slowly varying intensities.** It is natural to suppose intensity to be slowly varying treating the seasonal demand. As above we consider a family of systems $S_{\varepsilon}$ with intensity depending on $\varepsilon$. However now we assume $\lambda_{\varepsilon}(t) = \lambda(\varepsilon t)$ with periodic $\lambda(t)$. Thus period of $\lambda_{\varepsilon}(t)$ is equal to $T\varepsilon^{-1}$ and grows, as $\varepsilon \to 0$.

The main role in cyclic distribution analysis plays now function $\mu(t)$ given by

$$ \mu(t) = \max_{0 \le y \le t} [\mu(0) + \pi(t), \pi(t) - \pi(y)]. $$
If $b\lambda(y) \le 1$ for all $y$, setting $\mu(0) = 0$ we get $\mu(y) \equiv 0$. Otherwise we move the origin to the point of absolute minimum of $\mu(y)$ on interval of length $T$, getting $\mu(0) = 0$ and $\mu(t) = \pi(t) - \min_{0 \le y \le t} \pi(y)$.

**Theorem 4.** Let $\lambda(t)$ be piecewise continuous, $\rho = \lambda b < 1$ and $b_2 < \infty$. If $t_0$ is a continuity point for $\lambda(t)$, moreover, $\lambda(t_0) > 0$ and $\mu(t) = 0$ in some neighborhood of $t_0$, then $\Pi_{\varepsilon}(\varepsilon^{-1}t_0, z) \to P_{\lambda(t_0)}(z)$, as $\varepsilon \to 0$. Here $\Pi_{\varepsilon}(t, z)$ is a probability generating function of customers number in system $S_{\varepsilon}$, whereas $P_{\lambda}(z)$ corresponds to the system $M|GI|1|\infty$ with constant input intensity $\lambda$ and has the form

$$P_{\lambda}(z) = (1 - \rho)(1 - z)[b^{*}(\lambda(1 - z)) - z]^{-1}.$$

Moreover,

$$Eq_{\varepsilon}(\varepsilon^{-1}t_0) \to m(\lambda(t_0), b, b_2), \quad \text{as } \varepsilon \to 0,$$

$$\text{where } m(\lambda, b, b_2) = [2(1 - \rho)]^{-1}[\lambda b_2 + 2b(1 - \rho)].$$

*Proof* is based on the relation obtained for $H_{\varepsilon}^{*}(t, s)$ in Theorem 2.

Thus under assumptions of Theorem 4 one can use the stationary distribution for the system with input intensity $\lambda(t_0)$ as approximation for the cyclic distribution. The error of such approximation depends on the convergence rate to stationary distribution in system $M|GI|1|\infty$ with input intensity $\lambda(t_0)$.

Assume $\lambda(t)$ to be slowly varying step-function and $\lambda^* > b^{-1}$. In practice that means $b \ll T$. If $\lambda(t_0) > b^{-1}$ then $q(t_0)$ is well approximated by $\mu(t_0)/b$ according to the following

**Theorem 5.** If $\rho = \lambda b < 1$ then $\hat{q}_{\varepsilon}(t) = \varepsilon q_{\varepsilon}(\varepsilon^{-1}t)$ converges in probability to $\mu(t)$, as $\varepsilon \to 0$, for any $t \in [0, T]$.

*Proof* can be found in [1].

# 4 Application of Limit Theorems to Objective Function

Recall that under MTO strategy objective function is equal to

$$c_s \sum_{i=1}^{N} Eq_i = c_s T^{-1} \sum_{i=1}^{N} \int_{0}^{T} m(\alpha_i \lambda(t), b^{(i)}, b_2^{(i)}) dt.$$

The choice of approximation for $m(\cdot, \cdot, \cdot)$, as we have already seen, depends on values of $\rho_i = \alpha_i \lambda b^{(i)}$ and behavior of $\lambda(t)$ on interval $T$. Theorem 1 enables us to estimate the fluctuations of $m(\alpha_i \lambda(t), b^{(i)}, b_2^{(i)})$ and decide whether it can be replaced by the stationary mean $m(\alpha_i \lambda, b^{(i)}, b_2^{(i)})$ where $\lambda = T^{-1} \int_0^T \lambda(t) dt$ and expression for $m(\lambda, b, b_2)$ is given in Theorem 4.

If the answer is negative, it is necessary to use Theorem 2 for $\rho \approx 1$, that is, take $Eq_i = \lambda b_2^{(i)} (b^{(i)})^{-2} [2(1 - \alpha_i \lambda b^{(i)})]^{-1}$. For $\alpha_i$ such that $\rho \approx 0$, according to Theorem 3, we set $Eq_i = \alpha_i T^{-1} \int_0^\infty \lambda(t - y) [1 - B_i(y)] dy$. For $\rho_i$ sufficiently far
from 0 and 1 it is possible to apply Theorems 4 and 5. Thus, if $\alpha_i\lambda(t)b^{(i)} < 1$ for all $t \in [0, T]$ and $i = 1, \dots, N$ the objective function has the form

$$W(\alpha_1, \dots, \alpha_N) = c_s T^{-1} \sum_{i=1}^{N} \int_{0}^{T} \frac{\alpha_i \lambda(t) b_2^{(i)} + 2b^{(i)}(1 - \alpha_i \lambda(t)b^{(i)})}{2(1 - \alpha_i \lambda(t)b^{(i)})} dt.$$

If $\alpha_i\lambda(t)b^{(i)} > 1$ for some $i$ and $t$, then in a neighborhood of this point it is necessary to use $\mu(t)/b^{(i)}$ instead of corresponding $m(\alpha_i\lambda, b^{(i)}, b_2^{(i)})$.

Thus asymptotical analysis of systems behavior gives the possibility to propose an effective algorithm of $W(\alpha_1, \dots, \alpha_N)$ calculation and to find numerically the optimal policy of orders placement under seasonal demand. We can also use the reliability approach and choose parameters $\alpha_i$, $i = 1, \dots, N$, maximizing the service level of customers.

## 4.1 Example

Assume the input intensity to be a constant $\lambda$. Each supplier processes the orders as $M|GI|1$ system.

Objective function $\sum_{i=1}^{N} [\alpha_i \rho_i - \alpha_i^2 \rho_i^2 + (\alpha_i^2 \lambda^2 b_2^{(i)}/2)](1 - \alpha_i \rho_i)^{-1}$ is convex in $\alpha_1, \dots, \alpha_N$. We minimize it introducing the Lagrange multiplier $\beta$. Thus, we have to solve equations

$$\rho_i - \beta + \frac{\lambda^2 \alpha_i b_2^{(i)} (2 - \alpha_i \rho_i)}{2(1 - \alpha_i \rho_i)^2} = 0, \quad i = 1, \dots, N,$$

under assumptions $0 \le \alpha_i \le 1$, $\sum_{i=1}^{N} \alpha_i = 1$, $\alpha_i \rho_i < 1$, $i = 1, \dots, N$.

**Proposition 1.** Optimal probabilities $\alpha_i^*$, $i = 1, \dots, N$, have the form

$$\alpha_i^* = \rho_i^{-1} \left( 1 - \lambda \sqrt{b_2^{(i)}} \left( 2\rho_i(\beta - \rho_i) + \lambda^2 b_2^{(i)} \right)^{-1/2} \right).$$

Here $\beta$ is the solution of the following equation

$$\sum_{i=1}^{N} \rho_i^{-1} \left( 1 - \lambda \sqrt{b_2^{(i)}} \left( 2\rho_i(\beta - \rho_i) + \lambda^2 b_2^{(i)} \right)^{-1/2} \right) = 1$$

satisfying additionally $\rho_i \le \beta \le \rho_i + \lambda^2 b_2^{(i)} (2 - \rho_i)/2(1 - \rho_i)^2$, if $0 \le \rho_i < 1$, and $\beta \ge \rho_i$, if $\rho_i \ge 1$.

**Proposition 2.** The optimal $\alpha_i^*$, $i = 1, \dots, N$, can be obtained solving the Lagrange problem iff $\sum_{i=1}^{N} \rho_i^{-1} > 1$.

Using Maple 12 it is possible to calculate the optimal parameters.
Some numerical results are given below for $N = 2$, $\lambda = 2$, $b_2^{(1)} = b_2^{(2)} = 0.1$.

<table><thead><tr><th colspan="5">Case $\rho_1 < 1 < \rho_2$</th></tr><tr><th>$\rho_1$</th><th>$\rho_2$</th><th>$\beta$</th><th>$\alpha_1^*$</th><th>$\alpha_2^*$</th></tr></thead><tbody><tr><td>0,4</td><td>1,2</td><td>1,217325411</td><td>0,9597951872</td><td>0,04020481317</td></tr><tr><td>0,6</td><td>1,1</td><td>1,298578981</td><td>0,7194120625</td><td>0,2805879374</td></tr><tr><td>0,6</td><td>1,4</td><td>1,523859202</td><td>0,8084692052</td><td>0,1915307944</td></tr><tr><td>0,6</td><td>2</td><td>2,028167863</td><td>0,9416522358</td><td>0,05834776480</td></tr><tr><td>0,9</td><td>1,1</td><td>1,591743594</td><td>0,5632300892</td><td>0,4367699110</td></tr></tbody></table>

<table><thead><tr><th colspan="5">Case $\rho_1, \rho_2 < 1$</th></tr><tr><th>$\rho_1$</th><th>$\rho_2$</th><th>$\beta$</th><th>$\alpha_1^*$</th><th>$\alpha_2^*$</th></tr></thead><tbody><tr><td>0,2</td><td>0,4</td><td>0,5487697400</td><td>0,6947230170</td><td>0,3052769835</td></tr><tr><td>0,2</td><td>0,7</td><td>0,7222398081</td><td>0,9474492070</td><td>0,05255079300</td></tr><tr><td>0,4</td><td>0,8</td><td>0,9326329628</td><td>0,7603895062</td><td>0,2396104938</td></tr><tr><td>0,6</td><td>0,8</td><td>1,086054229</td><td>0,6036417640</td><td>0,3963582371</td></tr><tr><td>0,7</td><td>0,8</td><td>1,164767915</td><td>0,5471213116</td><td>0,4528786881</td></tr></tbody></table>

<table><thead><tr><th colspan="5">Case $\rho_2 > \rho_1 > 1$</th></tr><tr><th>$\rho_1$</th><th>$\rho_2$</th><th>$\beta$</th><th>$\alpha_1^*$</th><th>$\alpha_2^*$</th></tr></thead><tbody><tr><td>1,2</td><td>1,3</td><td>2,223372271</td><td>0,5214712735</td><td>0,4785287270</td></tr><tr><td>1,2</td><td>2,0</td><td>3,603736818</td><td>0,6211347510</td><td>0,3788652492</td></tr><tr><td>1,2</td><td>4,0</td><td>24,07589991</td><td>0,7624608414</td><td>0,2375391590</td></tr></tbody></table>

Other examples as well as results concerning make-to-stock strategy, omitted due to lack of space, will be published in the next paper.

## References

1. Afanas'eva, L.G.: On Periodic Distribution of Waiting Time Process. Lecture Notes in Math. Stability Problems for Stochastic Models, pp. 1–20 (1985)

2. Arda, Y., Hennet, J.-C.: Inventory Control in a Multi-supplier System. Int. J. Prod. Econom. 104, pp. 249–259 (2006)

3. Bashtova, E.: The Queue with Doubly Stochastic Poisson Input in the Light Traffic Situation. Trans. XXV Int. Sem. Stability Probl. Stoch Models. Maiori, Salerno, Italy, pp. 32–37 (2005)

4. Buchan, J., Koenigsberg, E.: Scientific Inventory Management. Prentice-Hall, Inc. Englewood Cliffs, New York (1963)

5. Bulinskaya, E.V.: Inventory Control in Case of Unknown Demand Distribution. Eng. Costs and Prod. Econom. 19, pp. 301–306 (1990)

6. Clarke, A.B.: A Waiting Time Process of Markov Type. Ann. Math. Stat. 27, 452–459 (1956)

7. Gupta, D., Weerawat, W.: Supplier-Manufacturer Coordination in Capacitated Two-Stage Supply Chains. Eur. J. Oper. Res. 175, pp. 67–89 (2006)

8. Harrison, J.M., Lemoine, A.: Limit Theorems for Periodic Queues. J. Appl. Probab. 14, 566–576 (1977)

9. Takacs, L.: Combinatorial Methods in the Theory of Stochastic Processes. Wiley, New York (1967)