=Paper= {{Paper |id=Vol-2278/paper04 |storemode=property |title=Classical Regeneration Based Approaches for Output Analysis of Multiserver Systems Simulation |pdfUrl=https://ceur-ws.org/Vol-2278/paper04.pdf |volume=Vol-2278 |authors=Irina Peshkova }} ==Classical Regeneration Based Approaches for Output Analysis of Multiserver Systems Simulation== https://ceur-ws.org/Vol-2278/paper04.pdf
     Classical Regeneration Based Approaches
    for Output Analysis of Multiserver Systems
                    Simulation

                                Irina Peshkova

                          Petrozavodsk State University,
                        Lenina Pr. 33, Petrozavodsk, Russia
                               iaminova@petrsu.ru
                              http://www.petrsu.ru



      Abstract. We present and discuss three types of regeneration that are
      based on classical regeneration for steady-state simulation and output
      analysis of multiserver systems. First method is artificial regeneration
      based on exponential splitting property of some distributions. Regen-
      erative envelopes allow to construct upper and lower bound for QoS
      characteristic of the system. Resampling of New Better than Used inter-
      arrivals or New Worse than Used service times accelerates simulation. All
      this techniques are reduced to classical regeneration and allow to receive
      shorter regeneration cycles.

      Keywords: Multiserver Systems · Regenerative Simulation· Artificial
      Regeneration · Splitting · Regenerative Envelopes · NBU distributions.


1   Introduction

Modern computer systems, as a rule, have a multiserver or multicore architec-
ture, which is caused both by the physical features of the modern technological
process of microchip production and by the need to create high-performance
computers suitable for solving computationally complex problems. Data trans-
fer to the users in modern networks also requires coordinated performance of
several devices, as well as work with several servers, for solving problems of re-
ducing delays and optimizing throughput. By this reason, reliable estimation
of Quality of Service (QoS) parameters of multiserver systems is an important
problem to study.
    Regeneration of stochastic process is one of the powerful and efficient meth-
ods of steady-state simulation and output analysis of complex queueing sys-
tems [3,5,8,7,6,11].
    Standard classical regeneration method uses visits of the process (describ-
ing the systems behavior) to a fixed return state as regeneration points. The
trajectory between such points can be divided into independent (for classical
regeneration) and identically distributed groups (in general, they can be weakly
dependent for wide-sense regeneration). It allows to apply well-developed classi-
cal statistical methods based on the special form of the Central Limit Theorem
(CLT) for confidence estimation of QoS parameters.
    However, in multiserver systems standard regeneration is not guaranteed by
stability condition of the system and regeneration points may be too rare to
be useful in practice or even they might not exist. We study an alternative
constructive method of obtaining the regeneration points.
    The so-called artificial regeneration uses the splitting property [2] and can
be applied for estimation of multiserver systems performance in two ways. First
method is based on constructing a new, equivalent to the original, system by
using the splitting property, which allows to change the internal structure of
continuous valued random variables in such a way to obtain a sequence of times
at which the distribution has some specific properties, and these times occur
with positive probability (possibly by enriching the probability space). Second
technique uses coupling method and the epochs when the process has a given
splitting distribution as regenerative times [1]. If the splitting is based on expo-
nential distribution, then the system at this particular points has memoryless
property, and thus, regenerates in classical sense in both cases.
    The regenerative envelopes method is reduced to replacing the original system
based on the coupling method [9] and the monotonicity properties of the studied
processes with a pair of minorant and majorant systems that regenerate in the
standard classical sense [10]. Then we can use upper (lower) bounds of confidence
intervals for QoS characteristic of majorant (minorant) system for the confidence
interval construction of the original system.
    We discuss the possibility of accelerating the simulation for multiserver sys-
tems in which service times and iterarrival times have New Worse than Used
distribution or New Better than Used distribution. If service times have New
Worse than Used distribution we resample service time each time when the re-
maining service time exceeds a fixed constant a. If interarrival times has New
Better than Used distribution we resample input each time when the remaining
interarrival time is less than a fixed constant b. In both cases we obtain more
frequent classical regeneration that allows to accelerate simulation.


2   Artificial regeneration based on splitting property in
    multiserver systems

The density f of continuous random variable (r.v.) T can be split if there exists
some 0 < p < 1 and density f0 , such that:

                                     f > pf0 .                                  (1)

We define a new Bernoulli random variable I called splitting indicator such that
P(I = 1) = p. We construct a new splitted r.v. T 0 (stochastically equivalent to
the r.v. T ) as follows:
                            T 0 = IT0 + (1 − I)T1 ,                          (2)
where T0 has density f0 , and T1 has density
                                            f − pf0
                                    f1 =            .                             (3)
                                             1−p
Now we can replace the original r.v. T with a triple [T0 , T1 , I].
    The particular case of splitting, the exponential splitting, could be efficiently
applied in simulation due to the memoryless property of exponential distribution.
    We say that a positive valued r.v. T with density f is exponentially splitted
if there exist constants λ > 0, τ0 > 0 and 0 < p < 1 such that
                          f (x) > pλe−λ(x−τ0 ) ,        x > τ0 .                  (4)
    That is, splitting representation T 0 of exponentially splitted r.v. T is a triple
(T0 , T1 , I) where T0 is left truncated exponentially distributed r.v. with trunca-
tion point τ0 , rate λ, density
                                      
                                        0,            x 6 τ0 ,
                             f0 (x) =
                                        λe−λ(x−τ0 ) , x > τ0 ,
and T1 has density
                                  ( f (x)
                                    1−p ,                x 6 τ0 ,
                       f1 (x) =     f (x)−pλe−λ(x−τ0 )
                                           1−p         , x > τ0 .

    Now we describe the application of exponential splitting technique to discrete
event simulation of multiserver systems.
    Consider a stochastic process Θ = {X(t), T(t)}t>0 , with discrete (vector)
component X = (X1 , . . . , Xn ) and continuous component T = (T1 , . . . , Tm ) > 0
(hereafter we denote vector-valued variables in bold and assume the dimension
is clear from context, if not sub-indexed explicitly). Moreover, we assume that T
decreases linearly with time, that is, T(t + δ) = T(t) − δ1, for δ > 0, and define
an event as a time epoch ti such that Tj ∗ (ti ) = 0 for some j ∗ ∈ {1, . . . , m} (in
this case we say that the ith event occurring at ti is of j ∗ th type). We assume
that the discrete component X is changed only at event epochs, e.g. by some
recurrent relation
                              X(ti +) = Gj ∗ (X(ti −)) ,                           (5)
where Gj ∗ is the recursion related to type j ∗ event, or, in general, assume the
change is governed by some Markovian probability
                   Pj ∗ (x, x0 ) = Pj ∗ {X(ti +) = x0 |X(ti −) = x}.
The component X(t) does not change if T(t) > 0, that formally means
                       X(t) = X(ti +), i = max{k : tk < t}.
We also assume that in the vector T the zeroed at ti − component Tj ∗ is initial-
ized at ti + from some distribution with density fj ∗ , that is,
           fj ∗ (u, x, x0 ) = P(Tj ∗ (ti +) ∈ du|X(ti −) = x, X(ti +) = x0 ).     (6)
    Now we assume, that if X(t) = x∗ for some specific x∗ , then the process
Θ allows to perform exponential splitting of the continuous component T, that
is for each continuous component Tj there exist constants τ0 (j), λ(j), p(j) such
that (4) holds for the density fj (u, x, x∗ ), with corresponding density f0,j (u) =
λj e−λj (u−τ0 (j)) and f1,j defined appropriately. At that, we introduce a split pro-
cess Θ0 = {X(t), T0 (t), Z(t)}t>0 enhanced with auxiliary phase variable Z ∈
{0, 1, 2}m denoting the phases of components of T defined componentwise in a
recursive manner as follows
                                         if X(ti +) = x∗ and Ij ∗ (ti +) = 1,
                        
                          1,
          Zj ∗ (ti +) =                                                           (7)
                          2, otherwise,
         Zk (ti +) = Zk (ti −), k 6= j ∗ ,                                       (8)

where ti is the epoch of type j ∗ event, and Ij ∗ is the splitting indicator corre-
sponding to the component Tj0∗ with success probability p(j ∗ ). However, starting
from Zj ∗ (ti +) = 1 at time ti , the component Zj ∗ makes a step between the
events:
                                     1, ti < t < ti + τ0 (j ∗ )
                                   
                        Zj ∗ (t) =                                              (9)
                                     0, t > ti + τ0 (j ∗ ),
that is, the phase Zj ∗ is changed from one to zero when not less than τ0 (j ∗ )
time passes since the epoch ti corresponding to the last event of j ∗ th type (we
call this time epoch a pseudo-event). In case Zj 6= 1, the component Zj does not
change during inter-event time.
    We define the artificial regeneration epochs as the (increasing) sequence of
pseudo-event epochs {βk }k>1

                  βk+1 = min {t > βk : X(t+) = x∗ , Z(t+) = 0} .                (10)

    This method of artificial regeneration could be applied to performance es-
timation of m-server system with energy efficiency control and single waiting
room. The discrete and continuous components of the simulated process Θ we
define in the following way. The discrete components τ0 (t) is the queue size at
time t, and Xi (t) is the mode of ith server, i = 1, . . . , m, X1 (t) = 1 if the
server is active, and X1 (t) = 0 if the server is in the sleep mode. Continuous
components: T0 (t) is the time before the next arrival epoch, and Ti (t) is the
remaining times of current activity (service/sleep) of ith server, i = 1, . . . , m.
For the regeneration epoch, we set x∗ = {k, 1, . . . , 1} for some k > 0 and recall
that artificial regeneration occurs at such pseudo-event epoch t that X(t+) = x∗
and all the exponentially split components T0 , . . . , Tm are in exponential phases,
that is, Z(t+) = 0.
    The second approach of defining the regeneration times is to construct stochas-
tically equivalent system with the same input, but sampling from either f0 or f1
at every random variable initialization time epoch. Then the original r.v. T is re-
placed by a stochastically equivalent splitting representation T 0 = IT0 +(1−I)T1 .
Coupling method allows to construct new artificial system, which is equivalent
to the original system, and by this reason has equivalent QoS parameters.
     To construct the artificial regeneration points we suggest the following method.
Sample the sequence of splitting indicators Ii , i > 1 as a Bernoulli 0-1 r.v. with
success probability p. Sample service times T̃j either from density f0 , or from
f1 , using relation (2). Fix integer ν0 (fixed number of customers in system) and
define
                          n
              βn+1 = inf k > βn : νk = ν0 6 m,
                                                            o
                           Ii = 1, i ∈ Mk , Mβn ∩ Mk = ∅ , n > 0,               (11)

where Mn = {i : ti 6 tn < tdi } is the set of the customers presented in new
system at arrival instant tn and tdi is the departure instant of the ith customer.
Note that at the constructed regeneration epoch, the remaining work possesses
memoryless property, since all the tasks being served have exponentially dis-
tributed service times. It is easy to see that Θn := {νn , Ij , j ∈ Mβn }, n > 1
is a Markov process, and that distribution of Θβn is independent of n > 1 and
pre-history {Θk , k < βn }.

3    Regenerative envelopes for multiserver systems
     simulation
Now we consider an infinite buffer FCFS m-server GI/G/m queueing system Σ,
with the renewal input with instants tn , the iid interarrival times Ln = tn+1 − tn
and the iid service times Sn , n > 1. Denote νn the number of customers at
instant t−                                                                   −
         n , Qn the number of customers waiting in the queue at instant tn , so
Qn = max(0, νn − m).
    The main idea of the regenerative envelopes method is the construction of
two new systems: a majorant system Σ and a minorant system Σ with the same
input as in the original system Σ, but (stochastically) smaller (in minorant
system) or larger (in majorant system) service times.
    We describe the construction of m-server majorant system Σ̃. The corre-
sponding variables in the system Σ̃ we endow with tildes. Assume that the
service time distributions in both systems are ordered as FS (x) 6 FS̃ (x), x > 0,
that is S >st S̃ (stochastically). Moreover, by construction, L =st L̃. Then the
coupling method allows to take t̃n = tn and Sn > S̃n w.p.1.
    Define the set Mn = {i : ti 6 tn < z i } of the customers, which are being
served in the system Σ at instant tn , S i (n) – the remaining service time of
customer i at instant t−n in Σ.
    Fix arbitrary integer ν0 > 0 (number of customers in the system), constants
0 6 a 6 b < ∞ and define

               n                                            \               o
    β n+1 = inf k > β n : ν k = ν0 , S i (k) ∈ (a, b), Mβ n   Mk = ∅, i ∈ Mk .
                                                                               (12)
    Then, at the (discrete) instant k = β n+1 , we replace all the non-zero remain-
ing times S i (k), i ∈ Mk , by the upper bound b.
     It is easy to see that Θn := {ν n , S i (n), i ∈ Mn }, n > 1, is a Markov pro-
cess and the distribution of Θβ n is independent of n and pre-history {Θk , k <
β n }, n > 1.
     The process {Θn , n > 1} classically regenerates at the instants {β k } and
has the i.i.d regeneration cycles {Θk , β n 6 k < β n+1 } with i.i.d cycle lengths
αn = β n+1 − β n , n > 1.
     In a minorant system Σ we define the regenerative moments β n , n > 0 in
a similar way (with possibly other values a, b, ν0 ). At each such instant β n we
replace the remaining times S i (k), i ∈ Mk , by the lower bound a.
     The systems Σ, Σ, Σ have zero initial state, and moreover (in an evident
notation) Ln = Ln = Ln , S n > Sn > Sn , n > 1. Then it follows from coupling
property that
                       ν n > νn > ν n , Qn > Qn > Qn , n > 1,                  (13)

and it allows us to use regenerative simulation of Σ and Σ to obtain an upper
and lower bounds of the mean stationary queue size (number of customers) in
original system.

                                  EQ 6 EQ 6 EQ.


4   Accelerated simulation of multiserver networks with
    New Better than Used distributions

In this section we will discuss how the properties of New Better (Worse) than
Used distributions of input or service times could be applied to accelerate regen-
eration simulation.
    We say that the distribution F of r.v. T has increasing (decreasing) failure
rate IFR (DFR) if

                                                          F (t + x) − F (t)
          Ft (x) = P (T 6 t + x|T > t) = P (T 6 x) =                                     (14)
                                                                F (t)

is increasing (decreasing) in t. Remark, that

                                Ft (x)   f (t)
                            lim        =       = r(t).                                   (15)
                            x→0   x      F (t)

                                                                                  Rt
                                                                              −        r(x)dx
   Function r(t) is called failure rate and it is easy to see that F (t) = e o .
Weibull and Gamma distributions with shape parameter α > 1 are IFR, with
α < 1 are DFR.
   It is known [4], that if ET = µ1 and F is IFR (DFR), then
                                          −t
                            F (t) > (6) e µ1 , t < µ1 .                                  (16)
    We say that the distribution F of r.v. T has increasing (decreasing) failure
rate in average IFRA (DFRA) if
                              Zt
                          1
                                   r(u)du increases in t.                   (17)
                          t
                              0

    Distribution F of r.v. T is called New Better (Worse) than Used NBU (NWU)
if F t 6 F or
                               F (t + x) 6 F (t) · F (x).                 (18)
The following relation between IFR, IFRA and NBU holds.

                              IFR ⇒ IFRA ⇒ NBU.

    Let S be the (generic) service time (with distribution B), L – (generic) in-
terarrival time (with distribution A) in m-server system Σ. Define r.v. St -
remaining service time at t− and B t (x) = P (S > t + x|S > t) = P (St > x), r.v.
Lt - remaining interarival time at t− and At = P (L > t + x|L > t) = P (Lt > x).
    If service times S have NBU distribution or interarrival times L have NWU
distribution in the original system Σ then we can apply coupling method to
construct a new system Σ̃ using property (14) of such distributions.
    If S has NWU (NBU) distribution, then B t > (6) B or

                          P (St > x) > (6) P (S > x).

Then we construct new system Σ̃ with the same input, but each time when
St > a = const we resample S with d.f. B. Using classical regeneration approach
for new system Σ̃ we obtain lower (upper) bound for QoS parameter of the
original system Σ (for example, queue size).
    If L has NBU (NWU) distribution, then At 6 (>) A or

                          P (Lt > x) 6 (>) P (L > x).

Then we construct new system Σ̃ with the same service times, but each time
when Lt 6 a = const we resample L with d.f. A. In this case we can use classical
regeneration approach for new system Σ̃ to obtain lower (upper) bound for QoS
parameter of the original system Σ.
    In the case of NWU service times or NBU interarrival times we obtain more
frequent classical regenerations that allows to accelerate simulation.


5   Conclusion
We present three techniques which allow to construct and speed-up classic regen-
erations in multiserver systems: artificial regeneration based on splitting prop-
erty, regenerative envelopes and accelerated simulation for systems with NWU
service times of NBU input Our future work will focus on the comparison of the
simulation results efficiency.
ACKNOWLEDGEMENTS

The publication was carried out under state order to the Karelian Research
Centre of the Russian Academy of Sciences (Institute of Applied Mathematical
Research KRC RAS). This research is partially supported by RFBR project 18-
07-00147 and by the Ministry of Education and Science of Russian Federation
(project no. 14.580.21.0009, unique identifier RFMEFI58017X0009).


References
 1. Andradóttir, S., Calvin, J.M., Glynn, P.W.: Accelerated Regeneration for Markov
    Chain Simulations. Probability in the Engineering and Informational Sciences
    9(04), 497–523 (1995). https://doi.org/10.1017/S0269964800004022
 2. Andronov, A.: Artificial regeneration points for stochastic simulation of complex
    systems. In: Simulation Technology: Science and Art. 10th European Simulation
    Symposium ESS’98. pp. 34–40. SCS, Delft, The Netherlands (1998)
 3. Asmussen, S.: Applied probability and queues. Springer, New York (2003)
 4. Barlow, R.E., Proschan, F.: Mathematical theory of reliability. With contributions
    by L. C. Hunter. The SIAM Series in Applied Mathematics, John Wiley and Sons
    (1965)
 5. Crane, M.A., Lemoine, A.J.: An Introduction to the Regenerative Method for
    Simulation Analys, Lecture Notes in Control and Information Sciences, vol. 4.
    Springer, Berlin, Heidelberg (1977). https://doi.org/10.1007/BFb0007339
 6. Foss, S.G., Kalashnikov, V.V.: Regeneration and renovation in queues. Queue-
    ing Systems 8(1), 211–223 (1991), http://link.springer.com/article/10.1007/
    BF02412251
 7. Glynn, P.: Wide-sense regeneration for Harris recurrent Markov pro-
    cesses: an open problem. Queueing Systems 68(3-4), 305–311 (2011).
    https://doi.org/10.1007/s11134-011-9238-x
 8. Glynn, P.W., Iglehart, D.L.: Simulation methods for queues: An overview. Queue-
    ing Systems 3(3), 221–255 (1988), http://link.springer.com/article/10.1007/
    BF01161216
 9. Morozov, E.: Coupling and stochastic monotonicity of queueing processes. PetrSU,
    Petrozavodsk (2013)
10. Morozov, E., Rumyantsev, A., Peshkova, I.: Monotonicity and stochastic bounds
    for simultaneous service multiserver systems. In: Ultra Modern Telecommuni-
    cations and Control Systems and Workshops (ICUMT), 2016 8th International
    Congress on. pp. 294–297. IEEE (2016), http://ieeexplore.ieee.org/abstract/
    document/7765374/
11. Sigman, K., Wolff, R.W.: A review of regenerative processes. SIAM Review 35(2),
    269–288 (1993)