{
"cells": [
{
"cell_type": "markdown",
"id": "reflected-joining",
"metadata": {},
"source": [
"# Homework assignment 9 (Lightweight Eastern edition)\n",
"\n",
"This homework assignment considers only the disease spreading models and their numerical solutions \n",
"already discussed in the interactive lecture from Wednesday \n",
"\n",
"**Deadline: Monday, 15th of April, 12:00.**"
]
},
{
"cell_type": "code",
"execution_count": 1,
"id": "crude-function",
"metadata": {},
"outputs": [],
"source": [
"%matplotlib widget\n",
"\n",
"import ipywidgets as widgets\n",
"from ipywidgets import interact, fixed\n",
"import numpy as np\n",
"from numpy import pi\n",
"from numpy.linalg import solve, norm \n",
"import matplotlib.pyplot as plt\n",
"\n",
"# Use a funny plotting style\n",
"plt.xkcd()\n",
"\n",
"newparams = {'figure.figsize': (6.0, 6.0),\n",
" 'axes.grid': True,\n",
" 'lines.markersize': 8, \n",
" 'lines.linewidth': 2,\n",
" 'font.size': 14}\n",
"plt.rcParams.update(newparams)"
]
},
{
"cell_type": "markdown",
"id": "southeast-grace",
"metadata": {},
"source": [
"## Problem 1 SIHR model\n",
"\n",
"**a)** Write down the final ODE system for the SIHR model developed during the first breakout session.\n",
" Denote the rates for the transition $I \\overset{\\gamma_h}{\\to} H$, $I \\overset{\\gamma_r}{\\to} R$,\n",
" and $H \\overset{\\delta}{\\to} R$\n",
" by $\\gamma_h$ and $\\gamma_r$, $\\delta$, respectively. \n",
" \n",
"**b)** As the SIHR model is a refinement of the SIR model where we\n",
" want to distinguish between infectious and hospitialized individuals, we assume the same total $\\gamma$\n",
" as for the simpler SIR model; that is,\n",
" $$\n",
" \\gamma_r + \\gamma_h =: \\gamma = 1/18\n",
" $$\n",
" Assume that\n",
" * that 3.5% of all infected individuals will be hospitalized which means that $\\gamma_h = 0.035\\gamma$\n",
" * hospitalized individuals stay 14 days in the hospital on average, that is $\\delta = 1/14$\n",
" * St. Olav's hospital has roughly 1000 beds with roughly 80% of them being occupied\n",
" \n",
" Use the same initial condition as in the Trondheim scenario as in the lecture, choose a tolerance $\\mathrm{tol}$ such that your numerical solutions are not off by more than 10 persons.\n",
" Now compute and plot the solution $S, I, H, R$ using both the adaptive Euler-Heun and Fehlberg's method.\n",
" Record and state the number of steps you needed with Euler-Heun vs Fehlberg.\n",
" What is (approximately) the largest basic reproduction number $R_0$ for which we will not exceed the maximal number of \n",
" available beds?\n",
" \n",
" *Hint* Before you start, it might be helpful to review the [[SimpleCovidModel.ipynb notebook]](https://wiki.math.ntnu.no/lib/exe/fetch.php?tok=0f251d&media=https%3A%2F%2Fwww.math.ntnu.no%2Femner%2FTMA4125%2F2021v%2Flectures%2FSimpleCovidModel.ipynb). Some of the code snippets provided below might be helpful as well."
]
},
{
"cell_type": "code",
"execution_count": 2,
"id": "exotic-craps",
"metadata": {},
"outputs": [],
"source": [
"# define SIHR class similar to SIR before\n",
"class SIHR:\n",
" # Create model with given parameters\n",
" def __init__(self, beta, gamma_r, gamma_h, delta):\n",
" self.beta = beta # infectional rate\n",
" self.gamma_r = gamma_r # removal rate\n",
" self.gamma_h = gamma_h # removal rate\n",
" self.delta = delta\n",
" \n",
" def __call__(self, t, y):\n",
" # TODO: Implement rhs of SIHR model\n",
" return ...\n",
" "
]
},
{
"cell_type": "code",
"execution_count": 3,
"id": "adjusted-confidentiality",
"metadata": {},
"outputs": [],
"source": [
"from rkm import EmbeddedExplicitRungeKutta\n",
"\n",
"#Heun-Euler\n",
"# TODO: Fill in missing steps indicated by ...\n",
"a = ...\n",
"bhat = ...\n",
"b = ...\n",
"c = ...\n",
"order = ...\n",
"\n",
"euler_heun = EmbeddedExplicitRungeKutta(a, b, c, bhat, order)\n",
"\n",
"#Fehlberg method combines a 4th and 5th order method\n",
"a = ...\n",
"bhat = ...\n",
"b = ...\n",
"c = ...\n",
"order = 4\n",
"fehlberg = EmbeddedExplicitRungeKutta(a, b, c, bhat, order)\n",
"\n",
"# max number of steps\n",
"Nmax = 10000\n",
"tol = ..."
]
},
{
"cell_type": "markdown",
"id": "horizontal-proxy",
"metadata": {},
"source": [
"## Problem 2 SIHRt model\n",
"\n",
"Redo problem 1, but this time develop and use an extension the SIHR model to account for time-limited immunity, assuming 1 year of immunity for each recovered person. Consider a time-period of 5 years and\n",
"find out how many \"infection waves\" will occur where the maximum capacity of beds are exceeded. "
]
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3",
"language": "python",
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.9.2"
}
},
"nbformat": 4,
"nbformat_minor": 5
}