Subscribe Now Subscribe Today
Research Article

Optimizing Multi-Product Multi-Constraint Inventory Control Systems with Stochastic Replenishments

Ata Allah Taleizadeh, Mir-Bahador Aryanezhad and Seyed Taghi Akhavan Niaki
Facebook Twitter Digg Reddit Linkedin StumbleUpon E-mail

Multi-periodic inventory control problems are mainly studied employing two assumptions. The first is the continuous review, where depending on the inventory level orders can happen at any time and the other is the periodic review, where orders can only happen at the beginning of each period. In this study, we relax these assumptions and assume that the periodic replenishments are stochastic in nature. Furthermore, we assume that the periods between two replenishments are independent and identically random variables. For the problem at hand, the decision variables are of integer-type and there are two kinds of space and service level constraints for each product. We develop a model of the problem in which a combination of back-order and lost-sales are considered for the shortages. Then, we show that the model is of an integer-nonlinear-programming type and in order to solve it, a search algorithm can be utilized. We employ a simulated annealing approach and provide a numerical example to demonstrate the applicability of the proposed methodology.

Related Articles in ASCI
Similar Articles in this Journal
Search in Google Scholar
View Citation
Report Citation

  How to cite this article:

Ata Allah Taleizadeh, Mir-Bahador Aryanezhad and Seyed Taghi Akhavan Niaki , 2008. Optimizing Multi-Product Multi-Constraint Inventory Control Systems with Stochastic Replenishments. Journal of Applied Sciences, 8: 1228-1234.

DOI: 10.3923/jas.2008.1228.1234



In multi-periodic inventory control models, the continuous review and the periodic review are the major vastly used policies. However, the underlying assumptions of the proposed models restrict their correct usage and utilization in real-world environments. In continuous review policy, the user has the freedom to act at anytime and replenish orders based upon the available inventory level. While in the periodic review policy, the user is allowed to replenish the orders only in specific and predetermined times.

The multi-periodic inventory control problems have been investigated in depth in different research. Chiang (2003) considered a periodic review model in which the period is partly long. Where discount has been considered, the costs associated with his modeling were the purchasing, holding and fixed order costs. The important aspect of his study was to introduce emergency orders to prevent shortages. The emergency orders were replenished with a higher cost compared to the normal ordering costs. He employed a dynamic programming approach to model the problem. Mohebbi and Posner (2002) investigated an inventory system with periodic review, multiple replenishment and multi-level delivery. They assumed that the demand was a Poisson random variable, shortages were allowed and that the lost sale policy could be employed. Feng and Rao (2007) considered a (R,nT) model in which in the first level a stochastic demand entered the system and the total unsatisfied demand were back-ordered at the second level. Ouyang and Chuang (2000) investigated a (R,T) model in which the period-length and lead-time were the decision variables, demand was a stochastic variable and the service level was a constraint. Chiang (2006) analyzed a periodic review problem in two cases of back-order and lost sales and employed the (R, T) policy. Qu et al. (1999) investigated a transportation model integrated with an inventory model with a periodic review policy. They designed a multi-periodic model with several providers in which the demand was stochastic and the planning horizon was finite. Eynan and Kropp (2007) have propounded the assumption of stochastic demand and variant warehousing costs on a periodic review system; while assuming nonzero lead-time and safety stock. Bylka (2005) investigated a model with constraints on the amounts of orders and back-order shortages in which the lead-time was constant and demand was stochastic. By analyzing the changes in the lead-time and ordering cost, they tried to define the optimized ordering time. For a continuous review inventory model, assuming stochastic lead-time, Mohebbi (2004) considered demand a compound Poisson random variable.

In summary, while there has been separate emphasis on the stochastic nature of the demands and lead-time, the real-world constraints of the systems have not been completely investigated. Specifically, there is no research in which both the demands and lead-time are considered to be probabilistic. Furthermore, some constraints have been partially studied, the decision variable has been considered integer and constraints such as budget and space have not been investigated.

Three main specifications of the proposed model of this research that have led to its novelty are the stochastic period length, the allowance of both multi-products and multi-constraints and the fact that the decision variables are integer. By deploying these conditions simultaneously, the created model is different from the other models in the periodic review literature.


