<a href="https://colab.research.google.com/github/john-s-butler-dit/Numerical-Analysis-Python/blob/master/Chapter%2003%20-%20Runge%20Kutta/3_ProblemSheet/301_Problem_Sheet.ipynb" target="_parent"></a>

# Application of 2nd order Runge Kutta to Populations Equations

This notebook implements the 2nd Order Runge Kutta method for three different population intial value problems.

# 2nd Order Runge Kutta
The general 2nd Order Runge Kutta method for to the first order differential equation
\begin{equation} y^{'} = f(t,y), \end{equation}
numerical approximates $y$ the at time point $t_i$ as $w_i$
with the  formula:
\begin{equation} w_{i+1}=w_i+\frac{h}{2}\big[k_1+k_2],\end{equation}
for $i=0,...,N-1$, where 
\begin{equation}k_1=f(t_i,w_i)\end{equation}
and
\begin{equation}k_2=f(t_i+h,w_i+hk_1)\end{equation}
and $h$ is the stepsize.

To illustrate the method we will apply it to three intial value problems:
## 1. Linear 
Consider the linear population Differential Equation
\begin{equation} y^{'}=0.1y, \ \ (2000 \leq t \leq 2020), \end{equation}
with the initial condition,
\begin{equation}y(2000)=6.\end{equation}

## 2. Non-Linear Population Equation 
Consider the non-linear population Differential Equation
\begin{equation} y^{'}=0.2y-0.01y^2, \ \ (2000 \leq t \leq 2020), \end{equation}
with the initial condition,
\begin{equation}y(2000)=6.\end{equation}

## 3. Non-Linear Population Equation with an oscillation 
Consider the non-linear population Differential Equation with an oscillation 
\begin{equation} y^{'}=0.2y-0.01y^2+\sin(2\pi t), \ \ (2000 \leq t \leq 2020), \end{equation}
with the initial condition,
\begin{equation}y(2000)=6.\end{equation}

#### Setting up Libraries


```python
## Library
import numpy as np
import math 
import pandas as pd
%matplotlib inline
import matplotlib.pyplot as plt # side-stepping mpl backend
import matplotlib.gridspec as gridspec # subplots
import warnings

warnings.filterwarnings("ignore")


```

## Discrete Interval
The continuous time $a\leq t \leq b $ is discretised into $N$ points seperated by a constant stepsize
\begin{equation} h=\frac{b-a}{N}.\end{equation}
Here the interval is $2000\leq t \leq 2020,$ 
\begin{equation} h=\frac{2020-2000}{200}=0.1.\end{equation}
This gives the 201 discrete points:
\begin{equation} t_0=2000, \ t_1=2000.1, \ ... t_{200}=2020. \end{equation}
This is generalised to 
\begin{equation} t_i=2000+i0.1, \ \ \ i=0,1,...,200.\end{equation}
The plot below shows the discrete time steps:


```python
N=200
t_end=2020.0
t_start=2000.0
h=((t_end-t_start)/N)
t=np.arange(t_start,t_end+h/2,h)
fig = plt.figure(figsize=(10,4))
plt.plot(t,0*t,'o:',color='red')
plt.title('Illustration of discrete time points for h=%s'%(h))
plt.show()
```

# 1. Linear Population Equation
## Exact Solution 
The linear population equation
\begin{equation}y^{'}=0.1y, \ \ (2000 \leq t \leq 2020), \end{equation}
with the initial condition,
\begin{equation}y(2000)=6.\end{equation}
has a known exact (analytic) solution
\begin{equation} y=6e^{0.1(t-2000)}. \end{equation}

## Specific 2nd Order Runge Kutta 
To write the specific 2nd Order Runge Kutta method for the linear population equation we need 
\begin{equation}f(t,y)=0.1y.\end{equation}


```python
def linfun(t,w):
    ftw=0.1*w
    return ftw
```

this gives
\begin{equation}k_1=f(t_i,w_i)=0.lw_i,\end{equation}
\begin{equation}k_2=f(t_i+h,w_i+hk_1)=0.1(w_i+hk_1),\end{equation}
and the difference equation
\begin{equation}w_{i+1}=w_{i}+\frac{h}{2}(k_1+k_2)\end{equation}
for $i=0,...,199$, where $w_i$ is the numerical approximation of $y$ at time $t_i$, with step size $h$ and the initial condition
\begin{equation}w_0=6.\end{equation}


```python
w=np.zeros(N+1)
w[0]=6.0
## 2nd Order Runge Kutta
for k in range (0,N):
    k1=linfun(t[k],w[k])
    k2=linfun(t[k]+h,w[k]+h*k1)
    w[k+1]=w[k]+h/2*(k1+k2)
```

## Plotting Results


```python
y=6*np.exp(0.1*(t-2000))
fig = plt.figure(figsize=(8,4))
plt.plot(t,w,'o:',color='purple',label='Runge Kutta')
plt.plot(t,y,'s:',color='black',label='Exact')
plt.legend(loc='best')
plt.show()
```

## Table
The table below shows the time, the Runge Kutta numerical approximation, $w$,  the exact solution, $y$, and the exact error $|y(t_i)-w_i|$ for the linear population equation:


```python

d = {'time t_i': t[0:10],    'Runge Kutta':w[0:10],'Exact (y)':y[0:10],'Exact Error':np.abs(np.round(y[0:10]-w[0:10],10))}
df = pd.DataFrame(data=d)
df
```




<div>
<style scoped>
    .dataframe tbody tr th:only-of-type {
        vertical-align: middle;
    }

    .dataframe tbody tr th {
        vertical-align: top;
    }

    .dataframe thead th {
        text-align: right;
    }
</style>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th></th>
      <th>time t_i</th>
      <th>Runge Kutta</th>
      <th>Exact (y)</th>
      <th>Exact Error</th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>0</th>
      <td>2000.0</td>
      <td>6.000000</td>
      <td>6.000000</td>
      <td>0.000000</td>
    </tr>
    <tr>
      <th>1</th>
      <td>2000.1</td>
      <td>6.060300</td>
      <td>6.060301</td>
      <td>0.000001</td>
    </tr>
    <tr>
      <th>2</th>
      <td>2000.2</td>
      <td>6.121206</td>
      <td>6.121208</td>
      <td>0.000002</td>
    </tr>
    <tr>
      <th>3</th>
      <td>2000.3</td>
      <td>6.182724</td>
      <td>6.182727</td>
      <td>0.000003</td>
    </tr>
    <tr>
      <th>4</th>
      <td>2000.4</td>
      <td>6.244861</td>
      <td>6.244865</td>
      <td>0.000004</td>
    </tr>
    <tr>
      <th>5</th>
      <td>2000.5</td>
      <td>6.307621</td>
      <td>6.307627</td>
      <td>0.000005</td>
    </tr>
    <tr>
      <th>6</th>
      <td>2000.6</td>
      <td>6.371013</td>
      <td>6.371019</td>
      <td>0.000006</td>
    </tr>
    <tr>
      <th>7</th>
      <td>2000.7</td>
      <td>6.435042</td>
      <td>6.435049</td>
      <td>0.000007</td>
    </tr>
    <tr>
      <th>8</th>
      <td>2000.8</td>
      <td>6.499714</td>
      <td>6.499722</td>
      <td>0.000009</td>
    </tr>
    <tr>
      <th>9</th>
      <td>2000.9</td>
      <td>6.565036</td>
      <td>6.565046</td>
      <td>0.000010</td>
    </tr>
  </tbody>
</table>
</div>



## 2. Non-Linear Population Equation 
\begin{equation} y^{'}=0.2y-0.01y^2, \ \ (2000 \leq t \leq 2020), \end{equation}
with the initial condition,
\begin{equation}y(2000)=6.\end{equation}
## Specific 2nd Order Runge Kutta for the Non-Linear Population Equation
To write the specific 2nd Order Runge Kutta method we need
\begin{equation}f(t,y)=0.2y-0.01y^2,\end{equation}
this gives
\begin{equation}k_1=f(t_i,w_i)=0.2w_i-0.01w_i^2,\end{equation}
\begin{equation}k_2=f(t_i+h,w_i+hk_1)=0.2(w_i+hk_1)-0.01(w_i+hk_1)^2,\end{equation}
and the difference equation
\begin{equation}w_{i+1}=w_{i}+\frac{h}{2}(k_1+k_2)\end{equation}
for $i=0,...,199$, where $w_i$ is the numerical approximation of $y$ at time $t_i$, with step size $h$ and the initial condition
\begin{equation}w_0=6.\end{equation}


```python
def nonlinfun(t,w):
    ftw=0.2*w-0.01*w*w
    return ftw
```


```python
w=np.zeros(N+1)
w[0]=6.0
## 2nd Order Runge Kutta
for k in range (0,N):
    k1=nonlinfun(t[k],w[k])
    k2=nonlinfun(t[k]+h,w[k]+h*k1)
    w[k+1]=w[k]+h/2*(k1+k2)
```

## Results
The plot below shows the Runge Kutta numerical approximation, $w$ (circles) for the non-linear population equation:


```python
fig = plt.figure(figsize=(8,4))
plt.plot(t,w,'o:',color='purple',label='Runge Kutta')
plt.legend(loc='best')
plt.show()
```

## Table
The table below shows the time and the Runge Kutta numerical approximation, $w$,  for the non-linear population equation:


```python
d = {'time t_i': t[0:10], 
     'Runge Kutta':w[0:10]}
df = pd.DataFrame(data=d)
df
```




<div>
<style scoped>
    .dataframe tbody tr th:only-of-type {
        vertical-align: middle;
    }

    .dataframe tbody tr th {
        vertical-align: top;
    }

    .dataframe thead th {
        text-align: right;
    }
</style>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th></th>
      <th>time t_i</th>
      <th>Runge Kutta</th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>0</th>
      <td>2000.0</td>
      <td>6.000000</td>
    </tr>
    <tr>
      <th>1</th>
      <td>2000.1</td>
      <td>6.084332</td>
    </tr>
    <tr>
      <th>2</th>
      <td>2000.2</td>
      <td>6.169328</td>
    </tr>
    <tr>
      <th>3</th>
      <td>2000.3</td>
      <td>6.254977</td>
    </tr>
    <tr>
      <th>4</th>
      <td>2000.4</td>
      <td>6.341270</td>
    </tr>
    <tr>
      <th>5</th>
      <td>2000.5</td>
      <td>6.428197</td>
    </tr>
    <tr>
      <th>6</th>
      <td>2000.6</td>
      <td>6.515747</td>
    </tr>
    <tr>
      <th>7</th>
      <td>2000.7</td>
      <td>6.603909</td>
    </tr>
    <tr>
      <th>8</th>
      <td>2000.8</td>
      <td>6.692672</td>
    </tr>
    <tr>
      <th>9</th>
      <td>2000.9</td>
      <td>6.782025</td>
    </tr>
  </tbody>
</table>
</div>



## 3. Non-Linear Population Equation with an oscilation 
\begin{equation} y^{'}=0.2y-0.01y^2+\sin(2\pi t), \ \ (2000 \leq t \leq 2020), \end{equation}
with the initial condition,
\begin{equation}y(2000)=6.\end{equation}

## Specific 2nd Order Runge Kutta for the Non-Linear Population Equation with an oscilation
To write the specific 2nd Order Runge Kutta difference equation for the intial value problem we need 
\begin{equation}f(t,y)=0.2y-0.01y^2+\sin(2\pi t),\end{equation}
which gives
\begin{equation}k_1=f(t_i,w_i)=0.2w_i-0.01w_i^2+\sin(2\pi t_i),\end{equation}
\begin{equation}k_2=f(t_i+h,w_i+hk_1)=0.2(w_i+hk_1)-0.01(w_i+hk_1)^2+\sin(2\pi (t_i+h)),\end{equation}
and the difference equation
\begin{equation}w_{i+1}=w_{i}+\frac{h}{2}(k_1+k_2)\end{equation}
for $i=0,...,199$, where $w_i$ is the numerical approximation of $y$ at time $t_i$, with step size $h$ and the initial condition
\begin{equation}w_0=6.\end{equation}


```python
def nonlin_oscfun(t,w):
    ftw=0.2*w-0.01*w*w+np.sin(2*np.math.pi*t)
    return ftw
```


```python
w=np.zeros(N+1)
w[0]=6.0
## 2nd Order Runge Kutta
for k in range (0,N):
    k1=nonlin_oscfun(t[k],w[k])
    k2=nonlin_oscfun(t[k]+h,w[k]+h*k1)
    w[k+1]=w[k]+h/2*(k1+k2)
```

## Results
The plot below shows the 2nd order Runge Kutta numerical approximation, $w$ (circles) for the non-linear population equation:


```python
fig = plt.figure(figsize=(8,4))
plt.plot(t,w,'o:',color='purple',label='Taylor')
plt.legend(loc='best')
plt.show()
```

## Table
The table below shows the time and the 2nd order Runge Kutta numerical approximation, $w$,  for the non-linear population equation:


```python
d = {'time t_i': t[0:10], 
     'Runge Kutta':w[0:10]}
df = pd.DataFrame(data=d)
df
```




<div>
<style scoped>
    .dataframe tbody tr th:only-of-type {
        vertical-align: middle;
    }

    .dataframe tbody tr th {
        vertical-align: top;
    }

    .dataframe thead th {
        text-align: right;
    }
</style>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th></th>
      <th>time t_i</th>
      <th>Runge Kutta</th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>0</th>
      <td>2000.0</td>
      <td>6.000000</td>
    </tr>
    <tr>
      <th>1</th>
      <td>2000.1</td>
      <td>6.113722</td>
    </tr>
    <tr>
      <th>2</th>
      <td>2000.2</td>
      <td>6.276109</td>
    </tr>
    <tr>
      <th>3</th>
      <td>2000.3</td>
      <td>6.458005</td>
    </tr>
    <tr>
      <th>4</th>
      <td>2000.4</td>
      <td>6.623032</td>
    </tr>
    <tr>
      <th>5</th>
      <td>2000.5</td>
      <td>6.741504</td>
    </tr>
    <tr>
      <th>6</th>
      <td>2000.6</td>
      <td>6.801784</td>
    </tr>
    <tr>
      <th>7</th>
      <td>2000.7</td>
      <td>6.814712</td>
    </tr>
    <tr>
      <th>8</th>
      <td>2000.8</td>
      <td>6.809444</td>
    </tr>
    <tr>
      <th>9</th>
      <td>2000.9</td>
      <td>6.822305</td>
    </tr>
  </tbody>
</table>
</div>




```python

```

8eb6739fc1948da9a12503e5c8a17c1e26eeb5f0
ipynb
Jupyter Notebook
Chapter 03 - Runge Kutta/Supplementary/301_Problem_Sheet.ipynb
jjcrofts77/Numerical-Analysis-Python
97e4b9274397f969810581ff95f4026f361a56a2
MIT
2019-09-05T21:39:12.000Z
2022-03-26T14:00:25.000Z
Chapter 03 - Runge Kutta/Supplementary/301_Problem_Sheet.ipynb
jjcrofts77/Numerical-Analysis-Python
97e4b9274397f969810581ff95f4026f361a56a2
MIT
Chapter 03 - Runge Kutta/Supplementary/301_Problem_Sheet.ipynb
jjcrofts77/Numerical-Analysis-Python
97e4b9274397f969810581ff95f4026f361a56a2
MIT
2021-06-17T15:34:04.000Z
2022-01-14T14:53:43.000Z
Qwen/Qwen-72B
1. YES
2. YES