Quick Start

The figure below illustrates the Dymodetron user workflow. There are two loops oriented around the Build Model activity.

%%{init: { 'theme': 'neutral' } }%% graph TB subgraph t1 [" "] direction TB A(Build Model) end subgraph s1 [" "] direction TB C(Generate Model Diagrams) F(View Model Diagrams) A(Build Model) -- model source --> C(Generate Model Diagrams) C(Generate Model Diagrams) -- model diagrams --> F(View Model Diagrams) F(View Model Diagrams) -- Iterate On Model Definition --> A(Build Model) end subgraph u1 [" "] direction TB B(Generate Model Simulation Code) D(Run Simulations) E(Analyze Results) A(Build Model) -- model source --> B(Generate Model Simulation Code) B(Generate Model Simulation Code) -- generated sim --> D(Run Simulations) D(Run Simulations) -- sim results --> E(Analyze Results) E(Analyze Results) -- Iterate On Model Definition --> A(Build Model) E(Analyze Results) -- Iterate On Model Parameters --> D(Run Simulations) end classDef class1 fill: white, stroke: black classDef class2 fill: white, stroke: white class t1 class2 class s1,u1 class1 style A fill:#aca,stroke:#333,stroke-width:4px

In this section, we’ll take a pre-built example model, and do the following:

Install Dymodetron

$ pip install git+https://github.com/InstituteforDiseaseModeling/dymodetron.git#egg=dymodetron

If you’ve already installed an earlier version, you can upgrade to a more recent version by specifying the version number at the end, as follows:

$ pip install git+https://github.com/InstituteforDiseaseModeling/dymodetron.git#egg=dymodetron==0.2.11

Generate Model Diagrams

This bash command generates diagrams for the model below.

$ python -m dymodetron.generators.state_machine_diagrams --model_definition_file=examples/mm1_queue.py --overwrite_existing=1
state_machine_diagrams
Generating state machine diagram for state machine [mm1_queue_arrivals]
Outputting to [generated/state_machine_diagrams/mm1_queue_arrivals.html]
Generating state machine diagram for state machine [mm1_queue_service]
Outputting to [generated/state_machine_diagrams/mm1_queue_service.html]

View Model Diagrams

The command above listed out the locations of the generated diagram files. Take note of the folder that the files are generated into. Open the listed .html files in a browser.

This example model has two state machines, each one gets its own diagram file:

mm1_queue_arrivals.html:

%%{init: { 'theme': 'base' } }%% stateDiagram-v2 arriving idle [*] --> idle : initial_event idle --> arriving : arrival_event arriving --> idle : done_arriving_event

mm1_queue_service.html:

%%{init: { 'theme': 'base' } }%% stateDiagram-v2 idle servicing [*] --> idle : initial_event idle --> servicing : service_event servicing --> idle : done_servicing_event

The diagrams provide a good picture of the “bones” of the model, but they don’t provide all the details. For that you have to look at the model definition. It’s often useful to have the model diagrams open side-by-side with the model definition, and look at both together.

Generate Model Simulation Code

This bash command generates simulation code for the same model.

$ python -m dymodetron.generators.python_arrayified --model_definition_file=examples/mm1_queue.py --overwrite_existing=1
python_arrayified
Generating code for entity type [MM1_Queue]
Outputting to [generated/python_arrayified/mm1_queue.py]

The command above listed the path of the generated simulation code. Take note of that path, that is the code you will use next to run the simulation.

Run Simulation

The following command illustrates how to run the simulation with two simulated entities for 100 time ‘ticks’.

Note

Dymodetron measures time in units of ‘ticks’. Time ‘ticks’ are whatever units of time you want them to be. Just be sure to be consistent in your use of time units throughout your model.

$ LOGLEVEL=INFO python generated/python_arrayified/mm1_queue.py --enable_logging=1 --end_time_ticks=100 --progress_interval_ticks=10 --num_entities=2

The simulation produces the output below. Some highlighted lines are described in the next section.

Sim Arguments:
{
  'enable_logging': 1,
  'end_time_ticks': 100.0,
  'num_entities': 2,
  'progress_interval_ticks': 10.0,
  'rand_seed': 1,
}

INFO:__main__:t = 0.0
INFO:__main__:t = 10.016675542743428
INFO:__main__:t = 20.027835278247828
INFO:__main__:t = 30.055466614166935
INFO:__main__:t = 40.07540843143107
INFO:__main__:t = 50.12466144061003
INFO:__main__:t = 60.12856024231663
INFO:__main__:t = 70.14824380249689
INFO:__main__:t = 80.15701700954823
INFO:__main__:t = 90.25112781277566

Runtime: [1.460378399999172 s]

-----------------------------------
State Tracking
-----------------------------------
instance(Entities_MM1_Queue_state_tracking):
  mm1_queue_arrivals__arriving: StateTracking(track_cumulative=None),
  mm1_queue_arrivals__idle: StateTracking(track_cumulative=None),
  mm1_queue_service__idle: StateTracking(track_cumulative=None),
  mm1_queue_service__servicing: StateTracking(track_cumulative=None)

Writing entity attribute history to [entity_attribute_history.nc]

-----------------------------------
Entity Attribute History Summary over ['time']
-----------------------------------

-------------------- mean --------------------
<xarray.Dataset>
Dimensions:      (entity_id: 2)
Coordinates:
  * entity_id    (entity_id) int64 0 1
Data variables:
    num_waiting  (entity_id) float64 0.3952 0.3241
----------------------------------------------------------

-------------------- std --------------------
<xarray.Dataset>
Dimensions:      (entity_id: 2)
Coordinates:
  * entity_id    (entity_id) int64 0 1
Data variables:
    num_waiting  (entity_id) float64 0.7688 0.6215
----------------------------------------------------------

-------------------- min --------------------
<xarray.Dataset>
Dimensions:      (entity_id: 2)
Coordinates:
  * entity_id    (entity_id) int64 0 1
    quantile     float64 0.0
Data variables:
    num_waiting  (entity_id) float64 0.0 0.0
