{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"\n",
"# Validating Runge Kutta Butcher tables using Truncated Taylor Series\n",
"## Authors: Zach Etienne & Brandon Clark\n",
"\n",
"\n",
"## This tutorial notebook is designed to validate the Butcher tables contained within the Butcher dictionary constructed in the [RK Butcher Table Dictionary](Tutorial-RK_Butcher_Table_Dictionary.ipynb) NRPy+ module. \n",
"\n",
"### NRPy+ Source Code for this module: \n",
"* [MoLtimestepping/RK_Butcher_Table_Validation.py](../edit/MoLtimestepping/RK_Butcher_Table_Validation.py) Stores the `Validate` function for validating convergence orders for Runge Kutta methods\n",
"* [MoLtimestepping/RK_Butcher_Table_Dictionary.py](../edit/MoLtimestepping/RK_Butcher_Table_Dictionary.py) [\\[**tutorial**\\]](Tutorial-RK_Butcher_Table_Dictionary.ipynb) Accesses the Butcher table dictionary `Butcher_dict` for known explicit Runge Kutta methods\n",
"\n",
"## Introduction:\n",
"\n",
"Starting with the ODE (ordinary differential equation) initial value problem:\n",
"$$\n",
"y'(t) = f(y,t)\\ \\ \\ y\\left(t=0\\right)=y_0,\n",
"$$\n",
"for various choices of $f(y,t)$, this module validates the Runge Kutta (RK) methods coded in [RK_Butcher_Table_Dictionary.py](../edit/MoLtimestepping/RK_Butcher_Table_Dictionary.py) [**tutorial notebook**](Tutorial-RK_Butcher_Table_Dictionary.ipynb) as follows:\n",
"\n",
"Given $y_0$, and a smooth $f(y,t)$, all explicit RK methods provide an estimate for $y_1 = y\\left(\\Delta t\\right)$, with an error term that is proportional to $\\left(\\Delta t\\right)^m$, where $m$ is an integer typically greater than zero. This error term corresponds to the *local* truncation error. For RK4, for example, while the *total accumulated truncation error* (i.e., the accumulated error at a fixed final time $t_f$) is proportional to $\\left(\\Delta t\\right)^4$, the *local* truncation error (i.e., the error after one arbitrarily chosen timestep $\\Delta t$) is proportional to $\\left(\\Delta t\\right)^5$.\n",
"\n",
"If the exact solution $y(t)$ is known as a closed-form expression, then $y\\left(\\Delta t\\right)$ can be *separately* written as a Taylor expansion about $y(t=0)$:\n",
"\n",
"$$\n",
"y\\left(\\Delta t\\right) = \\sum_{n=0}^\\infty \\frac{y^{(n)}(t=0)}{n!} \\left(\\Delta t\\right)^n,\n",
"$$\n",
"where $y^{(n)}(t=0)$ is the $n$th derivative of $y(t)$ evaluated at $t=0$.\n",
"\n",
"The above expression will be known exactly. Further if one chooses a numerical value for $y_0$ *and leaves $\\Delta t$ unspecified*, any explicit RK method will provide an estimate for $y\\left(\\Delta t\\right)$ of the form\n",
"\n",
"$$\n",
"y\\left(\\Delta t\\right) = \\sum_{n=0}^\\infty a_n \\left(\\Delta t\\right)^n,\n",
"$$\n",
"where $a_n$ *must* match the Taylor expansion of the *exact* solution at least up to and including terms proportional to $\\left(\\Delta t\\right)^m$, where $m$ is the order of the local truncation error. If this is *not* the case, then the Butcher table was almost certainly *not* typed correctly.\n",
"\n",
"Therefore, comparing the numerical result with unspecified $\\Delta t$ against the exact Taylor series provides a convenient (though not perfectly robust) means to verify that the Butcher table for a given RK method was typed correctly. Multiple typos in the Butcher tables were found using this approach.\n",
"\n",
"**Example from Z. Etienne's MATH 521 (Numerical Analysis) lecture notes:**\n",
"\n",
"Consider the ODE\n",
"$$\n",
"y' = y - 2 t e^{-2t},\\quad y(0)=y(t_0)=0.\n",
"$$\n",
"\n",
"* Solve this ODE exactly, then Taylor expand the solution about $t=0$ to\n",
"approximate the solution at $y(t=\\Delta t)$ to fifth order in $\\Delta\n",
"t$.\n",
"* Next solve this ODE using Heun's method (second order in total accumulated truncation error, third order in local truncation error) {\\it by hand} with a step size of\n",
"$\\Delta t$ to find $y(\\Delta t)$. Confirm that the solution obtained\n",
"when using Heun's method has an error term that is at worst\n",
"$\\mathcal{O}\\left((\\Delta t)^3\\right)$. If the dominant error is\n",
"proportional to a higher power of $\\Delta t$, explain the discrepancy.\n",
"\n",
"* Finally solve this ODE using the Ralston method {\\it by hand}\n",
" with a step size of $\\Delta t$ to find $y(\\Delta t)$. Is the\n",
" coefficient on the dominant error term closer to the exact solution\n",
" than Heun's method?\n",
"\n",
"We can solve this equation via the method of integrating factors,\n",
"which states that ODEs of the form:\n",
"$$\n",
"y'(t) + p(t) y(t) = g(t)\n",
"$$\n",
"are solved via \n",
"$$\n",
"y(t) = \\frac{1}{\\mu(t)} \\left[ \\int \\mu(s) g(s) ds + c \\right],\n",
"$$\n",
"where the integrating factor $\\mu(t)$ is given by\n",
"$$\n",
"\\mu(t) = \\exp\\left(\\int p(t) dt\\right)\n",
"$$\n",
"\n",
"Here, $p(t)=-1$ and $g(t) = - 2 t e^{-2t}$. Then\n",
"\n",
"\\begin{equation}\n",
"\\mu(t) = \\exp\\left(-\\int dt\\right) = e^{-t+c} = k e^{-t}\n",
"\\end{equation}\n",
"and\n",
"\n",
"\\begin{align}\n",
"y(t) &= e^t/k \\left[ \\int k e^{-s} (- 2 s e^{-2s}) ds + c \\right] = -2 e^t \\left[ \\int s e^{-3s} ds + c' \\right] \\\\\n",
"&= -2 e^t \\left[ e^{-3 t} \\left(-\\frac{t}{3}-\\frac{1}{9}\\right) + c' \\right] = -2 e^{-2t} \\left(-\\frac{t}{3}-\\frac{1}{9}\\right) -2 c' e^t \\\\\n",
"&= e^{-2t} \\left(2\\frac{t}{3}+\\frac{2}{9}\\right) + c'' e^t \\\\\n",
"\\end{align}\n",
"\n",
"If $y(0)=0$ then we can compute the integration constant $c''$, and\n",
"$y(t)$ becomes\n",
"$$\n",
"y(t) = \\frac{2}{9} e^{-2 t} \\left(3 t + 1 - e^{3 t}\\right).\n",
"$$\n",
"\n",
"The Taylor Series expansion of the exact solution about $t=0$\n",
"evaluated at $y(\\Delta t)$ yields\n",
"$$\n",
"y(\\Delta t) = -(\\Delta t)^2+(\\Delta t)^3-\\frac{3 (\\Delta t)^4}{4}+\\frac{23 (\\Delta\n",
" t)^5}{60}-\\frac{19 (\\Delta t)^6}{120}+O\\left((\\Delta t)^7\\right).\n",
"$$\n",
"\n",
"Next we evaluate $y(\\Delta t)$ using Heun's method. We know $y(0)=y_0=0$ and\n",
"$f(y,t)=y - 2 t e^{-2t}$, so\n",
"\\begin{align}\n",
"k_1 &= \\Delta t f(y(0),0) \\\\\n",
" &= \\Delta t \\times 0 \\\\\n",
" &= 0 \\\\\n",
"k_2 &= \\Delta t f(y(0)+k_1,0+\\Delta t) \\\\\n",
" &= \\Delta t f(y(0)+0,0+\\Delta t) \\\\\n",
" &= \\Delta t (-2 \\Delta t e^{-2\\Delta t}) \\\\\n",
" &= -2 (\\Delta t)^2 e^{-2\\Delta t} \\\\\n",
"y(\\Delta t) &= y_0 + \\frac{1}{2} (k_1 + k_2) + \\mathcal{O}\\left((\\Delta t)^3\\right) \\\\\n",
"&= 0 - (\\Delta t)^2 e^{-2\\Delta t} \\\\\n",
"&= - (\\Delta t)^2 ( 1 - 2 \\Delta t + 2 (\\Delta t)^2 + ...) \\\\\n",
"&= - (\\Delta t)^2 + 2 (\\Delta t)^3 + \\mathcal{O}\\left((\\Delta t)^4\\right).\n",
"\\end{align}\n",
"\n",
"Thus the coefficient on the $(\\Delta t)^3$ term is wrong, but\n",
"this is completely consistent with the fact that our stepping\n",
"scheme is only third-order accurate in $\\Delta t$.\n",
"\n",
"In the below approach, the RK result is subtracted from the exact Taylor series result, as a check to determine whether the RK Butcher table was coded correctly; if it was not, then the odds are good that the RK results will not match to the expected local truncation error order. Multiple $f(y,t)$ are coded below to improve the robustness of this test."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"# Table of Contents\n",
"$$\\label{toc}$$\n",
"\n",
"This notebook is organized as follows\n",
"\n",
"1. [Step 1](#initializenrpy): Initialize needed Python/NRPy+ modules\n",
"1. [Step 2](#table_validate) Validate Convergence Order of Butcher Tables\n",
" 1. [Step 2.a](#rhs): Defining the right-hand side of the ODE\n",
" 1. [Step 2.b](#validfunc): Defining a Validation Function\n",
" 1. [Step 2.c](#rkvalid): Validating RK Methods against ODEs\n",
"1. [Step 3](#latex_pdf_output): Output this notebook to $\\LaTeX$-formatted PDF file"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"# Step 1: Initialize needed Python/NRPy+ modules [Back to [top](#toc)\\]\n",
"$$\\label{initializenrpy}$$\n",
"\n",
"Let's start by importing all the needed modules from Python/NRPy+:"
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {
"execution": {
"iopub.execute_input": "2021-03-07T17:16:35.516909Z",
"iopub.status.busy": "2021-03-07T17:16:35.516113Z",
"iopub.status.idle": "2021-03-07T17:16:36.042355Z",
"shell.execute_reply": "2021-03-07T17:16:36.042843Z"
}
},
"outputs": [],
"source": [
"import sympy as sp # SymPy: The Python computer algebra package upon which NRPy+ depends\n",
"import numpy as np # NumPy: A numerical methods module for Python\n",
"from MoLtimestepping.RK_Butcher_Table_Dictionary import Butcher_dict"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"# Step 2: Validate Convergence Order of Butcher Tables [Back to [top](#toc)\\]\n",
"$$\\label{table_validate}$$\n",
"\n",
"\n",
"Each Butcher table/Runge Kutta method is tested by solving an ODE. Comparing the Taylor series expansions of the exact solution and the numerical solution as discussed in the **Introduction** above will confirm whether the method converges to the appropriate order. "
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"## Step 2.a: Defining the right-hand side of the ODE [Back to [top](#toc)\\]\n",
"$$\\label{rhs}$$\n",
"\n",
"Consider the form of ODE $y'=f(y,t)$. The following begins to construct a dictionary `rhs_dict` of right-hand side functions for us to validate explicit Runge Kutta methods. The most up-to-date catalog of functions stored in `rhs_dict` can be found in the [RK_Butcher_Table_Validation.py](../edit/MoLtimestepping/RK_Butcher_Table_Validation.py) module. "
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {
"execution": {
"iopub.execute_input": "2021-03-07T17:16:36.049940Z",
"iopub.status.busy": "2021-03-07T17:16:36.049261Z",
"iopub.status.idle": "2021-03-07T17:16:36.051301Z",
"shell.execute_reply": "2021-03-07T17:16:36.051767Z"
}
},
"outputs": [],
"source": [
"def fypt(y,t): # Yields expected convergence order for all cases\n",
" # except DP6 which converge to higher order (7, respectively)\n",
" return y+t\n",
"\n",
"def fy(y,t): # Yields expected convergence order for all cases\n",
" return y\n",
"\n",
"def feypt(y,t): # Yields expected convergence order for all cases\n",
" return sp.exp(1.0*(y+t))\n",
"\n",
"def ftpoly6(y,t): # Yields expected convergence order for all cases, L6 has 0 error\n",
" return 2*t**6-389*t**5+15*t**4-22*t**3+81*t**2-t+42\n",
"rhs_dict = {'ypt':fypt, 'y':fy, 'eypt':feypt, 'tpoly6':ftpoly6}"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"## Step 2.b: Defining a Validation Function [Back to [top](#toc)\\]\n",
"$$\\label{validfunc}$$\n",
"\n",
"To validate each Butcher table we compare the exact solutions to ODEs with the numerical solutions using the Runge Kutta scheme built into each Butcher table. The following is a function that\n",
"\n",
"1. Solves the ODE exactly,\n",
"2. Solves the ODE numerically for a given Butcher table, and\n",
"3. Compares the two solutions and checks for the order of convergence by returning their difference.\n",
"\n",
"The `Validate()` function inputs a specified `Butcher_key`, the starting guess solution and time `y_n`, `t_n` and the right-hand side of the ODE corresponding to a specified initial value problem, `rhs_key`.\n",
"\n",
"\n"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {
"execution": {
"iopub.execute_input": "2021-03-07T17:16:36.063763Z",
"iopub.status.busy": "2021-03-07T17:16:36.062979Z",
"iopub.status.idle": "2021-03-07T17:16:36.066476Z",
"shell.execute_reply": "2021-03-07T17:16:36.065801Z"
}
},
"outputs": [],
"source": [
"from MoLtimestepping.RK_Butcher_Table_Dictionary import Butcher_dict\n",
"def Validate(Butcher_key, yn, tn, rhs_key):\n",
" # 1. First we solve the ODE exactly\n",
" y = sp.Function('y')\n",
" sol = sp.dsolve(sp.Eq(y(t).diff(t), rhs_dict[rhs_key](y(t), t)), y(t)).rhs\n",
" constants = sp.solve([sol.subs(t,tn)-yn])\n",
" exact = sol.subs(constants)\n",
"\n",
" # 2. Now we solve the ODE numerically using specified Butcher table\n",
"\n",
" # Access the requested Butcher table\n",
" Butcher = Butcher_dict[Butcher_key][0]\n",
" # Determine number of predictor-corrector steps\n",
" L = len(Butcher)-1\n",
" # Set a temporary array for update values\n",
" k = np.zeros(L, dtype=object)\n",
" # Initialize intermediate variable\n",
" yhat = 0\n",
" # Initialize the updated solution\n",
" ynp1 = 0\n",
" for i in range(L):\n",
" #Initialize and approximate update for solution\n",
" yhat = yn\n",
" for j in range(i):\n",
" # Update yhat for solution using a_ij Butcher table coefficients\n",
" yhat += Butcher[i][j+1]*k[j]\n",
" if Butcher_key == \"DP8\" or Butcher_key == \"L6\":\n",
" yhat = 1.0*sp.N(yhat,20) # Otherwise the adding of fractions kills performance.\n",
" # Determine the next corrector variable k_i using c_i Butcher table coefficients\n",
" k[i] = dt*rhs_dict[rhs_key](yhat, tn + Butcher[i][0]*dt)\n",
" # Update the solution at the next iteration ynp1 using Butcher table coefficients\n",
" ynp1 += Butcher[L][i+1]*k[i]\n",
" # Finish determining the solution for the next iteration\n",
" ynp1 += yn\n",
"\n",
" # Determine the order of the RK method\n",
" order = Butcher_dict[Butcher_key][1]+2\n",
" # Produces Taylor series of exact solution at t=tn about t = 0 with the specified order\n",
" exact_series = sp.series(exact.subs(t, dt),dt, 0, order)\n",
" num_series = sp.series(ynp1, dt, 0, order)\n",
" diff = exact_series-num_series\n",
" return diff"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"## Step 2.c: Validating RK Methods against ODEs [Back to [top](#toc)\\]\n",
"$$\\label{rkvalid}$$\n",
"\n",
"The following makes use of the `Validate()` function above to demonstrate that each method within the Butcher table dictionary converges to the expected order for the given right-hand side expression. "
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {
"execution": {
"iopub.execute_input": "2021-03-07T17:16:36.111715Z",
"iopub.status.busy": "2021-03-07T17:16:36.080988Z",
"iopub.status.idle": "2021-03-07T17:16:58.312416Z",
"shell.execute_reply": "2021-03-07T17:16:58.311759Z"
},
"scrolled": false
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"RK method: \"Euler\".\n",
" When solving y'(t) = t + y(t), y(0)=1,\n",
" the first nonzero term should have local truncation error proportional to O(dt^2) or a higher power of dt.\n",
"Subtracting the numerical result from the exact Taylor expansion, we find a local truncation error of:\n",
" 2 ⎛ 3⎞\n",
"dt + O⎝dt ⎠\n",
" (Coefficients of order 1e-15 or less may generally be ignored, as these are at roundoff error.)\n",
"\n",
"RK method: \"RK2 Heun\".\n",
" When solving y'(t) = t + y(t), y(0)=1,\n",
" the first nonzero term should have local truncation error proportional to O(dt^3) or a higher power of dt.\n",
"Subtracting the numerical result from the exact Taylor expansion, we find a local truncation error of:\n",
" 3 \n",
"dt ⎛ 4⎞\n",
"─── + O⎝dt ⎠\n",
" 3 \n",
" (Coefficients of order 1e-15 or less may generally be ignored, as these are at roundoff error.)\n",
"\n",
"RK method: \"RK2 MP\".\n",
" When solving y'(t) = t + y(t), y(0)=1,\n",
" the first nonzero term should have local truncation error proportional to O(dt^3) or a higher power of dt.\n",
"Subtracting the numerical result from the exact Taylor expansion, we find a local truncation error of:\n",
" 3 \n",
"dt ⎛ 4⎞\n",
"─── + O⎝dt ⎠\n",
" 3 \n",
" (Coefficients of order 1e-15 or less may generally be ignored, as these are at roundoff error.)\n",
"\n",
"RK method: \"RK2 Ralston\".\n",
" When solving y'(t) = t + y(t), y(0)=1,\n",
" the first nonzero term should have local truncation error proportional to O(dt^3) or a higher power of dt.\n",
"Subtracting the numerical result from the exact Taylor expansion, we find a local truncation error of:\n",
" 3 \n",
"dt ⎛ 4⎞\n",
"─── + O⎝dt ⎠\n",
" 3 \n",
" (Coefficients of order 1e-15 or less may generally be ignored, as these are at roundoff error.)\n",
"\n",
"RK method: \"RK3\".\n",
" When solving y'(t) = t + y(t), y(0)=1,\n",
" the first nonzero term should have local truncation error proportional to O(dt^4) or a higher power of dt.\n",
"Subtracting the numerical result from the exact Taylor expansion, we find a local truncation error of:\n",
" 4 \n",
"dt ⎛ 5⎞\n",
"─── + O⎝dt ⎠\n",
" 12 \n",
" (Coefficients of order 1e-15 or less may generally be ignored, as these are at roundoff error.)\n",
"\n",
"RK method: \"RK3 Heun\".\n",
" When solving y'(t) = t + y(t), y(0)=1,\n",
" the first nonzero term should have local truncation error proportional to O(dt^4) or a higher power of dt.\n",
"Subtracting the numerical result from the exact Taylor expansion, we find a local truncation error of:\n",
" 4 \n",
"dt ⎛ 5⎞\n",
"─── + O⎝dt ⎠\n",
" 12 \n",
" (Coefficients of order 1e-15 or less may generally be ignored, as these are at roundoff error.)\n",
"\n",
"RK method: \"RK3 Ralston\".\n",
" When solving y'(t) = t + y(t), y(0)=1,\n",
" the first nonzero term should have local truncation error proportional to O(dt^4) or a higher power of dt.\n",
"Subtracting the numerical result from the exact Taylor expansion, we find a local truncation error of:\n",
" 4 \n",
"dt ⎛ 5⎞\n",
"─── + O⎝dt ⎠\n",
" 12 \n",
" (Coefficients of order 1e-15 or less may generally be ignored, as these are at roundoff error.)\n",
"\n",
"RK method: \"SSPRK3\".\n",
" When solving y'(t) = t + y(t), y(0)=1,\n",
" the first nonzero term should have local truncation error proportional to O(dt^4) or a higher power of dt.\n",
"Subtracting the numerical result from the exact Taylor expansion, we find a local truncation error of:\n",
" 4 \n",
"dt ⎛ 5⎞\n",
"─── + O⎝dt ⎠\n",
" 12 \n",
" (Coefficients of order 1e-15 or less may generally be ignored, as these are at roundoff error.)\n",
"\n",
"RK method: \"RK4\".\n",
" When solving y'(t) = t + y(t), y(0)=1,\n",
" the first nonzero term should have local truncation error proportional to O(dt^5) or a higher power of dt.\n",
"Subtracting the numerical result from the exact Taylor expansion, we find a local truncation error of:\n",
" 5 \n",
"dt ⎛ 6⎞\n",
"─── + O⎝dt ⎠\n",
" 60 \n",
" (Coefficients of order 1e-15 or less may generally be ignored, as these are at roundoff error.)\n",
"\n",
"RK method: \"DP5\".\n",
" When solving y'(t) = t + y(t), y(0)=1,\n",
" the first nonzero term should have local truncation error proportional to O(dt^6) or a higher power of dt.\n",
"Subtracting the numerical result from the exact Taylor expansion, we find a local truncation error of:\n",
" 6 \n",
" dt ⎛ 7⎞\n",
"- ──── + O⎝dt ⎠\n",
" 1800 \n",
" (Coefficients of order 1e-15 or less may generally be ignored, as these are at roundoff error.)\n",
"\n",
"RK method: \"DP5alt\".\n",
" When solving y'(t) = t + y(t), y(0)=1,\n",
" the first nonzero term should have local truncation error proportional to O(dt^6) or a higher power of dt.\n",
"Subtracting the numerical result from the exact Taylor expansion, we find a local truncation error of:\n",
" 6 \n",
"13⋅dt ⎛ 7⎞\n",
"────── + O⎝dt ⎠\n",
"231000 \n",
" (Coefficients of order 1e-15 or less may generally be ignored, as these are at roundoff error.)\n",
"\n",
"RK method: \"CK5\".\n",
" When solving y'(t) = t + y(t), y(0)=1,\n",
" the first nonzero term should have local truncation error proportional to O(dt^6) or a higher power of dt.\n",
"Subtracting the numerical result from the exact Taylor expansion, we find a local truncation error of:\n",
" 6 \n",
"dt ⎛ 7⎞\n",
"──── + O⎝dt ⎠\n",
"3600 \n",
" (Coefficients of order 1e-15 or less may generally be ignored, as these are at roundoff error.)\n",
"\n",
"RK method: \"DP6\".\n",
" When solving y'(t) = t + y(t), y(0)=1,\n",
" the first nonzero term should have local truncation error proportional to O(dt^7) or a higher power of dt.\n",
"Subtracting the numerical result from the exact Taylor expansion, we find a local truncation error of:\n",
" ⎛ 8⎞\n",
"O⎝dt ⎠\n",
" (Coefficients of order 1e-15 or less may generally be ignored, as these are at roundoff error.)\n",
"\n",
"RK method: \"L6\".\n",
" When solving y'(t) = t + y(t), y(0)=1,\n",
" the first nonzero term should have local truncation error proportional to O(dt^7) or a higher power of dt.\n",
"Subtracting the numerical result from the exact Taylor expansion, we find a local truncation error of:\n",
" 5 6 \n",
"- 1.0587911840678754238e-22⋅dt - 2.6469779601696885596e-23⋅dt + 0.0013227513\n",
"\n",
" 7 ⎛ 8⎞\n",
"227513227513⋅dt + O⎝dt ⎠\n",
" (Coefficients of order 1e-15 or less may generally be ignored, as these are at roundoff error.)\n",
"\n",
"RK method: \"DP8\".\n",
" When solving y'(t) = t + y(t), y(0)=1,\n",
" the first nonzero term should have local truncation error proportional to O(dt^9) or a higher power of dt.\n",
"Subtracting the numerical result from the exact Taylor expansion, we find a local truncation error of:\n",
" 2 \n",
"3.6854403535034607753e-18⋅dt + 5.8394451383711465375e-18⋅dt + 3.7764963953332\n",
"\n",
" 3 4 5 \n",
"980617e-18⋅dt + 9.542884942003761195e-19⋅dt + 1.2718729098615353529e-19⋅dt \n",
"\n",
" 6 7 \n",
"+ 3.9082629581905451582e-20⋅dt + 4.8075737201581968464e-21⋅dt + 5.1688448526\n",
"\n",
" 8 9 ⎛ 10⎞\n",
"907316834e-22⋅dt + 7.2078645877627939543e-9⋅dt + O⎝dt ⎠\n",
" (Coefficients of order 1e-15 or less may generally be ignored, as these are at roundoff error.)\n",
"\n"
]
}
],
"source": [
"t, dt = sp.symbols('t dt')\n",
"# Set initial conditions\n",
"t0 = 0\n",
"y0 = 1\n",
"# Set RHS of ODE\n",
"function = 'ypt'# This can be changed, just be careful that the initial conditions are satisfied\n",
"for key,value in Butcher_dict.items():\n",
" print(\"RK method: \\\"\"+str(key)+\"\\\".\")\n",
" y = sp.Function('y')\n",
" print(\" When solving y'(t) = \"+str(rhs_dict[function](y(t),t))+\", y(\"+str(t0)+\")=\"+str(y0)+\",\")\n",
" local_truncation_order = list(value)[1]+1\n",
" print(\" the first nonzero term should have local truncation error proportional to O(dt^\"+str(local_truncation_order)+\") or a higher power of dt.\")\n",
" print(\"Subtracting the numerical result from the exact Taylor expansion, we find a local truncation error of:\")\n",
" sp.pretty_print(Validate(key, y0, t0, function))\n",
"# print(\"\\n\")\n",
" print(\" (Coefficients of order 1e-15 or less may generally be ignored, as these are at roundoff error.)\\n\")"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"# Step 3: Output this notebook to $\\LaTeX$-formatted PDF file \\[Back to [top](#toc)\\]\n",
"$$\\label{latex_pdf_output}$$\n",
"\n",
"The following code cell converts this Jupyter notebook into a proper, clickable $\\LaTeX$-formatted PDF file. After the cell is successfully run, the generated PDF may be found in the root NRPy+ tutorial directory, with filename\n",
"[Tutorial-RK_Butcher_Table_Validation.pdf](Tutorial-RK_Butcher_Table_Validation.pdf) (Note that clicking on this link may not work; you may need to open the PDF file through another means.)"
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {
"execution": {
"iopub.execute_input": "2021-03-07T17:16:58.316815Z",
"iopub.status.busy": "2021-03-07T17:16:58.316182Z",
"iopub.status.idle": "2021-03-07T17:17:02.223642Z",
"shell.execute_reply": "2021-03-07T17:17:02.222688Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Created Tutorial-RK_Butcher_Table_Validation.tex, and compiled LaTeX file\n",
" to PDF file Tutorial-RK_Butcher_Table_Validation.pdf\n"
]
}
],
"source": [
"import cmdline_helper as cmd # NRPy+: Multi-platform Python command-line interface\n",
"cmd.output_Jupyter_notebook_to_LaTeXed_PDF(\"Tutorial-RK_Butcher_Table_Validation\")"
]
}
],
"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.8.2"
}
},
"nbformat": 4,
"nbformat_minor": 2
}