A Simulation Engine for Stochastic Timed Petri Nets and Application to Emergency Healthcare Systems

I. Introduction

BUSINESS process management (BPM) has been well accepted by industry as a discipline in operations management that studies the discovery, modeling, analysis, measurement, improvement, optimization, and automation of business processes. A business process is composed of a set of tasks that are performed by human stakeholders, software tools, machinery and/or equipment to achieve a concrete goal. BPM manages these tasks in a service company. Business processes are often repeated routinely with many instances running concurrently. Business processes span industries. Examples of service processes are financing, healthcare, travel, and banking processes.

Over the years, the emphasis of BPM has been shifting from the enactment of business processes towards the improvement of business processes [1]. Execution of a business process to all stakeholders’ satisfaction is not always easy. There are many challenges in reality. One such challenge arises from improper resource provisioning. In many service delivery systems, business processes are constrained by various available resources. The lack of sufficient resources can cause resource contention and delay the accomplishment of service goal. This is particularly true in an emergency response system where large quantities of resources, including emergency responders, ambulances, medical care personnel, medications, rooms, beds, etc., are required. Often potential delays can be avoided or reduced by using resource analysis to identify ways in which tasks can be executed in parallel and in efficient ways. A workflow with resource usage defined can help keep track of resource availability, disable the paths that are not executable, and present all executable paths, thus allowing the automatic selection of a feasible path for system execution [2].

Petri nets have been a popular formalism for modeling systems exhibiting behaviors of competition and concurrency for almost a half century [3], [4]. Petri nets are particularly powerful in the modeling and analysis of resources sharing systems [5], [6]. Stochastic timed Petri nets, an extension to regular Petri nets, are a powerful tool for system performance evaluation [7]. Some subclasses of STPNs can be analyzed numerically. For example, stochastic Petri nets (SPNs) [8], [9], in which each transition is associated with exponentially distributed firing times, can be mapped into continuous-time Markov chains. Generalized stochastic Petri nets (GSPN) [10], where each transition is associated with either exponentially distributed firing times, or a zero firing time, can be mapped into stochastic point processes. Extended stochastic Petri nets [11], [12] and deterministic and stochastic Petri nets (DSPNs) [13] can have exponential or non-exponential transition firing times, but non-exponential firing times can only be associated with exclusive or competitive transitions, in this way, these two types of stochastic timed Petri nets (STPNs) can be solved with semi-Markovian processes. When an arbitrary transition in an STPN is associated with general stochastic transition firing times, however, it would be very hard to solve it analytically due to inherited complexity with general stochastic processes. Unfortunately, this is the common case for real world service systems. For example, consider the emergency service in a hospital, where the processing time of many tasks, such as urine test and IV infusion, is close to a uniform distribution over an interval or normal distribution. Such kind of service systems also exhibit intensive concurrency, because multiple medical people work with multiple patients simultaneously. Thus, such a system cannot be specified by any of the above mentioned special STPN models. There are efforts reported that attempted to tackle systems with arbitrary event times [14]–[16], but because of the inherent computation complexity, those efforts didn’t generate practical applications. Therefore, we have to turn to software simulation tools [17].

A long list of Petri net software tools have been developed over the past three decades [18]. The majority of STPN tools are for specific classes of STPNs. For example, GreatSPN [19], which was initially developed at the University of Turin in the late eighties of last century, is mainly for the analysis of GSPNs. Later, a DSPN module was also developed and included in GreatSPN [20]. TimeNET [21] is also a tool for analyzing and simulating DSPNs. Many tools focus on other Petri net classes like SPNs, such as PIPE [22] and SPNP [23]; or like queuing Petri nets (QPNs), such as QPME [24]. Most software tools have their own graphical user interfaces (GUI) which support the design of Petri nets. However, most of them cannot be driven in the batch mode (only via the GUI), which does not ease the generation of multitudinous results.

In [25], a tool called TiPeNeSS (Timed Petri Net Simulator Software) is presented, which supports the simulation of timed Petri nets containing transitions with generally distributed firing delays. However, this tool does not take care of multi-server transitions in general case, i.e., not with exponential firing times.

This paper proposes a simulation framework that supports the analysis of optimal provisioning of resources that ensures timely care service and meanwhile incurs no excessive cost from unnecessary resource acquisitions. A tool based on this framework has been developed and tested. Notable features of this tool include:

1) Allow generic firing times for any transitions;

2) Support single-server and multi-server (infinite-server) transition models;

3) Support pre-selection and race transition firing policies;

4) Support step-by-step GUI-driven simulation mode and batch execution mode that facilitates large-scale simulations for system performance assessment;

5) Support resource provisioning, which is critical for the systems resource impact analysis.

6) STPN models are saved as JSON files.