----------------------------------------------------------

-------------------- 0.25 quantile --------------------
<xarray.Dataset>
Dimensions:      (entity_id: 2)
Coordinates:
  * entity_id    (entity_id) int64 0 1
    quantile     float64 0.25
Data variables:
    num_waiting  (entity_id) float64 0.0 0.0
----------------------------------------------------------

-------------------- 0.50 quantile --------------------
<xarray.Dataset>
Dimensions:      (entity_id: 2)
Coordinates:
  * entity_id    (entity_id) int64 0 1
    quantile     float64 0.5
Data variables:
    num_waiting  (entity_id) float64 0.0 0.0
----------------------------------------------------------

-------------------- 0.75 quantile --------------------
<xarray.Dataset>
Dimensions:      (entity_id: 2)
Coordinates:
  * entity_id    (entity_id) int64 0 1
    quantile     float64 0.75
Data variables:
    num_waiting  (entity_id) float64 1.0 0.5
----------------------------------------------------------

-------------------- max --------------------
<xarray.Dataset>
Dimensions:      (entity_id: 2)
Coordinates:
  * entity_id    (entity_id) int64 0 1
    quantile     float64 1.0
Data variables:
    num_waiting  (entity_id) float64 5.0 4.0
----------------------------------------------------------


-----------------------------------
Entity Attribute History Summary over ['time', 'entity_id']
-----------------------------------

-------------------- mean --------------------
<xarray.Dataset>
Dimensions:      ()
Data variables:
    num_waiting  float64 0.3597
----------------------------------------------------------

-------------------- std --------------------
<xarray.Dataset>
Dimensions:      ()
Data variables:
    num_waiting  float64 0.6999
----------------------------------------------------------

-------------------- min --------------------
<xarray.Dataset>
Dimensions:      ()
Coordinates:
    quantile     float64 0.0
Data variables:
    num_waiting  float64 0.0
----------------------------------------------------------

-------------------- 0.25 quantile --------------------
<xarray.Dataset>
Dimensions:      ()
Coordinates:
    quantile     float64 0.25
Data variables:
    num_waiting  float64 0.0
----------------------------------------------------------

-------------------- 0.50 quantile --------------------
<xarray.Dataset>
Dimensions:      ()
Coordinates:
    quantile     float64 0.5
Data variables:
    num_waiting  float64 0.0
----------------------------------------------------------

-------------------- 0.75 quantile --------------------
<xarray.Dataset>
Dimensions:      ()
Coordinates:
    quantile     float64 0.75
Data variables:
    num_waiting  float64 1.0
----------------------------------------------------------

-------------------- max --------------------
<xarray.Dataset>
Dimensions:      ()
Coordinates:
    quantile     float64 1.0
Data variables:
    num_waiting  float64 5.0
----------------------------------------------------------

Analyze Results

Entity Attribute History Summary

The model defines a single attribute named num_waiting, the number of items waiting in the queue at any one time. The model defines this attribute with probe=True, which causes it to be recorded at each event time while the simulation runs.

The highlighted lines in the mean and std blocks above illustrate the time-averaged mean and standard deviation of num_waiting for each entity in the simulation. We’ve run the simulation with two entities, so there are two entries for entity_id and num_waiting in each block.

Our model is configured with an interarrival rate of 3 and a service rate of 12. M/M/1 Queue theory says that we should expect the mean number of waiting queue requests to be 0.33, with a standard deviation of 0.66. You can see in the entity attribute history summary listed above that the simulated results match the theory well.

There are two Entity Attribute History Summary sections. One displays the results for each entity, aggregated over all simulated time. The other aggregates over all time and all entities.

Another highlighted output indicates that entity attribute data is written to the file entity_attribute_history.nc This is a NetCDF file containing the time-series history of the num_waiting entity attribute.

The section labeled State Tracking has a number of items listed as None, because this model doesn’t have state tracking enabled.

Entity Attribute History Plot

The simulation also generates a file named num_waiting_vs_entities_time.png. This is a 2-d plot of the value of the num_waiting attribute. The y axis is entity ID. The x axis is time. We ran this example with 2 entities, their ids are 0 and 1 in the plot below.

_images/mm1_queue_example_num_waiting_vs_entities_time.png

Output Data File

The simulation generates an output file named entity_attribute_history.nc. This is a NetCDF file. There are a lot of tools and libraries for accessing and using this data.

One example is the python library netCDF4. Below is an example of programmatically examining the file generated by this example.

$ pip install netCDF4
$ python -c 'import netCDF4; d = netCDF4.Dataset("entity_attribute_history.nc", "r"); print(d); print(d.variables["num_waiting"]);'

<class 'netCDF4._netCDF4.Dataset'>
root group (NETCDF4 data model, file format HDF5):
    dimensions(sizes): time(3039), entity_id(2)
    variables(dimensions): float64 time(time), int64 entity_id(entity_id), int64 num_waiting(time, entity_id)
    groups:
<class 'netCDF4._netCDF4.Variable'>
int64 num_waiting(time, entity_id)
unlimited dimensions:
current shape = (3039, 2)
filling on, default _FillValue of -9223372036854775806 used

The highlighted lines above tell you that this data has two dimensions, time and entity_id, a variable named num_waiting along those dimensions, that there are 30258 time coordinates, and 2 entities.

You can use netCDF4 or one of the many other libraries in various languages listed here to access the data for analysis.

The list of tools also has a number of GUI applications that you can use to explore and plot the simulation data in the NetCDF4 file. An especially useful one is Panoply.

Iterate on Model Parameters

We skipped this before, but you can ask the simulation to tell you what command line arguments it accepts by passing --help:

$ python generated/python_arrayified/mm1_queue.py --help

usage: mm1_queue.py [-h] [--num_entities NUM_ENTITIES] [--enable_logging ENABLE_LOGGING] [--rand_seed RAND_SEED] [--end_time_ticks END_TIME_TICKS] [--progress_interval_ticks PROGRESS_INTERVAL_TICKS]
                    [--PARAM_interarrival_time_distribution_rate PARAM_INTERARRIVAL_TIME_DISTRIBUTION_RATE] [--PARAM_service_time_distribution_rate PARAM_SERVICE_TIME_DISTRIBUTION_RATE]

