[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