Consider a periodic inventory control model for one provider in which the times required to order each of several available products are stochastic in nature. Let the time-periods between two product-replenishments be identical and independent random variables. The demands of the products are constant and distinct and in case of shortage, a fraction is considered back-order and a fraction as lost-sale. The costs associated with the inventory control system are holding, back-order, lost-sales and purchase costs. There are several products and the warehouse space and service level for each product are considered constraints of the problem. Furthermore, the lead-time is assumed zero and the decision variables are integer digits. We need to identify the inventory levels in each cycle such that the expected profit is maximized. In short, the assumptions involved in the problem are:

Time-periods between replenishments are independent and identical random variables
There are several products
There are several constraints
The decision variables are integer
There are holding, shortage and purchase costs
A fraction of a shortage is back-ordered
Only one provider exists
The demands are constant
All of the purchased product will be sold


For the problem at hand, since the time-periods between two replenishments are independent random variables, in order to maximize the expected profit of the planning horizon we need to consider only one period. Furthermore, since we assumed that the costs associated with the inventory control system are holding and shortage (back-order and lost-sale), we need to calculate the expected inventory level and the expected required storage space in each period. Before doing this, let us define the parameters and the variables of the model.

The parameters and the variables of the model: For, i = 1, 2,....n, let us define the parameters and the variables of the model as

Ri : : The inventory level of the ith product
Ti : A random variable denoting the time-period between two replenishments (cycle length) of the ith product
fTi(ti) : The Probability density function of Ti
hi : The holding cost per unit inventory of the ith product in each period
πi : The back-order cost per unit demand of the ith product
Wi : The purchasing cost per unit of the ith product
Pi : The sale price per unit of the ith product
Di : The constant demand rate of the ith product
tDi : The time at which the inventory level of the ith product reaches zero
βi : The percentage of unsatisfied demands of the ith product that is back-ordered
Ii : The expected amount of the ith product inventory per cycle
Li : The expected amount of the ith product lost-sale in each cycle
Bi : The expected amount of the ith product back-order in each cycle
Qi : The expected amount of the ith product order in each cycle
SLi : The lower limit of the service level for the ith product
fi : The required warehouse space per unit of the ith product
F : Total available warehouse space
Ch : The expected holding cost per cycle
Cb : The expected shortage cost in back-order state
Cl : The expected shortage cost in lost-sale state
Cp : The expected purchase cost
r : The expected revenue obtained from sales
Z : The expected profit obtained in each cycle

Fig. 1: Presenting the inventory cycle when Tmin≤T≤tD

Fig. 2: Presenting the inventory cycle when tD≤T≤Tmax

Inventory diagram: According to Ertogal and Rahim (2005) and considering the fact that the time-periods between replenishments are stochastic variables, two cases may occur. In the first case the time-period between replenishments is less than the amount of time required for the inventory level to reach zero (Fig. 1) and in the second case, it is greater (Fig. 2). Figure 3 shows the shortages in both cases.

Single product model–back ordered and lost sales cases
Calculating the costs and the profit: In order to calculate the expected profit in each cycle, we need to evaluate all of the terms in Eq. 1 (Ertogal and Rahim, 2005).


Based on Fig. 3, L, B, I and Q are evaluated by the following equations:




Fig. 3: Presenting shortages in two cases of compact back order and lost sales


Presenting the constraints: As the total available warehouse space is F, the space required for each unit of product is f and the upper limit for inventory is R, the space constraint will be:


Since the shortages only occur when the cycle time is more than tD and that the lower limit for the service level is SL, then;


In short, the complete mathematical model of the single product inventory is:



As introduced later, in the extending phase of the single-product model to the multiple product model, in addition to two cases of back-order and lost-sales, a combination of back-order and lost-sales is also needed.

Multi-product model-back order and lost sales case: The single-product inventory model can be easily extended to a multiple-product model as follows:



In what follows, we consider two probability density functions for Ti and hence we develop two models.

Ti follows a uniform distribution: In this case the probability density function of Ti is . Accordingly, (9) will change to:



Ti follows an exponential distribution: If follows an exponential distribution with parameter λi, then the probability density function of Ti will be fTi (ti) = λieGλiTi. In this case, the model is shown below:



Further, we will introduce a heuristic algorithm to solve the problem.