optional arguments:
  -h, --help            show this help message and exit
  --num_entities NUM_ENTITIES
                        Number of entity instances.
  --enable_logging ENABLE_LOGGING
                        Enable logging, or not. Set environment variable LOGLEVEL=DEBUG or INFO.
  --rand_seed RAND_SEED
                        Random number seed.
  --end_time_ticks END_TIME_TICKS
                        Simulation end time (ticks).
  --progress_interval_ticks PROGRESS_INTERVAL_TICKS
                        Interval between reporting sim progress (ticks).
  --PARAM_interarrival_time_distribution_rate PARAM_INTERARRIVAL_TIME_DISTRIBUTION_RATE
  --PARAM_service_time_distribution_rate PARAM_SERVICE_TIME_DISTRIBUTION_RATE

You will see that a couple of the arguments above are prefixed with PARAM. These arguments are used to configure parameters of the model that are defined in the model definition. (The rest of the arguments are for configuring other aspects of the simulation of the model.)

Below, we run the simulation again with different model parameters by providing values for the PARAM_... command line arguments. All the other arguments are specified the same as the previous example above.

The model parameters we’re setting here configure the random distributions used to determine the times of queue arrival and service events.

$LOGLEVEL=INFO python generated/python_arrayified/mm1_queue.py --enable_logging=1 --end_time_ticks=1000 --num_entities=2 --PARAM_interarrival_time_distribution_rate 5 --PARAM_service_time_distribution_rate 10

The highlighted lines below indicate what has changed in the simulation when running with the new PARAM_... arguments. Some of the output has been truncated.

The first highlighted lines are reading back to us the parameter values that we’ve provided. Lower down, you’ll see that the num_waiting values have changed relative to the previous example, which was using the default values for these parameters, which were different from the values that we’ve specified for this run.

If you check the M/M/1 Queue theory you’ll see that the highlighted results match what we’d expect given the new PARAM_... values.

Sim Arguments:
{
  'PARAM_interarrival_time_distribution_rate': 5.0,
  'PARAM_service_time_distribution_rate': 10.0,
  'enable_logging': 1,
  'end_time_ticks': 100.0,
  'num_entities': 2,
  'progress_interval_ticks': 10.0,
  'rand_seed': 1,
}

INFO:__main__:t = 0.0
INFO:__main__:t = 10.01873127619964
INFO:__main__:t = 20.084890768112547
INFO:__main__:t = 30.094167821482863
INFO:__main__:t = 40.108556547423525
INFO:__main__:t = 50.11398855853502
INFO:__main__:t = 60.149593728731936
INFO:__main__:t = 70.17381362676747
INFO:__main__:t = 80.18713329960623
INFO:__main__:t = 90.21419988309019

Runtime: [1.4623036999983015 s]

-----------------------------------
State Tracking
-----------------------------------
instance(Entities_MM1_Queue_state_tracking):
  mm1_queue_arrivals__arriving: StateTracking(track_cumulative=None),
  mm1_queue_arrivals__idle: StateTracking(track_cumulative=None),
  mm1_queue_service__idle: StateTracking(track_cumulative=None),
  mm1_queue_service__servicing: StateTracking(track_cumulative=None)

Writing entity attribute history to [entity_attribute_history.nc]

-----------------------------------
Entity Attribute History Summary over ['time']
-----------------------------------

-------------------- mean --------------------
<xarray.Dataset>
Dimensions:      (entity_id: 2)
Coordinates:
  * entity_id    (entity_id) int64 0 1
Data variables:
    num_waiting  (entity_id) float64 1.206 0.9456
----------------------------------------------------------

-------------------- std --------------------
<xarray.Dataset>
Dimensions:      (entity_id: 2)
Coordinates:
  * entity_id    (entity_id) int64 0 1
Data variables:
    num_waiting  (entity_id) float64 1.626 1.321
----------------------------------------------------------

Peruse Example Model Definition

The model below implements an M/M/1 Queue. This is the model fed to the diagram generator and the simulation code generator in the previous sections.

Note

