# HG changeset patch
# User Madhusudan.C.S <madhusudancs@gmail.com>
# Date 1270295071 -19800
# Node ID 4182c6c7e1c6c3f41f6d3b7a992d85477729ca6e
# Parent  a63d7dcba725df99c163c3c5fb0f6a729692158c
Added script on ODEs, session 6 of our workshop.

diff -r a63d7dcba725 -r 4182c6c7e1c6 odes.org
--- a/odes.org	Fri Apr 02 01:25:51 2010 +0530
+++ b/odes.org	Sat Apr 03 17:14:31 2010 +0530
@@ -7,4 +7,117 @@
 ********* working knowledge of arrays
 ********* working knowledge of functions
 *** Script
+    Welcome. 
+    
+    In this tutorial we shall look at solving Ordinary Differential Equations
+    using odeints in Python.
+
+    Let's consider a classical problem of the spread of epidemic in a
+    population.
+    This is given by dy/dt = ky(L-y) where L is the total population.
+    For our problem Let us use L=25000, k=0.00003.
+    Let the boundary condition be y(0)=250.
+
+    First of all run the magic command to import odeint to our program.
+
+    In []: from scipy.integrate import odeint
+
+
+    For now just remember this as a command that does some magic to obtain
+    the function odeint in to our program.
+    We will come back to the details of this command in subsequent sessions.
+
+    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.
+    Our function looks like this:
+    (Showing the slide should be sufficient)
+
+    In []: def epid(y, t):
+      ....     k = 0.00003
+      ....     L = 25000
+      ....     return k*y*(L-y)
+
+
+    Independent variable t can have 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 odeint of the ode we have already defined is as simple as
+    calling the Python's odeint function which we imported:
+    
+    In []: y = odeint(epid, 250, t)
+
+    We can plot the the values of y against t to get a graphical picture our ODE.
+
+
+    Let us move on to solving a system of two ordinary differential equations.
+    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 0 (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 is a 2-tuple containing the two
+    dependent variables in the system, theta and omega.
+    The second argument is the independent variable t.
+
+    In the function we assign theta and omega to first and second values of the
+    initial argument respectively.
+    Acceleration due to gravity, as we know is 9.8 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.
+    So first let us import odeint function into our program using the magic
+    import command
+
+    In []: from scipy.integrate import odeint
+
+    We can call ode_int as:
+
+    In []: pend_sol = odeint(pend_int, initial,t)
+
+    Plotting theta against t and omega against t we obtain the plots as shown
+    in the slide.
+
+    Thus we come to the end of this session on solving ordinary differential
+    equations in Python. Thanks for listening to this tutorial.
+
 *** Notes