Since the models in (10) and (11) are integer-nonlinear in nature, reaching an analytical solution (if any) to the problem is difficult (Gen, 1997). As a result, in this section we will try to solve the problem by the stochastic search algorithm of simulated annealing.

Simulated annealing: Simulated Annealing (SA) is a local search method that finds its inspiration in the physical annealing process studied in statistical mechanics (Aarts and Korst, 1989). It is an effective and efficient method that produces good suboptimal solutions and has diffused rapidly into many combinatorial optimization areas (Kirkpatrick et al., 1994). An SA algorithm repeats an iterative neighbor generation procedure and follows search directions that improve the objective function value. While exploring solution space, the SA method offers the possibility of accepting worse neighbor solutions in a controlled manner in order to escape from local minima. More precisely, in each iteration, for a current solution x characterized by an objective function value f (x), a neighbor x' is selected from the neighborhood of x denoted N (x) and defined as the set of all its immediate neighbors. For each move, the objective difference Δ = f (x')!f (x) is evaluated. For minimization problems, x' replaces x whenever Δ≤0. Otherwise, x' could also be accepted with a probability . The acceptance probability is compared to a number RN 0 Uni (0, 1) generated randomly and x' is accepted whenever P>RN.

The factors that influence the acceptance probability are the degree of objective function value Δ degradation (smaller degradations induce greater acceptance probabilities) and the parameter T called temperature (higher values of T give higher acceptance probability). The temperature can be controlled by a cooling scheme specifying how it should be progressively reduced to make the procedure more selective as the search progresses to neighborhoods of good solutions. There exist theoretical schedules, which require infinite computing time guaranteeing asymptotic convergence toward the optimal solution. In practice, much simpler and finite computing time schedules are preferred even if they do not guarantee an optimal solution.

A typical finite time implementation of SA consists of decreasing the temperature T in S steps, starting from an initial value T0 and using an attenuation factor α, (0<α<1). The initial temperature T0 is supposed to be high enough to allow acceptance of any new neighbor proposed in the first step. In each step s, the procedure generates a fixed number of neighbor solutions Nsol and evaluates them using the current temperature value TS = αST0. The whole process is commonly called cooling chain or also markov chain. Adaptation of SA to an optimization problem consists in defining its specific components: A neighbor generation procedure, objective function evaluation, a method for assigning the initial temperature, a procedure to change the temperature and a cooling scheme including stopping criteria. These adaptation steps for our new adaptations of SA for the inventory model of (10) and (11) are described below.

A neighbor generation procedure: The neighbor generation is an important component of SA. It has to be designed to allow easy neighbor generation and fast calculation of Δ and guarantee accessibility for the entire solution space.

In this study, the initial solutions are generated in two groups. In the first group, the procedure randomly chooses the initial solution from the possible solution space. For the other groups of initial solutions the procedure employs some solutions of a Genetic Algorithm.

Evaluating the objective function: After creating each solution, the objective function should be evaluated to compare the jth and the ith solutions. In a maximization problem, if the objective function of the new solution (j) is bigger than the previous solution (i), then (j) will be accepted. Otherwise, by generating a random number the better solution will be selected.

Assigning the initial temperature: Temperature is one of the important parameters of the SA algorithm that affects the acceptance or non-acceptance of the objective function changes. The initial temperature should be chosen such that at the first level a big number of inappropriate solutions are accepted. In this way, more changes are possible and hence more solutions are explored. Furthermore, the time required for the process repetition along with the annealing process depends on the initial temperature. In this research, the large values of 1000, 2000 and 5000 are chosen for initial temperatures.

Changing the temperature: One of the primary aspects of the annealing process in a SA algorithm is the range of temperature change. The temperature plays an important role in the possibility of selecting a bad solution. On one hand, when the temperature assumes a high value, a big number of bad solutions are accepted, leading to selecting a local optimum point. On the other hand, when the temperature is low, the probability of the solution to be a local optimum is high. In this paper, we change the temperature of the SA algorithm based on a geometric function shown in Eq. 12 with α = 0.9, 0.95 and 0.99.


Cooling scheme and stopping criterion: In a specific temperature of a SA algorithm, it is necessary to analyze the equilibrium state after a couple of renitence to see if the annealing process will be continued in that temperature or it will be stopped and moved to the next temperature. An epoch is the number changes that should occur in a specific temperature.