A stochastic timed Petri net based healthcare workflow and resource modeling technique is presented. Healthcare systems are one of areas where Petri nets find successful applications [26]–[29].

The rest of this paper is arranged as follows: Section II introduces stochastic timed Petri nets, their firing policies, and in particular, multiple enabledness. Section III discusses the proposed simulation engine, including its architecture and algorithms. Section IV presents a case study that applies the simulation tool in the resource requirement analysis of an emergency healthcare system. The last section concludes the paper.

II. Timed Petri Nets

A Petri net is a graphical and mathematical modeling tool [3]. It consists of places, transitions, and arcs that connect places and transitions. Input arcs connect places with transitions, while output arcs start at a transition and end at a place. Places can contain tokens; the current state of the modeled system (the marking) is given by the number (and type if the tokens are distinguishable) of tokens in each place. Transitions are active components. They model activities which can occur (the transition fires), thus changing the state of the system (the marking of the Petri net). Transitions are only allowed to fire if they are enabled, which means that all the preconditions for the activity must be fulfilled (there are enough tokens available in the input places). When the transition fires, it removes tokens from its input places and adds some at all of its output places. The number of tokens removed/added depends on the cardinality of each arc.

A Petri net is mathematically defined as a 5-tuple

where P = {p₁, p₂, …, p_m} is a finite set of places; T = {t₁, t₂, …, t_n} is a finite set of transitions; I: P × T → N is an input function that defines the directed arcs from places to transitions, where N is the set of all nonnegative integers; O: T × P → N is an output function that defines the directed arcs from places to transitions; M: P → N is the initial marking.

Petri nets are a promising tool for describing and studying systems that are characterized as being concurrent, asynchronous, distributed, parallel, nondeterministic, and/or stochastic. As a graphical tool, Petri nets can be used as a visual-communication aid similar to flow charts, block diagrams, and networks. In addition, tokens are used in these nets to simulate the dynamic and concurrent activities of systems. As a mathematical tool, it is possible to set up state equations, algebraic equations, and other mathematical models governing the behavior of systems [3], [30]. The characteristics exhibited by the activities of information processing systems such as concurrency, decision making, and synchronization are modeled very effectively with Petri nets. When extended with inhibitors, Petri nets can also model priorities with events [3].

To study performance issues of systems, it is necessary to include a timing concept into the model. There are several possibilities to do this for a Petri net; however, the most common way is to associate a firing delay with each transition [7]. This delay specifies the time that the transition has to be enabled, before it can actually fire. If the delay is a random distribution function, the resulting net class is called STPNs [8], [9].

STPNs as a formalism developed in the field of computer science have gained wide applications for modeling system performance. The most common practice is to assume all firing delays are exponentially distributed. Moreover, there are a lot of software tools developed to implement and solve models defined by using the STPN formalism.

A Modeling Power of Regular Petri Nets

Sequential execution, concurrency, competition and synchronization are typical characteristics exhibited in service systems. They can all be modeled by Petri nets nicely.

A sequence of events can be modeled with a Petri net construction shown in Fig. 1, in which transitions are connected with paces sequentially to each other. When the place p₁ gets a token, the three transitions will fire one after another starting from the transition t₁.

Figure  1.  Sequential execution.

DownLoad: Full-Size Img PowerPoint

In Fig. 2, t₁, t₂, t₃ are three parallel transitions. After t₀ fires, one token in p₀ is removed and p₁, p₂, p₃ each is added with one token, and thus all t₁, t₂, t₃ are enabled and can fire independently.

Figure  2.  Concurrent execution.

DownLoad: Full-Size Img PowerPoint

In Fig. 3, all t₁, t₂, and t₃ are enabled. However, they compete for the same resource (token) in p₀. When one transition fires, the other two will be disabled. Each transition can be selected to fire, though, which is the decidability of Petri nets. When transitions are associated with random firing times, the transition that is sampled with the minimum firing time will be selected to fire.

Figure  3.  Competition.

DownLoad: Full-Size Img PowerPoint

In Fig. 4, t₃ models synchronization. t₃ is enabled only when p₃ and p₄ each has at least one token.

Figure  4.  Synchronization

DownLoad: Full-Size Img PowerPoint

B STPN Modeling

STPNs are Petri nets extended with random event times [7]. Event times can be associated with arcs, places or transitions. Most STPN models introduce times to transitions, called transition firing times. The time can be deterministic or stochastic. Stochastic transition firing time model is widely adopted because it covers general cases. STPN-based system performance evaluation starts with the modeling of the system under consideration with STPN. The STPN model has an underling process that governs the system behavior. This stochastic process is then analyzed using known analytical techniques or through software simulation.

Formally, an STPN is a tuple

