# Introduction to simmer

#### Bart Smeets, Iñaki Ucar

#### 2019-03-10

## Basic usage

First, load the package and instantiate a new simulation environment.

    library(simmer)

    set.seed(42)

    env <- simmer("SuperDuperSim")
    env
    #> simmer environment: SuperDuperSim | now: 0 | next:
    #> { Monitor: in memory }

Set-up a simple trajectory. Let’s say we want to simulate an ambulatory consultation where a patient is first seen by a nurse for an intake, next by a doctor for the consultation and finally by administrative staff to schedule a follow-up appointment.

    patient <- trajectory("patients' path") %>%
      ## add an intake activity
      seize("nurse", 1) %>%
      timeout(function() rnorm(1, 15)) %>%
      release("nurse", 1) %>%
      ## add a consultation activity
      seize("doctor", 1) %>%
      timeout(function() rnorm(1, 20)) %>%
      release("doctor", 1) %>%
      ## add a planning activity
      seize("administration", 1) %>%
      timeout(function() rnorm(1, 5)) %>%
      release("administration", 1)

In this case, the argument of the `timeout` activity is a function, which is evaluated dynamically to produce an stochastic waiting time, but it could be a constant too. Apart from that, this function may be as complex as you need and may do whatever you want: interact with entities in your simulation model, get resources’ status, make decisions according to the latter…

Once the trajectory is known, you may attach arrivals to it and define the resources needed. In the example below, three types of resources are added: the nurse and administration resources, each one with a capacity of 1, and the doctor resource, with a capacity of 2. The last method adds a generator of arrivals (patients) following the trajectory `patient`. The time between patients is about 10 minutes (a Gaussian of `mean=10` and `sd=2`). (Note: returning a negative interarrival time at some point would stop the generator).

    env %>%
      add_resource("nurse", 1) %>%
      add_resource("doctor", 2) %>%
      add_resource("administration", 1) %>%
      add_generator("patient", patient, function() rnorm(1, 10, 2))
    #> simmer environment: SuperDuperSim | now: 0 | next: 0
    #> { Monitor: in memory }
    #> { Resource: nurse | monitored: TRUE | server status: 0(1) | queue status: 0(Inf) }
    #> { Resource: doctor | monitored: TRUE | server status: 0(2) | queue status: 0(Inf) }
    #> { Resource: administration | monitored: TRUE | server status: 0(1) | queue status: 0(Inf) }
    #> { Source: patient | monitored: 1 | n_generated: 0 }

The simulation is now ready for a test run; just let it simmer for a bit. Below, we specify that we want to limit the runtime to 80 time units using the `until` argument. After that, we verify the current simulation time (`now`) and when will be the next 3 events (`peek`).

    env %>%
      run(80) %>%
      now()
    #> [1] 80
    env %>% peek(3)
    #> [1] 81.08101 81.08101 81.37616

It is possible to run the simulation step by step, and such a method is chainable too.

    env %>%
      stepn() %>% # 1 step
      print() %>%
      stepn(3)    # 3 steps
    #> simmer environment: SuperDuperSim | now: 81.0810054927203 | next: 81.0810054927203
    #> { Monitor: in memory }
    #> { Resource: nurse | monitored: TRUE | server status: 1(1) | queue status: 2(Inf) }
    #> { Resource: doctor | monitored: TRUE | server status: 2(2) | queue status: 0(Inf) }
    #> { Resource: administration | monitored: TRUE | server status: 0(1) | queue status: 0(Inf) }
    #> { Source: patient | monitored: 1 | n_generated: 7 }
    #> simmer environment: SuperDuperSim | now: 81.3761633133839 | next: 81.3761633133839
    #> { Monitor: in memory }
    #> { Resource: nurse | monitored: TRUE | server status: 1(1) | queue status: 2(Inf) }
    #> { Resource: doctor | monitored: TRUE | server status: 1(2) | queue status: 0(Inf) }
    #> { Resource: administration | monitored: TRUE | server status: 0(1) | queue status: 0(Inf) }
    #> { Source: patient | monitored: 1 | n_generated: 8 }
    env %>% peek(Inf, verbose=TRUE)
    #>       time  process
    #> 1 81.37616 patient2
    #> 2 85.04420 patient4
    #> 3 86.20007 patient7
    #> 4 86.20007  patient
    #> 5 92.41640 patient3

