<!DOCTYPE article PUBLIC "-//NLM//DTD JATS (Z39.96) Journal Archiving and Interchange DTD v1.0 20120330//EN" "JATS-archivearticle1.dtd">
<article xmlns:xlink="http://www.w3.org/1999/xlink">
  <front>
    <journal-meta />
    <article-meta>
      <title-group>
        <article-title>Parallel Discrete Event Simulation as a Paradigm for Large Scale Modeling Experiments</article-title>
      </title-group>
      <contrib-group>
        <contrib contrib-type="author">
          <string-name>© Lev Shchur</string-name>
          <email>shchur@chg.ru</email>
          <xref ref-type="aff" rid="aff0">0</xref>
          <xref ref-type="aff" rid="aff1">1</xref>
        </contrib>
        <aff id="aff0">
          <label>0</label>
          <institution>Proceedings of the XVII International Conference «Data Analytics and Management in Data Intensive Domains» (DAMDID/RCDL'2015)</institution>
          ,
          <addr-line>Obninsk</addr-line>
          ,
          <country country="RU">Russia</country>
        </aff>
        <aff id="aff1">
          <label>1</label>
          <institution>Liudmila Shchur Science Center in Chernogolovka</institution>
          ,
          <addr-line>142432</addr-line>
          <country country="RU">Russia</country>
        </aff>
      </contrib-group>
      <fpage>107</fpage>
      <lpage>113</lpage>
      <abstract>
        <p>Parallel Discrete Event Simulation (PDES) is a method introduced to conduct simulation of a complex system, which consists of a large number of objects. The PDES has been known for about 30 years and we discuss why further analysis of the method is important nowadays. We briefly introduce the main features of the method and discuss the current state of PDES. We present a survey of the method as well as promising applications.</p>
      </abstract>
    </article-meta>
  </front>
  <body>
    <sec id="sec-1">
      <title>1 Introduction</title>
      <p>Methods of parallel algorithms developed in the
middle of the last century are weakly applicable on the
modern supercomputers. During the last two decades
the high-performance computing evolved from the big
mainframes with a complicated central processing unit
(CPU) into the parallel processing of large numbers of
CPUs and GPGPUs (general-purpose graphics
processing units). Running hundreds of thousands or
even millions of CPUs in a single task is one of the
challenging problems of high-performance computing
(HPC). We will focus our discussion on how the parallel
discrete event simulation (PDES) method may be used
in large-scale simulations.</p>
      <p>
        There are two ways to parallelize a simulation model:
time-parallel and space-parallel [
        <xref ref-type="bibr" rid="ref1">1</xref>
        ]. The time-parallel
approach partitions the simulated time axis into intervals
[T1,T2], [T2,T3], …, [Ti,Ti+1], … and assigns each
interval to a separate process. A process i computes a
sample path for the interval [Ti,Ti+1]. The final state of
the system in the sample path [Ti-1,Ti] must be the
initial state in the sample path for [Ti,Ti+1]. It is clear
that the time-parallel simulations are restricted to special
models that fulfil that requirement.
      </p>
      <p>The space-parallel approach partitions the system
being modelled into a collection of subsystems, and
assigns a logical process (LP) to each subsystem.
Therefore, each LP is just a sequential event-driven
simulator. Sending an event to a LP is equivalent to
scheduling that event in the receiving list of LP pending
events.</p>
      <p>
        The PDES is a space-parallel approach with an
execution of a single discrete event simulation program
on a parallel computer or on a cluster of computers [
        <xref ref-type="bibr" rid="ref1">1</xref>
        ].
It is widely used in economics and logistics, and
sometimes applied in physics and computer science. A
system that should be simulated is divided into disjoint
subsystems, which are mapped onto the programming
objects [
        <xref ref-type="bibr" rid="ref2">2</xref>
        ]. Subsystems are not isolated during the
simulations, and dependencies between objects should
be resolved properly. There are three main essentials of
PDES. First, it is assumed that changes of possible, or
potential, dependencies occur at some particular
moments of time. It is supposed that these moments of
time are spread on a scale that is large enough compared
to an elementary unit of simulation time. Therefore,
changes are considered to be discrete (although random)
in time, and are called discrete events. The objects
evolve independently in time, as soon as there are no
dependencies generated. We need to have some protocol
in order to process dependencies correctly, in other
words we would like to keep causality. Using general
language we can say that we need some protocol for
synchronization of discrete events.
      </p>
      <p>
        The main idea of PDES is to use an asynchronous
mechanism and to implement the concept of Virtual
Time [
        <xref ref-type="bibr" rid="ref2">2</xref>
        ]. A local virtual time variable i is associated
with each object i. As soon as a new time event is
generated by an object i, that event is stamped with the
local virtual time (LVT) i. A message containing
information on the event is generated and is sent to
other objects. A message sending mechanism is the
second essential feature of PDES. The third essential
feature of PDES - there is no information exchange
between subsystems through common variables and
there is no access of objects to a shared memory.
Instead, all the information is distributed via messages.
The ensemble of LVT i generates a profile of local
virtual times. A minimal value of the profile defines
time horizon, and an evaluation of this value defines
Global Virtual Time (GVT) of simulations. An evolution
of LVT profile and of GVT in time depends on a
particular scheme of managing causality, and the
analysis of time profile properties (which is related to
the problem of synchronization!) is one of subjects of
the present paper.
      </p>
      <p>
        Originally, the PDES was designed to deal with
complex and heavy problems. Suppose, one has to
simulate a complex system with the condition that
simulation is resource consuming (huge processor time,
big RAM, and big data storage, large number of
sensors) and could not be fit in a single computer. An
appropriate simulation method should be able to handle
any problem. The main requirement was not in the
effective use of CPU time, or storage, or other hardware
loading level. The main requirement was just to have
some possibility to simulate. It is not a usual
requirement of scalability, which is traditionally
associated with the effective use of hardware. Contrary,
scalability in the PDES is just a possibility to evolve a
system regardless of how large and complex this system
is. During the last years, the PDES method attracts
interest of a computational physics community, which is
looking for methods to simulate huge systems that can
be simulated on computers with millions of nodes
[
        <xref ref-type="bibr" rid="ref4 ref5 ref6">4,5,6</xref>
        ]. A discussion of that issue is another subject of
the paper.
      </p>
      <p>
        The third subject of the paper is a discussion of the
possibility to use the PDES method on a lower level, on
a level of system software that will link together billions
of nodes (computational units, CPUs, cores, sensors,
etc.) within one task. This was mentioned in the White
Paper “Performance Technologies for Peta-Scale
Systems” prepared by a group of researchers from US
National Laboratories and associated Universities [
        <xref ref-type="bibr" rid="ref7">7</xref>
        ].
      </p>
      <p>The paper is organized as follows. In Section 2 we
describe briefly the PDES method. In Section 3 we
review the research on the evolution of LVT time
profiles. In Section 4 we discuss an application of PDES
for in physics in informatics. In Section 5 we discuss a
possible application of the PDES method in system
software. In Section 6 we summarize our discussion and
point out possible future research.</p>
    </sec>
    <sec id="sec-2">
      <title>2 Time evolutions in Parallel Discrete Event</title>
    </sec>
    <sec id="sec-3">
      <title>Simulation method</title>
      <p>
        We discuss an evolution of the local time profile and
of the GVT in PDES algorithms [
        <xref ref-type="bibr" rid="ref8 ref9">8,9</xref>
        ] using simplified
models of that evolution. In a sense, it is a modelling of
the possible PDES algorithms. As we already
emphasized, the evolution of the LVT and the GVT is
associated with synchronization problems. The
importance of the model-of-model approach becomes
clear if we take into account that using that approach we
may prove, for example, that deadlocks are never
happens in conservative algorithms. In terms of LVT
profile time evolution, it means that the speed of the
time profile is always positive because there exists at
least one object with the LVT lower than LVTs of other
objects, therefore that object may proceed in time. In the
case of local algorithms it is possible to prove existence
of a lower bound of that speed [
        <xref ref-type="bibr" rid="ref10">10</xref>
        ].
      </p>
      <p>
        We use an extension [
        <xref ref-type="bibr" rid="ref11">11</xref>
        ] of Virtual Time concept on
the current multi-core and multi-thread architecture of
hardware and software. In our discussion, the LVT is
associated with logical processes (LP) and not with
processing elements, as it was the case in “before-core
era” of hardware [
        <xref ref-type="bibr" rid="ref9">9</xref>
        ]. A processing element (PE) is
associated with hardware, and a LP is more flexible
because it is the virtual. Nowadays this extension is
even more important with GPGPU in which threads
may run LP in parallel and concurrently.
      </p>
      <p>
        Let us define an abstract model of simulation in
PDES. The whole physical system is decomposed into
N subsystems. Each physical subsystem is associated
with a logical process. Each logical process develops in
time, changing its internal state at some moments of
local simulation time. Time of simulation is measured in
the Global Virtual Time [
        <xref ref-type="bibr" rid="ref2">2</xref>
        ], which is a minimal time of
the LVT profile. We call a change in the internal state of
simulated subsystem i as an event that happens at logical
process i, it is denoted LP(i). We suppose that intervals
of time between events are distributed according to a
Poisson law. A logical process generates a message with
information on the changes in the internal state, marked
up (stamped) by the value of LVT. Sequence of LVT
times i is not decreasing. Each logical process LP(i)
receives a message and may include the received
information in simulation, or ignore information if it is
not necessary for simulation. In the case when LPs
perform simulations in parallel, the following problem
may occur. It may happen that LVT j of LP(j) is higher
than i in the received message from LP(i), and LP(j)
does depend on the information contained in the
message generated by LP(i). Therefore, causality is
violated! Fujimoto classified PDES algorithms into two
groups [
        <xref ref-type="bibr" rid="ref1 ref9">1,9</xref>
        ]. They are called optimistic algorithms and
conservative algorithms reflecting the way in which
synchronization is performed. FaS algorithm is
introduced and described in the paper [
        <xref ref-type="bibr" rid="ref10">10</xref>
        ], it may be
viewed as representatives of the third group of PDES.
      </p>
    </sec>
    <sec id="sec-4">
      <title>2.1 Conservative PDES algorithm</title>
      <p>
        In the case of a conservative algorithm, a subsystem
waits for all events to happen in every subsystem from
each it depends on, before proceed in time. Clearly it is
possible only for the simulations in which we know in
advance all dependencies between subsystems in the
physical system (named as arcs on graphs of LPs).
Therefore, in a conservative algorithm LVT i of LP(i)
should not be greater than LVT k of all LPs from which
LP(i) depends on [
        <xref ref-type="bibr" rid="ref1">1</xref>
        ].
      </p>
      <p>In realizations, one has to define matrix of
dependencies M(i,j) with elements equal to 0 or 1. If an
element M(i,j) equals to 1, than LP(j) does may depend
on events generated by LP(i). This matrix may be
viewed as a matrix of sender-receiver for message
sending. For each LP(i) the corresponding column of a
matrix is defined, and each time event happens in
simulations of LP(i), it sends messages to all LP(j) for
which element of matrix M(i,j) is equal to one.
Accordingly, LP(j) will never go forward in simulation
with LVT j larger than minimal k in its queue of
messages stamped with the LVTs.</p>
      <p>In subsection 3.1 we will show that a conservative
algorithm never leads to a deadlock and will
demonstrate how simulation speed of GVT behaves for
some simple cases</p>
      <sec id="sec-4-1">
        <title>2.2 Optimistic PDES algorithm</title>
        <p>
          In an optimistic algorithm, all LPs are evaluated in
time with assumption that causality is fulfilled. Clearly,
LVT of subsystems increased independently, and it may
happen that LP(j) will receive message with time stamp
i smaller than l in the last processed message.
Causality is broken! LP(j) proceeds using the wrong
information, while ignoring an event which happened
earlier in simulation time, i &lt; l. An optimistic
algorithm introduces a protocol of rollback in time that
resolves the problem of causality. This is done again
using a message sending procedure. The most well
known optimistic protocol is Time Warp paradigm [
          <xref ref-type="bibr" rid="ref2">2</xref>
          ].
The event causing rollback is called a straggler. A
recovery is organized by undoing effects of all events
that have been processed prematurely by the LP
receiving the straggler.
        </p>
        <p>In addition to the Input Message Queue and Output
Message Queue, each LP keeps State Queue, which is a
list of the recently processed events. It gives a
possibility to restore an old state vector on rollback. LP
sends anti-message which annihilates the original one
when it reaches its destination LP. So, rollback
generates an anti-message (negative), which annihilates
while meeting the corresponding primary (positive)
message. This happens like annihilation in particle
physics, with an electron and a positron. If LP(i)
receives an anti-message that corresponds to a positive
message, which it has already processed, then the
process must be rolled back itself to undo the effect of
processing that positive message. In the case if both an
anti-message and a positive message are in LP(i) queue
they just annihilate within queue and do not cause any
rollback by LP(i). Recursive repetition of this procedure
allows cancelling all erroneous calculations. This
recursive procedure may lead to the avalanche in the
rollback, but the length of the avalanche is limited by
GVT. No event with time-stamp smaller than GVT will
ever be rolled back, so LPs with GVT may discard that
event from the State Queue.</p>
        <p>In subsection 3.2 we will describe a simple model of
evolution of LVT profile, and demonstrate an analogy
with a process known in physics as an unrestricted
surface growth.</p>
      </sec>
      <sec id="sec-4-2">
        <title>2.3 Freeze-and-Shift (FaS) algorithm</title>
        <p>
          To make the classification of the PDES algorithms
possible, two steps are performed. First, concept of
Virtual Time [
          <xref ref-type="bibr" rid="ref2">2</xref>
          ] was generalized, as LVT is associated
with the logical processes (LP) rather than with the
processing elements (PE). Second, it is assumed that it
is possible to handle several LPs within one PE, and that
communication within PE may be efficient using
conservative algorithm. In that case we can use
correspondence of time profile evolution in the
conservative PDES [
          <xref ref-type="bibr" rid="ref10">10</xref>
          ] with the evolution of surface
growth on substrate [
          <xref ref-type="bibr" rid="ref12">12</xref>
          ], described by
Kardar-ParisiZhang partial differential equation, named KPZ
equation. There is one-to-one correspondence between
classification of boundary conditions in KPZ equation
and classification of PDES algorithms [
          <xref ref-type="bibr" rid="ref11">11</xref>
          ]. In KPZ
equation, we may divide space into domains, and
classify possible boundary conditions (BC) between
domains. There are three possibilities of BC:
continuous, free, and fixed. Mapping of KPZ onto
PDES is as follows. Domain corresponds to PE, and BC
between domains corresponds to the communication
protocol. Therefore, with many LPs running within one
PE and updated using conservative algorithm, there are
three possible algorithms for communication protocol
between PEs (associated with corresponding BC):
conservative, optimistic, and FaS. FaS is abbreviation of
the Freeze-and-Shift algorithm.
        </p>
        <p>Let us consider the simplest parallel simulation task,
in which LPs are placed along the horizontal lines (see
Figure 1), and send time-stamped messages only to the
left and right LP. Let us assume we distribute l logical
processes on each of N processor elements.</p>
        <p>Fig. 1. Possible realization of the first two “steps” of the
time profile evolution in the conservative algorithm.
Number of LPs is equal to the number of PEs. Each LP
sends message only to the neighboring LP at the left and
right side.</p>
        <p>Figure 1 shows time profile of simulation with 12 LPs
on the 12 PEs. We start with all LVT i = 0. The system
evaluates to the second raw of black LPs. At that
“moment” GVT is defined by the process with the
minimal LVT time, it is process number 3. At the
“next” step of simulation only LPs with the minimal
local virtual time (number 1, 3, 6, 8, 10, and 12) are
allowed to proceed to the possible time marked by pink,
according with the conservative algorithm. The third
“step” consists in possible evaluation of those processes,
which are sitting in local minima of LVT profile; it is
processes with number 4, 7, 9, and 11.</p>
      </sec>
    </sec>
    <sec id="sec-5">
      <title>3 Evolution of time profile in PDES</title>
      <p>In this section we will discuss evolution of LVT time
profile in PDES. We do not discuss case studies in
which particular realization of PDES are investigated
empirically. Instead, we will discuss properties of
models which mimics essential features of PDES and
which can be investigated using methods of applied
mathematics and mathematical physics. First, we
discuss results of the research done in series of papers
by group of Korniss, mainly on the conservative
algorithm. Next, we discuss results of our research with
Mark Novotny, mainly for the optimistic algorithm.</p>
      <sec id="sec-5-1">
        <title>3.1 Conservative PDES and Kardar-Parisi-Zhang equation</title>
        <p>Let us consider the simplest case of application of
conservative algorithm to the chain of LPs
communicating with neighbors only as already
considered in Subsection 2.3 and depicted schematically
in the Figure 1. LVT i is associated with each LP(i). We
start with all i = 0 (LVT time profile is flat at t = 0,
where t is artificial time, and the discrete one for
simplicity1). Evolution of LVT profile may be written
as iterative process: i (t+1) = i (t)+ i(t) if i (t) min
{ i-1 (t), i+1 (t)}, and i (t+1) = i (t), otherwise. Here i(t)
are random variables. When LP(i) advanced in time it
sends messages to the right LP(i+1) and left LP(i-1)
stamped with the LVT i (t+1). In this process causality
is never violated.</p>
        <p>
          This algorithm is free of deadlock. This can be seen
qualitatively from the fact that statistically there are four
possible local configurations for LVT of LP(i), and only
one of them with the local minima. Therefore, average
speed of GVT (time horizon) should be 1/4. In fact,
numerical simulation of Korniss et al model gives
average speed of the time horizon is 0.246 410(7) [
          <xref ref-type="bibr" rid="ref10">10</xref>
          ]
for the Poisson distribution of i(t), 0.267(4) [
          <xref ref-type="bibr" rid="ref11">11</xref>
          ] for the
uniformly distributed i(t), and 0.258(5) [
          <xref ref-type="bibr" rid="ref11">11</xref>
          ] for the
Gaussian noise. So, average utilization of CPU time in
this simplified local version of conservative algorithm is
close to 25 per cent.
        </p>
        <p>
          Korniss et al argued that in the continuum limit (the
increment of the artificial time t is infinitesimal)
iteration process for the LVT time profile obeys the
equation introduced by Kardar-Parisi-Zhang (KPZ) for
the surface growth on the solid substrate [
          <xref ref-type="bibr" rid="ref10 ref11">10,11</xref>
          ]. This is
very important observation because it provides mapping
of conservative PDES algorithm with local message
exchange onto KPZ problem. In turn, KPZ problem is
well investigated [
          <xref ref-type="bibr" rid="ref12">12</xref>
          ], and solution demonstrates
universal properties of the profile fluctuations
depending on the time and system size.
        </p>
        <p>First, from the fluctuations of surface profile it is
possible to estimate the width growth of the profile.
Mapping this result from KPZ onto the properties of
LVT time profile in PDES shows us that width of LVT
grows with time as t3/2. This means that LPs are
unsynchronized with simulation time growth. This is a
bad news, but from KPZ analysis it follows, that width
growth is saturated at some moment of time. Saturation
level of LVT width is proportional to the square root
from the number of LPs. This is a good news.</p>
      </sec>
      <sec id="sec-5-2">
        <title>3.2 Optimistic PDES and unrestricted surface</title>
        <p>1 One should not confuse artificial time t we introduce for the sake of
simplicity with the time at which events happens and time-stamped
messages are generated.
growth</p>
        <p>
          Model for the evolution of time profile, which mimics
essential features of optimistic algorithm, was
introduced first in [
          <xref ref-type="bibr" rid="ref13">13</xref>
          ]. Let us consider again for the
sake of simplicity the linear chain of LPs. The first step
is the optimistic unrestricted evolution of LVT at which
simulation evolves forward in time and the second is the
rollback (or backward) algorithm of sending
antimessages. For this purpose one may introduce two
parameters, F and B associated with two steps of
optimistic scheme. Namely, at the first step one
evaluates the time horizon by updating time of the
randomly chosen LPs with value F following the
Poisson distribution, and every LP can be chosen with
the non-zero probability. Then one have to relax
(rollback!) LVT profile B times: for all LP(i) value of
LVT time changed to the value of LVT time of the
nearest left LP(i-1) or right LP(i+1). Value of B satisfies
the Poisson distribution. The average value of the speed
profile u(t) evolves proportional to (q-qc)a, with
q=1/(1+B), with a value of a=1.74 close to the critical
exponent of directed percolation and describing
roughening transition in one-dimensional unrestricted
growth process [
          <xref ref-type="bibr" rid="ref14">14</xref>
          ]. The value of qc=0.23(3)
corresponds to the value of the parameter B=3.3. It
should be mentioned that the width of the LVT profile
acts better for values q far away from the critical value
qc. When close to qc the system practically does not
evolve in time. The growth of width provides an
analogy with the roughening transition, which is in the
same universality class as a directed percolation [
          <xref ref-type="bibr" rid="ref14">14</xref>
          ].
Practically this means that if the length of avalanche
(see subsection 2.2) is larger than B=3.3 in average,
optimistic algorithm will not proceed in time but will be
rather stock. It is quite important to perform case studies
of optimistic algorithm in order to understand that
prediction in details.
        </p>
      </sec>
    </sec>
    <sec id="sec-6">
      <title>4 Using PDES in physics and informatics</title>
      <p>The most important feature of PDES for using in
physics is that PDES algorithms are scaled very
naturally with the physical system size (number of
objects, or logical processes) and with the hardware size
(number of nodes, cores, threads) with the only
requirement that the time for message distribution is
smaller than computation time. This is valid for the
simulation in which time to simulate object is much
larger than time to send message. It is clear that more
complex and hard simulations the more gain can be
reached.</p>
      <p>
        Last years a number of papers with the analysis of
implementation of PDES ideas in the large-scale
computing of models in physics have been published.
Most of discussions were done for the analysis of
application of PDES algorithm to the Ising model,
which is the simplest model of ferromagnetism, and
which is standard model to test algorithms in Monte
Carlo simulations [
        <xref ref-type="bibr" rid="ref15">15</xref>
        ]. Here we review some of the
papers that may be of general interest.
      </p>
      <p>
        SPPARKS is the project running by Sandia National
Laboratory and developing kinetic Monte Carlo
simulator [
        <xref ref-type="bibr" rid="ref16">16</xref>
        ]. SPPARKS runs on single processors or
in parallel using message-passing techniques and a
spatial-decomposition of the simulation domain. The
code is designed to be easy to modify or extend with
new functionality. Main idea is in partitioning of
simulation space into the computational domains, and
performing simulations in parallel within some time
window. In some sense, this idea is similar to FaS
algorithm. It implements several KMC solvers whose
serial computational complexity ranges from O(N) to
O(NlogN) to O(1) in the number of events N owned by a
processor. In a generic sense the solvers are catalog a
list of "events", each with an associated probability,
choose a single event to perform, and advance time by
the correct amount. Events may be chosen individually
at random, or by sweeping over sites in a more ordered
fashion. Problems they addressed are magnetic spin
models, surface growth on substrate, etc. One of the
simulations is connected with the microstructural
evolution during sintering [
        <xref ref-type="bibr" rid="ref17">17</xref>
        ]. The model is a grain
growth with physical effects of particle diffusion,
vacancy diffusion, and vacancy annihilation. This
sintering model is able to capture all the necessary
mechanism to simulate simple solid-state sintering
correctly. It was demonstrated by comparing it to the
three--dimensional, in-situ images taken in a
highenergy synchrotron during the sintering of Cu particles.
Similar project is running by Lawrence Livermore
National Laboratory, also for the parallel kinetic Monte
Carlo [
        <xref ref-type="bibr" rid="ref18 ref19">18,19</xref>
        ]. It is based on ideas of virtual time and
random order of simulating computational domains
within time window. It is tested on the model of billion
atoms. The approach was tested with the Ising model
(the simplest ferromagnetic model) and demonstrated
that scalability with the number of nodes is close to the
ideal one.
      </p>
      <p>
        A number of projects connected with simulation of
communication networks are using PDES. One of such
projects is running by INRIA, and it is connected with
the simulation of Border Gateway Protocol using
optimistic PDES [
        <xref ref-type="bibr" rid="ref20">20</xref>
        ] with number of nodes from 10 to
100 thousands. In the last case simulation lasts up to 36
hours on 4-core Intel Xeon W5580 3.2Ghz with 64GB
of RAM running 64-bit Fedora 12 on a Linux kernel
2.6.32.
      </p>
    </sec>
    <sec id="sec-7">
      <title>5 PDES and system software</title>
      <p>
        The goal of high end computing (HEC) initiative in
US and Europe is to deploy large-scale computing
platforms with hundreds thousands of nodes. In the
White Paper “Performance Technologies for Peta-Scale
Systems” the group of US researchers emphasized
importance of the research on the program tools and
software facilities, and in particular in system
simulation [
        <xref ref-type="bibr" rid="ref7">7</xref>
        ]. They see that as “an open-source
architectural simulation framework and API that enables
plug-and-play between separately-developed simulators
for different architectural features… and would also
enable zoom-out and zoom-in between
statisticallybased and cycle-accurate simulation techniques”.
Solution seen by them as “the simulations will be
decomposed into logical processes, and will be
synchronized by either conservative or optimistic
methods … as developed in PDES community”. The
reason is that by such approach the high degree of
computational parallelism in the simulation will
perfectly match the high degree of real parallelism of
HEC systems.
      </p>
      <p>There are several simulators available as open source
software.</p>
      <p>
        The most famous realization of Time Warp
synchronization algorithm is Rensselaer's Optimistic
Simulation System, named as ROSS simulator [
        <xref ref-type="bibr" rid="ref21">21</xref>
        ]. It is
public domain software, which can be installed on the
cluster. Recently group of researchers from Lawrence
Livermore National Laboratory, Rensselaer Polytechnic
Institute and University of Illinois at
UrbanaChampaign tested ROSS on the Blue Gene architecture
with successful use of almost 2 billions of cores [
        <xref ref-type="bibr" rid="ref22">22</xref>
        ]
yielding speed of 504 billion events per second.
Modification of optimistic algorithm named
shockresistant Time WARP (SRTW) designed by group in
Westminster University and analysis was performed
using the Grid’5000 platform [
        <xref ref-type="bibr" rid="ref23">23</xref>
        ].
      </p>
      <p>
        Group from Bologna University developed PDES
software named ErlangTW based on Time Warp
algorithm. Work was presented at 1st ACM SIGPLAN
workshop on Functional high-performance computing
(FHPC '12) [
        <xref ref-type="bibr" rid="ref24">24</xref>
        ], where some preliminary performance
results on multicore and distributed architectures using
the PHOLD benchmark.
      </p>
      <p>
        There are at least two attempts to implement PDES in
Cloud architecture. Fujimoto group introduced
TWSMIP (the Time Warp Straggler Message Identification
Protocol) protocol resigned for environment with shared
hardware resources for which it is known that traditional
Time Warp (TW) algorithm shows poor performance
[
        <xref ref-type="bibr" rid="ref25">25</xref>
        ]. TW-SMIP protocol defines dynamic
synchronization points for individual LPs based on
straggler messages (a new event with time stamp
smaller than other events it has already processed),
which improves efficiency of simulation.
      </p>
      <p>Second example is realization of PDES simulator on
Cloud/Virtual machine platforms developed in Oak
Ridge National Laboratory [26]. The scalability of
virtual time scheduler has been tested on 128 virtual
machines multiplexed on 32 cores, showing
improvement in the runtime relative to the default
Cloud/VM scheduler. Authors believe that algorithmic
design, observations, and results of their work are
timely for emerging cloud/VM-based installations,
highlighting the need for PDES-specific support in high
performance PDES on the cloud/VM platforms.</p>
    </sec>
    <sec id="sec-8">
      <title>6 Summary and discussion</title>
      <p>Parallel Discrete Event Simulation paradigm is one of
the paradigms for the large-scale high performance
computing. It is under extensive development by many
groups of researches. PDES was successfully tested on
the new architectures, from Blue Gene and conventional
clusters and CPUs to the Virtual Machine/Cloud
approach. Efficiency of PDES does depend on the
problem and it is promising to investigate properties of
synchronization, scalability and efficiency.</p>
    </sec>
    <sec id="sec-9">
      <title>7 Acknowledgments</title>
      <p>Russian Science Foundation supported this work
under the grant number 14-21-00158.</p>
    </sec>
  </body>
  <back>
    <ref-list>
      <ref id="ref1">
        <mixed-citation>
          [1]
          <string-name>
            <given-names>R.M.</given-names>
            <surname>Fujimoto</surname>
          </string-name>
          .
          <source>Parallel Discrete Event Simulation. Comm. of the ACM</source>
          ,
          <volume>33</volume>
          (
          <issue>10</issue>
          ), p.
          <fpage>30</fpage>
          -
          <lpage>53</lpage>
          ,
          <year>1990</year>
          .
        </mixed-citation>
      </ref>
      <ref id="ref2">
        <mixed-citation>
          [2]
          <string-name>
            <given-names>A.M.</given-names>
            <surname>Law</surname>
          </string-name>
          and
          <string-name>
            <given-names>W.D.</given-names>
            <surname>Kelton</surname>
          </string-name>
          .
          <article-title>Simulation modeling and analysis</article-title>
          .
          <source>McGraw-Hill. Third edition</source>
          .
          <year>2000</year>
          .
        </mixed-citation>
      </ref>
      <ref id="ref3">
        <mixed-citation>
          [3]
          <string-name>
            <given-names>D.R.</given-names>
            <surname>Jefferson</surname>
          </string-name>
          . Virtual Time.
          <source>ACM Trans. Programming Languages and Systems</source>
          ,
          <volume>7</volume>
          (
          <issue>3</issue>
          ), p.
          <fpage>404</fpage>
          -
          <lpage>425</lpage>
          ,
          <year>1985</year>
          .
        </mixed-citation>
      </ref>
      <ref id="ref4">
        <mixed-citation>
          [4]
          <string-name>
            <given-names>S.</given-names>
            <surname>Miller</surname>
          </string-name>
          and
          <string-name>
            <given-names>S.</given-names>
            <surname>Luding</surname>
          </string-name>
          .
          <article-title>Event-Driven Molecular Dynamics in Parallel</article-title>
          .
          <source>J. Comp. Phys</source>
          .
          <volume>193</volume>
          (
          <issue>1</issue>
          ), p.
          <fpage>306</fpage>
          -
          <lpage>316</lpage>
          ,
          <year>2004</year>
          .
        </mixed-citation>
      </ref>
      <ref id="ref5">
        <mixed-citation>
          [5]
          <string-name>
            <given-names>D.C.</given-names>
            <surname>Richardson</surname>
          </string-name>
          ,
          <string-name>
            <given-names>K.J.</given-names>
            <surname>Walsh</surname>
          </string-name>
          ,
          <string-name>
            <given-names>N.</given-names>
            <surname>Murdoch</surname>
          </string-name>
          , and
          <string-name>
            <surname>P. Michel P</surname>
          </string-name>
          .
          <article-title>Numerical simulations of granular dynamics: I. Hard-sphere discrete element method and tests</article-title>
          .
          <source>Icarus</source>
          <volume>212</volume>
          , p.
          <fpage>427</fpage>
          -
          <lpage>437</lpage>
          ,
          <year>2011</year>
          .
        </mixed-citation>
      </ref>
      <ref id="ref6">
        <mixed-citation>
          [6]
          <string-name>
            <given-names>J.P.</given-names>
            <surname>Nilmeier</surname>
          </string-name>
          ,
          <string-name>
            <given-names>J. Marian J. A Rigorous</given-names>
            <surname>Sequential</surname>
          </string-name>
          <article-title>Update Strategy for Parallel Kinetic Monte Carlo Simulation</article-title>
          . (accepted to Communications in Computer Physics.
          <year>2014</year>
          ) arXiv link: http://arxiv.org/pdf/1403.4674.pdf
        </mixed-citation>
      </ref>
      <ref id="ref7">
        <mixed-citation>
          [7]
          <string-name>
            <given-names>D.H.</given-names>
            <surname>Bailey</surname>
          </string-name>
          et al.
          <article-title>Performance Technologies for Peta-Scale Systems: A White Paper Prepared by the Performance Evaluation Research Center and Collaborators</article-title>
          . May 14
          <year>2003</year>
          .
        </mixed-citation>
      </ref>
      <ref id="ref8">
        <mixed-citation>
          [8]
          <string-name>
            <given-names>R.</given-names>
            <surname>Fujimoto</surname>
          </string-name>
          .
          <source>Parallel and Distributed Simulation Systems</source>
          . Wiley Interscience,
          <year>2000</year>
          .
        </mixed-citation>
      </ref>
      <ref id="ref9">
        <mixed-citation>
          [9]
          <string-name>
            <given-names>G.</given-names>
            <surname>Korniss</surname>
          </string-name>
          ,
          <string-name>
            <given-names>Z.</given-names>
            <surname>Toroczkai</surname>
          </string-name>
          ,
          <string-name>
            <given-names>M.A.</given-names>
            <surname>Novotny</surname>
          </string-name>
          , and
          <string-name>
            <given-names>P.A.</given-names>
            <surname>Rikvold</surname>
          </string-name>
          .
          <article-title>From Massively Parallel Algorithms and Fluctuating Time Horizons to Nonequilibrium Surface Growth</article-title>
          .
          <source>Phys. Rev. Lett.</source>
          ,
          <volume>84</volume>
          , p.
          <fpage>1351</fpage>
          -
          <lpage>1354</lpage>
          ,
          <year>2000</year>
          .
        </mixed-citation>
      </ref>
      <ref id="ref10">
        <mixed-citation>
          [10]
          <string-name>
            <given-names>L.N.</given-names>
            <surname>Shchur</surname>
          </string-name>
          and
          <string-name>
            <given-names>M.A.</given-names>
            <surname>Novotny</surname>
          </string-name>
          .
          <article-title>Evolution of time horizons in parallel and grid simulations</article-title>
          .
          <source>Phys. Rev. E</source>
          ,
          <volume>70</volume>
          ,
          <fpage>026703</fpage>
          ,
          <year>2004</year>
          .
        </mixed-citation>
      </ref>
      <ref id="ref11">
        <mixed-citation>
          [11]
          <string-name>
            <given-names>M.</given-names>
            <surname>Kardar</surname>
          </string-name>
          , G. Parisi, and
          <string-name>
            <given-names>Y.C.</given-names>
            <surname>Zhang</surname>
          </string-name>
          .
          <article-title>Dynamic Scaling of Growing Interfaces</article-title>
          .
          <source>Phys. Rev. Lett.</source>
          ,
          <volume>56</volume>
          , p.
          <fpage>889</fpage>
          -
          <lpage>892</lpage>
          ,
          <year>1986</year>
          .
        </mixed-citation>
      </ref>
      <ref id="ref12">
        <mixed-citation>
          [12]
          <string-name>
            <given-names>L.N.</given-names>
            <surname>Shchur</surname>
          </string-name>
          ,
          <string-name>
            <given-names>M.A.</given-names>
            <surname>Novotny</surname>
          </string-name>
          .
          <article-title>Evolution of time horizon in parallel discrete event simulations</article-title>
          .
          <source>Proceedings of Conference “Scientific services in Internet: multicore computer world”, AbrauDyurso</source>
          , p.
          <fpage>197</fpage>
          -
          <lpage>200</lpage>
          ,
          <year>2007</year>
          (in Russian).
        </mixed-citation>
      </ref>
      <ref id="ref13">
        <mixed-citation>
          [13]
          <string-name>
            <given-names>U.</given-names>
            <surname>Alon</surname>
          </string-name>
          ,
          <string-name>
            <given-names>M.R.</given-names>
            <surname>Evans</surname>
          </string-name>
          ,
          <string-name>
            <given-names>H.</given-names>
            <surname>Hinrinchsen</surname>
          </string-name>
          , and
          <string-name>
            <given-names>D.</given-names>
            <surname>Mukamel</surname>
          </string-name>
          .
          <article-title>Roughening Transition in a OneDimensional Growth Process</article-title>
          .
          <source>Phys. Rev. Lett.</source>
          ,
          <volume>76</volume>
          , p.
          <fpage>2746</fpage>
          -
          <lpage>2749</lpage>
          ,
          <year>1996</year>
          .
        </mixed-citation>
      </ref>
      <ref id="ref14">
        <mixed-citation>
          [14]
          <string-name>
            <given-names>D.P.</given-names>
            <surname>Landau</surname>
          </string-name>
          and
          <string-name>
            <given-names>K.</given-names>
            <surname>Binder</surname>
          </string-name>
          . A Guide to Monte Carlo Simulations in
          <source>Statistical Physics. 3d edition</source>
          . Cambridge: Cambridge University Press.
          <year>2013</year>
          .
        </mixed-citation>
      </ref>
      <ref id="ref15">
        <mixed-citation>
          [15]
          <string-name>
            <given-names>SPPARKS</given-names>
            <surname>Kinetic Monte Carlo</surname>
          </string-name>
          Simulator
        </mixed-citation>
      </ref>
      <ref id="ref16">
        <mixed-citation>[16] . Sandia National Laboratory. http://spparks.sandia.gov/index.html</mixed-citation>
      </ref>
      <ref id="ref17">
        <mixed-citation>
          [17]
          <string-name>
            <given-names>C. G.</given-names>
            <surname>Cardona</surname>
          </string-name>
          ,
          <string-name>
            <given-names>V.</given-names>
            <surname>Tikare</surname>
          </string-name>
          ,
          <string-name>
            <given-names>S. J.</given-names>
            <surname>Plimpton</surname>
          </string-name>
          .
          <source>Parallel Simulation of 3D Sintering. Int J Comp Materials Science and Surface Engineering</source>
          ,
          <volume>4</volume>
          , p.
          <fpage>37</fpage>
          -
          <lpage>54</lpage>
          ,
          <year>2011</year>
          .
        </mixed-citation>
      </ref>
      <ref id="ref18">
        <mixed-citation>
          [18]
          <string-name>
            <given-names>J.P.</given-names>
            <surname>Nimeier</surname>
          </string-name>
          and
          <string-name>
            <given-names>J.</given-names>
            <surname>Marian</surname>
          </string-name>
          .
          <article-title>A rigorous sequential update strategy for parallel kinetic Monte Carlo simulation</article-title>
          .
          <source>Computer Physics Communications</source>
          ,
          <volume>185</volume>
          , p.
          <fpage>2479</fpage>
          -
          <lpage>2486</lpage>
          ,
          <year>2014</year>
          .
        </mixed-citation>
      </ref>
      <ref id="ref19">
        <mixed-citation>
          [19]
          <string-name>
            <given-names>D.</given-names>
            <surname>Coudert</surname>
          </string-name>
          ,
          <string-name>
            <given-names>L.</given-names>
            <surname>Hogie</surname>
          </string-name>
          ,
          <string-name>
            <given-names>A.</given-names>
            <surname>Lancin</surname>
          </string-name>
          ,
          <string-name>
            <given-names>D.</given-names>
            <surname>Papadimitrou</surname>
          </string-name>
          ,
          <string-name>
            <given-names>S.</given-names>
            <surname>Perennes</surname>
          </string-name>
          , and
          <string-name>
            <given-names>Isaam</given-names>
            <surname>Tahiri</surname>
          </string-name>
          .
          <article-title>Feasibility study on distributed simulation of BGP</article-title>
          . Preprint arxiv:
          <volume>1304</volume>
          .
          <fpage>4750</fpage>
          .
        </mixed-citation>
      </ref>
      <ref id="ref20">
        <mixed-citation>
          [20]
          <string-name>
            <surname>ROSS -</surname>
          </string-name>
          <article-title>Rensselaer's Optimistic Simulation System</article-title>
          . http://www.cs.rpi.edu/~chrisc/ross.html
        </mixed-citation>
      </ref>
      <ref id="ref21">
        <mixed-citation>
          [21]
          <string-name>
            <surname>P.D. Barnes</surname>
            Jr,
            <given-names>D.R.</given-names>
          </string-name>
          <string-name>
            <surname>Jefferson</surname>
            ,
            <given-names>M.</given-names>
          </string-name>
          <string-name>
            <surname>Schordan</surname>
            ,
            <given-names>D.</given-names>
          </string-name>
          <string-name>
            <surname>Quinlan</surname>
            ,
            <given-names>C.D.</given-names>
          </string-name>
          <string-name>
            <surname>Carothers</surname>
            ,
            <given-names>L.V.</given-names>
          </string-name>
          <string-name>
            <surname>Kale</surname>
          </string-name>
          .
          <article-title>Extreme scale optimistic parallel discrete event simulation with dynamic load balancing</article-title>
          .
          <source>Proceedings of the 2014 Winter Simulation Conference</source>
          .
          <year>2014</year>
          .
        </mixed-citation>
      </ref>
      <ref id="ref22">
        <mixed-citation>
          [22]
          <string-name>
            <given-names>G.</given-names>
            <surname>Krafft</surname>
          </string-name>
          and
          <string-name>
            <given-names>V.</given-names>
            <surname>Getov</surname>
          </string-name>
          .
          <article-title>Transaction-oriented simulation in Ad Hoc Grids: design and experience</article-title>
          .
          <source>Procedings of HPCS</source>
          , pp.
          <fpage>38</fpage>
          -
          <lpage>44</lpage>
          ,
          <year>2008</year>
          .
        </mixed-citation>
      </ref>
      <ref id="ref23">
        <mixed-citation>
          [23]
          <string-name>
            <given-names>L.</given-names>
            <surname>Toscano</surname>
          </string-name>
          ,
          <string-name>
            <surname>G. D'Angelo</surname>
            , and
            <given-names>M.</given-names>
          </string-name>
          <string-name>
            <surname>Marzolla</surname>
          </string-name>
          .
          <article-title>Parallel discrete event simulation with Erlang</article-title>
          .
          <source>FHPC '12 Proceedings of the 1st ACM SIGPLAN workshop on Functional high-performance computing</source>
          . p.
          <fpage>83</fpage>
          -
          <lpage>92</lpage>
          ,
          <year>2012</year>
          .
        </mixed-citation>
      </ref>
      <ref id="ref24">
        <mixed-citation>
          [24]
          <string-name>
            <given-names>A.W.</given-names>
            <surname>Malik</surname>
          </string-name>
          ,
          <string-name>
            <given-names>A.J.</given-names>
            <surname>Park</surname>
          </string-name>
          ,
          <string-name>
            <given-names>R.M.</given-names>
            <surname>Fujimoto</surname>
          </string-name>
          .
          <article-title>An Optimistic Parallel Simulation Protocoal for Cloud Computing Environments</article-title>
          .
          <string-name>
            <surname>SCS M&amp;S Magazin</surname>
          </string-name>
          ,
          <volume>4</volume>
          , p.
          <fpage>1</fpage>
          -
          <lpage>9</lpage>
          ,
          <year>2010</year>
          .
        </mixed-citation>
      </ref>
      <ref id="ref25">
        <mixed-citation>
          [25]
          <string-name>
            <given-names>S.B.</given-names>
            <surname>Yoginati</surname>
          </string-name>
          and
          <string-name>
            <given-names>K.S.</given-names>
            <surname>Permulla</surname>
          </string-name>
          .
          <source>Efficient Parallel Discrete Event Simulation on Cloud/Virtual Machine platform. ACM Transactions on Modeling and Computer Simulation</source>
          .
          <year>2014</year>
          .
        </mixed-citation>
      </ref>
    </ref-list>
  </back>
</article>