where (P, T, I, O, M₀) is a regular Petri net; Φ = (Φ₁, Φ₂, …, Φ_n}, where Φ_i is the random distribution of the firing time of the transition t_i. We denote by τ_i the firing time of the transition t_i. Then

Transitions that take zero time to fire are called immediate transitions. Normally, immediate transitions are used to model events whose durations can be ignored comparing to other events in a system.

For a regular Petri nets, the set of all possible transition firing sequences starting from initial marking M₀ is denoted by L(M₀). The timed execution of an STPN can be specified by a sequence in L(M₀) augmented with a non-decreasing sequence of real nonnegative values representing the instants of transition firings, such that consecutive transitions t_i and t_i+1 correspond to the ordered firing epochs τ_i and τ_i+1, τ_i ≤ τ_i+1. That is, if $\sigma \in$ L(M₀) and $\sigma=t_{1}, t_{2}, \ldots, t_i , t_{i+1}, \ldots$.

Then the augmented sequence is

Assume the marking of the STPN after firing transition t_i is M_i. The time intervals [τ_i, τ_i+1) between consecutive epochs represent the periods in which the STPN sojourns in the marking M_i. A history of an STPN up to the kth epoch τ_k, or the sequence of σ’ up to τ_k, is denoted by Z(k). We assume that for all k, Z, and M, the following distribution functions can be uniquely determined

where the random variable d represents the time elapsed from the STPN’s entry to M up to the next transition firing epoch, i.e., the time difference τ_i+1 –τ_i.

In modeling customer arrival processes, a common assumption is that they follow Poisson distribution. Also, the customer service time is normally treated as exponentially distributed. The biggest advantage of the two assumptions is the memoryless property associated with Poisson process and exponential distribution. If in a marking all enabled transitions are associated with the exponentially distributed transition firing time, the above distribution function only depends on the current marking. That is

C Policy for Next Transition to be Fired

When multiple transitions are enabled in a marking, a policy is needed to decide which transition to fire next. This is particularly important when we investigate the performance of an STPN model through simulation. We consider two policies [13]:

Preselection: With this policy, transitions that are enabled in a marking are assigned fixed probabilities to fire. Let E(M) be the set of all transitions that are enabled in the marking M. When the preselection policy is applied, each enabled transition t_i is assigned a firing probability p_i(M, Z), and

Notice that such a probability distribution is defined at the marking level. In practice, it is difficult to implement, because one has to generate the entire state space of an STPN and then find proper markings to define the distribution. In this research, we propose to define the preselection policy at the net structure level. We make the following assumptions:

1) Transitions in a preselection set will be either all enabled or all disabled in any marking.

2) All transitions in preselection set are immediate transitions, that is, they all have 0 firing delay.

Assumption 1) makes sense because in a practical system, we do not expect a randomly formed group of events that need to enforce the preselection policy. Assumption 2) puts no restriction of the application of the preselection policy, because for any timed transition, we can expand it to two sequential transitions: one immediate transition and one timed transition. With these two assumptions

It significantly simplifies the practical use of the preselection policy.

At runtime, the selection of the next firing transition is performed based on some simple rule of sampling. For example, consider the three transitions in Fig. 3. Assume the firing probabilities of t₁, t₂, t₃ are 0.3, 0.5 and 0.2, respectively. When the place p₀ gets marked, they are all enabled. Our simulation tool will generate a random number between 0 and 1 to determine which transition to fire. A simple rule is, if we get a number in the range [0, 0.3), we fire t₁; If the number drops in [0.3, 0.8), we fire t₂; Otherwise, we fire t₃.

Race: Let E(M) be the set of all enabled transitions at M. In the race policy, for each transition t_i in the set E(M), a random sample τ is extracted from the joint distribution

The minimum of these samples determines which transition will actually fire, and the sojourn time in M. We only consider the case in which all random samples are independent. In this way, (6) is uniquely determined by the marginal distributions

D Multiple Enabledness and Multi-Servers

In an STPN, concurrency is not only exhibited when different transitions are enabled and can run in parallel, but also may be exposed by a single transition that can have multiple instances. For example, we can use a single transition to model three machines of the same type. If the transition is triple-enabled, i.e., the transition can fire three times continuously, then the transition actually works like three transitions. We call this transition a multi-server. On the other hand, if a transition only represents a single machine or processor or any other task performer, then the transition is called a single-server [31], [32].

In the single-server case, each transition can only have one instance and thus it can perform only one operation at one time. Fig. 5 illustrates the case that a single-server processes three inputs. There are three tokens in p₁. Assume t₁ is a single server. Then t₁ has to fire three times, each time it will take the assigned firing time to finish. That is, if the firing time is a fixed value, say 3 time units, and the transition starts firing at time 0, then the first token moves from p₁ to p₂ in three time units. The second token and third token move to p₂ in 6 time unites and 9 time units, respectively.