This section is focused on how to use and iterate on a model, but we’re not describing all the elements of the model here. The concepts and constructs used to build this model are described in the Overview and Modeling Constructs sections.

  1import numpy
  2from dymodetron import \
  3    EntityType, \
  4    EntityAttribute, \
  5    StateMachine, \
  6    ModelDescription, \
  7    Params, \
  8    State, \
  9    Event, \
 10    Transition, \
 11    Transitions, \
 12    dymaction, \
 13    initial_state, \
 14    initial_event, \
 15    types, \
 16    random as dyrandom, \
 17    action
 18
 19# This is the model description. The name of the class is the name of the model.
 20# The docstring is a place to summarize the model.
 21class mm1_queue(ModelDescription):
 22    """An mm1 queue model. See https://en.wikipedia.org/wiki/M/M/1_queue"""
 23
 24# This model has a single entity type. The name of this class is the name of the entity type.
 25class MM1_Queue(EntityType):
 26
 27    # We'll track how many items are waiting in the queue to be serviced.
 28    num_waiting = EntityAttribute(attribute_type=types.scalar_int, initial_value=0.0, probe=True)
 29
 30# Model parameters go here.
 31class ModelParameters(Params):
 32
 33    # We'll sample from this distribution for interarrival times. Note that the distribution parameter
 34    # is in terms of rate (1/time), rather than time directly. The samples will be inter-arrival times, though.
 35    interarrival_time_distribution = dyrandom.ExponentialDistribution(rate=3.0)
 36
 37    # We'll sample from this distribution for time spent waiting to be serviced by the queue. Same comment
 38    # as above about rate vs. time.
 39    service_time_distribution = dyrandom.ExponentialDistribution(rate=12.0)
 40
 41# These events are used in the state machines below.
 42class arrival_event(Event):
 43    pass
 44
 45class done_arriving_event(Event):
 46    pass
 47
 48class service_event(Event):
 49    pass
 50
 51class done_servicing_event(Event):
 52    pass
 53
 54# This is one of two state machines in this model.
 55class mm1_queue_arrivals(StateMachine):
 56    """Handles arrivals for mm1 queue."""
 57
 58    # The model has to specify what entity type is associated with this state machine.
 59    # This state machine models a portion of the behavior of an 'MM1_Queue'.
 60    entity_type = MM1_Queue()
 61
 62    # This is how states are defined within a state machine.
 63    class idle(State):
 64
 65        # These statements are executed when entities enter this state.
 66        @dymaction
 67        def entry_action(queue: MM1_Queue):
 68            """Kick off the next arrival."""
 69
 70            # Declare an entity variable to store the next arrival time. There's one value in this
 71            # variable for each entity in the entity set 'queue'.
 72            t_next_arrival = action.declare(
 73                entity_set=queue,
 74                var_type=numpy.float64
 75            )
 76
 77            # Sample from the random distribution, and use those samples to populate the values of the
 78            # entity variable 't_next_arrival'. One value is sampled and stored for each entity in the
 79            # entity set 'queue'.
 80            action.assign(
 81                entity_set=queue,
 82                entity_var=t_next_arrival,
 83                dist=ModelParameters.interarrival_time_distribution
 84            )
 85
 86            # Generate 'arrival_event' on each entity in the entity set 'queue'. The times for the events
 87            # are taken from the entity variable that we declared and assigned to above.
 88            action.generate_event_rel(
 89                entity_type=MM1_Queue(),
 90                entity_set=queue,
 91                event=arrival_event(),
 92                time_ticks_rel=t_next_arrival.value()
 93            )
 94
 95    class arriving(State):
 96
 97        @dymaction
 98        def entry_action(queue: MM1_Queue):
 99            """Count the arrival, and kick off transition back to mm1_queue_arrivals.idle state."""
100
101            # This is how we track how many are waiting in the queue to be serviced.
102            queue.num_waiting = queue.num_waiting + 1
103
104            # This will kick the state machine back to state 'idle' (see transition table below),
105            # for each entity in the entity set 'queue'.
106            action.generate_event_rel(
107                entity_type=MM1_Queue(),
108                entity_set=queue,
109                event=done_arriving_event(),
110                time_ticks_rel=0
111            )
112
113    # The transition table defines how the state machine moves between states when it receives events.
114    # The information in this table is also visualized in the generated state machine diagram,
115    # which is arguably easier to interpret.
116    transitions = Transitions([
117        Transition(event_type=initial_event(),   source_state=initial_state(), target_state=idle()),
118        Transition(           arrival_event(),                idle(),                       arriving()),
119        Transition(           done_arriving_event(),          arriving(),                   idle()),
120    ])
121
122# This is the other state machine for this model. A model can have as many state machines
123# as you want.
124class mm1_queue_service(StateMachine):
125    """Handles servicing mm1 queue."""
126
127    entity_type = MM1_Queue()
128
129    class idle(State):
130
131        @dymaction
132        def entry_action(queue: MM1_Queue):
133            """Kick off the next service."""
134
135            # Declare entity variable to store the next service time.
136            t_next_service = action.declare(
137                entity_set=queue,
138                var_type=numpy.float64
139            )
140
141            # Sample from the random distribution, and use those samples to populate the values
142            # of 't_next_service'. One value per entity in entity set 'queue'.
143            action.assign(
144                entity_set=queue,
145                entity_var=t_next_service,
146                dist=ModelParameters.service_time_distribution
147            )
148
149            # Generate 'service_event' on each entity in the entity set 'queue', using the times we
150            # sample above.
151            action.generate_event_rel(
152                entity_type=MM1_Queue(),
153                entity_set=queue,
154                event=service_event(),
155                time_ticks_rel=t_next_service.value()
156            )
157
158    class servicing(State):
159
160        @dymaction
161        def entry_action(queue: MM1_Queue):
162            """Count the service, and kick off transition back to mm1_queue_service.idle state."""
163
164            # If the queue has more than zero waiting, we'll decrement.
165            with action.entity_match(
166                entity_type=MM1_Queue(),
167                entity_set=queue,
168                criteria=lambda q: q.num_waiting > 0
169            ) as non_empty_queue:
170
171                # Update how many are waiting in the queue to be serviced.
172                non_empty_queue.num_waiting = non_empty_queue.num_waiting - 1
173
174            # This will kick the state machine back to state 'idle' (see transition table below),
175            # for each entity in the entity set 'queue'.
176            action.generate_event_rel(
177                entity_type=MM1_Queue(),
178                entity_set=queue,
179                event=done_servicing_event(),
180                time_ticks_rel=0
181            )
182
183    # The transition table defines how the state machine moves between states when it receives events.
184    # The information in this table is also visualized in the generated state machine diagram,
185    # which is arguably easier to interpret.
186    transitions = Transitions([
187        Transition(event_type=initial_event(),   source_state=initial_state(), target_state=idle()),
188        Transition(           service_event(),                idle(),                       servicing()),
189        Transition(           done_servicing_event(),         servicing(),                  idle()),
190    ])
191

Iterate On Model Definition

Now we’ll change the model to capture some more real-life effects. We’ll leave the ‘arrivals’ part of the model alone, and modify the ‘service’ portion to model queue processors being on and off shift, i.e. services are being handled or not. Here’s the new state machine diagram:

shifty_mm1_queue_service.html:

%%{init: { 'theme': 'base' } }%% stateDiagram-v2 mm1_queue_service off_shift [*] --> off_shift : initial_event off_shift --> mm1_queue_service : shift_start_event mm1_queue_service --> off_shift : shift_end_event state mm1_queue_service { idle servicing [*] --> idle : initial_event idle --> servicing : service_event servicing --> idle : done_servicing_event }

The changes in the model definition, relative to the previous example model, are highlighted below.

