Discrete Models for Intermediate Storages under Stochastic Operational Conditions

In this paper we investigate the operation of an intermediate storage assuming discrete stochastic operational conditions. The storage operates as a buffer. The material is collected into it and is withdrawn from the buffer for production. Deterministic constant withdrawal rate is supposed and the filling process is described by discrete random variables. The main question addressed here is the probability of material shortage as a function of the initial amount of material. For this purpose, an auxiliary function is defined and a difference equation is set for it. In a special case we present the solution of the difference equation, we compare the exact and simulated results, and we propose a method for the determination of unknown constants. We apply the results for expressing the initial amount of material necessary to a given reliability level. Finally, we compare the results of the discrete model with those of the continuous one.


Introduction
Intermediate storages are frequently used in operational systems.They are able to balance out the differences between the operations of different subsystems.Storages are suitable for providing spare material in case of faults or maintenance problems and so on.They serve as a buffer in chemical, pharmaceutical or food industry, in water supply but they can be used in informatics as well (see for example [1,2,3]).The investigation of their operation is an important task in process engineering.
It can be easily realized that the operation of such systems is often stochastic.Random faults are supposed in [4] and a general stochastic model is set up in [5].In stochastic models an important question is the reliability of the system [4,5,6].Stochastic systems can be studied by simulation or by analytical methods.The analytical methods require plenty of mathematical tools in probability theory and analysis.Sometimes it is worth combining them with simulations to get useful solutions.In this paper we present a method which applies both for providing the solution of the model.
Moreover, the practical problems are rather discrete than continuous, but mathematical models are frequently continuous.The reason might be that the theory of differential and integral equations is well-developed, while the theory of difference equations is less elaborated.Some discrete models are also presented, see for example [3,7].In this paper we turn our attention to a model defined by discrete random variables.
In this presentation we first introduce the model investigated and we define an auxiliary function which can handle the reliability and the expected time of material shortage at the same time.Then we set up a difference equation for it and we provide its solution in the case of special discrete distributions.We compare this solution with the results given by the Monte-Carlo simulation supporting their correctness.Applying the form of the exact solution we suggest a numerical method for providing an approximate solution and we compute the necessary initial amount of material for given reliability levels.Finally, we present an example to compare the results of the discrete model and the corresponding continuous one.

The investigated model
Let us consider an intermediate storage which collects the material and stores it until processing.For processing, the material is withdrawn from the buffer and it is used as raw material.The input process is supposed to be a batch process.Moreover, the batches filled into the storage are assumed to have random amounts and the fillings arrive at random time as well.The output process is assumed to be linear as a function of time.The process begins at time T 0 =0 and the time interval between the (i-1) th and i th filling is denoted by T i , i = 1,2,... .Consequently, the i th filling happens at the time T s s i = ∑ 1 , i = 1,2,... .We suppose that T s , s = 1,2,..., are independent, identically distributed nonnegative random variables.In this paper we suppose that they are discrete random variables with integer values 0,1,2,..., j,... and P(T s = j) = f (j), 0 ≤ f (j) and . The inequality 0 < f (0) stands for the case when more than one batch can be filled into the storage at the same time with positive probability.N(t) provides the number of fillings in the interval [0,t], that is The amount of material filled into the buffer during the i th filling (Y i , i = 1,2,...) is described by independent, identically distributed discrete random variables with nonnegative integer values 0,1,2,... k,... and P(Y i = k) = g(k).The inequalities 0 ≤ g(k) and the equality . The withdrawn material is equal to m, consequently the amount of material in the intermediate storage at time m is equal to We emphasise that it contains the material filled at time m and withdrawn at the same time as well.The amount of material in the storage before the last filling is equal to with U(0−) = n We define shortage as lack of material in the storage.This means that at a given time the function describing the amount of material in the storage is less than or equal to zero.This is the event {U(m−)≤ 0 for some 0 ≤ m} and the time of shortage is "the first" time when the event {U(m−)≤ 0} occurs.By mathematical formula, Note that it is a function of the amount of initial material n.The definition of shortage given by (3) expresses that the amount of material in the storage falls to zero before the actual filling happens.This time is the first time point when the storage becomes empty.Other possibilities for this definition are: a definitely negative value and/or taking the actual fillings into account.These variations can not essentially change neither the method nor the results of this presentation.The model given by continuous random variables is investigated in [8].In that case these variations are not important.For many points of view discrete models are more complicated than the continuous models, although they are more realistic.
The amount of material in the storage as a function of time is a point process and it can be seen in Fig. 1.The initial amount of material equals n=5.The function can be interpreted as follows: the process starts from 5, and the material decreases one by one until the first filling.It happens at T 1 =3.Before the first filling the material in the storage is 2 units, the filled material is Y 1 =5 units, hence after the first filling there are 7 units T n m U m m if there exists m for which U m U ( ) min : ( ) , , ( ) T U (5) of material in the storage.This decreases until the next filling, which happens at T 1 +T 2 =9.The amount of the filled material is Y 2 =3 units, consequently the function U(m) jumps upwards 3 units.Then the amount of material is decreasing again and it reaches level 0 at m=13, consequently there is material shortage in the buffer.The time point when it is realized first is m=13, consequently T U (5)=13 in this case.As T i and Y i are random, in another example the level zero might be reached at another time point, or there might not be a shortage at all.The question is the probability and the expected time of the shortage during such a stochastic process.
To determine these quantities we define the following auxiliary function: where E denotes the expectation and It can be seen that and Consequently the function x z (n) is suitable for investigating the probability and the expected time of the shortage together.Moreover, right hand side of Eq. (5) shows that x z (n) is the z-transform of the distribution of the finite time of shortage.The z-transform is the discrete version of the Laplace transform [9].