Figure  5.  t₁ is a single-server.

DownLoad: Full-Size Img PowerPoint

In the multiple-server case, each transition represents k servers, and thus it can perform k operations at one time. For example, assume the transition t₁ is a triple-server in Fig. 6 and there are three tokens in p₁ at time 0 in the initial marking. t₁ can fire and consume the three tokens at the same time. That is, all the three tokens will move from p₁ to p₂ in three time units. The triple-server transition here is equivalent to three single-server transitions.

Figure  6.  t₁ is a multi-server.

DownLoad: Full-Size Img PowerPoint

If a multi-server transition can be enabled up to k times simultaneously (the case we use a single transition to model up to k identical machines, servers, processors, etc.), there will be up to k instances of the transition created in a simulation tool. Each instance fires independently of other instances and the firing time of each instance follows the same random distribution as defined in (1). When all k instances are busy in firing, the transition cannot be further enabled. But, when one instance completes its firing, the transition can be enabled again. In the implementation of multi-server transitions, we put all enabled instances of each transition in a queue. We will explain that in the next section.

III. The Simulation Engine

In this section, we give an overview of the simulation engine for resource oriented timed Petri nets. As shown in Fig. 7, the simulation system presented in this paper has seven components:

Figure  7.  Architecture of the proposed simulation engine

DownLoad: Full-Size Img PowerPoint

1) A simulation engine. It is a stochastic timed Petri net based state transition system and supports non-Markovian transition firing times, race and preselection transition firing policies, and both single-server and multi-server transitions.

2) STPN editor. It is a graphic tool for users to construct STPN models.

3) An STPN bank. It stores all STPN models built with the tool in the JSON file format.

4) A resource configurator. It specifies and adjusts, manually or automatically, the resource allocation of a service system modeled by an STPN.

5) Algorithms. It provides different algorithms to the simulation engine. Each algorithm represents a simulation strategy. Currently there are single-server and multi-server two cases.

6) Event data. It stores event names and all timing parameters. Event data can be edited independently of the STPN structure.

7) Results. Simulation results are stored in the MS Excel spreadsheets.

A Class Diagram

The overall design of the proposed simulation engine for multiple server STPN models is illustrated in the class diagram in Fig. 8. Ten major classes and their functions are explained as follows.

Figure  8.  Class diagram of simulation engine of STPNs with multiple-server transitions

DownLoad: Full-Size Img PowerPoint

The PNModel class is used to import external files of STPN models and store the models in internal data structures.

The PNElement class defines the elements in a Petri net model. The Arc, Place and Transition inherit from the PNElement class for common attributes and methods. In the Arc class, we define two types of arcs, one is from place to transition, the other from transition to place. In the Transition class, we define a priority queue for each transition to save firing times of its instances. In the Place class, we set the number of tokens for each place.

The ResourceMgr class specifies resource type and resource quantity.

The RandomNumGenerator class generates the random number from different stochastic distributions. For an exponential distribution, we only need the average time to get a random number. For a uniform distribution, we need two numbers to define the minimum and maximum values,r espectively.

The TranFiringTimeMgr manages transition firing times in the priority queue for each transition.

The Petrinet class is to compute the fireable transition list for the Simulation class. In the Petrinet class, we define a global arraylist to save fireable transitions and update them in the Simulation class.

The Simulation class manages the operation of STPN models and saves the simulation results.

B Implementing Multi-Servers

Like in any typical discrete event simulation system, we use a global clock in our tool to record time. The clock time increases at each simulation step. The simulation step length can be fixed or dynamic. In the fixed step case, the step length is set to 1. In the dynamic step case, the step length is determined by the gap between the current clock time and the next transition’s absolute firing time. We implemented both options.

When a transition becomes enabled, the transition firing time is sampled from its distribution function. If it is a single-server, then at each step of the simulation, the global clock increases and the transition’s keeps firing time decrements by the length of the simulation step. When the firing time decreases to 0, the transition fires and the token distribution changes. Of course, not all enabled transitions will fire. If an enabled transition becomes disabled before its firing time decreases to 0, it cannot fire. This is dictated by the race firing policy. If before the transition fires, more tokens come to its input places that are sufficient to make the transition double or triple enabled, there will be no change to the transition’s firing behavior. After the transition is fired, it will be resampled for firing time at the new marking.

If the newly enabled transition is a multi-server, we will follow the three-phase procedure to handle its firing:

Phase 1: At the time the transition is enabled, tokens are removed from the incoming places. The transition firing time is sampled and inserted into a queue.

Phase 2: the firing time decreases when the global clock ticks until it goes down to 0,

Phase 3: tokens are added to the outgoing places.

