SEE HERE FOR MORE INFORMATION: https://en.wikipedia.org/wiki/Damping#Example:_mass.E2.80.93spring.E2.80.93damper
In [1]:
%matplotlib inline
import numpy as np
from scipy.integrate import odeint
from pylab import *
from bokeh.plotting import *
from bokeh.models import ColumnDataSource
#from bokeh.io import output_notebook
# DEFINE A FUNCTION THAT OUTPUTS POSITION & VELOCITY TIME SERIES,
#GIVEN AN ARRAY OF TIME POINTS (x) AND SOME CONSTANTS (A,B)
#(RELATED TO SPRNG AND DAMPING CONSTANTS)
def g(x,A,B):
#A = 1
#B = 3
#DEFINE OUR SYSTEM OF DIFFERENTIAL EQUATIONS dZ/dt = f(Z, t)
def f(Z,t):
X = Z[0]
Y = Z[1]
#the model equations
f0 = Y
f1 = -A*X + B*Y
return [f0 , f1]
#initial conditions
X0 = 1.1
Y0 = 3
Z0 = [ X0 , Y0 ]
#t = np.linspace(0, 50, 1000) # time grid
#Solve the DEs
return odeint( f , Z0 , x)
#GENERATE THE FIRST SOLUTION
x = np.linspace(0, 50, 1000)
soln = g(x , 3 , -0.3)
y = soln[:,0]
z = soln[:,1]
In [2]:
#LOAD THE NOTEBOOK
output_notebook()
#GET ALL DATA TOGETHER TO DEFINE THE 1ST PLOTS
source = ColumnDataSource( data = dict(x=x,y=y,z=z))
#output_file("lines1.html")
p1 = figure( width=350, height=350, title='Chemical species Y')
p1.line('x','y', line_width = 4 , source = source )
p1.xaxis.axis_label = 'time'
p1.yaxis.axis_label = 'concentration'
p2 = figure(width=350, height=350, title='Chemical species Z')
p2.line('x','z', line_width = 4 , source = source , color = 'red')
p2.xaxis.axis_label = 'time'
p2.yaxis.axis_label = 'concentration'
p3 = figure( width=350, height=350, title='Phase plane')
p3.line('y','z', line_width = 2 , source = source , color = 'green')
p3.yaxis.axis_label = '[Z]'
p3.xaxis.axis_label = '[Y]'
p = gridplot([[p1 , p2 ], [p3]])
#THIS IS HOW WE UPDATE THE PLOTS
def update(A=1 , B =3):
source.data['y'] = g(x , A, B)[:,0]
source.data['z'] = g(x , A, B)[:,1]
source.push_notebook()
#LETS LOOK AT THE CONCENTRATION PROFILES AND PHASE PLANE NOW
show(p)
In [ ]:
#HERE WE GENERATE SOME SLIDERS THAT ALLOW US TO CHANGE THE REACTION RATES A,B ON THE FLY
from IPython.html.widgets import interact
interact( update , A = (-5,5,0.1), B = (-10,10,0.1))
#output_file("lines1.html")
In [ ]: