5 ******* How are we going to do? |
5 ******* How are we going to do? |
6 ******* Arsenal Required |
6 ******* Arsenal Required |
7 ********* working knowledge of arrays |
7 ********* working knowledge of arrays |
8 ********* working knowledge of functions |
8 ********* working knowledge of functions |
9 *** Script |
9 *** Script |
10 Welcome. |
10 Welcome friends. |
11 |
11 |
12 In this tutorial we shall look at solving Ordinary Differential Equations |
12 In this tutorial we shall look at solving Ordinary Differential Equations, |
13 using odeints in Python. |
13 ODE henceforth using odeint in Python. |
14 |
14 |
15 Let's consider a classical problem of the spread of epidemic in a |
15 Let's consider the classic problem of the spread of an epidemic in a |
16 population. |
16 population. |
17 This is given by dy/dt = ky(L-y) where L is the total population. |
17 This is given by the ordinary differential equation dy/dt = ky(L-y) |
18 For our problem Let us use L=25000, k=0.00003. |
18 where L is the total population and k is an arbitrary constant. For our |
|
19 problem Let us use L=25000, k=0.00003. |
19 Let the boundary condition be y(0)=250. |
20 Let the boundary condition be y(0)=250. |
20 |
21 |
21 Lets start ipython -pylab interpreter. |
22 Lets fire up IPython by typing ipython -pylab interpreter. |
22 |
23 |
23 As we saw in one of earlier session, sometime pylab wont 'import' all |
24 As we saw in one of earlier session, sometimes pylab wont 'import' all |
24 packages. For solving 'ordinary differential equations' also we shall |
25 packages. For solving 'ordinary differential equations' also we shall |
25 import 'odeint' function which is part SciPy package. So we run the |
26 have to import 'odeint' function which is a part of the SciPy package. |
26 magic command: |
27 So we run the magic command: |
27 |
28 |
28 In []: from scipy.integrate import odeint |
29 In []: from scipy.integrate import odeint |
29 |
30 |
30 # For now just remember this as a command that does some magic to obtain |
31 For now just remember this as a command that does some magic to obtain |
31 # the function odeint in to our program. |
32 the function odeint in to our program. |
32 We will cover more details regarding 'import' in subsequent sessions. |
33 The details regarding `import' shall be covered in a subsequent tutorial. |
33 |
34 |
34 We can represent the given ODE as a Python function. |
35 We can represent the given ODE as a Python function. |
35 This function takes the dependent variable y and the independent variable t |
36 This function takes the dependent variable y and the independent variable t |
36 as arguments and returns the ODE. |
37 as arguments and returns the ODE. |
37 Our function looks like this: |
38 Our function looks like this: |
38 (Showing the slide should be sufficient) |
39 (Showing the slide should be sufficient) |
39 |
40 |
|
41 Let us now define our function. |
|
42 |
40 In []: def epid(y, t): |
43 In []: def epid(y, t): |
41 .... k = 0.00003 |
44 .... k = 0.00003 |
42 .... L = 25000 |
45 .... L = 25000 |
43 .... return k*y*(L-y) |
46 .... return k*y*(L-y) |
44 |
47 |
45 |
48 Independent variable t can be assigned the values in the interval of |
46 Independent variable t can have be assigned the values in the interval of |
|
47 0 and 12 with 61 points using linspace: |
49 0 and 12 with 61 points using linspace: |
48 |
50 |
49 In []: t = linspace(0, 12, 61) |
51 In []: t = linspace(0, 12, 61) |
50 |
52 |
51 Now obtaining the odeint of the ode we have already defined is as simple as |
53 Now obtaining the solution of the ode we defined, is as simple as |
52 calling the Python's odeint function which we imported: |
54 calling the Python's odeint function which we just imported: |
53 |
55 |
54 In []: y = odeint(epid, 250, t) |
56 In []: y = odeint(epid, 250, t) |
55 |
57 |
56 We can plot the the values of y against t to get a graphical picture our ODE. |
58 We can plot the the values of y against t to get a graphical picture our ODE. |
57 |
59 |
58 plot(y, t) |
60 plot(y, t) |
59 Lets close this plot and move on to solving ordinary differential equation of |
61 Lets now close this plot and move on to solving ordinary differential equation of |
60 second order. |
62 second order. |
61 Here we shall take the example ODEs of a simple pendulum. |
63 Here we shall take the example ODEs of a simple pendulum. |
62 |
64 |
63 The equations can be written as a system of two first order ODEs |
65 The equations can be written as a system of two first order ODEs |
64 |
66 |
89 |
91 |
90 It takes two arguments. The first argument itself containing two |
92 It takes two arguments. The first argument itself containing two |
91 dependent variables in the system, theta and omega. |
93 dependent variables in the system, theta and omega. |
92 The second argument is the independent variable t. |
94 The second argument is the independent variable t. |
93 |
95 |
94 In the function we assign theta and omega to first and second values of the |
96 In the function we assign the first and second values of the |
95 initial argument respectively. |
97 initial argument to theta and omega respectively. |
96 Acceleration due to gravity, as we know is 9.8 meter per second sqaure. |
98 Acceleration due to gravity, as we know is 9.81 meter per second sqaure. |
97 Let the length of the the pendulum be 0.2 meter. |
99 Let the length of the the pendulum be 0.2 meter. |
98 |
100 |
99 We create a list, f, of two equations which corresponds to our two ODEs, |
101 We create a list, f, of two equations which corresponds to our two ODEs, |
100 that is d(theta)/dt = omega and d(omega)/dt = - g/L sin(theta). |
102 that is d(theta)/dt = omega and d(omega)/dt = - g/L sin(theta). |
101 We return this list of equations f. |
103 We return this list of equations f. |