If after Phase 1 and before Phase 3, the transition becomes enabled again, we follow the same three-phase procedure to handle the new instance. That is, we will maintain a queue that stores the remaining firing time of each instance of a multi-server transition.

C Performance Evaluation

Fig. 9 shows the external view of typical service systems. The transition A is a source transition, which models arrival of customers or, in general, service requests. The token in the place i represents a customer. Notice that A and i constitute a self-loop, modeling consecutive arrival of customers. The number of tokens in the sink place o counts the number of customers being served by the system. Fig. 10 models a single shot execution of the system. Transitions D₁, D₂, …, D_k are sink transitions, which model different ways that customers depart the service system. The STPN box models the service workflow. To validate the correctness (no deadlock, no dangling transitions, etc.), we can add a transition to the model in Fig. 10 and connect o to it as an input and also connect i to it as output, and then check to see if the new STPN processes all desired properties we want.

Figure  9.  Service system input and output.

DownLoad: Full-Size Img PowerPoint

Figure  10.  A single shot execution.

DownLoad: Full-Size Img PowerPoint

It is a widely accepted assumption that the customer arrival process is a Poisson process, and the interval time between customer arrivals is exponentially distributed. That is, the firing time of the transition A is exponentially distributed. Therefore, we can use a single parameter, the average rate of customer arrival, to specify the customer arrival behavior.

Assume we sample the exponential distribution of transition i N times, and the sequence of the samples we get is x₁, x₂, …, x_N, where x_i is the arrival time of the ith customer relative to the arrival time of the (i –1)th customer. Assume the system starts to run at time 0. Let a_i be the arrival time of the i th customer according to the global clock. Then

Let b_i be the departure time of the ith customer. The ith customer entering the system and the ith customer leaving the system may or may not be the same one. These times are recorded when any of the transitions D₁, D₂, …, D_k fires. The average sojourn time, s, of customers in the system can be calculated by

A customer’s sojourn time in a system includes service time and waiting time. The average service time can be obtained as the statistics of the sojourn time of a custom entering into the system in which no resources are occupied. In other words, we can run the model in Fig. 10 N times independently, starting at time 0, and record the time, c_i, that the place o is marked. By this way, the average service time is

Then the average waiting time of the system is

IV. Algorithms

In this section, we introduce the major algorithms used in our simulation engine for STPNs with multi-server transitions (We omit the algorithms for STPNs with only single-server STPNs). We will also use a toy example to show the difference between simulation results when different servers are used, which further indicates the importance to have multi-server implementation in a simulation tool.

A Algorithms

We present four major algorithm here. Algorithm 1 depicts the main control flow of the simulation engine. It also shows how a queue of the firing time of each instance of a transition is managed, which is the most critical part of the simulation engine. Algorithm 2 specifies how the list of enabled transitions is maintained and updated. Algorithm 3 is an abstract of the firing behaviors of the transition that follow the pre-selection policy. Algorithm 4 is about the data collection for system performance evaluation.

Selecting the length of simulation step (clock step) is tricky. A bigger simulation step would mean a faster simulation but a bigger simulation error. Simulation error is incurred due to the discretization of continuous transition firing time delays. As a rule of thumb, choosing a small simulation step length will allow a fine-grained sampling time for transition firing delays and reduce the error.

Algorithm 1: Multiple server engine

Input: A timed Petri net model

Output: Delay time for each customer

Initialization:

　transitions and timing parameters

　places and tokens

　arcs and weights

　a global clock (CLK)

　an enabled transition list (ETL)

　a transition instance queue (TIQ)

Main Flow:

　update ETL;

　while (ETL is not empty) {

　　traverse ETL to find the transitions-to-fire;

　　peek TIQ of transitions-to-fire to find the firing instance;

　　if (simulation step is fixed at 1)

　　　step ← 1;

　　else if variable step length

　　　step ← the firing time of transitions-to-fire;

　　add step to CLK;

　　for each transitions in ETL

　　　decrease by step the firing time of each instance in TIQ;

　　while (TIQ of transitions-to-fire not empty && head value is 0)

　　{

　　　for all outgoing arcs

　　　　update number of tokens of connected places;

　　　fire immediate transitions;

　　　fire pre-selection transitions.

　　　remove entries of 0’s from TIQ of transitions-to-fire;

　　　if transitions-to-fire is source/sink transitions

　　　　handle source/sink transitions;

　　　remove transitions with empty PIQ from ETL;

　　}

　　update ETL;

　}

　display results;

Algorithm 2: Update enabled transition list (ETL)

Input: the timed Petri net model

Output: an ETL with all enabled transitions

Initialization: an array to save all enabled transitions

Main Flow:

　for each transition t in the Petri net

