Abstract

Translate

In this paper we present an exact numerical model for the evaluation of a three-echelon supply chain with multiple retailers. Poisson demand, exponentially distributed transportation times and lost sales at the retailers are assumed. The system is modeled as a continuous time Markov chain, and the analysis is based on matrix analytic methods. We analyze the infinitesimal generator matrix of the process and develop an algorithm for its construction. Performance measures for the system are calculated algorithmically from the stationary probabilities vector. The algorithm is used for an extensive numerical investigation of the system so that conclusions of managerial importance may be drawn.

Full text

Turn on search term navigation

Translate

1. Introduction

Multi-echelon inventory systems with multiple retailers have practical importance, but their modeling can become highly complex, sometimes even to the point of mathematical intractability. In general, to tackle complexity either simplifying assumptions must be made (deterministic parameters, no stock-outs, nested policies) or approximate methods have to be employed. Both approaches have their respective drawbacks, and there is an ongoing search for more realistic models that would allow us to assess longer and more complex supply chains.

In this paper we propose an algorithm for the exact numerical evaluation of performance measures for a three-echelon system with multiple retailers. Our modeling approach is based on the analysis of the infinitesimal generator matrix, and it allows us to use more realistic assumptions such as stochastic demand, stochastic transportation times and lost sales in case of stock-out. The resulting model is used for an extensive numerical research so that conclusions of managerial importance can be drawn.

2. Literature Review

A great part of the literature on divergent systems is concerned with two-echelon systems with one wholesaler and multiple retailers (OWMRs). An introduction to the evaluation of such systems is given in [1]. A systematic review of systems with periodic review inventory control policies can be found in [2]. An overview of some of the OWMR-related literature is given in Table 1.

Analyses of divergent systems with more than two echelons are less common. In [3] the authors analyze a multi-echelon system with general arborescent topology. They assume backordering, demand with Poisson characteristics and zero lead times. They define the system as a Markov process and characterize the optimal cost, applying dynamic programming that provides numerical approximations for global optimization.

Table 1

Related Literature.

Paper/Review Policy	Demand	Transportation Time	Shortages at the Retailer	Method
[4]/Mixed continuous–periodic review policy	Poisson	Deterministic	Backorders	Approximation method
[5]/Single cycle	Constant	Instantaneous	No shortages	Closed form
[6]/Stationary nested policies	Deterministic	Negligible	No shortages	Heuristic
[7]/Stationary nested policies	Constant	Instantaneous	No shortages	Search algorithm based on the characteristics of the optimal cost curve
[8]/Periodic	Random	Constant	Backorders	Closed form
[9]/Single cycle	Deterministic	Negligible	No shortages	Order interval division and recursive tightening
[10]/Integer-ratio policies	Deterministic	Negligible	No shortages	Iterative procedure based on the balance between thereplenishment and the inventory holding costs
[11]/Periodic–finite planning horizon	Stochastic (Poisson/normal)	Zero	Lost sales	Approximate dynamic programming/stochastic dynamic programming
[12]/Periodic	Correlated/bivariate Poisson	Constant (1 period)	Lost sales	Restricted observation Markov decision process
[13]/Periodic	Correlated with inventory and price	Constant	No shortages	Genetic algorithm/fuzzy simulation
[14]/Periodic base stock	Deterministic	Constant	Lost sales	Mathematical programming/genetic algorithms
[15]/Continuous review (R,Q) policies	Poisson	Constant	Backorders	Approximate optimal solutions
[16]/Finite horizon policy	Deterministic dynamic	Negligible	No shortages	Integer programming
[17]/Periodic (s,nq) for warehouse/discrete base stock policy for retailers	Stationary i.i.d.random	Deterministic	Backorders	Decomposition
[18]/Periodic, integer ratio, echelon stock, order up to policies	Poisson	Constant	Backorders	Exact solution based on the regenerative cycle
[19]/Periodic finite horizon policy	Deterministic dynamic	Deterministic	Backorders/lost sales	Decomposition
[20]/Cross-docking warehouse/(1,T) policy for retailers	Poisson	Constant	Lost sales	Heuristic
[21]/Periodic review base stock, finite horizon	Deterministic dynamic	Constant	Backorders	Mixed-integer linear programming
[22]/Continuous review (S−1,S) for the retailers, no replenishment policy for the warehouse	Stochastic	Stochastic	Backorders	Simulation
[23]/Continuous review	Deterministic	Normal	Backorders	Numerical iterative algorithm
[24]/VMI policy	Deterministic	Deterministic	No shortages	Integer non-linear programming/imperialist competitive algorithm
[25]/Periodic review base stock	Deterministic	Constant	Backorders	Mixed-integer linear programming/heuristic
[26]/Shipments at regular intervals	Deterministic	Deterministic	No shortages	Iterative procedure based on the concavity property of total profit
[27]/Periodic, finite horizon, option for delayed distribution	Stochastic	Negligible	Lost sales	Recursive solution algorithm based on multi-stage stochastic programming
[28]/Continuous review	Compound Poisson	Constant	Backorders	Exact analysis of the system based on the probability mass function (pmf) of each retailer’s inventory level
[29]/Periodic with an integer-ratio ordering schedule	Poisson	Constant	Backorders	Exact solution based on regenerative cycles
[30]/Periodic	Poisson/Gaussian	Constant	Lost sales/partial or complete backordering	Deep reinforcement learning
[31]/Periodic	Seasonal/stochastic	No lead time	Backorders	Deep reinforcement learning and multi-stage stochastic programming

In [32] the authors investigate a system with multiple stocking echelons and multiple retailers. They formulate a model for a mixed produce-to-order and produce-in-advance inventory system and seek to determine the optimal inventory at each installation on a single-period basis. They analyze the system for uniform and normal external demand with allowed transshipments between the retailers, and they conclude that in both cases the problem can be solved as a constrained optimization problem.

In [33] a model for a three-tier system with one supplier, one manufacturer and multiple retailers is developed. The retailers’ demands are random variables with generic probability density functions and lost sales are allowed. The model is based on a single cycle basis and a procedure for its solution is proposed. The authors focus on the comparison of traditional policies, where every node acts independently, to consignment policies where the retailers and the supplier are subordinates of the manufacturer.

The author in [34] develops approximation algorithms for optimal solutions for k-echelon divergent systems in a deterministic setting, with or without backorders. His approach is based on the decomposition technique.

In [35,36] an integrated three-level supply chain with a distribution structure, where one production plant serves multiple warehouses and each warehouse may serve many retailers, is examined. Both works assume dynamic and known demand over a discrete and finite planning horizon. In [35] the authors assume un-capacitated shipments, and they expand on the OWMR solutions. The authors compare different mixed-integer programming formulations for the optimal lot sizing, scheduling, transportation and warehousing decisions. In [36] limited production and transportation capacity, retailers that may change supplier over the planning horizon, as well as the possibility that the demand from a specific period may be satisfied by deliveries over multiple periods is assumed. The authors propose two heuristic algorithms for the optimal solution of the problem as well as an exact brunch and cut algorithm.

In [37] the authors employ mixed-integer programming techniques. They investigate three-level systems with multiple retailers, each one of which is supplied by a predefined warehouse, while the warehouses are supplied by a single production plant. The retailers face deterministic dynamic demand over a finite planning horizon, while it is assumed there is no restriction on the amount that can be produced or transported in a given period. The authors propose approaches for more effective and efficient optimal solutions with regard to the total cost incurred.

In [38] a multi-product integrated four-level supply chain consisting of a supplier, a producer, a wholesaler and multiple retailers is considered. Demand and lead times are deterministic, the review policy is periodic, and shortages occur at the retailers. The authors assume that the order quantity of products in each one of the levels has a normal distribution, while all the levels of the supply chain orders have the same number of products during the same period length. They formulate the problem as a non-linear programming model, and they apply two different approximation algorithms (sequential quadratic programming and interior point with super-linear convergence rates) in order to optimize the number and volume of the stockpiles.

The authors in [39] investigate both divergent and general structure systems that are centrally controlled. The authors assume periodic review policies with backorders and take into consideration both lead time and demand uncertainties. They model the problem as a Markov decision process and then they apply deep reinforcement learning and the proximal policy optimization algorithm for an approximate solution that minimizes the holding and backorder costs.

A dynamic supply chain member selection algorithm based on conditional generative adversarial networks (CGANs) is presented in [40]. For the analysis and the prediction of purchase and inventory links in the supply chain, machine learning is also used. They also examine the vehicle scheduling module where the path is reasonably planned to improve the operation efficiency. Finally, they use the SSH framework for the integrated implementation of the SCM system.

Inventory policies for Lindley systems with possibly unbounded costs, using as an objective the minimization of the expected discounted total cost by ordering (production) strategies are examined in [41]. The authors also show the existence of a subsequence of minimizers of the value iteration functions that converge to an optimal inventory system policy.

A systematic review of the existing state-of-the-art literature on machine learning (ML) in logistics and supply chain management (LSCM) is given in [42]. A wide collection of eight databases from 1994 to 2019 are explored. In total, 110 articles are analyzed showing that only nine literature reviews have been published in this area. The most important key findings show that 53.8% of publications were closely clustered on transportation and manufacturing industries and 54.7% were centered on mathematical models and simulations.