What we’ve changed:

  • New model name and summary description (lines 21-22).

  • New random distributions for modeling on/off shift times (lines 39-43).

  • New events for transitioning on/off shift (lines 58-62).

  • New top-level state machine shifty_mm1_queue_service to define how queues are serviced (starting on line 129).

  • If you look closely, you’ll see that the original mm1_queue_service state machine is now a sub-state-machine within the top-level shifty_mm1_queue_service. The top-level state machine toggles back and forth between the old behavior happening (on shift), vs. nothing happening (off shift).

  • We’ve also modified the mm1_queue_service state machine (which is now a sub-state-machine) to have an entry action that schedules the next shift_end_event to take us back to state off_shift.

  • The off_shift state entry action schedules the next time the processor is servicing requests again.

  1import numpy
  2from dymodetron import \
  3    EntityType, \
  4    EntityAttribute, \
  5    StateMachine, \
  6    SubStateMachine, \
  7    ModelDescription, \
  8    Params, \
  9    Param, \
 10    State, \
 11    Event, \
 12    Transition, \
 13    Transitions, \
 14    dymaction, \
 15    initial_state, \
 16    initial_event, \
 17    types, \
 18    random as dyrandom, \
 19    action
 20
 21class mm1_queue_with_shifts(ModelDescription):
 22    """The example mm1 queue model, modified for queue processing happening only during certain periods of time."""
 23
 24class MM1_Queue(EntityType):
 25
 26    # We'll track how many items are waiting in the queue to be serviced.
 27    num_waiting = EntityAttribute(attribute_type=types.scalar_int, initial_value=0.0, probe=True)
 28
 29class ModelParameters(Params):
 30
 31    # We'll sample from this distribution for interarrival times. Note that the distribution parameter
 32    # is in terms of rate (1/time), rather than time directly. The samples will be inter-arrival times, though.
 33    interarrival_time_distribution = dyrandom.ExponentialDistribution(rate=3.0)
 34
 35    # We'll sample from this distribution for time spent waiting to be serviced by the queue. Same comment
 36    # as above about rate vs. time.
 37    service_time_distribution = dyrandom.ExponentialDistribution(rate=12.0)
 38
 39    # How many hours a queue processor stays on shift.
 40    shift_on_distribution = dyrandom.NormalDistribution(mean=8.0, std=0.1)
 41
 42    # How many hours a queue processor stays off shift.
 43    shift_off_distribution = dyrandom.NormalDistribution(mean=16.0, std=0.5)
 44
 45# These events are used in the state machines below.
 46class arrival_event(Event):
 47    pass
 48
 49class done_arriving_event(Event):
 50    pass
 51
 52class service_event(Event):
 53    pass
 54
 55class done_servicing_event(Event):
 56    pass
 57
 58class shift_start_event(Event):
 59    pass
 60
 61class shift_end_event(Event):
 62    pass
 63
 64
 65class mm1_queue_arrivals(StateMachine):
 66    """Handles arrivals for mm1 queue."""
 67
 68    entity_type = MM1_Queue()
 69
 70    class idle(State):
 71
 72        @dymaction
 73        def entry_action(queue: MM1_Queue):
 74            """Kick off the next arrival."""
 75
 76            # Declare an entity variable to store the next arrival time. There's one value in this
 77            # variable for each entity in the entity set 'queue'.
 78            t_next_arrival = action.declare(
 79                entity_set=queue,
 80                var_type=numpy.float64
 81            )
 82
 83            # Sample from the random distribution, and use those samples to populate the values of the
 84            # entity variable 't_next_arrival'. One value is sampled and stored for each entity in the
 85            # entity set 'queue'.
 86            action.assign(
 87                entity_set=queue,
 88                entity_var=t_next_arrival,
 89                dist=ModelParameters.interarrival_time_distribution
 90            )
 91
 92            # Generate 'arrival_event' on each entity in the entity set 'queue'. The times for the events
 93            # are taken from the entity variable that we declared and assigned to above.
 94            action.generate_event_rel(
 95                entity_type=MM1_Queue(),
 96                entity_set=queue,
 97                event=arrival_event(),
 98                time_ticks_rel=t_next_arrival.value()
 99            )
100
101    class arriving(State):
102
103        @dymaction
104        def entry_action(queue: MM1_Queue):
105            """Count the arrival, and kick off transition back to mm1_queue_arrivals.idle state."""
106
107            # This is how we track how many are waiting in the queue to be serviced.
108            queue.num_waiting = queue.num_waiting + 1
109
110            # This will kick the state machine back to state 'idle' (see transition table below),
111            # for each entity in the entity set 'queue'.
112            action.generate_event_rel(
113                entity_type=MM1_Queue(),
114                entity_set=queue,
115                event=done_arriving_event(),
116                time_ticks_rel=0
117            )
118
119    # The transition table defines how the state machine moves between states when it receives events.
120    # The information in this table is also visualized in the generated state machine diagram,
121    # which is arguably easier to interpret.
122    transitions = Transitions([
123        Transition(event_type=initial_event(),   source_state=initial_state(), target_state=idle()),
124        Transition(           arrival_event(),                idle(),                       arriving()),
125        Transition(           done_arriving_event(),          arriving(),                   idle()),
126    ])
127
128
129class shifty_mm1_queue_service(StateMachine):
130    """Handles servicing mm1 queue when on shift, otherwise idle."""
131
132    entity_type = MM1_Queue()
133
134    class off_shift(State):
135        @dymaction
136        def entry_action(queue: MM1_Queue):
137            """Schedule the next 'on' shift."""
138
139            t_next_shift = action.declare(
140                entity_set=queue,
141                var_type=numpy.float64
142            )
143
144            action.assign(
145                entity_set=queue,
146                entity_var=t_next_shift,
147                dist=ModelParameters.shift_off_distribution
148            )
149
150            action.generate_event_rel(
151                entity_type=MM1_Queue(),
152                entity_set=queue,
153                event=shift_start_event(),
154                time_ticks_rel=t_next_shift.value()
155            )
156
157    class mm1_queue_service(SubStateMachine):
158        """Handles servicing mm1 queue."""
159
160        @dymaction
161        def entry_action(queue: MM1_Queue):
162            """Schedule the next 'off' shift."""
163
164            t_shift_ends = action.declare(
165                entity_set=queue,
166                var_type=numpy.float64
167            )
168
169            action.assign(
170                entity_set=queue,
171                entity_var=t_shift_ends,
172                dist=ModelParameters.shift_on_distribution
173            )
174
175            action.generate_event_rel(
176                entity_type=MM1_Queue(),
177                entity_set=queue,
178                event=shift_end_event(),
179                time_ticks_rel=t_shift_ends.value()
180            )
181
182
183        class idle(State):
184
185            @dymaction
186            def entry_action(queue: MM1_Queue):
187                """Kick off the next service."""
188
189                # Declare entity variable to store the next service time.
190                t_next_service = action.declare(
191                    entity_set=queue,
192                    var_type=numpy.float64
193                )
194
195                # Sample from the random distribution, and use those samples to populate the values
196                # of 't_next_service'. One value per entity in entity set 'queue'.
197                action.assign(
198                    entity_set=queue,
199                    entity_var=t_next_service,
200                    dist=ModelParameters.service_time_distribution
201                )
202
203                # Generate 'service_event' on each entity in the entity set 'queue', using the times we
204                # sample above.
205                action.generate_event_rel(
206                    entity_type=MM1_Queue(),
207                    entity_set=queue,
208                    event=service_event(),
209                    time_ticks_rel=t_next_service.value()
210                )
211
212        class servicing(State):
213
214            @dymaction
215            def entry_action(queue: MM1_Queue):
216                """Count the service, and kick off transition back to mm1_queue_service.idle state."""
217
218                # If the queue has more than zero waiting, we'll decrement.
219                with action.entity_match(
220                    entity_type=MM1_Queue(),
221                    entity_set=queue,
222                    criteria=lambda q: q.num_waiting > 0
223                ) as non_empty_queue:
224
225                    # Update how many are waiting in the queue to be serviced.
226                    non_empty_queue.num_waiting = non_empty_queue.num_waiting - 1
227
228                # This will kick the state machine back to state 'idle' (see transition table below),
229                # for each entity in the entity set 'queue'.
230                action.generate_event_rel(
231                    entity_type=MM1_Queue(),
232                    entity_set=queue,
233                    event=done_servicing_event(),
234                    time_ticks_rel=0
235                )
236
237        # The transition table defines how the state machine moves between states when it receives events.
238        # The information in this table is also visualized in the generated state machine diagram,
239        # which is arguably easier to interpret.
240        transitions = Transitions([
241            Transition(event_type=initial_event(),   source_state=initial_state(), target_state=idle()),
242            Transition(           service_event(),                idle(),                       servicing()),
243            Transition(           done_servicing_event(),         servicing(),                  idle()),
244        ])
245
246    # This is the transition table for the top-level state machine 'shifty_mm1_queue_service'.
247    transitions = Transitions([
248        Transition(event_type=initial_event(),   source_state=initial_state(), target_state=off_shift()),
249        Transition(           shift_start_event(),            off_shift(),                  mm1_queue_service()),
250        Transition(           shift_end_event(),              mm1_queue_service(),          off_shift()),
251    ])
252

We can generate the code for this as follows:

$ python -m dymodetron.generators.python_arrayified --model_definition_file=examples/mm1_queue_with_shifts.py --overwrite_existing=1
python_arrayified
Generating code for entity type [MM1_Queue]
Outputting to [generated/python_arrayified/mm1_queue_with_shifts.py]

If we ask the generated code for --help, it will tell us about the command line options to configure the new parameters:

$ python generated/python_arrayified/mm1_queue_with_shifts.py --help
usage: mm1_queue_with_shifts.py [-h] [--num_entities NUM_ENTITIES] [--enable_logging ENABLE_LOGGING] [--rand_seed RAND_SEED] [--end_time_ticks END_TIME_TICKS] [--progress_interval_ticks PROGRESS_INTERVAL_TICKS]
                                [--PARAM_interarrival_time_distribution_rate PARAM_INTERARRIVAL_TIME_DISTRIBUTION_RATE] [--PARAM_service_time_distribution_rate PARAM_SERVICE_TIME_DISTRIBUTION_RATE]
                                [--PARAM_shift_on_distribution_mean PARAM_SHIFT_ON_DISTRIBUTION_MEAN] [--PARAM_shift_on_distribution_std PARAM_SHIFT_ON_DISTRIBUTION_STD]
                                [--PARAM_shift_off_distribution_mean PARAM_SHIFT_OFF_DISTRIBUTION_MEAN] [--PARAM_shift_off_distribution_std PARAM_SHIFT_OFF_DISTRIBUTION_STD]

optional arguments:
  -h, --help            show this help message and exit
  --num_entities NUM_ENTITIES
                        Number of entity instances.
  --enable_logging ENABLE_LOGGING
                        Enable logging, or not. Set environment variable LOGLEVEL=DEBUG or INFO.
  --rand_seed RAND_SEED
                        Random number seed.
  --end_time_ticks END_TIME_TICKS
                        Simulation end time (ticks).
  --progress_interval_ticks PROGRESS_INTERVAL_TICKS
                        Interval between reporting sim progress (ticks).
  --PARAM_interarrival_time_distribution_rate PARAM_INTERARRIVAL_TIME_DISTRIBUTION_RATE
  --PARAM_service_time_distribution_rate PARAM_SERVICE_TIME_DISTRIBUTION_RATE
  --PARAM_shift_on_distribution_mean PARAM_SHIFT_ON_DISTRIBUTION_MEAN
  --PARAM_shift_on_distribution_std PARAM_SHIFT_ON_DISTRIBUTION_STD
  --PARAM_shift_off_distribution_mean PARAM_SHIFT_OFF_DISTRIBUTION_MEAN
  --PARAM_shift_off_distribution_std PARAM_SHIFT_OFF_DISTRIBUTION_STD