　　if (t can fire && t not a pre-selection transition

　　　&& t not an immediate transition)

　　　add t into the array;

　for each transition t in the array {

　　if (t is still enabled)

　　　sample firing time for t;

　　else

　　　set it as null;

　}

　t ← transition with smallest firing time;

　for incoming arcs of t

　　remove adequate number of tokens from connected places;

　if (t is not in ETL)

　add t into ETL;

B Example

Consider a small canteen. Customers arrive at random times that follow the Poisson distribution. There are three kinds of foods served in the canteen: pizza, salad and sandwich, each served by one person. The food service time is uniformly distributed in the range of [2, 8] (minute) for all three types of foods. There are three check-out service stations. The check-out time is exponentially distributed with the average being 5 minutes for each station.

Algorithm 3: Fire pre-selection transitions

Input: pre-selection transitions and current marking

Output: new marking after firing

Initialization: firing probability

　for each pair of pre-selection trans

　　if they can fire

　　　p ← random number in [0, 1];

　　　if p smaller than preset probability of one transition

　　　　fire the transition;

　　　else

　　　　fire the other transition;

Algorithm 4: Handle source/sink transitions

Input: CLK, transitions_to_fire, number of inCustomer, number of outCustomer

Output: updated quantity for in and out customers

Main Flow:

　if (transitions-to-fire is a source transition)

　　inTime ← inTime + CLK;

　　inCustomer ++;

　if (transition-to-fire is a sink transition)

　　outTime ← outTime + CLK;

　　outCustomer ++;

The STPN model of this service system is illustrated in Fig. 11. We performed simulation on two options:

Figure  11.  STPN model of a canteen service.

DownLoad: Full-Size Img PowerPoint

1) The transition checkout is implemented as a single-server.

2) The transition checkout is implemented as a multi-server.

We want to show the difference between the two options.

The firing time of the transition arrive is exponentially distributed. Its average firing time, τ, which specifies the average arrival interval of customers, is a variable in the simulation. We simulated the cases when τ = 2, 3, 5, 7, 9, 12, 15, 20 min. For each case in each option, we simulated the service process for 5 hours. The simulation results are presented in Fig. 12.

Figure  12.  Simulation results of canteen service.

DownLoad: Full-Size Img PowerPoint

Customers may have to queue for food to be served and for checkout. For both single-server and multi-server implementations, with the arrival interval of customers increases, the total waiting time of a customer decreases. When the arrival interval is big enough such that customers arrive at a slower rate than food service and checkout combined together, the average waiting time approaches 0. This trend is clearly shown in Fig. 12.

More importantly, Fig. 12 clearly shows the difference between the simulation result when the transition checkout is implemented as a single-server and that when it is implemented as a multi-server. That is, the single-server implementation results in a longer waiting time, because the single-server implementation does not reflect the concurrency among the three checkout stations.

Notice that with the increase of customer arrival intervals, the system becomes less congested and the difference between the two implementations become less dramatic. When the customer arrival interval is sufficiently large, there will be no difference between using one cashier and using more cashiers, and thus the simulation results of the two implementations will be statistically the same.

V. Case Study

In this section, we present a case study on healthcare emergency services. A health emergency is a situation that poses an immediate risk to somebody’s health and even life, and requires an urgent medical intervention to prevent things from going worse. This case study aims to find optimal resource provisioning in an emergency department (ED) regarding the treatment of patients who suffer from strokes. By optimal we mean the quantity of a resource, say doctors, is adequate such that it is sufficient to ensure timely service to patient flow and meanwhile it is not over allocated to incur any unnecessary cost to the ED.

A The Workflow

Heart attack, stroke and arrhythmia are common heart and blood vessel problems. Although each problem has a unique clinical workflow, they do share some common tasks, such as CT (computerized tomography) scan, X-Ray test, blood test, etc. As an example, we model the clinic workflow of stroke treatment with STPN and study the related resource impact. The workflow handling a patient who is experiencing symptoms of a stroke is as follows: As soon as the patient arrives at the ED, he and she is taken to the point-of-care desk where he or she gets registered by a registrar. A registered nurse asks some preliminary questions, takes his or her vitals, assigns the patient an ESI, or emergency severity index, and orders CT scan. The patient takes a seat in the waiting room. Then a CT tech wheels him off for diagnostic test. After the test, the patient is returned the waiting room and waits for initial assessment by a doctor.

A doctor asks the patient questions about his previous history, signs and symptoms, medications, emotional state, drug abuse screen, etc. Then according to all information he has collected, including the CT scan result, the doctor assesses the situation of patient and decides if the patient needs more care at the ED. If it is not the case, the patient is disposed and told to follow up with his regular physician. Otherwise, he is assigned a bed and over there is gowned, hooked up with EKG (electrocardiogram) and has an IV (infusion intravenous) in place under the help of a nurse. A phlebotomy tech takes his blood sample and then conducts lab test. Sometime later, the doctor comes with the lab results and EKG readings and decides if the patient needs more emergency care. If not, the patient is discharged. Otherwise, the patient has to undergo another round of treatment with an adjusted medication.

All major tasks involved in the clinical workflow described above are listed in Table I. The column Res.Required lists resources that are required to execute the corresponding task, while the column Res.Released lists resources that are released when the task execution is finished. This table clearly shows the shared resources among tasks. For example, doctors as a resource are shared among several tasks, namely tasks 9, 19, 20, and 22. It also shows that an examination bed is occupied after executing task 10 but won’t be released until the execution of task 20.

Table  I.  Task Specification

	Tasks	Abbr.	Res. Required	Res. Released	Time Distr.
1	Registration	RG	Registrar	Registrar	Exp. 10
2	Take Vital	TV	Nursing Assistant	Nursing Assistant	Uni. 5–15
3	Initial Assessment	IA	Physician Assistant	Physician Assistant	Exp. 15
4	Chest X-ray	CXR	X-ray Technician, X-ray Machine	X-ray Technician, X-ray Machine	Exp. 25
5	CT Scan	CT	CT Tech, CT Machine	CT Tech, CT Machine	Exp. 25
6	Triage	TR	Triage Nurse	Triage Nurse	Exp. 17
7	Take Urine	TU	Urine Test Nurse	Urine Test Nurse	Exp. 8
8	Urine Test	UR	Specimen Cup, Urine Test Technician, Urine Test Machine	Urine Test Technician, Urine Test Machine	Uni. 25–35
9	Reassessment	RA	Doctor	Doctor	Exp. 20
10	Bed Assignment	BA	Nurse, bed	Nurse	Exp. 18
11	Gown	GN	Nurse, gown	Nurse	Exp. 12
12	EKG Hookup	EH	Nurse, EKG Set	Nurse	Exp. 28.5
13	Infusion Intravenous	IV	Nurse, IV Bag	Nurse	Uni. 140–160
14	IV Infusion Start	IS	Nurse	Nurse	Uni. 1–5
15	IV Infusion End	IE	Nurse	Nurse	Uni. 1–5
16	Draw Blood	DB	Nurse, Needle, 3 Tubes	Nurse	Exp. 28
17	EKG Test	EKG	Nurse	Nurse	Exp. 12
18	Blood Test	BT	Phlebotomy Tech	Phlebotomy Tech	Uni. 50–70
19	Assessment	AM	Doctor	Doctor	Exp. 20
20	Disposition	DP	Doctor	Doctor, EKG Set, Grown, Bed	Exp. 10
21	Adjust Drug	AD	Nurse, Drug	Nurse	Exp. 32
22	Early Disposition	ED	Doctor	Doctor	Exp. 10

 | Show Table

DownLoad: CSV

The last column of the table shows the random distribution of each task’s duration. For example, the duration of the task Registration follows an exponential distribution and its average time is 25 minutes. The task IV Infusion takes a time which is uniformly distributed in the range of [140, 160].

B STPN Modeling

The STPN model of the clinic workflow is illustrated in Fig. 13. Each task in Table I has a corresponding transition in this model. In addition, there are also transitions that are used to model different assessment results. For example, the transitions d₁₁ and d₁₂ are a pair. d₁₁ represents that the doctor’s assessment result indicates the patient can have an early departure, while d₁₂ represents that the doctor’s assessment result indicates the patient needs to have an examination bed assigned for further treatment. d₁₁ and d₁₂ are immediate transitions and they follow the pre-selection firing policies. So are d₂₁ and d₂₂. The unnamed place and transition following EH are introduced solely for readability purpose. They do not alter the logic of the workflow.

Figure  13.  Stochastic timed Petri net model of the workflow of stroke care.

DownLoad: Full-Size Img PowerPoint

Tasks RG, TR, TV and IA are executed sequentially for each patient, and thus in the STPN model they form a sequential structure. Tasks CT, CXR and TU can be executed in an arbitrary order. Therefore, these three tasks plus UR form a random execution structure. But different from a decision structure in which only one transition will be executed, here all these transitions have to be executed. Tasks BA, GN and EH are also executed sequentially. After that the patient receives IV and undertakes other tests. This process may be repeated upon a doctor’s decision. Therefore, there is a loop structure in the model, marked by transitions d₂₂ and AD on the feedback path.

Transition firing times are specified in the right-most column. Two types of firing time distributions are adopted, exponential distribution and uniform distribution, just to show our tool can handle non-Markovian STPNs.

Resources are modeled by the places P_doctor, P_register, P_bed, and P_nurse. All transitions with one or more of these places as input must be implemented as multi-servers. Otherwise, the concurrency among multiple processors which utilize the same resource will be lost.

C Resource Provisioning

In this case study, we only consider reusable resources: doctors, nurses, registrars and beds. Non-reusable resources can only be consumed by a single transition and are not shared, and thus not studied here. For each set of resource provisioning, we run the STPN model with various patient arrival rates and collect the average patient waiting time in the system. We stop the simulation when further increasing the quantity of a resource does not reduce the patient waiting time. According to the cost of resources, we follow the order doctor→nurse→registrar→bed to determine a proper resource quantity. In the simulation, we first set the quantities of nurses, registrars and beds to be unlimited, and slightly increase the number of doctors. Then, we need to compare the average waiting time of patients with different provision for doctors to determine the most reasonable doctor number. Second, after we get the reasonable number of doctors, we fix the number for doctors and set the other two resources unlimited to determine the adequate number of nurses. At last, we can get all the resources provisioned properly.

D Results

Fig. 14 shows the average waiting time of patients vs. different numbers of available nurses when all transitions are implemented as single-servers. As a contrast, Fig. 15 shows the average waiting time of patients vs. different numbers of available nurses when all related transitions are implemented as multi-servers as they should be. One can observe that in any cases in terms of the number of available nurses and average arrival interval, patients have longer waiting times in the single-server implementation than in the multi-server implementation. This difference is more dramatic when the patient arrival rate is high, because it causes more congestion in a single-server implementation.

Figure  14.  Patient waiting time vs. # of available nurses when all transitions are implemented as single-servers.

DownLoad: Full-Size Img PowerPoint

Figure  15.  Patient waiting time vs. # of available nurses when all transitions are implemented as multi-servers.

DownLoad: Full-Size Img PowerPoint

It is worth pointing out that in the single-server implementation, when we increase the number of nurses, the average patient waiting time is still going down. This is because there are several concurrent events in the system performed by nurses.

Fig. 16 shows that after the number of doctors is increased over 5, providing more doctors does not result in any obvious improvement on the patient waiting time. Therefore, we can reasonably say 5 doctors are enough for this service settings. The same approach is applied to decide the optimal provisioning of nurses, registrars and beds. At last, we get the relatively optimal numbers of resources as 5 doctors, 12 nurses, 20 beds and 5 registrars. Simulation results for nurses, registrars and beds are illustrated in Figs. 15, 17 and 18.

Figure  16.  Patient waiting time vs. # of available doctors.

DownLoad: Full-Size Img PowerPoint

Figure  17.  Patient waiting time vs. # of available beds.

DownLoad: Full-Size Img PowerPoint

Figure  18.  Patient waiting time vs. # of available registrars.

DownLoad: Full-Size Img PowerPoint

E Discussion

This case study is only to show an application of the proposed STPN simulation tool. The order of resource provisioning that we used in this case study is just one option. There are many different simulation strategies that could help find proper quantities of resources. For example, one can first find the best match among the quantities of these resources, then use the simulation to decide the quantity of each. In fact, this is an optimization problem with multi-variables.

Of course, we can also model the availability or more types of resources in the emergency care systems. Examples are the capacity of the waiting room, number of stretchers, number of CT machines, etc.

VI. Conclusion

Resource availability plays critical roles in service providing systems. This paper proposed a simulation framework that supports the analysis of optimal provisioning of resources that ensures timely service and meanwhile incurs no excessive cost from unnecessary resource acquisitions. The contribution of the work is two-fold: First, it implemented a simulation system that is based on stochastic timed Petri nets and supports non-Markovian transition firing times, race and preselection transition firing policies, and both single-server and multi-server transitions. An important feature of the simulation system is that it allows a systematic system performance vs. resource provisioning analysis. Second, a meaningful case study on an emergency healthcare system is developed to show the power of the simulation system.

In the current implementation of the STPN tool, we assume that all transitions work in the enabling memory mode. The distributions depend on the history Z through the concept of an enabling age variable associated with each transition. The age variable accounts for the work done by the activity corresponding to the transition from the last instant at which it has become enabled up to the considered epoch. This is the most applicable memory mode in real world systems. There are other two memory modes though, resampling and age memory. There are occasions that these two modes need to be considered [7]. In our future work, we will implement them in our tool. From the application perspective, right now we have to manually provide the quantity of resources and check out the corresponding system performance. In our future work, we will consider the automated optimal resource provisioning based on user input criteria and all constraints that a system is subject to. Of course, we will also try to extend applications of the simulation tool to more time-critical systems, such as those mentioned in [33].

It will also be interesting to explore the use of statistical model checkers to validate the desired properties of the STPN model of a service delivery system of interest, especially when many details of the system are considered in the model [34].