The difference equation for x z (n) and its solution
First we note that the mathematical background of this paper has not been published yet, but the mathematical details are out of the frame of this paper.
On the basis of ( 5), it can be seen that |x z (n)|≤z −n .With the help of renewal theory it can be proved that, for any values of 1≤z, x z (n) satisfies the following difference equation: Eq.( 8) is not a recursion.Its solution can not be computed step by step because the formula concerning x z (n) contains the members of the sequence with indices greater than n.This equation has a unique bounded solution and it was partially investigated in [10] for geometrically distributed random times and amount of material.Now we extend the set of distributions to the case when the amount of material has an arbitrary discrete distribution and the time distribution satisfies a linear difference equation with constant coefficients.
We introduce the following notations for a sequence b = b(.), Let us assume that f(.) satisfies the difference equation p(D) (f)≡0.Then, solution of Eq.( 8) has the following form: where μ i (z) are the roots of the characteristic equation of the difference equation ( 8) with multiplicities k i , are constant values depending on the initial conditions on x z (n), which originate from the difference equation satisfied by f.One can prove that (8) has L roots with the property |μ i (z)|<1 for any 1 < z and also in the case z . This inequality expresses that the expectation of the filled amount of material during a unit time interval is more than the withdrawn material during that time.The characteristic equation of the difference equation can be written as where the polynomial q(y) is of the form q y a q y Consequently, based on (9), Taking the derivative with respect to z we obtain

Determination of coefficients and exponents on the basis of simulation results
Distributions which satisfy a linear difference equation form a wide class of discrete distributions, therefore (11) and (12) can be frequently applied.But they require the determination of μ i (z), c i,j (z) and also the derivatives at the point z=1.It is a hard task even for simple polynomials p.The first problem can be the computation of the initial values of (8).We note that if A j =0 for j = 0,1,2,...,L−2, then Eq.( 10) contains an infinite series.Applying complex analysis (5) x n x n k j f j g k z f j z it can be proved that (10) has exactly L roots (with multiplicities) with the property |μ i (z)|<1.Unfortunately, the computation of these roots is complicated, and so are the derivatives, as well.Nevertheless, the form of the solution is very useful, because it gives us the possibility to get an approximation for the probability and expectation.First we present an example in which we compare the results arising from the exact formula (9) and the probability arising from Monte-Carlo simulation.Let p(y) = (s + y) 2 with some 0 < s < 1 and suppose that p(D)(f) ≡ 0, with A = 0 and A 1 = s 2 .One can prove that f(n) = ns 2 (1 − s) n − 1 , n = 0,1,... .In case of g(1) = 0.5 , g(2) = 0.5, Eq.( 10) has the form If z = 1, then μ(1) = 1 is a root of ( 13) and two further roots can be computed by solving a quadratic equation.Denoting them by μ 1 (1) and μ 2 (1), the constant values are c 1 2 If ).The process is investigated up to a fixed (large) value of time.The probability of shortage is estimated by the relative frequency of shortages.The expectation to be estimated is , and we estimate it by averaging the realized shortage times and zeros.The zero values come up when there is no shortage during the realization.The simulations were carried out by MATLAB.
The simulated, the exact and the approximated probabilities can be seen in Fig. 2.These probabilities are very close to each other.The differences can hardly be seen in Fig. 2, to present them we plotted Fig. 3.We can see that the maximal difference between the exact and the approximate values of the probabilities is less than 0.004 and the maximal difference between the exact and simulated results is less than 0.01.As the number of simulations was N=10000, the above mentioned error is less than the error of the Monte-Carlo simulation.The approximate values of the derivatives were determined numerically by the simulated values of x z (n), n=1,2,3,4,5 z =1, z = 1.005, z = 1.01.Numerical values for them are Substituting these values into (12), the approximate values of the expected time of shortage and their simulated values are plotted in Fig. 4. The approximate values follow quite well the simulated values even in the case of large initial amount of material, which is important from practical point of view.
Finally, if we want to find the necessary initial amount of material for the reliability level 1 − α , we have to find the smallest n for which x 1 (n) ≤ α.In our example this means that we have to find the smallest value for which x 1 app (n) = 1.5764 • 0.6718 n − 0.5764 • 0.1024 n ≤ α. (13) This can be executed numerically, and the necessary amounts of initial material n for different values of are the values contained in Table 1.