In this paper we propose an algorithm for the exact numerical evaluation of a three- echelon system consisting of a distribution center, a wholesaler and multiple retailers. We assume continuous review inventory control policies with each installation following an independent (s,Q) policy. Both demand and transportation times are assumed to be stochastic, while demand that cannot be met from inventory on hand at the retailers is assumed to be lost. The system is modeled as a continuous time discrete state Markov process, and the analysis is based on the infinitesimal generator matrix. In comparison with existing models, we investigate a longer network under more realistic stochastic conditions, and we offer an exact solution. Moreover, whereas most of the existing literature focuses on the optimal solution, our investigation is concerned with the general behavior of the system. In practice, optimal policies are not always easy to find or follow, and we consider it important to understand the effect of small changes in structural and operational parameters on the overall system performance, even when operating in sub-optimal conditions. The proposed algorithm could be used as the evaluative tool in the context of an optimization algorithm.

3. System Description

We investigate a pull system of three tiers with multiple retailers. A distribution center (DC) orders from a plant and supplies a wholesaler. In its turn, the wholesaler supplies n independent retailers (Figure 1). Each member holds an inventory and follows an independent continuous review inventory control policy (s,Q). The external demand has pure Poisson characteristics (exponentially distributed inter-arrival times and unitary demand per customer), while the external demand that cannot be met from the inventory on hand at each retailer is lost. In case of a stock-out at the wholesaler, the highest indexed retailer always has priority (retailer i has priority over retailer i − 1). Transportation times are exponentially distributed and independent of the replenishment order quantity, while both the DC and the wholesaler may send partial orders. The transportation processes are modeled as virtual independent stations. On replenishment order initiation, the respective inventory is subtracted from the upstream node and remains “in transit” until its delivery to the downstream node. The plant is assumed to be saturated, always sending complete orders of Q_d units to the DC. It is also assumed that at any given time at most one order can be in transit to any given node (one outstanding order assumption). Such an assumption is common in analytic models, and it is necessary in order to maintain a tractable level of complexity [43].

For our analysis we denote as the decision variables the number of retailers and the parameters of the inventory control policies at each node:

n: the number of retailers

s_d: the reorder point at the distribution center

Q_d: the replenishment order quantity at the distribution center

s_w: the reorder point at the wholesaler

Q_w: the replenishment order quantity at the wholesaler

s_i: the reorder point at retailer i

Q_i: the replenishment order quantity at retailer i

The other parameters that are necessary to completely define the system are as follows:

μ_d: the transportation rate for orders from the plant to the distribution center

μ_w: the transportation rate for orders from the distribution center to the wholesaler

μ_i: the transportation rate for orders from the wholesaler to retailer i

λ_i: the arrival rate of external customers at retailer i

4. Modeling Approach

Our modeling approach is based on the analysis of the infinitesimal generator matrix of the process. The system is modeled as a continuous time, discreet state Markov chain with a finite and multi-dimensional state space. We start by defining the state space of the process and by prescribing rules for ordering the permissible states. The resulting infinitesimal generator matrix G is partitioned into recurring blocks of states (sub-matrices), with different kinds of blocks corresponding to different kinds of transitions. In general, matrix G has a three-tier structure. Blocks on the diagonal describe transitions in the retailers for a given state of the upstream part of the system, with each basic level corresponding to a different level of inventory on hand at the distribution center. Blocks above the diagonal describe the arrival of a replenishment order at the distribution center (“birth” transitions), while blocks below the diagonal correspond to changes at the wholesaler. Transitions between non-adjacent basic levels are generally allowed.

Given the generator matrix, it is easy to construct a linear system of balance equations and compute a vector of stationary probabilities. For the results presented here, we have used LU factorization with partial pivoting from the Matlab© toolbox. System performance measures can be computed algorithmically from the stationary probabilities, taking advantage of the predefined ordering of states.

This section may be divided by subheadings. It should provide a concise and precise description of the experimental results, their interpretation, as well as the experimental conclusions that can be drawn.

4.1. States Definition

The system is modeled as a (2n + 3) dimensional continuous time Markov chain. At any given time t ≥ 0, the state of the system can be defined by a 2n + 3 dimensional vector $({I_{d}^{t}, T_{w,}^{t} I}_{w}^{t}, T_{n}^{t}, I_{n}^{t}, T_{n - 1}^{t}, I_{n - 1}^{t} {, \dots, T}_{1}^{t}, I_{1}^{t})$ , where

$I_{d}^{t}$ : the inventory on hand at the distribution center at time t with 0 ≤ $I_{d}^{t}$ ≤ s_d + Q_d

$T_{w}^{t}$ : the inventory in transit to the wholesaler at time t with 0 ≤ $T_{w}^{t}$ ≤ Q_w

$I_{w}^{t}$ : the inventory on hand at the wholesaler at time t with 0 ≤ $I_{w}^{t}$ ≤ s_w + Q_w

$T_{i}^{t}$ : the inventory in transit to retailer i at time t with 0 ≤ $T_{i}^{t}$ ≤ Q_i

$I_{i}^{t}$ : the inventory on hand at retailer i at time t with 0 ≤ $I_{i}^{t}$ ≤ s_i + Q_i

The Markov process has a finite state space and its dimension, which depends on the decision variables, can be calculated algorithmically. In general, for computational reasons, transient states are not taken into consideration, except where it simplifies the algorithmic analysis. The states are ordered using a lexicographical order. The subset of all states corresponding to the fixed inventory at the distribution center (I_d) is taken as the basic level, and the basic levels are ordered from lower to higher. Within each basic level, the states are ordered according to the inventory in transit to the wholesaler (T_w), then according to the inventory at the wholesaler (I_w) and finally according to T_i and I_i_, with higher priority retailers preceding lower priority ones and lower values preceding higher values.

The four kinds of events that may instantaneously alter the state of the system are as follows:

The arrival of an outstanding order from the Plant to the Distribution Centre $(I_{d}^{t + d t} = I_{d}^{t} + Q_{d})$

The arrival of an outstanding order from the DC to the wholesaler $(I_{w}^{t + d t} = I_{w}^{t} + T_{w}^{t})$

The arrival of an outstanding order at retailer i $(I_{i}^{t + d t} = I_{i}^{t} + T_{i}^{t})$

The occurrence of external demand at retailer i $(I_{i}^{t + d t} = I_{i}^{t} - 1)$

For example, for the simple system with two retailers and parameters s_d = 0, Q_d = 2, s_w = 0, Q_w = 2, s₁ = 2, Q₁ = 2, s₂ = 0, Q₂ = 1, the system is modeled as a seven-dimensional continuous time Markov chain ${I_{d}^{t}, T_{w}^{t}, I_{w}^{t}, T_{2}^{t}, I_{2}^{t}, T_{1}^{t}, I_{1}^{t}, t \geq 0}$ with 163 possible states. Assuming that we start at state (0,0,2,0,1,1,1), the possible transitions and associated events are

-. To state (0,0,0,0,1,2,2) with transition rate μ₁—Arrival of a replenishment order at retailer 1 and triggering of a new replenishment order from the wholesaler to retailer 1;

-. To state (0,0,1,1,0,1,1) with transition rate λ₂—Occurrence of external demand at retailer 2 and triggering of a replenishment order from the wholesaler to retailer 2.

-. To state (0,0,2,0,1,1,0) with transition rate λ₁—occurrence of external demand at retailer 1.

-. To state (2,0,2,0,1,1,1) with transition rate μ_d—Arrival of a replenishment order at the DC.

4.2. The Infinitesimal Generator Matrix

Diagonal blocks describe transitions where only the state of the retailers may change (inventory at the DC, inventory at the wholesaler and inventory in transit to the wholesaler do not change). The structure of these blocks can be defined recursively. We use as the “seed” the block describing transitions for the lowest priority retailer (retailer 1) and its associated transport station. The block for each successive retailer is constructed using as the “building block” the block for the previous retailer and according to rules that hold for all retailers. The specific structure for each such block depends on whether it corresponds to $I_{w}^{t} = 0$ , or $I_{w}^{t} > 0$ .

As an example we analyze the structure of the diagonal blocks when $I_{w}^{t} > 0$ . We define the following:

$b s d$ : the greatest common divisor of $(Q_{d}, Q_{w})$

$b s w$ : The greatest common divisor of $(b s d, Q_{1}, Q_{2}, \dots, Q_{n})$

$q_{i}$ : the integer part of $(\frac{Q_{i}}{b s w})$

Block C₁ for transitions where only the state of retailer 1 is changed ( $I_{1}^{t}$ or $T_{1}^{t}$ ) is a square matrix of $C l_{1} = Q_{1} + q_{1} \cdot (s_{1} + 1)$ dimension (Figure 2). C₁ can be further reduced into smaller blocks. For example, block $C_{1}^{z}$ is a $Q_{1} \times Q_{1}$ matrix for $T_{1}^{t} = 0$ , $I_{1}^{t} > s_{1}$ :

Block C₂ includes transitions where only the state of retailer 1 or 2 is changed (Figure 3). For every state of retailer 2, there are $C l_{1}$ possible states of retailer 1. C₂ is a square matrix with a dimension of $C l_{2} = Q_{2} \cdot C l_{1} + q_{2} \cdot (s_{2} + 1) \cdot C l_{1}$ . Correspondingly, $C_{2}^{z}$ for $T_{2}^{t} = 0$ , $I_{2}^{t} > s_{2}$ is a $Q_{2} \cdot C l_{1} \times Q_{2} \cdot C l_{1}$ matrix. If I₁ is the identity matrix of Cl₁ dimension