Also, it is possible to resume the automatic execution simply by specifying a longer runtime. Below, we continue the execution until 120 time units.

    env %>%
      run(120) %>%
      now()
    #> [1] 120

You can also reset the simulation, flush all results, resources and generators, and restart from the beginning.

    env %>%
      reset() %>%
      run(80) %>%
      now()
    #> [1] 80

## Replication

It is very easy to replicate a simulation multiple times using standard R functions.

    envs <- lapply(1:100, function(i) {
      simmer("SuperDuperSim") %>%
        add_resource("nurse", 1) %>%
        add_resource("doctor", 2) %>%
        add_resource("administration", 1) %>%
        add_generator("patient", patient, function() rnorm(1, 10, 2)) %>%
        run(80)
    })

The advantage of the latter approach is that, if the individual replicas are heavy, it is straightforward to parallelise their execution (for instance, in the next example we use the function `mclapply` from the parallel) package. However, the external pointers to the C++ simmer core are no longer valid when the parallelised execution ends. Thus, it is necessary to extract the results for each thread at the end of the execution. This can be done with the helper function `wrap` as follows.

    library(parallel)

    envs <- mclapply(1:100, function(i) {
      simmer("SuperDuperSim") %>%
        add_resource("nurse", 1) %>%
        add_resource("doctor", 2) %>%
        add_resource("administration", 1) %>%
        add_generator("patient", patient, function() rnorm(1, 10, 2)) %>%
        run(80) %>%
        wrap()
    })

This helper function brings the simulation data back to R and makes it accessible through the same methods that would ordinarily be used for a `simmer` environment.

    envs[[1]] %>% get_n_generated("patient")
    #> [1] 7
    envs[[1]] %>% get_queue_count("doctor")
    #> [1] 0
    envs[[1]] %>% get_queue_size("doctor")
    #> [1] Inf
    envs %>%
      get_mon_resources() %>%
      head()
    #>   resource     time server queue capacity queue_size system limit
    #> 1    nurse 12.46836      1     0        1        Inf      1   Inf
    #> 2    nurse 24.86070      1     1        1        Inf      2   Inf
    #> 3    nurse 27.42167      1     0        1        Inf      1   Inf
    #> 4   doctor 27.42167      1     0        2        Inf      1   Inf
    #> 5    nurse 37.14626      1     1        1        Inf      2   Inf
    #> 6    nurse 42.88991      1     0        1        Inf      1   Inf
    #>   replication
    #> 1           1
    #> 2           1
    #> 3           1
    #> 4           1
    #> 5           1
    #> 6           1
    envs %>%
      get_mon_arrivals() %>%
      head()
    #>       name start_time end_time activity_time finished replication
    #> 1 patient0  12.468365 53.13489      40.66653     TRUE           1
    #> 2 patient1  24.860695 67.66651      40.24485     TRUE           1
    #> 3 patient0  11.152832 51.88635      40.73352     TRUE           2
    #> 4 patient1  19.444382 67.43141      40.14704     TRUE           2
    #> 5 patient0   9.470693 47.60155      38.13085     TRUE           3
    #> 6 patient1  25.735504 64.36031      38.46786     TRUE           3

Unfortunately, as the C++ simulation cores are destroyed, the downside of this kind of parallelization is that one cannot resume execution of the replicas.

## Basic visualisation tools

You may want to try the `simmer.plot` package, a plugin for `simmer` that provides some basic visualisation tools to help you take a quick glance at your simulation results or debug a trajectory object:

- Plotting simmer statistics
- Plotting simmer trajectories.