Different types of stopping criteria have been presented for different SA algorithms in the literature. Some of them are:

A lower bound on the final temperature
An upper bound on the total number of states
Reaching a solidification point based on the value of a function evaluated on the value of the objective function and the temperature
An upper bound on the total number of exchanges accepted along the annealing process
An upper bound on the total number of refused changes

In this research, we select the first criterion and stop the algorithm after it reaches the final temperature Ti. Furthermore, we use 100, 200 and 300 for different values of N (t).

In short, the steps involved in the SA algorithm used in this research are:

(1) Choosing an initial solution i from the group of the feasible solutions S.
(2) Choosing the initial temperature T0>0.
(3) Selecting the number of iterations N (t) at each temperature.
(4) Selecting the final temperature TF.
(5) Determining the process of the temperature reduction until it reaches TF.
(6) Setting the temperature exchange counter n to zero for each temperature. (Balancing process).
(7) Creating the j solution at the neighborhood of the i solution.
(8) Evaluating the objective function f = Z at any temperature.
(9) Calculating Δ = f (j)!f (i).
(10) Accepting the solution j, if Δ<0. Otherwise, generating a random number RN-Uniform. If then selecting the j solution.
(11) Setting n = n+1. If n is equal to N (t) then go to 12. Otherwise, go to 7.
(12) Reducing the temperature. If it reaches TF then stop. Otherwise, go to 6.

In order to demonstrate the proposed SA algorithm and evaluate its performances, in the next section we bring a numerical example used in Ertogal and Rahim (2005). In this example, two cases of the uniform and the exponential distributions for the time-period between two replenishments are investigated.


Consider a multi-product inventory control problem with eight products and general data shown in Table 1. Table 2 and 3 show the parameters of the uniform and the exponential distributions used for the time-period between two replenishments, respectively. The total available warehouse space is 18000.

Table 1 : General data

Table 2 : Data for uniform distribution

Table 3 : Data for exponential distribution

Table 4 : The parameters of the SA algorithm

Table 4 shows different values of the parameters of the SA method. In this research all the possible combinations of the parameters in SA (N (t), T0 and α) method are employed and using the max(max) criterion the best combination of the parameters has been selected. Table 5 and 6 shows the best result and best combination, respectively. Furthermore, the convergence paths of the objective function values of the SA algorithm in uniform and exponential distributions are shown in Fig. 4 and 5.

Table 5: The best result for Ri

Table 6: The best combination of the SA parameters

Fig. 4: The convergence paths of the best result in Uniform example of the SA

Fig. 5: The convergence paths of the best result in Exponential example of SA


In this study, a stochastic replenishment multi-product inventory model was investigated. Two mathematical modeling for two cases of uniform and exponential distribution of the time between two replenishments have been developed and shown to be integer-nonlinear programming problems. Then, a meta-heuristic solution algorithm of SA has been proposed to solve the integer non-linear problems. Finally, two numerical examples were given to demonstrate the applicability of the proposed methodology.

Some recommendation on future works follow:

As we may consider the model parameters fuzzy, the analysis may be performed using fuzzy systems.
In addition to the SA algorithm, some other meta-heuristic algorithms like Genetic Algorithm, Tabu-Search and Ant-Colony optimization may be employed to solve the integer non-linear problems.
Some other probability density functions may be considered for the time between replenishments.
The discount factor may be augmented to the problem.
1:  Allameh, A., A.R. Safamehr, S.A. Mirhadi, M. Shivazad and M. Razzaghi-Abyaneh, 2005. Evaluation of biochemical and production parameters of broiler chicks fed ammonia treated aflatoxin contaminated maize grains. Anim. Feed. Sci. Technol., 122: 289-301.
CrossRef  |  Direct Link  |  

2:  Aravind, K.L., V.S. Patil, G. Devegowda, B. Umakantha and S.P. Ganpule, 2003. Efficacy of esterified glucomannan to counteract mycotoxicosis in naturally contaminated feed on performance and serum biochemical and hematological parameters in broilers. Poult. Sci., 82: 571-576.
Direct Link  |  