In general, block C_i (i > 2) includes transitions where only the state of retailers 1 to i may be changed (Figure 4). For every state of retailer i, there are $C l_{i - 1}$ possible states of the retailers with lower priority. C_i is a square matrix with dimension:

$C l_{i} = Q_{i} \cdot C l_{i - 1} + q_{i} (s_{i} + 1) \cdot C l_{i - 1}, where C l_{i - 1} = Q_{i - 1} \cdot C l_{i - 2} + q_{i - 1} (s_{i - 1} + 1) \cdot C l_{i - 2}$

Accordingly, $C_{i}^{z}$ for $T_{i}^{t} = 0$ , $I_{i}^{t} > s_{i}$ is a $Q_{i} \cdot C l_{i - 1} \times Q_{i} \cdot C l_{i - 1}$ matrix, where I_i−1 is the identity matrix of the Cl_i−1 dimension:

Above the diagonal there are diagonal blocks of μ_d describing the arrival of a replenishment order at the distribution center (DC). An outstanding replenishment order from the plant to the DC occurs as long as $I_{d}^{t} \leq s_{d}$ , and according to our assumptions, Q_d units are always delivered. When there is no outstanding demand from the wholesaler, the incoming order increases the available inventory on hand at the DC $(I_{d}^{t + d t} = I_{d}^{t} + Q_{d})$ . In the case of DC stock-out ( $I_{d}^{t} = 0$ , $T_{w}^{t} = 0$ , $I_{w}^{t} \leq s_{w}$ ), part or all of the incoming replenishment order is immediately forwarded to the wholesaler, and the μ_d blocks are moved appropriately to the left.

The blocks below the diagonal correspond either to the arrival of a replenishment order at the wholesaler or the triggering of a replenishment order to the retailers. In the first case, there are blocks of μ_w which occur as long as $I_{w}^{t} \leq s_{w}$ . A full or a partial order may be delivered, and the exact position of each μ_w element in the infinitesimal generator matrix depends on the values of $I_{d}^{t}$ and $I_{w}^{t}$ . When $I_{w}^{t} = 0$ and there are pending orders from the retailers, part or all of the incoming replenishment order is forwarded to the retailers, according to the priority of each retailer. When $I_{d}^{t} > 0,$ a new replenishment order will be triggered if $I_{w}^{t + d t} \leq s_{w}$ .

The triggering of a replenishment order to a retailer can occur when (A) external demand occurs at retailer i, while $I_{i}^{t} = s_{i} + 1$ and (B) a replenishment order arrives at the retailer i, and the updated $I_{i}^{t + d t}$ is less than or equal to s_i. In some cases a replenishment order from the DC to the wholesaler can also be triggered by these events. The corresponding blocks consist of λ_i and μ_i and are constructed recursively starting from the block for the retailer with the lowest priority (1), in a way similar to that presented for the diagonal blocks. The exact structure and positioning of the blocks depend on whether there is more than enough inventory to meet the retailer’s demand $(I_{w}^{t + d t} > 0)$ , or there is just enough or not enough inventory to meet the retailer’s demand $(I_{w}^{t + d t} = 0)$ . The position of the blocks also depends on whether there is inventory on hand at the DC ( $I_{d}^{t} > 0$ or $I_{d}^{t} = 0$ ) and on whether there is already a replenishment order in transit to the wholesaler ( $T_{w}^{t} > 0$ or $T_{w}^{t} = 0$ ).

The solution algorithm can be summarized in the following steps:

Step 1: Calculate infinitesimal generator matrix parameters from system parameters (number of basic levels, number of possible inventory levels at each stage, maximum inventory levels where outstanding orders occur)

Step 2: Create diagonal blocks

Step 2.1: Create diagonal block B₁ for lowest priority retailer 1 and $I_{w}^{t} = 0$

Step 2.2: Through an iterative process, create diagonal block B_n for highest priority retailer n and $I_{w}^{t} = 0$ based on block B_n−1

Step 2.3: Create diagonal block C₁ for lowest priority retailer 1 and $I_{w}^{t} > 0$

Step 2.4: Through an iterative process, create diagonal block C_n for highest priority retailer n and $I_{w}^{t} > 0$ based on block C_n−1

Step 2.5: Create the diagonal tier of the infinitesimal generator matrix using blocks B_n and C_n

Step 3: Insert upper diagonal blocks

Step 3.1: Insert blocks for the cases when $I_{d}^{t} > 0$

Step 3.2: Insert blocks for the cases when $I_{d}^{t} = 0$ (DC stockout)

Step 4: Insert below the diagonal blocks

Step 4.1: Add blocks corresponding to the arrival of a replenishment order at the wholesaler

Step 4.2: Add blocks corresponding to triggering of a replenishment order to a retailer

Step 5: Calculate numerically the vector of stationary probabilities (steady state solution)

using the following algorithm (see Calculation of the steady-state probabilities below)

Step 6: Calculate performance measures algorithmically from the vector of stationary probabilities

The calculation of the stationary probabilities (Step 5) of the previous algorithm can be achieved using the following algorithm:

Calculation of the steady-state probabilities

Let n be the number of states of the system, A be the infinitesimal generator matrix with dimension (n × n), X be the matrix of the unknown steady-state probabilities with dimension (n × 1), and B be the zero matrix with dimension (n × 1)

Step 1: Calculate matrix A − I (I is the identity matrix)

Step 2: Calculate ${(A - I)}^{T}$ (Transposition of A − I)

Step 3: Take ${(A - I)}^{T}$ and replace the last row with a new row where all the elements of this new row are equal to 1. Call this new matrix P

Step 4: Solve the system $P \cdot X = {[0 0 0 \dots 1]}^{T}$

The solution of this system provides the steady-state probabilities.

4.3. Performance Measures

Our analysis is based on the steady-state solution of the system and the numerical calculation of the vector of stationary probabilities $\vec{p}$ , where the i-th element of the vector $(p (i))$ corresponds to the i-th state in the hierarchy of states defined according to the lexicographical ordering. Performance measures about the system are computed algorithmically using the stationary probabilities. Some examples are given below.

Average inventory at the Distribution Centre—WIP_d

$W I P_{d} = \sum_{L e v e l = 1}^{N L d} (L e v e l \cdot b s d \cdot \sum_{j = 1}^{L_{1}} p (L_{0} + (L e v e l - 1) \cdot L_{1} + j))$

where

NLd: the integer part of $(\frac{s_{d} + Q_{d}}{b s d})$ . NLd is the number of the different positive values of the inventory at the distribution centre (basic levels);

L₀: the dimension of basic levels for I_d = 0 (the number of all states where I_d = 0);

L₁: the dimension of basic levels for I_d > 0 (the number of all states for any different value of I_d > 0).

Utilization of resource for transportation towards the wholesaler—u_w

Utilization of the resource for transportation to the wholesaler is the percentage of time that there is a replenishment order in transit to the wholesaler. To calculate u_w we sum the stationary probabilities of states corresponding to $T_{w}^{t} > 0$

$u_{w} = \sum_{i = i_{0}}^{L_{0}} p (i) + \sum_{i = 1}^{N L d} \sum_{j = j_{0}}^{L_{1}} p (L_{0} + (i - 1) \cdot L_{1} + j)$

$i_{0} = B l_{n} + N L w \cdot C l_{n} + 1$

$j_{0} = (N L w - n s w) \cdot C l_{n} + 1$

NLw: the integer part of $(\frac{s_{w} + Q_{w}}{b s w})$ . NLw is the number of levels of inventory on hand at the Wholesaler for $I_{w}^{t} > 0$

nsw: the integer part of $(\frac{s_{w}}{b s w})$ . nsw is the greatest $I_{w}^{t}$ level where the wholesaler asks for a replenishment order from the DC

Bl_n: the dimension of the block of states comprising all possible states for the retailers for a given state of the rest of the system and $I_{w}^{t} = 0$ . The block describes transitions where only the state of the retailers is changed and its dimension can be computed algorithmically:

$B l_{i} = (s_{i} + Q_{i} + 1) \cdot B l_{i - 1} + n Q_{i} \cdot (s_{i} + 1) \cdot B l_{i - 1}, i \geq 1$

$B l_{1} = (s_{1} + Q_{1} + 1) + n Q_{1} \cdot (s_{1} + 1)$

nQ_i: the integer part of $(\frac{Q_{i}}{b s w})$ . nQ_i is the number of the permissible values for inventory in transit to retailer i when $T_{i}^{t} > 0$ .

Average inventory at the wholesaler—WIP_w

WIP_w is the average inventory on hand at the wholesaler. We define a vector W such that the i-th element W(i) is the probability of $I_{w}^{t} = i$ . The value of each probability is computed algorithmically by summing up the appropriate steady-state probabilities. To facilitate the analysis, we break the process into four distinct cases depending on $I_{d}^{t}$ and $T_{w}^{t}$ values. Defining b = Bl_n + nsw × Cl_n and assuming $L e v e l$ , $L e v e l_{T}$ , and $L e v e l_{D}$ positive integers,

${I_{d}^{t} = 0, T}_{w}^{t} = 0$

$W (L e v e l) = \sum_{i = 1}^{C l_{n}} p (B l_{n} + (L e v e l - 1) \cdot C l_{n} + i)$

$1 \leq L e v e l \leq N L w$

${I_{d}^{t} = 0, T}_{w}^{t} > 0$