Comparison of discrete and continuous models
Finally we compare the discrete and the continuous models in order to be able to replace one by the other if it is necessary.Our real world seems to be discrete rather than continuous in processing systems, but continuous models are better elaborated mathematically, as the theory of differential and integral equations is more widely discussed in the literature than that of the difference equations.
For the sake of simplicity let us consider the case when both the polynomial applied for defining the distribution of interarrival time and the distribution of the amount of material is p(y) = y + s.Now p(D)(f) ≡ 0 implies f(n + 1) + sf(n) = 0, n = 0,1,2,... under the condition s g with 0 < s g < 1 (s = s g ).After computations, the characteristic equation ( 10) can be written as Determining its roots at z = 1 we get Consequently, the results in [8] imply that φ(x,0) = exp(− c c µ λ − x).

Substituting c = 1 we get
To compare ( 15) and ( 17) numerically, we have to find the mappings between the parameters of the discrete and the continuous models.If the polynomials defining the distributions and probability density functions are equal, that is λ = s f and μ = s g , this does not imply the equality of expectations.As the behaviour of the model at z = 1, δ = 0 strongly depends on the expected inter-arrival time and the expectation of the filled amount of material, it is worth requiring the equality of these expectations.
In the continuous model, with . Supposing equal expectations in the continuous and discrete models, the following equalities have to be satisfied: and To compare the probability of shortage in the continuous and the discrete models we present Figs.5-9.Fig. 5 shows the distributions, Fig. 6 presents the probabilities of shortage as the function of the amount of the initial material.Fig. 7 shows the differences between the reliabilities in case of the discrete and continuous model, while Fig. 8 depicts the relative error.
The parameters are λ = 0.1, μ = 0.05,   The functions given by (20) and (21) are plotted in Fig. 9.The differences between the results of the discrete and continuous models represent that they can replace each other in case of appropriate values of parameters.

Summary
In this paper the operation of intermediate storages under stochastic operational conditions is investigated.Both the inter-arrival time distribution and the distribution of the filled amount of material are supposed to be discrete.The probability and the expected time of material shortage are handled with the help of an auxiliary function.A difference equation is set up for it and solved for a wide set of inter-arrival time distributions and a general class of the distribution of the filled amount.The results given by analytical solutions and Monte-Carlo simulation are compared and small differences are obtained.Approximate solutions are presented on the basis of the form of the

Fig. 1 .
Fig. 1.Change of the amount of material as a function of time (* U(m-), o U(m))

1 app
, for example s = 0.85, μ 1 (1) = 0.6760, μ 2 (1) = 0.0921, c 1 (1) = 1.5549 and c 2 (1) = -0.5549and x 1 (n) = 1.5549 • 0.6760 n − 0.5549 • 0.0921 n .On the other hand, we can determine the roots and constants by parameter fitting applying least squares method.Fitting a function of form (9) on the simulation points we get μ = -0.5764.Then the approximate formula for the probabilities is x 1 app (n) = 1.5764 • 0.6718 n − 0.5764 • 0.1024 n The shortage probabilities and the expectation of the shortage time can be investigated by Monte-Carlo simulation as well.The amount of material in the storage should be considered in those points when fillings happen.Generating random variables T s and Y s s=1,2,3,…, if the inequality 0

Fig. 4 .
Fig. 4. Approximate values of the expectation of shortage time (-) and the simulated values (*)

µ
To compare the expected time of shortage we have taken the derivatives with respect to δ and z.Omitting the details, in the continous model we finally obtain that while in the discrete one

Fig. 7 .Fig. 8 .Fig. 9 .Fig. 5 .Fig. 6 .
Fig. 7. Differences between the results of the continuous and the discrete models Necessary initial amount of material for fixed reliability levels