{
"cells": [
{
"cell_type": "markdown",
"id": "56889fcd-f08a-4475-9928-af817748f5c8",
"metadata": {
"deletable": false,
"editable": false,
"nbgrader": {
"cell_type": "markdown",
"checksum": "3cdba29419235462e9e49cbaa50af966",
"grade": false,
"grade_id": "cell-5e2b58ff19f0997c",
"locked": true,
"schema_version": 3,
"solution": false,
"task": false
}
},
"source": [
"## Numerical Analysis - Fall semester 2025\n",
"# Serie 11 - Romberg integration and ODEs"
]
},
{
"cell_type": "markdown",
"id": "5d0c899b-58da-49e8-ace0-828c8d5b60ca",
"metadata": {
"deletable": false,
"editable": false,
"nbgrader": {
"cell_type": "markdown",
"checksum": "10edd7f8d5f42ef44d458651f4e462e8",
"grade": false,
"grade_id": "cell-e56f5ef08590e711",
"locked": true,
"schema_version": 3,
"solution": false,
"task": false
}
},
"source": [
"First, we will need to import some of the usual packages. You will have to run this cell every time you restart your notebook."
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "7d640d3e-d688-47ad-95b6-96eed7152684",
"metadata": {
"deletable": false,
"editable": false,
"nbgrader": {
"cell_type": "code",
"checksum": "c421cd6444adfd593791bcaa50c264d3",
"grade": false,
"grade_id": "cell-6bbe893666203fcf",
"locked": true,
"schema_version": 3,
"solution": false,
"task": false
}
},
"outputs": [],
"source": [
"import matplotlib.pyplot as plt\n",
"import numpy as np\n",
"import scipy as sp"
]
},
{
"cell_type": "markdown",
"id": "f63dc15a-a707-43e2-9722-59d5ba8eb67f",
"metadata": {
"deletable": false,
"editable": false,
"nbgrader": {
"cell_type": "markdown",
"checksum": "21b3201c9b8ea4fa934e23e4f9db1540",
"grade": false,
"grade_id": "cell-cf5c15060376e2e7",
"locked": true,
"schema_version": 3,
"solution": false,
"task": false
}
},
"source": [
"
\n",
"\n",
"### Romberg integration method"
]
},
{
"cell_type": "markdown",
"id": "27fd5981-14a6-408f-9f16-00a774a2be5f",
"metadata": {
"deletable": false,
"editable": false,
"nbgrader": {
"cell_type": "markdown",
"checksum": "f74419a853831f3571a636e4dc33443d",
"grade": false,
"grade_id": "cell-95f85d5601bd8908",
"locked": true,
"schema_version": 3,
"solution": false,
"task": false
}
},
"source": [
"Given $a, b \\in \\mathbb{R}$, $a < b$, $n \\in \\mathbb{N} \\setminus \\{0, 1\\}$, and $f \\in \\mathrm{C}^6([a, b], \\mathbb{R})$, the Euler--Maclaurin expansion ensures the existence of $c_1, c_2 \\in \\mathbb{R}$ such that\n",
"$$\n",
"\\int_a^b f = T(n) + c_1 h^2 + c_2 h^4 + \\mathrm{O}(h^6),\n",
"$$\n",
"where $h := (b-a)/n$ and $T(n) := h\\left((f(a)+f(b))/2+\\sum_{i=1}^{n-1} f(a+hi)\\right)$ is the composite trapezium formula. "
]
},
{
"cell_type": "markdown",
"id": "4ceed02c-e5db-4c1d-9b9f-9be98c53aeb6",
"metadata": {
"deletable": false,
"editable": false,
"nbgrader": {
"cell_type": "markdown",
"checksum": "34a3a7383e5a96874b1698bf06e200bd",
"grade": false,
"grade_id": "cell-6378ec916260fd6f",
"locked": true,
"points": 0,
"schema_version": 3,
"solution": false,
"task": true
}
},
"source": [
"\n",
" \n",
"**Exercise 1 (Theoretical):** \n",
"Propose a formula $T_2(n)$ relying only on $T(n)$, $T(2n)$, and $T(4n)$ such that\n",
"$$\n",
"\\int_a^b f = T_2(n) + \\mathrm{O}(h^6).\n",
"$$\n",
"How many evaluations of $f$ does $T_2(n)$ require?\n",
"
\n"
]
},
{
"cell_type": "markdown",
"id": "c8c1936d-5795-41cf-845d-3fb1a490f02a",
"metadata": {
"deletable": false,
"editable": false,
"nbgrader": {
"cell_type": "markdown",
"checksum": "689417c985357a67c506c75709f083a5",
"grade": false,
"grade_id": "cell-d1011fe74ddbe047",
"locked": true,
"schema_version": 3,
"solution": false,
"task": false
}
},
"source": [
"
\n",
"\n",
"### Theoretical tools for ODEs"
]
},
{
"cell_type": "markdown",
"id": "a7651f50-0778-46d9-a531-c96c24ce8038",
"metadata": {
"deletable": false,
"editable": false,
"nbgrader": {
"cell_type": "markdown",
"checksum": "64ab56af452c024436e4ae92cb0a818b",
"grade": false,
"grade_id": "cell-3c49e73133deb1ac",
"locked": true,
"points": 0,
"schema_version": 3,
"solution": false,
"task": true
}
},
"source": [
"\n",
" \n",
"**Exercise 2 (Theoretical):** \n",
"Consider the autonomous vector field associated with the logistic growth model, namely $f : \\mathbb{R} \\to \\mathbb{R} : x \\mapsto x(1-x)$. For every $u_0 \\in [0, \\infty)$, compute the maximal integral curve $u$ of $f$ such that $u(0) = u_0$.\n",
"
\n"
]
},
{
"cell_type": "markdown",
"id": "4e019a05-0667-4229-8d66-77f1f60f3cda",
"metadata": {
"deletable": false,
"editable": false,
"nbgrader": {
"cell_type": "markdown",
"checksum": "cac9dfd27c3c8d4484e05e0427d63bf5",
"grade": false,
"grade_id": "cell-b6258d4954827230",
"locked": true,
"points": 0,
"schema_version": 3,
"solution": false,
"task": true
}
},
"source": [
"\n",
" \n",
"**Exercise 3 (Theoretical):** \n",
"Compute the solution to the IVP\n",
"\\begin{align*}\n",
"\\left\\{\\begin{array}{l}\n",
"u'(t) = Au(t) \\text{ for all } t \\in \\mathbb{R}\\\\\n",
"u(0) = u_0\n",
"\\end{array}\\right.&&\n",
"\\text{where}&&\n",
"A \\coloneq \\begin{bmatrix} 1 & -1\\\\ 1 & 1 \\end{bmatrix}&&\n",
"\\text{and}&&\n",
"u_0 \\in \\mathbb{R}^2.\n",
"\\end{align*}\n",
"
\n"
]
},
{
"cell_type": "markdown",
"id": "2a23c5ca-7434-4632-b3a5-1d37c1b24a69",
"metadata": {
"deletable": false,
"editable": false,
"nbgrader": {
"cell_type": "markdown",
"checksum": "24eb5096d130d8411241451578448596",
"grade": false,
"grade_id": "cell-4e32ab71461f3c03",
"locked": true,
"points": 0,
"schema_version": 3,
"solution": false,
"task": true
}
},
"source": [
"\n",
" \n",
"**Exercise 4 (Theoretical):** \n",
"Consider the autonomous vector field associated with a damped pendulum, namely $f : \\mathbb{R}^2 \\to \\mathbb{R}^2 : (x_1, x_2) \\mapsto (x_2, - \\sin x_1 - x_2)$. What are the equilibrium points of $f$? Study their stability using Lyapunov's first method.\n",
"
\n"
]
},
{
"cell_type": "markdown",
"id": "7f504b66-6e8b-454d-8a9b-7810e53e49b0",
"metadata": {
"deletable": false,
"editable": false,
"nbgrader": {
"cell_type": "markdown",
"checksum": "0e012d327e697c8da928e860c61eaa8b",
"grade": false,
"grade_id": "cell-2f76ac45431af2b3",
"locked": true,
"schema_version": 3,
"solution": false,
"task": false
}
},
"source": [
"
\n",
"\n",
"### Explicit methods for population model"
]
},
{
"cell_type": "markdown",
"id": "40813ae6-1d05-46ac-93e5-6a9cccdd71d9",
"metadata": {
"deletable": false,
"editable": false,
"nbgrader": {
"cell_type": "markdown",
"checksum": "ff678ca7823d57bae25b67484c84fb8a",
"grade": false,
"grade_id": "cell-c92451a8da9db056",
"locked": true,
"schema_version": 3,
"solution": false,
"task": false
}
},
"source": [
"Let $f : \\mathbb{R} \\times \\mathbb{R}^n \\to \\mathbb{R}^n$ and $(t_0, u_0) \\in \\operatorname{dom} f$. Let us consider an ordinary differential equation of the form\n",
"\n",
"$$\n",
"\\begin{cases}\n",
" u'(t) = f(t, u(t)), & t \\in (0, \\tau], \\\\\n",
" u(0)=u_0.\n",
"\\end{cases}\n",
"$$\n",
"\n",
"We want to approximate $u(t_i)$ at the discrete time-steps $t_i = i h$, $i = 0, 1, 2, \\ldots,m$, where $h = \\tau/m$ is the step size and $m$ the number of time-steps. "
]
},
{
"cell_type": "markdown",
"id": "4b008597-b7bf-4761-a266-5c6b696a90ea",
"metadata": {
"deletable": false,
"editable": false,
"nbgrader": {
"cell_type": "markdown",
"checksum": "4967e6b72cf33fcd597bf65a2e199d07",
"grade": false,
"grade_id": "cell-c21721b11e77180b",
"locked": true,
"schema_version": 3,
"solution": false,
"task": false
}
},
"source": [
"\n",
" \n",
"**Exercise 5:** Complete the function `forward_euler` which implements the forward Euler method to compute the approximations $u_i \\approx u(t_i)$ for $i = 0, 1, 2, \\dots, m$.\n",
"
"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "654958ce-11b9-4b26-a87c-a92d60668008",
"metadata": {
"deletable": false,
"nbgrader": {
"cell_type": "code",
"checksum": "af3c42090354b7abb26dd70d479e12a6",
"grade": false,
"grade_id": "cell-a41bad7f587b2cb7",
"locked": false,
"schema_version": 3,
"solution": true,
"task": false
}
},
"outputs": [],
"source": [
"def forward_euler(f, u_0, tau, m):\n",
" # YOUR CODE HERE\n",
" raise NotImplementedError()\n",
" return u"
]
},
{
"cell_type": "markdown",
"id": "79c04a9d-0932-4de2-9508-8c64be7f13b0",
"metadata": {
"deletable": false,
"editable": false,
"nbgrader": {
"cell_type": "markdown",
"checksum": "8330ef15d3abb3d8a2160e3ba83ada94",
"grade": false,
"grade_id": "cell-a9a51259c1344ec6",
"locked": true,
"schema_version": 3,
"solution": false,
"task": false
}
},
"source": [
"\n",
" \n",
"**Exercise 6:** Complete the function `heun` which implements the Heun method to compute the approximations $u_i \\approx u(t_i)$ for $i = 0, 1, 2, \\dots, m$.\n",
"
"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "c2c7cb2a-fd72-4ca2-9888-88b4da58951c",
"metadata": {
"deletable": false,
"nbgrader": {
"cell_type": "code",
"checksum": "205d83cb81a8722acbe6a5c3d3086c92",
"grade": false,
"grade_id": "cell-3e22f0ffd534b765",
"locked": false,
"schema_version": 3,
"solution": true,
"task": false
}
},
"outputs": [],
"source": [
"def heun(f, u_0, tau, m):\n",
" # YOUR CODE HERE\n",
" raise NotImplementedError()\n",
" return u"
]
},
{
"cell_type": "markdown",
"id": "96ca18ab-baf9-4424-9584-7b5b246ef134",
"metadata": {
"deletable": false,
"editable": false,
"nbgrader": {
"cell_type": "markdown",
"checksum": "f59f1016c1b31b73715a7f58e9b72221",
"grade": false,
"grade_id": "cell-546d0e879a54168d",
"locked": true,
"schema_version": 3,
"solution": false,
"task": false
}
},
"source": [
"Let us consider a population of $u$ individuals in an environment\n",
"where at most $u_{\\max}=1000$ individuals can coexist. We suppose that initially\n",
"the number of individuals is $u_0 = 100$ and that the growth factor equals\n",
"a constant $C=2/15$. The considered model for the evolution of the\n",
"population is the following:\n",
"\n",
"$$\n",
"\\left\\{\\begin{array}{ll}\n",
" u'(t)=Cu(t)\\left(1-\\frac{u(t)}{u_{\\max}}\\right), & t \\in (0,100], \\\\\n",
" u(0)=u_0. & \n",
"\\end{array}\n",
"\\right.\n",
"$$"
]
},
{
"cell_type": "markdown",
"id": "9783ac31-f06b-4a8f-b84b-4327fbf300d7",
"metadata": {
"deletable": false,
"editable": false,
"nbgrader": {
"cell_type": "markdown",
"checksum": "4f9bbba6711e41743c5ca71a1682d83e",
"grade": false,
"grade_id": "cell-be34d926e36a1dbd",
"locked": true,
"schema_version": 3,
"solution": false,
"task": false
}
},
"source": [
"\n",
" \n",
"**Exercise 7:** Use the two functions you have implemented above to compute the approximate solution $u_i$, $i=0,1,...,m$ for $m=20$ time-steps and plot on the same graph the obtained numerical solutions in terms of time. Compare the obtained approximations for both forward Euler and Heun with the exact solution given by\n",
" $$\n",
" u(t) = \\frac{u_0u_\\max}{u_0+(u_\\max-u_0)\\exp(-Ct)}.\n",
" $$\n",
" \n",
" Which method gives a better approximation?\n",
"
\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "d85d7412-63ff-4bb5-920d-f6e3fcb8aea5",
"metadata": {
"deletable": false,
"nbgrader": {
"cell_type": "code",
"checksum": "5a913a87f141b23f173a6c9c183eb1cd",
"grade": false,
"grade_id": "cell-ee50b5b4981d323b",
"locked": false,
"schema_version": 3,
"solution": true,
"task": false
}
},
"outputs": [],
"source": [
"# YOUR CODE HERE\n",
"raise NotImplementedError()"
]
},
{
"cell_type": "markdown",
"id": "e6b8708e-df6d-46d0-a55a-632402c287a8",
"metadata": {
"deletable": false,
"editable": false,
"nbgrader": {
"cell_type": "markdown",
"checksum": "3d44ec64ac92546e1e05be494ba97b3d",
"grade": false,
"grade_id": "cell-e07241bfa6f26b16",
"locked": true,
"schema_version": 3,
"solution": false,
"task": false
}
},
"source": [
"\n",
" \n",
"**Exercise 8:** Repeat the last exercise using $m=2000$ time-steps. Is the obtained approximation with $m=2000$ better than the one obtained with $m=20$?\n",
"
\n",
"\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "62cca3be-75ce-4b98-b278-a460c05e8394",
"metadata": {
"deletable": false,
"nbgrader": {
"cell_type": "code",
"checksum": "1e4e12352df40eaa039c062b21873de3",
"grade": false,
"grade_id": "cell-ecc3d186889f7cfa",
"locked": false,
"schema_version": 3,
"solution": true,
"task": false
}
},
"outputs": [],
"source": [
"# YOUR CODE HERE\n",
"raise NotImplementedError()"
]
},
{
"cell_type": "markdown",
"id": "3065e471-f22d-42dd-a7dd-7d6156139e36",
"metadata": {
"deletable": false,
"editable": false,
"nbgrader": {
"cell_type": "markdown",
"checksum": "f7436f9f26b5d47f34dc4a33ad1bd7c3",
"grade": false,
"grade_id": "cell-5b39b5b957674ceb",
"locked": true,
"schema_version": 3,
"solution": false,
"task": false
}
},
"source": [
"
\n",
"\n",
"## The end\n",
"\n",
"Great! You have reached the end of the eleventh exercise notebook."
]
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3 (ipykernel)",
"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.12.3"
}
},
"nbformat": 4,
"nbformat_minor": 5
}