$W (L e v e l) = W (L e v e l) + \sum_{i = 1}^{C l_{n}} p (B l_{n} + N L w \cdot C l_{n} + (L e v e l_{T} - 1) \cdot b + B l_{n} + (L e v e l - 1) \cdot C l_{n} + i)$

$1 \leq L e v e l_{T} \leq n Q_{w}$

$1 \leq L e v e l \leq n s w$

nQ_w: the integer part of $(\frac{Q_{w}}{b s d})$ . nQ_w is the number of levels for inventory in transit towards the wholesaler when $T_{w}^{t} > 0$ .

${I_{d}^{t} > 0, T}_{w}^{t} = 0$

$W (L e v e l + n s w) = W (L e v e l + n s w) + \sum_{i = 1}^{C l_{n}} p (L_{0} + (L e v e l_{D} - 1) \cdot L_{1} + (L e v e l - 1) \cdot C l_{n} + i)$

$1 \leq L e v e l_{D} \leq N L d$

$1 \leq L e v e l \leq N L w - n s w$

${I_{d}^{t} > 0, T}_{w}^{t} > 0$

$\begin{array}{l} W (L e v e l) = W (L e v e l) \\ + \sum_{i = 1}^{C l_{n}} p (L_{0} + (L e v e l_{D} - 1) \cdot L_{1} + (N L w - n s w) \cdot C l_{n} + (L e v e l_{T} - 1) \cdot b + B l_{n} + (L e v e l - 1) \\ \cdot C l_{n} + i) \end{array}$

$1 \leq L e v e l_{D} \leq N L d$

$1 \leq L e v e l_{T} \leq n Q_{w}$

$1 \leq L e v e l \leq n s w$

Having constructed vector W, WIP_w can be easily calculated as the sum:

$W I P_{w} = \sum_{i = 1}^{N L w} (i \cdot b s w \cdot W (i))$

In order to calculate the retailers performance measures we define the following:

B-type blocks (B_i) such that each B_i block comprises the possible states for retailers 1 to i for a given state of the rest of the system and $I_{w}^{t} = 0$ .

C-type blocks (C_i) such that each C_i block comprises the possible states for retailers 1 to i for a given state of the rest of the system and $I_{w}^{t} > 0$ .

To facilitate the analysis, the performance measures concerning the retailers are computed by evaluating each B_n and C_n block separately. We denote L₀: the dimension of basic levels for I_d = 0; L₁: the dimension of basic levels for I_d > 0; and b = Bl_n + nsw × Cl_n. If lp + 1 the number of the first state of the B_n or C_n block under consideration and, $L e v e l_{T}$ , $L e v e l_{D}$ , and $L e v e l_{W}$ positive integers:

B_n blocks

For ${I_{d}^{t} = 0, T}_{w}^{t} = 0$ :

$l p = 0$

For ${I_{d}^{t} = 0, T}_{w}^{t} > 0$ :

$l p = B l_{n} + N L w \cdot C l_{n} + (L e v e l_{T} - 1) \cdot b$

$1 \leq L e v e l_{T} \leq n Q_{w}$

For ${I_{d}^{t} > 0, T}_{w}^{t} > 0$ :

$l p = L_{0} + (L e v e l_{D} - 1) \cdot L_{1} + (N L w - n s w) \cdot C l_{n} + (L e v e l_{T} - 1) \cdot b$

$1 \leq L e v e l_{D} \leq N L d$

$1 \leq L e v e l_{T} \leq n Q_{w}$

C_n blocks

For ${I_{d}^{t} = 0, T}_{w}^{t} = 0$ :

$l p = B l_{n} + (L e v e l_{W} - 1) \cdot C l_{n}$

$1 \leq L e v e l_{W} \leq N L w$

For ${I_{d}^{t} = 0, T}_{w}^{t} > 0$ :

$l p = B l_{n} + N L w \cdot C l_{n} + (L e v e l_{T} - 1) \cdot b + B l_{n} + (L e v e l_{W} - 1) \cdot C l_{n}$

$1 \leq L e v e l_{T} \leq n Q_{w}$

$1 \leq L e v e l_{W} \leq n s w$

For ${I_{d}^{t} > 0, T}_{w}^{t} = 0$ :

$l p = L_{0} + (L e v e l_{D} - 1) \cdot L_{1} + (L e v e l - 1) \cdot C l_{n}$

$1 \leq L e v e l_{D} \leq N L d$

$1 \leq L e v e l_{W} \leq N L w - n s w$

For ${I_{d}^{t} > 0, T}_{w}^{t} > 0$ :

$\begin{matrix} l p = L_{0} + (L e v e l_{D} - 1) \cdot L_{1} + (N L w - n s w) \cdot C l_{n} + (L e v e l_{T} - 1) \cdot b + B l_{n} \\ + (L e v e l_{w} - 1) \cdot C l_{n} \end{matrix}$

$1 \leq L e v e l_{T} \leq n Q_{w}$

$1 \leq L e v e l_{W} \leq n s w$

Average inventory at the retailers—WIP_r

WIP_r is the average inventory on hand at retailer r, 1 ≤ r ≤ n.

B-type blocks ( $I_{w}^{t} = 0$ )

Retailer 1

In B₁ blocks, for the states where $T_{1}^{t} = 0$ , positive inventory at retailer 1 corresponds to s₁ + Q₁ states. For $T_{1}^{t} > 0$ , there are nQ₁ different levels of T₁, while in each level s₁ state, it corresponds to $I_{1}^{t} > 0$ . If $b_{1} = B l_{n} / B l_{1}$ the number of B₁ blocks in B_n and lp + 1, the number of the first state of the B_n block under consideration is as follows:

$\begin{array}{l} W I P_{1}^{B, l p} = \sum_{i = 1}^{b_{1}} (\sum_{j = 1}^{s_{1} + Q_{1}} j \cdot p (l p + (i - 1) \cdot B l_{1} + j + 1) \\ + \sum_{j = 1}^{n Q_{1}} \sum_{z = 1}^{s_{1}} z \\ \cdot p (l p + (i - 1) \cdot B l_{1} + s_{1} + Q_{1} + 1 + (j - 1) \cdot (s_{1} + 1) + z \\ + 1)) \end{array}$

Retailer r, 2 ≤ r ≤ n

With the same iterative approach, we can calculate the average inventory on hand for higher priority retailers. In B_r blocks, for each state of retailer r correspond Bl_r−1 states of the lower priority retailers. In B-type blocks, for the states where $T_{r}^{t} = 0$ , inventory at retailer r can take s_r + Q_r values. For $T_{r}^{t} > 0$ , there are nQ_r different levels of $T_{r}^{t}$ , while in each level s_r values correspond to $I_{r}^{t} > 0$ . If $b_{r} = B l_{n} / B l_{r}$ and lp + 1, the number of the first state of the B_n block under consideration is as follows:

$W I P_{r}^{B 1} = \sum_{i = 1}^{b_{r}} (\sum_{j = 1}^{s_{r} + Q_{r}} \sum_{z = 1}^{B l_{r - 1}} j \cdot p (l p + (i - 1) \cdot B l_{r} + B l_{r - 1} + (j - 1) \cdot B l_{r - 1} + z))$

$\begin{array}{l} W I P_{r}^{B 2} = \sum_{i = 1}^{b_{r}} \sum_{j = 1}^{n Q_{r}} \sum_{z = 1}^{s_{r}} \sum_{y = 1}^{B l_{r - 1}} z \\ \cdot p (l p + (i - 1) \cdot B l_{r} + (s_{r} + Q_{r} + 1) \cdot B l_{r - 1} + (j - 1) \cdot (s_{r} \\ + 1) \cdot B l_{r - 1} + B l_{r - 1} + (z - 1) \cdot B l_{r - 1} + y) \end{array}$

$W I P_{r}^{B, l p} = W I P_{r}^{B 1} + W I P_{r}^{B 2}$

C-type blocks ( $I_{w}^{t} > 0$ )

Retailer 1

When $T_{1}^{t} = 0$ , there are Q₁ different states where $I_{1}^{t} > 0$ . For $T_{1}^{t} > 0$ , there are nQ₁ different levels of $T_{1}^{t}$ , and in each level s₁ state, it corresponds to $I_{1}^{t} > 0$ . If $c_{1} = C l_{n} / C l_{1}$ the number of C₁ blocks in C_n and lp + 1, the number of the first state of the C_n block under consideration is as follows:

$W I P_{1}^{C 1} = \sum_{j = 1}^{c_{1}} \sum_{z = 1}^{Q_{1}} (z + s_{1}) \cdot p (l p + (j - 1) \cdot C l_{1} + z$

$W I P_{1}^{C 2} = \sum_{j = 1}^{c_{1}} \sum_{z = 1}^{n Q_{1}} \sum_{y = 1}^{s_{1}} y \cdot p (l p + (j - 1) \cdot C l_{1} + Q_{1} + (z - 1) \cdot (s_{1} + 1) + y + 1)$

$W I P_{1}^{C, l p} = W I P_{1}^{C 1} + W I P_{1}^{C 2}$

Retailer r, 2 ≤ r ≤ n

In C_r blocks each retailer r state corresponds to Cl_r−1 states of the lower priority retailers. When $T_{r}^{t} = 0$ , there are Q_r different values for $I_{r}^{t} > 0$ . For $T_{r}^{t} > 0$ , there are nQ_r different levels of $T_{r}^{t}$ , while in each level correspond s_r different levels of positive $I_{r}^{t}$ . If $c_{r} = \frac{C l_{n}}{C l_{r}},$ the number of C_r blocks in C_n is as follows:

$W I P_{r}^{C 1} = \sum_{j = 1}^{c_{r}} \sum_{z = 1}^{Q_{r}} \sum_{y = 1}^{C l_{r - 1}} (s_{r} + z) \cdot p (l p + (j - 1) \cdot C l_{r} + (z - 1) \cdot C l_{r - 1} + y)$

$\begin{array}{l} W I P_{r}^{C 2} = \sum_{j = 1}^{c_{r}} \sum_{z = 1}^{n Q_{r}} \sum_{y = 1}^{s_{r}} \sum_{x = 1}^{C l_{r - 1}} y \\ \cdot p (l p + (j - 1) \cdot C l_{r} + Q_{r} \cdot C l_{r - 1} + (z - 1) \cdot (s_{r} + 1) \cdot C l_{r - 1} \\ + (y - 1) \cdot C l_{r - 1} + C l_{r - 1} + x) \end{array}$

$W I P_{r}^{C, l p} = W I P_{r}^{C 1} + W I P_{r}^{C 2}$

The average inventory at retailer r (1 ≤ r ≤ n) will be the sum $W I P_{r}^{B, l p}$ for B_n blocks and $W I P_{r}^{C, l p}$ for C_n blocks for all possible lp.

$W I P_{r} = \sum_{l p, I_{w}^{t} = 0} W I P_{r}^{B, l p} + \sum_{l p, I_{w}^{t} > 0} W I P_{r}^{C, l p}$

Stock-out probability for Retailer r—SO_r

SO_r is the probability that the external demand at retailer r will become lost sales. Since external demand at the retailers is independent and uniformly distributed in time, SO_r will be same as the probability of inventory on hand at retailer r being zero.

B-type Blocks ( $I_{w}^{t} = 0$ )

Retailer 1

$I_{1}^{t}$ is zero in the first state of each B₁ block of states, where also $T_{1}^{t} = 0$ . For $T_{1}^{t} > 0$ , $I_{1}^{t}$ is zero in the first state of each (s₁ + 1)-dimension sub-block corresponding to a different $T_{1}^{t}$ value. If $b_{1} = B l_{n} / B l_{1},$

$\begin{array}{l} S O_{1}^{B, l p} = \sum_{i = 1}^{b_{1}} (p (l p + (i - 1) \cdot B l_{1} + 1) \\ + \sum_{j = 1}^{n Q_{1}} p (l p + (i - 1) \cdot B l_{1} + s_{1} + Q_{1} + 1 + (j - 1) \cdot (s_{1} + 1) \\ + 1)) \end{array}$

Retailer r, 2 ≤ r ≤ n

In B_r blocks, each retailer r state (r > 1) corresponds to Bl_r−1 states of the lower priority retailers. $I_{r}^{t}$ is zero in the first Bl_r−1 states of each B_r block of states, where $T_{r}^{t}$ is also zero. For $T_{r}^{t} > 0$ , $I_{r}^{t} = 0$ in the first Bl_r−1 states of each (s_r + 1)∙Bl_r−1 dimension sub-block corresponding to a particular $T_{r}^{t}$ value. If $b_{r} = B l_{n} / B l_{r}$ ,

$\begin{array}{l} S O_{r}^{B, l p} = \sum_{i = 1}^{b_{r}} (\sum_{j = 1}^{B l_{r - 1}} (p (l p + (i - 1) \cdot B l_{r} + j)) \\ + \sum_{j = 1}^{n Q_{r}} \sum_{z = 1}^{B l_{r - 1}} p (l p + (i - 1) \cdot B l_{r} + (s_{r} + Q_{r} + 1) \cdot B l_{r - 1} + (j \\ - 1) \cdot (s_{r} + 1) \cdot B l_{r - 1} + z)) \end{array}$

C-type blocks ( $I_{w}^{t} > 0$ )

Retailer 1

In C₁ blocks $I_{1}^{t} = 0$ in the first state of each of (s₁ + 1)-dimension sub-block corresponding to a different $T_{1}^{t}$ value and only when $T_{1}^{t} > 0$ . If $c_{1} = C l_{n} / C l_{1},$

$S O_{1}^{C, l p} = \sum_{j = 1}^{c_{1}} \sum_{z = 1}^{n Q_{1}} p (l p + (j - 1) \cdot C l_{1} + Q_{1} + (z - 1) \cdot (s_{1} + 1) + 1)$

Retailer r, 2 ≤ r ≤ n

In C_r blocks each retailer r state corresponds to Cl_r−1 states of the lower priority retailers. $I_{r}^{t} = 0$ in the first state of each (s_r + 1)∙Bl_r−1 dimension sub-block corresponding to a particular $T_{r}^{t}$ value and only when $T_{r}^{t} > 0$ . If $c_{r} = C l_{n} / C l_{r},$

$\begin{matrix} S O_{r}^{C, l p} = \sum_{j = 1}^{c_{r}} \sum_{z = 1}^{n Q_{r}} \sum_{y = 1}^{C l_{r - 1}} p (l p + (j - 1) \cdot C l_{r} + Q_{r} \cdot C l_{r - 1} + (z - 1) \cdot (s_{r} + 1) \\ \cdot C l_{r - 1} + y) \end{matrix}$

Stock-out probability for retailer r (1 ≤ r ≤ n) will be the sum $S O_{r}^{B, l p}$ for B_n blocks and $S O_{r}^{C, l p}$ for C_n blocks for all possible lp.

$S O_{r} = \sum_{l p, I_{w}^{t} = 0} S O_{r}^{B, l p} + \sum_{l p, I_{w}^{t} > 0} S O_{r}^{C, l p}$

Fill Rate of retailer r—FR_r

Fill rate is the percentage of external customers arriving at retailer r whose demand is met by inventory on hand at the retailer:

$F R_{r} = 1 - S O_{r}$

Throughput of retailer r—Thr_r

Throughput is the number of product units per time unit that flow through retailer r. Alternatively, Thr_r could be defined as the rate of sales at retailer r:

$T h r_{r} = λ_{r} \cdot F R_{r}$

4.4. Validation and Model Performance

The validity of the developed algorithm was tested with simulation. A total of 900 different scenarios were tested for supply chains with one to five retailers and for various parameter relations. In all cases and for all the tested performance measures, the difference between analytic results and simulation results was within the expected deviation attributed to the experimental nature of the simulation approach. Some examples are given in Figure 5.

With regard to performance, our model shares the common problem of Markovian models, namely the increasing number of states as the system under consideration becomes bigger or more complex. As a general trend, the dimension of the infinitesimal generator matrix increases with an increasing number of retailers and increasing values for the inventory policy parameters. Regarding the computational complexity limitations, theoretically the model has no structural (mathematical) limitations. However, in practical terms, the computational feasibility is constrained by hardware resources like RAM, CPU, etc. Larger systems remain solvable but at the cost of significantly increased runtime and resource demands.

Despite size limitations, the proposed algorithm still offers certain advantages. The exact algorithm is significantly faster than simulation, in some cases the difference in computation time being several orders of magnitude. Moreover, the exact solution poses no limits on precision in contrast to simulation or approximation methods. Practitioners should consider using our model over alternatives when

The impact of decision variables (quantity order, reorder point, etc.) on throughput and WIP must be explicitly modeled.

The system properties, behavior and its sensitivity can be identified because the model is analysed explicitly. This does not hold when approximation methods are used.

Validation against simulation is required for credibility and decision support.

System design requires precise performance estimates, validated via both analytical and simulation methods.

Decision-makers seek to evaluate various maintenance policies under uncertainty.

A hybrid approach can also be effective: use the exact Markovian solution for small-scale cases to calibrate, validate or benchmark approximate/simulation models for larger networks.

Finally, the proposed algorithm can be easily integrated with other components in the framework of a more generic model as, for example, in the context of an optimization algorithm.

5. Numerical Results

Two different scenarios have been examined in order to check the behavior of the system and the proposed algorithm under different conditions.

The first scenario examines the case where the system is balanced (i.e., the total supply rate is equal to the total demand rate). The second scenario examines the case where the system is supply constrained. In this case the total supply rate is less than the total demand rate (i.e., demand exceeds supply). For all cases examined first the transition matrix is generated, second the steady-state probabilities are calculated, and finally the performance measures are computed according to the procedures described in the previous sections.

5.1. The Effect of Policy Parameters—Balanced Systems

We investigate “balanced” systems where upstream and downstream transportation rates are balanced. The following numerical results refer to a system with two retailers (n = 2), where $μ_{1} = μ_{2} = λ_{1} = λ_{2} = \frac{μ_{w}}{2} = \frac{μ_{d}}{2}$ .

5.1.1. Distribution Center’s Policy (s_d, Q_d) for Balanced Systems

The performance of the retailers increases with increasing s_d but at the cost of an increase in the average total inventory in the system (average total inventory—WIP_total is the average inventory from the DC and downstream). Sometimes this increase is almost linear (Figure 6). The effect of Q_d depends on the value of the other variables, but in general both the fill rate at the retailers and total inventory tend to increase with increasing Q_d. In some cases jagged patterns are observed (Figure 7). The effects of the parameters on both retailers are similar.

DC policy is important for the retailers’ performance only for low s_d and Q_d values when the distribution center is the “bottleneck” of the system (Figure 8). The presence of some safety stock at the DC seems preferable, as for s_d = 0 the system is less stable and more difficult to predict. In some cases, with minor policy adjustments it is possible to achieve lower total inventory without any serious negative effect on the retailers’ performance.

5.1.2. Wholesaler Policy (s_w, Q_w) for Balanced Systems

By increasing s_w we can achieve better service levels at the retailers but again at the cost of an increase in the average total inventory (Figure 9). The effect of Q_w is less straightforward. In general, fill rates at the retailers and the average total inventory tend to increase with increasing Q_w, but jagged patterns may be observed (Figure 10). There is strong interplay between the parameters, and for certain scenarios it is actually possible to enhance retailers’ performance, while actually decreasing the average inventory in the system (Figure 11). The presence of some safety stock at the wholesaler (s_w > 0) can protect the retailers from significant deviations in fill rates (Figure 12). The effect of the wholesaler’s policy on both retailers is similar.

5.1.3. Retailers’ Policy (s_i, Q_i) for Balanced Systems

Increasing s_i causes an increase in the fill rate at retailer i but also an increase in the average inventory at retailer i. Although both average inventory at the DC and average inventory at the wholesaler is negatively correlated with s_i, the overall effect on WIP_total is a clear increase. Up to a point, increasing the reorder point at retailer i causes a transfer of available inventory downstream. This has a negative effect on the performance of the other retailers (Figure 13).

Similarly, an increase in Q_i increases both the fill rate and average inventory at retailer i. The effect on the wholesaler’s average inventory is more dynamic. In general, WIP wholesaler tends to decrease with an increasing Q_i, but saw-like patterns can also be observed (Figure 14). Average inventory at the distribution center decreases with Q_i. The average total inventory in the system is generally positively correlated with Q_i, but due to the WIP_wholesaler contribution, a jagged pattern may be observed. As was the case with s_i, the increase in Q_i causes a downstream transfer of available inventory, and this has a negative effect on the performance of the other retailers. The managerial implication of the dynamic system behavior is that there are good reasons for the fine-tuning of the system. Small changes in inventory policies may achieve an enhanced performance in terms of both customer satisfaction and total system inventory.

5.2. The Effect of Policy Parameters—Supply-Constrained Systems

We investigate supply-constrained systems where $λ_{1} = λ_{2} > μ_{1} = μ_{2} > μ_{w} > μ_{d}$ . For the numerical examples that follow, the parameters were λ₁ = λ₂ = 3, μ₁ = μ₂ = 2, μ_w = 1, and μ_d = 0.8.

5.2.1. Distribution Center’s Policy (s_d, Q_d) for Supply Constrained Systems

Compared to balanced systems, in supply-constrained systems the retailers are more sensitive to DC policy changes. In general the reorder quantity Q_d has a greater effect on the performance measures than s_d (Figure 15). Jagged patterns may also be observed, so attention should be given to fine-tuning the system.

5.2.2. Wholesaler’s Policy (s_w, Q_w) for Supply Constrained Systems

Again, supply-constrained systems are more sensitive to design variable changes. In general the effect of reorder quantity Q_w is more important than that of s_w, while jagged patterns in the performance measures may also occur. The wholesaler’s policy can be an effective way to enhance retailers’ performance. Higher s_w and Q_w values lead to higher availability of products at the wholesaler and a better service for the retailers. This causes an increase in average inventories not only at the wholesaler but also downstream at the retailers (Figure 16).

5.2.3. Retailers’ Policy (s_i, Q_i) for Supply Constrained Systems

As with balanced systems, increasing s_i or Q_i causes a transfer of available inventory downstream, and the fill rate at retailer i increases. In general, total inventory increases, but in some scenarios a small decrease in WIP_total may be initially observed (Figure 17). Compared to balanced systems, changing the policy of one retailer has a greater impact on the performance of the other retailers.

5.3. The Effect of Number of Retailers

Each retailer acts independently from the others, and no demand correlation is assumed. However, since all retailers are supplied by a finite-capacity wholesaler, some interaction amongst them occurs. We study a system with three retailers, and we focus on how the inventory policy of one retailer affects the performance of the others. We focus on the retailer with the highest (retailer 3) and the lowest (retailer 1) priority. For the following examples we assume for the balanced systems μ₁ = μ₂ = μ₃ = λ₁ = λ₂ = λ₃ = μ_w/3 = μ_d/3, and for the supply-constrained systems λ₁ = λ₂ = λ₃ = 3, μ₁ = μ₂ = μ₃ = 2, μ_w = 1, and μ_d = 0.8.

5.3.1. Interactions Between the Retailers

Changing the reorder point s_i causes small (sometimes insignificant) but consistent decreases in the performance of the other retailers (Figure 18 and Figure 19).

Parameter Q directly affects the availability of inventory at the wholesaler, and thus the ability of the wholesaler to respond to demand from the other retailers, so greater effects are observed (Figure 20 and Figure 21). For balanced systems, in some scenarios an increase in Q in one retailer caused small increases in fill rates in two or more retailers (Figure 20). Such cases indicate the coordination between the inventory policies of the various network members and are of interest from a managerial point of view.

The lowest priority retailers were found to be more sensitive to changes in the other retailers’ policies.

5.3.2. The Effect of Retailer Addition

We investigate the behavior of the system with an increasing number of retailers supplied by the same wholesaler. All retailers are assumed to follow the same inventory control policy, while replenishment times are also the same for all retailers.

Increasing the number of retailers causes the fill rates to decrease as the available inventory at the wholesaler decreases. Initially, the effect is almost the same for all retailers, irrespective of their designated priority. Beyond a point the antagonism between the retailers becomes more intense as the stock-out probability for the wholesaler becomes important and priorities start to have an effect on retailers’ performance (Figure 22).

Considering total system performance, by increasing the number of retailers we increase the total output, but at the same time total lost sales also increase. The ratio of lost sales to output increases with n. The changes in the ratio depend on the other parameters of the system as well, and they are more important for higher retailer transportation rates (Figure 23). The optimal number of retailers will depend on the specific costs parameters.

6. Conclusions

In this paper we presented a model based on continuous time discrete space Markov processes for the exact numerical analysis of a single-product, three-echelon inventory system with multiple retailers. A solution algorithm based on the properties of the infinitesimal generator matrix was developed, and the model was used to investigate the effect of the decision variables on the performance measures of the system. The algorithm was validated by comparing the results with simulation, and for all cases there was a very good agreement.

Conclusions were drawn based on extensive numerical research on different scenarios for balanced and supply-constrained systems. The results indicate a dynamic behavior. There is interdependence between the different members of the network and interplay between system parameters. More importantly from a managerial point of view, our analysis also indicates that in some cases it is possible to coordinate the system, increasing customer satisfaction (fill rates), while at the same time decreasing average total inventory in the network. For balanced systems some safety stock was found to be desirable at both the DC and the wholesaler, while the possibility of coordination of the supply network members could be explored by fine-tuning either the wholesaler’s or the retailers’ policies. In regard to supply-constrained systems, DC and wholesaler policies were found to be more important for overall system performance, while the effect of each retailer’s policy on the performance of the others was also stronger.

The novelty of our approach in comparison to previous three-echelon models is the following:

Our model provides an exact analytical solution for the examined system in contrast with other three-echelon models where approximation methods are used.

The exact Markov chain solution in this paper can reveal important properties and the general behavior of the system that cannot be identified when approximation techniques are used.

Beyond its theoretical contribution, the model serves as a benchmarking tool for simulation and heuristic methods, thereby providing both methodological rigor and actionable managerial guidance.

As further research, four directions are proposed. Regarding external demand, more general demand distributions can be investigated. Although the Poisson distribution is commonly used to model demand in inventory systems, it may not be appropriate to describe more complex and more realistic cases. A compound Poisson distribution combining Poisson arrivals for customers with an empirical distribution for individual demand would offer more modeling flexibility, and it would be relatively easy to integrate in the presented model. In a second direction, more general distributions can be used for replenishment times. The application of phase-type distributions (Erlang, Coxian) to model times would allow for more realistic modeling, but it must be kept in mind that the total number of states, and so the computing requirements, would increase. In a third direction, such a model may be used as a decomposition block for the analysis of larger systems using approximate methods like decomposition, aggregation techniques, etc. Finally, the development of a cost function would allow the model to be used for optimization purposes. This could be achieved either through an exhaustive enumeration of possible policies or by combining the evaluative algorithm with an optimization heuristic.

Author Contributions

Conceptualization, G.V., S.K., A.D. and E.I.; methodology, G.V., S.K., A.D. and E.I.; software, G.V., S.K., A.D. and E.I.; validation, G.V., S.K., A.D. and E.I.; formal analysis, G.V., S.K., A.D. and E.I.; investigation, G.V., S.K., A.D. and E.I.; resources, G.V., S.K., A.D. and E.I.; data curation, G.V., S.K., A.D. and E.I.; writing—original draft preparation, G.V., S.K., A.D. and E.I. All authors have read and agreed to the published version of the manuscript.

Data Availability Statement

The original contributions presented in this study are included in the article. Further inquiries can be directed to the corresponding author.

Conflicts of Interest

The authors declare no conflict of interest.

Footnotes

Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Figures and Table

Figure 1 System layout.

Figure 2 The submatrices of block C₁.

Figure 3 The submatrices of block C₂.

Figure 4 The submatrices of block C_i.

Figure 5 % deviation (100 × (analytic-simulation)/analytic) between algorithmic solution and simulation results for a system with five retailers, s_d = {2,4}, Q_d = 2, s_w = 0, Q_w = 2, s₁ = 0, Q₁ = 1, s₂ = 1, Q₂ = 1, s₃ = 0, Q₃ = {1,2}, s₄ = 0, Q₄ = {1,2}, 0 ≤ s₅ ≤ 3, 1 ≤ Q₅ ≤ 2 and parameters μ_d = 2.5, μ_w = 3.6, μ₁ = 1, μ₂ = 1.2, μ₃ = 1.4, μ₄ = 1.6, μ₅ = 1.8, λ₁ = 0.5, λ₂ = 0.7, λ₃ = 0.9, λ₄ = 1.2, λ₅ = 1.5. Simulation parameters: one replication of 2,000,000 time units with a warm up period of 10,000 time units.

Figure 6 The effect of s_d on system performance measures; Q_d = 2, s_w = 2, Q_w = 1, s₁ = 1, Q₁ = 2, s₂ = 1, Q₂ = 1.

Figure 7 The effect of Q_d on system performance measures: s_d = 2, s_w = 2, Q_w = 2, s₁ = 1, Q₁ = 2, s₂ = 1, Q₂ = 2.

Figure 8 Compound effect of DC policy for Balanced Systems: s_w = 2, Q_w = 2, s₁ = 2, Q₁ = 3, s₂ = 2, Q₂ = 1.

Figure 9 The effect of s_w on system performance measures: s_d = 2, Q_d = 2, Q_w = 1, s₁ = 1, Q₁ = 2, s₂ = 1, Q₂ = 1.

Figure 10 The effect of Q_w on system performance measures: s_d = 4, Q_d = 4, s_w = 2, s₁ = 1, Q₁ = 2, s₂ = 1, Q₂ = 2.

Figure 11 Compound effect of wholesaler’s policy: s_d = 2, Q_d = 6, s₁ = 1, Q₁ = 2, s₂ = 1, Q₂ = 1.

Figure 12 Compound effect of wholesaler’s policy: s_d = 3, Q_d = 5, s₁ = 1, Q₁ = 2, s₂ = 1, Q₂ = 2.

Figure 13 The effect of s₁ on system performance measures: s_d = 4, Q_d = 2, s_w = 2, Q_w = 4, Q₁ = 2, s₂ = 1, Q₂ = 2.

Figure 14 The effect of Q₁ on system performance measures: s_d = 2, Q_d = 2, s_w = 6, Q_w = 2, s₁ = 2, s₂ = 1, Q₂ = 2.

Figure 15 Compound effect of DC policy for Supply Constrained Systems: s_w = 2, Q_w = 2, s₁ = 2, Q₁ = 3, s₂ = 2, Q₂ = 1.

Figure 16 Compound effect of wholesaler’s policy: s_d = 2, Q_d = 6, s₁ = 1, Q₁ = 2, s₂ = 1, Q₂ = 2.

Figure 17 Compound effect of retailer 1’s policy: s_d = 2, Q_d = 6, s_w = 2, Q_w = 4, s₂ = 1, Q₂ = 2.

Figure 18 Balanced system. Effect of s, sd = 0, Qd = 4, sw = 2, Qw = 4, si = 1, Qi = 2, μ₁ = μ₂ = μ₃ = λ₁ = λ₂ = λ₃ = μ_w/3 = μ_d/3.

Figure 19 Supply-constrained system. Effect of s, s_d = 0, Q_d = 4, s_w = 2, Q_w = 4, Q₁ = 2, s_i = 1, Q_i = 2λ₁ = λ₂ = λ₃ > μ₁ = μ₂ = μ₃ > μ_w > μ_d.

Figure 20 Balanced system. Effect of Q, s_d = 0, Q_d = 4, s_w = 2, Q_w = 4, s_i = 1, Q_i = 2, μ₁ = μ₂ = μ₃ = λ₁ = λ₂ = λ₃ = μ_w/3 = μ_d/3.

Figure 21 Supply-constrained system. Effect of lowest priority retailer—Q, s_d = 0, Q_d = 4, s_w = 2, Q_w = 4, s_i = 1, Q_i = 2, λ₁ = λ₂ = λ₃ > μ₁ = μ₂ = μ₃ > μ_w > μ_d.

Figure 22 Order fill rate of the highest and lowest priority retailers as a function of the number of retailers: s_d = 0, Q_d = 2, s_w = 0, Q_w = 2, s_i = 0, Q_i = 1, i = [1,7]. λ = 1, μ_d = μ_w = 4.

Figure 23 Ratio of lost sales to output as a function of the number of retailers.

References

1. Axsäter, S. Inventory Control; International Series in Operations Research & Management Science 225; Springer International Publishing: Cham, Switzerland, 2015; pp. 191-222.

2. Agrawal, N.; Smith, S.A. Multi-location Inventory Models for Retail Supply Chain Management—A Review of Recent Research. Retail Supply Chain Management; Agrawal, N.; Smith, S.A. International Series in Operations Research & Management Science 223; Springer: Berlin/Heidelberg, Germany, 2009; pp. 319-347. Available online: https://link.springer.com/book/10.1007/978-0-387-78902-6 (accessed on 22 September 2025).

3. Rofman, E.; González, R.; Sagastizábal, C. Global Optimization of Arborescent Multilevel Inventory Systems. J. Glob. Optim.; 1995; 6, pp. 269-292. [DOI: https://dx.doi.org/10.1007/bf01099465]

4. Ahire, S.L.; Schmidt, C.P. A model for a Mixed Continuous-Periodic Review One-Warehouse, N-Retailer inventory system. Eur. J. Oper. Res.; 1996; 92, pp. 69-82. [DOI: https://dx.doi.org/10.1016/0377-2217(94)00309-2]

5. Yang, P.; Wee, H. An arborescent inventory model in a supply chain system. Prod. Plan. Control; 2001; 12, pp. 728-735. [DOI: https://dx.doi.org/10.1080/09537280010024063]

6. Abdul-Jalbar, B.; Gutiérrez, J.M.; Sicilia, J. Single cycle policies for the one-warehouse N-retailer inventory/distribution system. Omega; 2006; 34, pp. 196-208. [DOI: https://dx.doi.org/10.1016/j.omega.2004.10.003]

7. Yao, M.J.; Wang, Y. A new algorithm for one-warehouse multi-retailer systems under stationary nested policy. Optim. Methods Softw.; 2006; 21, pp. 41-56.

8. Lagodimos, A.; Koukoumialos, S. Service performance of two-echelon supply chains under linear rationing. Int. J. Prod. Econ.; 2008; 112, pp. 869-884. [DOI: https://dx.doi.org/10.1016/j.ijpe.2007.07.007]

9. Hsiao, Y.-C. Optimal single-cycle policies for the one-warehouse multi-retailer inventory/distribution system. Int. J. Prod. Econ.; 2008; 114, pp. 219-229. [DOI: https://dx.doi.org/10.1016/j.ijpe.2008.01.008]

10. Abdul-Jalbar, B.; Segerstedt, A.; Sicilia, J.; Nilsson, A. A new heuristic to solve the one-warehouse N-retailer problem. Comput. Oper. Res.; 2010; 37, pp. 265-272. [DOI: https://dx.doi.org/10.1016/j.cor.2009.04.012]

11. Geng, W.; Qiu, M.; Zhao, X. An inventory system with single distributor and multiple retailers: Operating scenarios and performance comparison. Int. J. Prod. Econ.; 2010; 128, pp. 434-444. [DOI: https://dx.doi.org/10.1016/j.ijpe.2010.08.002]

12. Helper, C.M.; Davis, L.B.; Wei, W. Impact of demand correlation and information sharing in a capacity constrained supply chain with multiple-retailers. Comput. Ind. Eng.; 2010; 59, pp. 552-560. [DOI: https://dx.doi.org/10.1016/j.cie.2010.06.014]

13. Panda, D.; Maiti, M.K.; Maiti, M. Two warehouse inventory models for single vendor multiple retailers with price and stock dependent demand. Appl. Math. Model.; 2010; 34, pp. 3571-3585. [DOI: https://dx.doi.org/10.1016/j.apm.2010.03.007]

14. Paul, B.; Rajendran, C. Rationing mechanisms and inventory control-policy parameters for a divergent supply chain operating with lost sales and costs of review. Comput. Oper. Res.; 2011; 38, pp. 1117-1130. [DOI: https://dx.doi.org/10.1016/j.cor.2010.11.002]

15. Guan, R.; Zhao, X. Pricing and inventory management in a system with multiple competing retailers under (r, Q) policies. Comput. Oper. Res.; 2011; 38, pp. 1294-1304. [DOI: https://dx.doi.org/10.1016/j.cor.2010.12.005]

16. Solyalı, O.; Süral, H. The one-warehouse multi-retailer problem: Reformulation, classification, and computational results. Ann. Oper. Res.; 2012; 196, pp. 517-541. [DOI: https://dx.doi.org/10.1007/s10479-011-1022-0]

17. Tempelmeier, H. A multi-level inventory system with a make-to-order supplier. Int. J. Prod. Res.; 2013; 51, pp. 6880-6890. [DOI: https://dx.doi.org/10.1080/00207543.2013.776190]

18. Wang, Q. A periodic-review inventory control policy for a two-level supply chain with multiple retailers and stochastic demand. Eur. J. Oper. Res.; 2013; 230, pp. 53-62. [DOI: https://dx.doi.org/10.1016/j.ejor.2013.04.004]

19. Gayon, J.-P.; Massonnet, G.; Rapine, C.; Stauffer, G. Constant approximation algorithms for the one warehouse multiple retailers problem with backlog or lost-sales. Eur. J. Oper. Res.; 2016; 250, pp. 155-163. [DOI: https://dx.doi.org/10.1016/j.ejor.2015.10.054]

20. Tayebi, H.; Haji, R.; Jeddi, B.G. Joint order (1, T) policy for a two-echelon, single-item, multi-retailer inventory system with Poisson demand. Comput. Ind. Eng.; 2018; 119, pp. 353-359. [DOI: https://dx.doi.org/10.1016/j.cie.2018.04.009]

21. John, K.; Rajendran, C.; Ziegler, H. A comparative study on allocation/rationing mechanisms operational with/without backorder clearing in divergent supply chains. Sadhana; 2019; 44, 231. [DOI: https://dx.doi.org/10.1007/s12046-019-1217-7]

22. Tsai, S.C.; Ho, I.-Y. Sample average approximation for a two-echelon inventory system with service-level constraints. J. Oper. Res. Soc.; 2019; 70, pp. 675-688. [DOI: https://dx.doi.org/10.1080/01605682.2018.1457479]

23. Hoque, M.A. An optimal solution policy to an integrated manufacturer-retailers problem with normal distribution of lead times of delivering equal and unequal-sized batches. OPSEARCH; 2021; 58, pp. 483-512. [DOI: https://dx.doi.org/10.1007/s12597-020-00485-2]

24. Najafnejhad, E.; Roodsari, M.T.; Sepahrom, S.; Jenabzadeh, M. A mathematical inventory model for a single-vendor multi-retailer supply chain based on the Vendor Management Inventory Policy. Int. J. Syst. Assur. Eng. Manag.; 2021; 12, pp. 579-586. [DOI: https://dx.doi.org/10.1007/s13198-021-01120-z]

25. Kurian, J.; Brijesh, P.; Chandrasekharan, R.; Ziegler, H. Priority fractional rationing (PFR) policy and a hybrid metaheuristic for managing stock in divergent supply chains. Sådhanå; 2022; 47, 254.

26. Lu, C.-J.; Gu, M.; Lee, T.-S.; Yang, C.-T. Integrated multistage supply chain inventory model of multiple retailers with imperfect production and inspection systems. Soft Comput.; 2022; 26, pp. 12057-12075. [DOI: https://dx.doi.org/10.1007/s00500-022-07490-1]

27. Smirnov, D.; Herer, Y.T.; Avrahami, A. The continuous delayed distribution problem. Comput. Oper. Res.; 2022; 148, 105976. [DOI: https://dx.doi.org/10.1016/j.cor.2022.105976]

28. Andersson, J.; Malmberg, F.; Marklund, J. Exact analysis of One-Warehouse-Multiple-Retailer inventory systems with quantity restricted deliveries. Eur. J. Oper. Res.; 2023; 309, pp. 1161-1172. [DOI: https://dx.doi.org/10.1016/j.ejor.2023.02.026]

29. Wang, Q.; Wan, G. Fixed-interval order-up-to policies and myopic optimal warehouse stock allocation for one-warehouse multiple-retailer systems. Eur. J. Oper. Res.; 2023; 309, pp. 1112-1124. [DOI: https://dx.doi.org/10.1016/j.ejor.2023.02.025]

30. Kaynov, I.; van Knippenberg, M.; Menkovski, V.; van Breemen, A.; van Jaarsveld, W. Deep Reinforcement Learning for One-Warehouse Multi-Retailer inventory management. Int. J. Prod. Econ.; 2024; 267, 109088. [DOI: https://dx.doi.org/10.1016/j.ijpe.2023.109088]

31. Stranieri, F.; Fadda, E.; Stella, F. Combining deep reinforcement learning and multi-stage stochastic programming to address the supply chain inventory management problem. Int. J. Prod. Econ.; 2024; 268, 109099. [DOI: https://dx.doi.org/10.1016/j.ijpe.2023.109099]

32. Wu, Y.; Cheng, T.C.E.; Zhang, J. A serial mixed produce-to-order and produce-in-advance inventory model with multiple retailers. Int. J. Prod. Econ.; 2012; 136, pp. 378-383. [DOI: https://dx.doi.org/10.1016/j.ijpe.2011.12.025]

33. Islam, S.S.; Hoque, A.; Hamzah, N. Single-supplier single-manufacturer multi-retailer consignment policy for retailers’ generalized demand distributions. Int. J. Prod. Econ.; 2017; 184, pp. 157-167. [DOI: https://dx.doi.org/10.1016/j.ijpe.2016.11.021]

34. Stauffer, G. Approximation algorithms for k-echelon extensions of the one warehouse multi-retailer problem. Math. Methods Oper. Res.; 2018; 88, pp. 445-473. [DOI: https://dx.doi.org/10.1007/s00186-018-0642-4]

35. Gruson, M.; Bazrafshan, M.; Cordeau, J.-F.; Jans, R. A comparison of formulations for a three-level lot sizing and replenishment problem with a distribution structure. Comput. Oper. Res.; 2019; 111, pp. 297-310. [DOI: https://dx.doi.org/10.1016/j.cor.2019.07.005]

36. Gruson, M.; Cordeau, J.-F.; Jans, R. Split demand and deliveries in an integrated three-level lot sizing and replenishment problem. Comput. Oper. Res.; 2024; 161, 106434. [DOI: https://dx.doi.org/10.1016/j.cor.2023.106434]

37. Cunha, J.O.; Melo, R.A. Valid inequalities, preprocessing, and an effective heuristic for the uncapacitated three-level lot-sizing and replenishment problem with a distribution structure. Eur. J. Oper. Res.; 2021; 295, pp. 874-892. [DOI: https://dx.doi.org/10.1016/j.ejor.2021.03.029]

38. Gharaei, A.; Amjadian, A.; Shavandi, A. An integrated reliable four-level supply chain with multi-stage products under shortage and stochastic constraints. Int. J. Syst. Sci. Oper. Logist.; 2023; 10, 1958023. [DOI: https://dx.doi.org/10.1080/23302674.2021.1958023]

39. Geevers, K.; van Hezewijk, L.; Mes, M.R.K. Multi-echelon inventory optimization using deep reinforcement learning. Central Eur. J. Oper. Res.; 2023; 32, pp. 653-683. [DOI: https://dx.doi.org/10.1007/s10100-023-00872-2]

40. Lin, H.; Lin, J.; Wang, F. An innovative machine learning model for supply chain management. J. Innov. Knowl.; 2022; 7, 100276. [DOI: https://dx.doi.org/10.1016/j.jik.2022.100276]

41. Blancas-Rivera, R.; Cruz-Suárez, H.; Portillo-Ramírez, G.; López-Ríos, R. $(s, S)$ Inventory policies for stochastic controlled system of Lindley-type with lost-sales. AIMS Math.; 2023; 8, pp. 19546-19565. [DOI: https://dx.doi.org/10.3934/math.2023997]

42. Akbari, M.; Do, T.N.A. A systematic review of machine learning in logistics and supply chain management: Current trends and future directions. Benchmarking Int. J.; 2021; 28, pp. 2977-3005. [DOI: https://dx.doi.org/10.1108/BIJ-10-2020-0514]

43. Bijvank, M.; Vis, I.F. Lost-sales inventory theory: A review. Eur. J. Oper. Res.; 2011; 215, pp. 1-13. [DOI: https://dx.doi.org/10.1016/j.ejor.2011.02.004]

© 2025 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https://creativecommons.org/licenses/by/4.0/). Notwithstanding the ProQuest Terms and Conditions, you may use this content in accordance with the terms of the License.

Analysis of a Three-Echelon Supply Chain System with Multiple Retailers, Stochastic Demand and Transportation Times

Content area

Abstract

Full text

1. Introduction

2. Literature Review

3. System Description

4. Modeling Approach

4.1. States Definition

4.2. The Infinitesimal Generator Matrix

4.3. Performance Measures

4.4. Validation and Model Performance

5. Numerical Results

5.1. The Effect of Policy Parameters—Balanced Systems

5.1.1. Distribution Center’s Policy (sd, Qd) for Balanced Systems

5.1.2. Wholesaler Policy (sw, Qw) for Balanced Systems

5.1.3. Retailers’ Policy (si, Qi) for Balanced Systems

5.2. The Effect of Policy Parameters—Supply-Constrained Systems

5.2.1. Distribution Center’s Policy (sd, Qd) for Supply Constrained Systems

5.2.2. Wholesaler’s Policy (sw, Qw) for Supply Constrained Systems

5.2.3. Retailers’ Policy (si, Qi) for Supply Constrained Systems

5.3. The Effect of Number of Retailers

5.3.1. Interactions Between the Retailers

5.3.2. The Effect of Retailer Addition

6. Conclusions

5.1.1. Distribution Center’s Policy (s_d, Q_d) for Balanced Systems

5.1.2. Wholesaler Policy (s_w, Q_w) for Balanced Systems

5.1.3. Retailers’ Policy (s_i, Q_i) for Balanced Systems

5.2.1. Distribution Center’s Policy (s_d, Q_d) for Supply Constrained Systems

5.2.2. Wholesaler’s Policy (s_w, Q_w) for Supply Constrained Systems

5.2.3. Retailers’ Policy (s_i, Q_i) for Supply Constrained Systems