Home Subscribe

3. Differentiation and Integration

3.1. Using Python to find derivatives and integrals.

Idea

Example 1

Use the Python sympy package to solve the following

(a) \(\frac{dy}{dx} = x^2\)
(b) \(\frac{dy}{dx} = cos(x)\)
(c) \(\frac{dy}{dx} = e^{-x^2}\)
(d) \(\frac{d^3y}{dx^3} = x^4\)
(e) \(\int cos(x)\)
(f) \(\int_0^\infty e^{-x}\)
(g) \(\int_{-\infty}^{\infty}\int_{-\infty}^{\infty} e^{- x^{2} - y^{2}}\, dx\, dy\)

Use the Python scipy package to solve and plot the following

(a) \(\frac{dy(t)}{dt} = -k \ y(t)\) \(k = 0.3; y_0 = 5\)
(b) \(\frac{dy(t)}{dt} = ke^t + y\)

Solution:

sympy

 1from sympy import *
 2# create a "symbol" called x
 3x = Symbol('x')
 4#Define function
 5f = x**2
 6#Calculating Derivative
 7derivative_f = f.diff(x)
 8print(derivative_f)
 9x,y,z = symbols('x y z')
10#sympy.init_printing()
11print(diff(cos(x), x))
12print(diff(exp(x**2), x))
13# third derivative of x^4
14print(diff(x**4, x, x, x))
15print(diff(x**4, x, 3))
16
17print(integrate(cos(x), x))
18print(integrate(exp(-x), (x, 0, oo)))
19print(integrate(exp(-x**2 - y**2), (x, -oo, oo), (y, -oo, oo)))
2*x
-sin(x)
2*x*exp(x**2)
24*x
24*x
sin(x)
1
pi

scipy

\[\begin{align*} \frac{dy(t)}{dt} &= -k \; y(t) \\ k &= 0.3 \\ y_0 &= 5 \\ \frac{dy(t)}{dt} &= ke^t + y \end{align*}\]
 1import numpy as np
 2import matplotlib.pyplot as plt
 3from scipy.integrate import odeint
 4
 5def model(y,t):
 6    k = 0.3
 7    dydt = -k * y 
 8    return dydt 
 9
10y01 = 5
11t1 = np.linspace(0,20)
12y1 = odeint(model, y01, t1)
13
14def model2(y, t, a):
15    dydt = a*np.exp(t) + y 
16    return dydt 
17
18y02 = 1 
19t2 = np.linspace(0, 5) 
20
21y2= odeint(model2, y02, t2, args=(13,)) 
22
23fig, ax = plt.subplots(1,2, figsize=(10,4))
24
25ax[0].plot(t1,y1)
26ax[0].set_xlabel('t')
27ax[0].set_ylabel('y')
28ax[0].set_title(r'$\frac{dy(t)}{dt} = -k \; y(t)$')
29
30ax[1].plot(t2,y2)
31ax[1].set_xlabel('t')
32ax[1].set_title(r'$\frac{dy(t)}{dt} = ke^t + y$')
33#plt.savefig(r'../images/scipy_ode.png')
34plt.show()
scipy-ode
Figure 2. Scipy ode calculation

3.2. Spread of virus using the SIR Model.

Idea

Example 2

A new infectious flu virus is discovered to have infected \(10\%\) of a community. The rate of infection is \(0.35\) and that of recovery is \(0.1\). Given that no one has recovered so far, use the SIR Model to simulate the spread of this virus.

Solution:
 1import numpy as np
 2import matplotlib.pyplot as plt 
 3from scipy.integrate import odeint 
 4
 5# model 
 6def SIR_model(y, t, betta, gamma):
 7    S, I, R = y
 8    dSdt = -beta*S*I 
 9    dIdt = beta*S*I - gamma*I 
10    dRdt = gamma*I 
11    return ([dSdt, dIdt, dRdt]) 
12
13# initial conditions 
14S0 = 0.9 # 90 % 
15I0 = 0.1 
16R0 = 0.0 
17beta = 0.35 
18gamma = 0.1 
19t = np.linspace(0, 100, 10000) # create the simul time 
20
21# result 
22
23sol = odeint(SIR_model, [S0, I0, R0], t, args=(beta, gamma)) 
24sol = np.array(sol) 
25
26plt.figure(figsize=[6,4]) 
27plt.plot(t, sol[:, 0], 'b', label=r"$\frac{dS}{dt}$") 
28plt.plot(t, sol[:, 1], 'r', label=r"$\frac{dI}{dt}$") 
29plt.plot(t, sol[:, 2], 'lime', label=r"$\frac{dR}{dt}$") 
30plt.xlabel("Unit Time"); plt.ylabel("Percentage")
31plt.legend()
32plt.savefig("../images/sir_model.png")
33plt.show()
sir-model
Figure 3. SIR model output.


Add Comment

* Required information
1000
Drag & drop images (max 3)
What is the opposite word of weak?

Comments

No comments yet. Be the first!