3:  Devegoda, G., B.I.R. Aravind and M.G. Morton, 1997. Biotechnology in the feed industry. Proceedings of Alltech's 13th Annual Symposium, April 30, 1997, Nottingham University Press, pp: 205-215.

4:  Edrington, T.S.L., F. Kubena, R.B. Harvey and G.E. Rottinghaus, 1997. Influence of a super activated charcoal on the toxic effects of aflatoxin or T-2 toxin in growing broilers. Poult. Sci., 76: 1205-1211.
PubMed  |  

5:  Edrington, T.S., A.B. Sarr, L.F. Kubena, R.B. Harvey and T.D. Phillips, 1996. Hydrated sodium calcium aluminosilicate (HSCAS), acidic HSCAS and activated charcoal reduce urinary excretion of aflatoxin M1 in turkey poults. Lack of effect by activated charcoal on aflatoxicosis. Toxicol. Lett., 89: 115-122.
CrossRef  |  

6:  Eraslan, G., B.C. Liman, B.K. Guclu, A. Atasever, A.N. Koc and L. Beyaz, 2004. Evaluation of aflatoxin toxicity in Japanese quails given various doses of hydrated sodium calcium aluminosilicate. Bull. Vet. Res. Inst. Pulawy, 48: 511-517.
Direct Link  |  

7:  Teleb, H.M., A.A. Hegazy and Y.A. Hussein, 2004. Efficiency of kaolin and activated charcoal to reduce the toxicity of low level of aflatoxin in broilers. Scient. J. King Faisal Univ. (Basic Applied Sci.), 5: 145-160.
Direct Link  |  

8:  Jindal, N., S.K. Mahipal and N.K. Mahajan, 1994. Toxicity of aflatoxin B1 in broiler chicks and its reduction by activated charcoal. Res. Vet. Sci., 56: 37-40.
CrossRef  |  PubMed  |  Direct Link  |  

9:  Kececi, T., H. Oguz, V. Kurtoglu and O. Demet, 1998. Effects of polyvinylpolypyrrolidone, synthetic zeolite and bentonite on serum biochemical and haematological character of broiler chickens during aflatoxicosis. Br. Poult. Sci., 39: 452-458.
PubMed  |  

10:  Kubena, L.F., R.B. Harvey, R.H. Bailey, S.A. Buckley and G.E. Rottinghaus, 1998. Effects of a hydrated sodium calcium aluminosilicate (T-Bind) on mycotoxicosis in young broiler chicks. Poult. Sci., 77: 1502-1509.
PubMed  |  

11:  Kubena, L.F., R.B. Harvey, T.D. Phillips and B.A. Clement, 1992. The use of sorbent compounds to modify the toxic expression of mycotoxins in poultry. Proceedings of the 19th World's Poultry Congress, September 20-24, 1992, Bangladesh, pp: 357-360.

12:  Kubena, L.F., R.B. Harvey, W.E. Huff, D.E. Corrier, T.D. Phillips and G.E. Rottinghaus, 1990. Efficacy of a hydrated sodium calcium aluminosilicate to reduce the toxicity of aflatoxin and T-2 toxin. Poult. Sci., 69: 1078-1086.
PubMed  |  

13:  Kubena, L.F., R.B. Harvey, W.E. Huff, M.H. Elissalde, A.G. Yersin, T.D. Phillips and G.E. Rottinghaus, 1993. Efficacy of a hydrated sodium calcium aluminosilicate to reduce the toxicity of aflatoxin and diacetoxyscirpenol. Poult. Sci., 72: 51-59.
PubMed  |  

14:  Ledoux, D.R., G.E. Rottinghaus, A.J. Bermudez and M. Alonso-Debolt, 1999. Efficacy of a hydrated sodium calcium aluminosilicate to ameliorate the toxic effects of aflatoxin in broiler chicks. Poult. Sci., 78: 204-210.
CrossRef  |  Direct Link  |  

15:  Legator, M., 1966. Biological effects of aflatoxin in cell culture. Bacterial Rev., 30: 471-477.
PubMed  |  

16:  Mgbodile, M.U.K., M. Holscher and R.A. Neal, 1975. Possible protective role toxicity: Effect of pretreatment of reduced glutathione in aflatoxin B1 rats with phenobarbitol and 3-methylchoranthrene on aflatoxin toxicity. Toxicol. Applied Pharmacol., 34: 128-142.

17:  Miazzo, R., C.A. Rosa, E.C. De Queiroz Carvalho, C. Magnoli and S.M. Chiacchiera et al., 2000. Efficacy of synthetic zeolite to reduce the toxicity of aflatoxin in broiler chicks. Poult. Sci., 79: 1-6.
PubMed  |  Direct Link  |  

18:  Miazzo, R.D., M.F. Peralta, C. Magnoli, M. Salvano and S. Ferrero et al., 2005. Efficacy of sodium bentonite as a detoxifier of broiler feed contaminated with aflatoxin and fumonisin. Poult. Sci., 84: 1-8.
CrossRef  |  Direct Link  |  

19:  NRC. (National Research Council), 1994. Nutrient Requirement of Poultry. 3rd Edn., National Academy Press, Washington, DC., USA., pp: 3.

20:  Oguz, H., T. Kececi, Y.O. Birdane, F. Onder and V. Kurtoglu, 2000. Effect of clinoptilolite on serum biochemical and haematological characters of broiler chickens during aflatoxicosis. Res. Vet. Sci., 69: 89-93.
CrossRef  |  Direct Link  |  

21:  Ortatatli, M., H. Oguz, F. Hatipoglu and M. Karaman, 2005. Evaluation of biochemical characters of broiler chickens during dietary aflatoxin (50-100 ppb) and clinoptilolite exposure. Res. Vet. Sci., 78: 61-68.

22:  Parlat, S.S., A.O.H. Yildiz and H. Oguz, 1999. Effect of clinoptilolite on fattening performance of Japanese quail (Coturinix coturnix japonica) during experimental aflatoxicosis. Br. Poult. Sci., 40: 495-500.

23:  Pasha, T.N., M.U. Farooq, F.M. Khattak and M.A. Jabbar, 2006. Effectiveness of sodium bentonite and two commercial products as aflatoxin absorbents in diets for broiler chickens. J. Anim. Feed Sci. Technol., 132: 103-110.
Direct Link  |  

24:  Raju, M.V. and G. Devegowda, 2000. Influence of esterified-glucomannan on performance and organ morphology, serum biochemistry and haematology in broilers exposed to individual and combined mycotoxicosis (aflatoxin, ochratoxin and T-2 toxin). Br. Poult. Sci., 41: 640-650.
PubMed  |  Direct Link  |  

25:  Rosa, C.A., R.D. Miazzo, C. Magnoli, M. Salvano and S.M. Chiacchiera et al., 2001. Evaluation of efficacy of bentonite from the South of Argentina to ameliorate the toxic effects of aflatoxin in broilers. Poult. Sci., 80: 139-144.
Direct Link  |  

26:  Shane, M.S., 2001. Mannan Oligosaccharides in Poultry Nutrition: Mechanism and Benefits. In: Science and Technology in the Feed Industry: Proceedings of Alltech's 17th Annual Symposium, Lyons, T.P. and K.A. Jacques (Eds.). Nottingham University Press, Nottingham, UK., pp: 65-77.

27:  Smith, J.W. and P.B. Hamilton, 1970. Aflatoxicosis in the broiler chicken. Poult. Sci., 49: 207-215.
PubMed  |  

28:  Stanley, V.G. and O. Jor, 1993. The use of Saccharomyces cerevisiae to suppress the effect of aflatoxicosis in broiler chicks. Poult. Sci., 72: 1867-1872.
Direct Link  |  

29:  Trucksess, M.W., M.E. Stack, S. Nesheim, R.H. Albert and T.R. Romer, 1994. Multifunctional column coupled with liquid chromatography for determination of aflatoxins B1, B2, G1 and G2 in corn, almonds, Brazil nuts, peanuts and pistachio nuts: Collaborative study. J. AOAC Int., 77: 1512-1521.
PubMed  |  

30:  Yildiz, A.O., S.S. Parlatand and I. Yildirim, 2004. Effect of dietary addition of live yeast (Saccharomyces cerevisiae) on some performance parameters of adult Japanese quail (Coturnix japonica) induced by aflatoxicosis. Revue. Med. Vet., 155: 38-41.
Direct Link  |  

©  2021 Science Alert. All Rights Reserved