#Nth Order Delay Demo

This is a fun demonstration designed to build intuition around the idea that balancing feedback loops with delays lead to oscillation. It uses a vensim model as the ‘system’ but provides a way for a user to interact with the simulation in realtime - essentially acting as the controller - a balancing feedback loop around the model output.

##About this Technique

This is a way to interact with the models in realtime using your keyboard.

Ingredients

The Game

The student is asked to use the ‘up’ and ‘down’ arrow keys to bring a blue line (the system output) to the value of the dashed red line (the target). However, the inputs from the keyboard go through a delay process (here using either the ‘first order delay’ model, or the ‘third order delay’ model).

When we run this cell, the student will have 60 seconds to bring the blue line to the level of the red line.

%pylab
import pysd
from matplotlib import animation
import numpy as np

Using matplotlib backend: MacOSX
Populating the interactive namespace from numpy and matplotlib

#import the model (need to import each time to reinitialize)
#choose one of the following lines:
#model = pysd.read_vensim('../../models/Basic_Structures/First_Order_Delay.mdl')
model = pysd.read_vensim('../../models/Basic_Structures/Third_Order_Delay.mdl')

#set the delay time in the model
model.set_components({'delay':5})

#set the animation parameters
fps=4
seconds=60
dt=1./fps

#set up the figure axes
fig, ax = plt.subplots()
ax.set_xlim(0,1)
ax.set_ylim(-10, 20)
ax.set_xticks([])
title = ax.set_title('Time %.02f'%0)

#draw the target line
ax.plot([0,1], [10,10], 'r--')

#draw the moving line, just for now. We'll change it later
line, = ax.plot([0,1], [0,0], lw=2)

#set up variables for simulation
input_val = 1
model.components.input = lambda: input_val

#capture keyboard input
def on_key_press(event):
    global input_val
    if event.key == 'up':
        input_val += .25
    elif event.key == 'down':
        input_val -= .25
    sys.stdout.flush()

fig.canvas.mpl_connect('key_press_event', on_key_press)

#make the animation
def animate(t):
    #run the simulation forward
    time = model.components.t+dt
    stocks = model.run(return_columns=['input', 'delay_buffer_1', 'delay_buffer_2', 'delay_buffer_3', 'output'],
                       return_timestamps=[time],
                       initial_condition='current', collect=True)

    #make changes to the display
    level = stocks['output']
    line.set_data([0,1], [level, level])
    title.set_text('Time %.02f'%time)

# call the animator.
anim = animation.FuncAnimation(fig, animate, repeat=False,
                               frames=seconds*fps, interval=1000./fps,
                               blit=False)

record = model.get_record()
record.head()

	input	delay_buffer_1	delay_buffer_2	delay_buffer_3	output
0.25	1	0.221199	0.026499	0.002161	0.002161
0.50	1	0.393469	0.090204	0.014388	0.014388
0.75	1	0.527633	0.173359	0.040505	0.040505
1.00	1	0.632121	0.264241	0.080301	0.080301
1.25	1	0.713495	0.355364	0.131532	0.131532

record.plot();

<matplotlib.axes._subplots.AxesSubplot at 0x10f08a250>

Display student input vs model output

To show how we did, we can plot the input and output over time. Here we start to see the oscillatory behavior (for higher order and longer delays)

plt.plot(x,input_collector, label='Your Input')
plt.plot(x,y, label='Model Response')
plt.legend(loc='lower right')
plt.xlabel('Time [Seconds]')
plt.ylabel('Value');

<matplotlib.text.Text at 0x108ab4bd0>

Display the value of each of the buffer stocks over time

If we plot the stock levels over time, we can see (especially for the third order case) how the delay works to smooth out the input values.

import pandas as pd
delay_stock_values = pd.DataFrame(stocks_collector)
delay_stock_values.plot()
plt.xlabel('Time [Seconds]')
plt.ylabel('Stock Level');

<matplotlib.text.Text at 0x108c7c590>