--- a/odes.org Mon Sep 13 18:35:56 2010 +0530
+++ /dev/null Thu Jan 01 00:00:00 1970 +0000
@@ -1,123 +0,0 @@
-* Solving ODEs
-*** Outline
-***** Introduction
-******* What are we going to do?
-******* How are we going to do?
-******* Arsenal Required
-********* working knowledge of arrays
-********* working knowledge of functions
-*** Script
- Welcome friends.
-
- In this tutorial we shall look at solving Ordinary Differential Equations,
- ODE henceforth using odeint in Python. In this tutorial we shall be using
- the concepts of arrays, functions and lists which we have covered in
- previous tutorials.
-
- Let's consider the classic problem of the spread of an epidemic in a
- population.
- This is given by the ordinary differential equation dy/dt = ky(L-y)
- where L is the total population and k is an arbitrary constant. For our
- problem Let us use L=250000, k=0.00003.
- Let the boundary condition be y(0)=250.
-
- Let's now fire up IPython by typing ipython -pylab interpreter.
-
- As we saw in one of the earlier sessions, sometimes we will need more than
- pylab to get our job done. For solving 'ordinary differential equations'
- we shall have to import 'odeint' function, which is a part of the SciPy package.
- So we run the magic command:
-
- In []: from scipy.integrate import odeint
-
- The details regarding `import' shall be covered in a subsequent tutorial.
-
- We can represent the given ODE as a Python function.
- This function takes the dependent variable y and the independent variable t
- as arguments and returns the ODE.
-
- Let us now define our function.
-
- In []: def epid(y, t):
- .... k = 0.00003
- .... L = 250000
- .... return k*y*(L-y)
-
- Independent variable t can be assigned the values in the interval of
- 0 and 12 with 61 points using linspace:
-
- In []: t = linspace(0, 12, 61)
-
- Now obtaining the solution of the ode we defined, is as simple as
- calling the Python's odeint function which we just imported:
-
- In []: y = odeint(epid, 250, t)
-
- We can plot the the values of y against t to get a graphical picture our ODE.
-
- plot(y, t)
- Lets now close this plot and move on to solving ordinary differential equation of
- second order.
- Here we shall take the example ODEs of a simple pendulum.
-
- The equations can be written as a system of two first order ODEs
-
- d(theta)/dt = omega
-
- and
-
- d(omega)/dt = - g/L sin(theta)
-
- Let us define the boundary conditions as: at t = 0,
- theta = theta naught = 10 degrees and
- omega = 0
-
- Let us first define our system of equations as a Python function, pend_int.
- As in the earlier case of single ODE we shall use odeint function of Python
- to solve this system of equations by passing pend_int to odeint.
-
- pend_int is defined as shown:
-
- In []: def pend_int(initial, t):
- .... theta = initial[0]
- .... omega = initial[1]
- .... g = 9.81
- .... L = 0.2
- .... f=[omega, -(g/L)*sin(theta)]
- .... return f
- ....
-
- It takes two arguments. The first argument itself containing two
- dependent variables in the system, theta and omega.
- The second argument is the independent variable t.
-
- In the function we assign the first and second values of the
- initial argument to theta and omega respectively.
- Acceleration due to gravity, as we know is 9.81 meter per second sqaure.
- Let the length of the the pendulum be 0.2 meter.
-
- We create a list, f, of two equations which corresponds to our two ODEs,
- that is d(theta)/dt = omega and d(omega)/dt = - g/L sin(theta).
- We return this list of equations f.
-
- Now we can create a set of values for our time variable t over which we need
- to integrate our system of ODEs. Let us say,
-
- In []: t = linspace(0, 20, 101)
-
- We shall assign the boundary conditions to the variable initial.
-
- In []: initial = [10*2*pi/360, 0]
-
- Now solving this system is just a matter of calling the odeint function with
- the correct arguments.
-
- In []: pend_sol = odeint(pend_int, initial,t)
-
- In []: plot(pend_sol)
- This gives us 2 plots. The green plot is omega vs t and the blue is theta vs t.
-
- Thus we come to the end of this tutorial. In this tutorial we have learnt how to solve ordinary
- differential equations of first and second order.
-
-*** Notes