Below we run with the arrival and service rates such that the processor is just keeping up with demand while on shift. Note that we’re running with 5 queues now. We’ve configured so that queues processors on average run for 8 hours with 16 hour breaks.

$ LOGLEVEL=DEBUG python generated/python_arrayified/mm1_queue_with_shifts.py --enable_logging=1 --end_time_ticks=100 --num_entities=5 --PARAM_interarrival_time_distribution_rate=10 --PARAM_service_time_distribution_rate=20 --PARAM_shift_on_distribution_mean=8 --PARAM_shift_off_distribution_mean=16 | grep -E "shift_(start|end)_event" | grep Processing

We’ve enabled LOGLEVEL=DEBUG above, and are applying some bash and grep magic to filter the output in order to look at shift start and end events:

DEBUG:__main__:Processing event: [ScheduledEvent(time_ticks=15.839870402347834, event_id=7, event_obj=<examples.mm1_queue_with_shifts.shift_start_event object at 0x7fb0a56f7a30>, entity_ids=array([2]))]
DEBUG:__main__:Processing event: [ScheduledEvent(time_ticks=16.411418394602425, event_id=8, event_obj=<examples.mm1_queue_with_shifts.shift_start_event object at 0x7fb0a56f7a30>, entity_ids=array([0]))]
DEBUG:__main__:Processing event: [ScheduledEvent(time_ticks=16.48581232444172, event_id=9, event_obj=<examples.mm1_queue_with_shifts.shift_start_event object at 0x7fb0a56f7a30>, entity_ids=array([4]))]
DEBUG:__main__:Processing event: [ScheduledEvent(time_ticks=16.516611992313333, event_id=10, event_obj=<examples.mm1_queue_with_shifts.shift_start_event object at 0x7fb0a56f7a30>, entity_ids=array([3]))]
DEBUG:__main__:Processing event: [ScheduledEvent(time_ticks=16.587737090994786, event_id=11, event_obj=<examples.mm1_queue_with_shifts.shift_start_event object at 0x7fb0a56f7a30>, entity_ids=array([1]))]
DEBUG:__main__:Processing event: [ScheduledEvent(time_ticks=24.00950271981257, event_id=1642, event_obj=<examples.mm1_queue_with_shifts.shift_end_event object at 0x7fb0a56f77f0>, entity_ids=array([2]))]
DEBUG:__main__:Processing event: [ScheduledEvent(time_ticks=24.34930968829747, event_id=1719, event_obj=<examples.mm1_queue_with_shifts.shift_end_event object at 0x7fb0a56f7b50>, entity_ids=array([0]))]
DEBUG:__main__:Processing event: [ScheduledEvent(time_ticks=24.389206084520275, event_id=1751, event_obj=<examples.mm1_queue_with_shifts.shift_end_event object at 0x7fb0a56f7af0>, entity_ids=array([3]))]
DEBUG:__main__:Processing event: [ScheduledEvent(time_ticks=24.50686468621639, event_id=1738, event_obj=<examples.mm1_queue_with_shifts.shift_end_event object at 0x7fb0a56f7ca0>, entity_ids=array([4]))]
DEBUG:__main__:Processing event: [ScheduledEvent(time_ticks=24.71094645767817, event_id=1784, event_obj=<examples.mm1_queue_with_shifts.shift_end_event object at 0x7fb0a56f7940>, entity_ids=array([1]))]
DEBUG:__main__:Processing event: [ScheduledEvent(time_ticks=39.39334645813642, event_id=4160, event_obj=<examples.mm1_queue_with_shifts.shift_start_event object at 0x7fb0a56f7a60>, entity_ids=array([4]))]
DEBUG:__main__:Processing event: [ScheduledEvent(time_ticks=39.58149606852187, event_id=4051, event_obj=<examples.mm1_queue_with_shifts.shift_start_event object at 0x7fb0a56f7b80>, entity_ids=array([2]))]
DEBUG:__main__:Processing event: [ScheduledEvent(time_ticks=40.290511055040334, event_id=4130, event_obj=<examples.mm1_queue_with_shifts.shift_start_event object at 0x7fb0a56f79a0>, entity_ids=array([0]))]
DEBUG:__main__:Processing event: [ScheduledEvent(time_ticks=40.43694464902083, event_id=4195, event_obj=<examples.mm1_queue_with_shifts.shift_start_event object at 0x7fb0a56f7fa0>, entity_ids=array([1]))]
DEBUG:__main__:Processing event: [ScheduledEvent(time_ticks=40.86773007523732, event_id=4143, event_obj=<examples.mm1_queue_with_shifts.shift_start_event object at 0x7fb0a56f7b20>, entity_ids=array([3]))]
DEBUG:__main__:Processing event: [ScheduledEvent(time_ticks=47.32726868032039, event_id=5670, event_obj=<examples.mm1_queue_with_shifts.shift_end_event object at 0x7fb0a56f7a00>, entity_ids=array([4]))]
DEBUG:__main__:Processing event: [ScheduledEvent(time_ticks=47.58026972021337, event_id=5691, event_obj=<examples.mm1_queue_with_shifts.shift_end_event object at 0x7fb0a56f7f10>, entity_ids=array([2]))]
DEBUG:__main__:Processing event: [ScheduledEvent(time_ticks=48.14047961449103, event_id=5830, event_obj=<examples.mm1_queue_with_shifts.shift_end_event object at 0x7fb0a56f77f0>, entity_ids=array([0]))]
DEBUG:__main__:Processing event: [ScheduledEvent(time_ticks=48.36089674700067, event_id=5859, event_obj=<examples.mm1_queue_with_shifts.shift_end_event object at 0x7fb0a56f7d00>, entity_ids=array([1]))]
DEBUG:__main__:Processing event: [ScheduledEvent(time_ticks=48.86332882350327, event_id=5970, event_obj=<examples.mm1_queue_with_shifts.shift_end_event object at 0x7fb0a56f7be0>, entity_ids=array([3]))]
DEBUG:__main__:Processing event: [ScheduledEvent(time_ticks=63.355878950906636, event_id=8018, event_obj=<examples.mm1_queue_with_shifts.shift_start_event object at 0x7fb0a56f7ee0>, entity_ids=array([2]))]
DEBUG:__main__:Processing event: [ScheduledEvent(time_ticks=63.54038311499454, event_id=7961, event_obj=<examples.mm1_queue_with_shifts.shift_start_event object at 0x7fb0a56f78e0>, entity_ids=array([4]))]
DEBUG:__main__:Processing event: [ScheduledEvent(time_ticks=63.636522730745384, event_id=8129, event_obj=<examples.mm1_queue_with_shifts.shift_start_event object at 0x7fb0a56f7f70>, entity_ids=array([0]))]
DEBUG:__main__:Processing event: [ScheduledEvent(time_ticks=64.25318121513186, event_id=8170, event_obj=<examples.mm1_queue_with_shifts.shift_start_event object at 0x7fb0a56f79a0>, entity_ids=array([1]))]
DEBUG:__main__:Processing event: [ScheduledEvent(time_ticks=64.8634887873948, event_id=8253, event_obj=<examples.mm1_queue_with_shifts.shift_start_event object at 0x7fb0a56f7bb0>, entity_ids=array([3]))]
DEBUG:__main__:Processing event: [ScheduledEvent(time_ticks=71.54804911256912, event_id=9702, event_obj=<examples.mm1_queue_with_shifts.shift_end_event object at 0x7fb0a56f7e20>, entity_ids=array([2]))]
DEBUG:__main__:Processing event: [ScheduledEvent(time_ticks=71.60075531927286, event_id=9745, event_obj=<examples.mm1_queue_with_shifts.shift_end_event object at 0x7fb0a56f7ee0>, entity_ids=array([4]))]
DEBUG:__main__:Processing event: [ScheduledEvent(time_ticks=71.62732639263827, event_id=9762, event_obj=<examples.mm1_queue_with_shifts.shift_end_event object at 0x7fb0a56f7b80>, entity_ids=array([0]))]
DEBUG:__main__:Processing event: [ScheduledEvent(time_ticks=72.14682560697014, event_id=9901, event_obj=<examples.mm1_queue_with_shifts.shift_end_event object at 0x7fb0a56f7f70>, entity_ids=array([1]))]
DEBUG:__main__:Processing event: [ScheduledEvent(time_ticks=72.83586559363486, event_id=10052, event_obj=<examples.mm1_queue_with_shifts.shift_end_event object at 0x7fb0a56f7f40>, entity_ids=array([3]))]
DEBUG:__main__:Processing event: [ScheduledEvent(time_ticks=87.01674377115683, event_id=12080, event_obj=<examples.mm1_queue_with_shifts.shift_start_event object at 0x7fb0a56f7850>, entity_ids=array([4]))]
DEBUG:__main__:Processing event: [ScheduledEvent(time_ticks=87.22284147193372, event_id=12067, event_obj=<examples.mm1_queue_with_shifts.shift_start_event object at 0x7fb0a56f7fd0>, entity_ids=array([2]))]
DEBUG:__main__:Processing event: [ScheduledEvent(time_ticks=88.04646050718927, event_id=12275, event_obj=<examples.mm1_queue_with_shifts.shift_start_event object at 0x7fb0a56f7e20>, entity_ids=array([3]))]
DEBUG:__main__:Processing event: [ScheduledEvent(time_ticks=88.41109984957293, event_id=12085, event_obj=<examples.mm1_queue_with_shifts.shift_start_event object at 0x7fb0a56f7bb0>, entity_ids=array([0]))]
DEBUG:__main__:Processing event: [ScheduledEvent(time_ticks=88.53995986997506, event_id=12176, event_obj=<examples.mm1_queue_with_shifts.shift_start_event object at 0x7fb0a56f7ca0>, entity_ids=array([1]))]
DEBUG:__main__:Processing event: [ScheduledEvent(time_ticks=94.99741638388537, event_id=13696, event_obj=<examples.mm1_queue_with_shifts.shift_end_event object at 0x7fb0a56f7f70>, entity_ids=array([4]))]
DEBUG:__main__:Processing event: [ScheduledEvent(time_ticks=95.21440310860663, event_id=13721, event_obj=<examples.mm1_queue_with_shifts.shift_end_event object at 0x7fb0a56f78e0>, entity_ids=array([2]))]
DEBUG:__main__:Processing event: [ScheduledEvent(time_ticks=95.9447754901954, event_id=13876, event_obj=<examples.mm1_queue_with_shifts.shift_end_event object at 0x7fb0a56f7a90>, entity_ids=array([3]))]
DEBUG:__main__:Processing event: [ScheduledEvent(time_ticks=96.4479386316862, event_id=13945, event_obj=<examples.mm1_queue_with_shifts.shift_end_event object at 0x7fb0a56f7f40>, entity_ids=array([0]))]
DEBUG:__main__:Processing event: [ScheduledEvent(time_ticks=96.56639366711914, event_id=13978, event_obj=<examples.mm1_queue_with_shifts.shift_end_event object at 0x7fb0a56f7850>, entity_ids=array([1]))]

The summary plot shows us that the queues are having trouble keeping up with demand:

_images/mm1_queue_with_shifts_overwhelmed.png

Now we’ll modify the shift times, so that queue processors run for 20 hours alternating with 4 hour breaks.

$ LOGLEVEL=DEBUG python generated/python_arrayified/mm1_queue_with_shifts.py --enable_logging=1 --end_time_ticks=100 --num_entities=5 --PARAM_interarrival_time_distribution_rate=10 --PARAM_service_time_distribution_rate=20 --PARAM_shift_on_distribution_mean=20 --PARAM_shift_off_distribution_mean=4 | grep -E "shift_(start|end)_event" | grep Processing

The summary plot shows that the queue processors are now able to manage a little better. Note that the color scale is not the same in this plot as the last one.

_images/mm1_queue_with_shifts_hanging_in_there.png