[Tutor] System of ODEs Question
Eric Lofgren
Eric.Lofgren at unc.edu
Thu Feb 3 08:05:13 CET 2011
So I'm in the process of learning Python, and have been working on a program to solve a very simple S-I-R model, used to study infectious diseases. Basically, it's just a smallish set of differential equations which need to be numerically integrated over time.
Working off of a template program, I came up with the following:
---
#Python implementation of continuous SIR model
#Import Necessary Modules
import numpy as np
import pylab as pl
import scipy.integrate as spi
#Parameter Values
S0 = 0.99999
I0 = 0.00001
R0 = 0.0
PopIn= (S0, I0, R0)
beta=1.4547
gamma=1/7.
t_end = 70
t_start = 1
t_step = 1
t_interval = np.arange(t_start, t_end, t_step)
def eq_system(PopIn,x):
'''Defining SIR System of Equations'''
#Creating an array of equations
Eqs= np.zeros((3))
Eqs[0]= -beta * PopIn[0]*PopIn[1]
Eqs[1]= beta * PopIn[0]*PopIn[1] - gamma*PopIn[1]
Eqs[2]= gamma*PopIn[1]
return Eqs
SIR = spi.odeint(eq_system, PopIn, t_interval)
print SIR
#Plot Everything
#This is really ugly, but works for now
pl.plot(SIR[:,0])
pl.plot(SIR[:,1])
pl.plot(SIR[:,2])
pl.show()
---
The part that is confusing me is defining the function. Namely, it seems to need two arguments - the first needs to be the initial states of the population being modeled, but the second...it doesn't seem to matter what it is. Originally, I didn't include it at all, and the program was very, very unhappy - looking back at the example, it had a second argument, so I put one in, and magically, it worked. But it can be x, it can be t, it can be anything.
Which raises the question: What is it *doing*? Honestly, I'm baffled - can anyone point me in the right direction?
Thanks,
Eric
More information about the Tutor
mailing list