{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"\n",
"# Explicit Runge Kutta methods and their Butcher tables\n",
"## Authors: Brandon Clark & Zach Etienne\n",
"\n",
"## This tutorial notebook stores known explicit Runge Kutta-like methods as Butcher tables in a Python dictionary format. \n",
"\n",
"**Notebook Status:** Validated \n",
"\n",
"**Validation Notes:** This tutorial notebook has been confirmed to be **self-consistent with its corresponding NRPy+ module**, as documented [below](#code_validation). In addition, each of these Butcher tables has been verified to yield an RK method to the expected local truncation error in a challenging battery of ODE tests, in the [RK Butcher Table Validation tutorial notebook](Tutorial-RK_Butcher_Table_Validation.ipynb).\n",
"\n",
"### NRPy+ Source Code for this module: [MoLtimestepping/RK_Butcher_Table_Dictionary.py](../edit/MoLtimestepping/RK_Butcher_Table_Dictionary.py)\n",
"\n",
"## Introduction:\n",
"\n",
"The family of explicit [Runge Kutta](https://en.wikipedia.org/w/index.php?title=Runge%E2%80%93Kutta_methods&oldid=898536315)-like methods are commonly used when numerically solving ordinary differential equation (ODE) initial value problems of the form\n",
"\n",
"$$ y'(t) = f(y,t),\\ \\ \\ y(t_0)=y_0.$$\n",
"\n",
"These methods can be extended to solve time-dependent partial differential equations (PDEs) via the [Method of Lines](https://en.wikipedia.org/w/index.php?title=Method_of_lines&oldid=855390257). In the Method of Lines, the above ODE can be generalized to $N$ coupled ODEs, all written as first-order-in-time PDEs of the form\n",
"\n",
"$$ \\partial_{t}\\mathbf{u}(t,x,y,u_1,u_2,u_3,...)=\\mathbf{f}(t,x,y,...,u_1,u_{1,x},...),$$\n",
"\n",
"where $\\mathbf{u}$ and $\\mathbf{f}$ are vectors. The spatial partial derivatives of components of $\\mathbf{u}$, e.g., $u_{1,x}$, may be computed using approximate numerical differentiation, like finite differences.\n",
"\n",
"As any explicit Runge-Kutta method has its own unique local truncation error, can in principle be used to solve time-dependent PDEs using the Method of Lines, and may be stable under different Courant-Friedrichs-Lewy (CFL) conditions, it is useful to have multiple methods at one's disposal. **This module provides a number of such methods.**\n",
"\n",
"More details about the Method of Lines is discussed further in the [Tutorial-RK_Butcher_Table_Generating_C_Code](Tutorial-RK_Butcher_Table_Generating_C_Code.ipynb) module where we generate the C code to implement the Method of Lines, and additional description can be found in the [Numerically Solving the Scalar Wave Equation: A Complete C Code](Tutorial-Start_to_Finish-ScalarWave.ipynb) NRPy+ tutorial notebook."
]
},
{
"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 modules\n",
"1. [Step 2](#introbutcher): The Family of Explicit Runge-Kutta-Like Schemes (Butcher Tables)\n",
" 1. [Step 2a](#codebutcher): Generating a Dictionary of Butcher Tables for Explicit Runge Kutta Techniques \n",
" 1. [Step 2.a.i](#euler): Euler's Method\n",
" 1. [Step 2.a.ii](#rktwoheun): RK2 Heun's Method\n",
" 1. [Step 2.a.iii](#rk2mp): RK2 Midpoint Method\n",
" 1. [Step 2.a.iv](#rk2ralston): RK2 Ralston's Method\n",
" 1. [Step 2.a.v](#rk3): Kutta's Third-order Method\n",
" 1. [Step 2.a.vi.](#rk3heun): RK3 Heun's Method\n",
" 1. [Step 2.a.vii](#rk3ralston): RK3 Ralston's Method\n",
" 1. [Step 2.a.viii](#ssprk3): Strong Stability Preserving Runge-Kutta (SSPRK3) Method\n",
" 1. [Step 2.a.ix](#rkfour): Classic RK4 Method\n",
" 1. [Step 2.a.x](#dp5): RK5 Dormand-Prince Method\n",
" 1. [Step 2.a.xi](#dp5alt): RK5 Dormand-Prince Method Alternative\n",
" 1. [Step 2.a.xii](#ck5): RK5 Cash-Karp Method\n",
" 1. [Step 2.a.xiii](#dp6): RK6 Dormand-Prince Method\n",
" 1. [Step 2.a.xiv](#l6): RK6 Luther Method\n",
" 1. [Step 2.a.xv](#dp8): RK8 Dormand-Prince Method\n",
"1. [Step 3](#code_validation): Code Validation against `MoLtimestepping.RK_Butcher_Table_Dictionary` NRPy+ module\n",
"1. [Step 4](#latex_pdf_output): Output this notebook to $\\LaTeX$-formatted PDF file"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"# Step 1: Initialize needed Python modules [Back to [top](#toc)\\]\n",
"$$\\label{initializenrpy}$$\n",
"\n",
"Let's start by importing all the needed modules from Python:"
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {
"execution": {
"iopub.execute_input": "2021-03-07T17:16:28.146682Z",
"iopub.status.busy": "2021-03-07T17:16:28.145902Z",
"iopub.status.idle": "2021-03-07T17:16:28.462938Z",
"shell.execute_reply": "2021-03-07T17:16:28.462271Z"
}
},
"outputs": [],
"source": [
"# Step 1: Initialize needed Python modules\n",
"import sympy as sp # SymPy: The Python computer algebra package upon which NRPy+ depends"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"# Step 2: The Family of Explicit Runge-Kutta-Like Schemes (Butcher Tables) [Back to [top](#toc)\\]\n",
"$$\\label{introbutcher}$$\n",
"\n",
"In general, a predictor-corrector method performs an estimate timestep from $n$ to $n+1$, using e.g., a Runge Kutta method, to get a prediction of the solution at timestep $n+1$. This is the \"predictor\" step. Then it uses this prediction to perform another, \"corrector\" step, designed to increase the accuracy of the solution.\n",
"\n",
"Let us focus on the ordinary differential equation (ODE)\n",
"\n",
"$$ y'(t) = f(y,t), $$\n",
"\n",
"which acts as an analogue for a generic PDE $\\partial_{t}u(t,x,y,...)=f(t,x,y,...,u,u_x,...)$.\n",
"\n",
"The general family of Runge Kutta \"explicit\" timestepping methods are implemented using the following scheme:\n",
"\n",
"$$y_{n+1} = y_n + \\sum_{i=1}^s b_ik_i $$\n",
"\n",
"where \n",
"\n",
"\\begin{align}\n",
"k_1 &= \\Delta tf(y_n, t_n) \\\\\n",
"k_2 &= \\Delta tf(y_n + [a_{21}k_1], t_n + c_2\\Delta t) \\\\\n",
"k_3 &= \\Delta tf(y_n +[a_{31}k_1 + a_{32}k_2], t_n + c_3\\Delta t) \\\\\n",
"& \\ \\ \\vdots \\\\\n",
"k_s &= \\Delta tf(y_n +[a_{s1}k_1 + a_{s2}k_2 + \\cdots + a_{s, s-1}k_{s-1}], t_n + c_s\\Delta t)\n",
"\\end{align}\n",
"\n",
"Note $s$ is the number of right-hand side evaluations necessary for any given method, i.e., for RK2 $s=2$ and for RK4 $s=4$, and for RK6 $s=7$. These schemes are often written in the form of a so-called \"Butcher tableau\". or \"Butcher table\":\n",
"\n",
"$$\\begin{array}{c|ccccc}\n",
" 0 & \\\\\n",
" c_2 & a_{21} & \\\\\n",
" c_3 & a_{31} & a_{32} & \\\\\n",
" \\vdots & \\vdots & & \\ddots \\\\\n",
" c_s & a_{s_1} & a_{s2} & \\cdots & a_{s,s-1} \\\\ \\hline \n",
" & b_1 & b_2 & \\cdots & b_{s-1} & b_s\n",
"\\end{array} $$\n",
"\n",
"As an example, the \"classic\" fourth-order Runge Kutta (RK4) method obtains the solution $y(t)$ to the single-variable ODE $y'(t) = f(y(t),t)$ at time $t_{n+1}$ from $t_n$ via:\n",
"\n",
"\\begin{align}\n",
"k_1 &= \\Delta tf(y_n, t_n), \\\\\n",
"k_2 &= \\Delta tf(y_n + \\frac{1}{2}k_1, t_n + \\frac{\\Delta t}{2}), \\\\\n",
"k_3 &= \\Delta tf(y_n + \\frac{1}{2}k_2, t_n + \\frac{\\Delta t}{2}), \\\\\n",
"k_4 &= \\Delta tf(y_n + k_3, t_n + \\Delta t), \\\\\n",
"y_{n+1} &= y_n + \\frac{1}{6}(k_1 + 2k_2 + 2k_3 + k_4) + \\mathcal{O}\\big((\\Delta t)^5\\big).\n",
"\\end{align}\n",
"\n",
"Its corresponding Butcher table is constructed as follows:\n",
"\n",
"$$\\begin{array}{c|cccc}\n",
" 0 & \\\\\n",
" 1/2 & 1/2 & \\\\ \n",
" 1/2 & 0 & 1/2 & \\\\\n",
" 1 & 0 & 0 & 1 & \\\\ \\hline\n",
" & 1/6 & 1/3 & 1/3 & 1/6\n",
"\\end{array} $$\n",
"\n",
"\n",
"This is one example of many explicit [Runge Kutta methods](https://en.wikipedia.org/w/index.php?title=List_of_Runge%E2%80%93Kutta_methods&oldid=896594269). Throughout the following sections we will highlight different Runge Kutta schemes and their Butcher tables from the first-order Euler's method up to and including an eighth-order method."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"## Step 2.a: Generating a Dictionary of Butcher Tables for Explicit Runge Kutta Techniques [Back to [top](#toc)\\]\n",
"$$\\label{codebutcher}$$\n",
"\n",
"We can store all of the Butcher tables in Python's **Dictionary** format using the curly brackets {} and 'key':value pairs. The 'key' will be the *name* of the Runge Kutta method and the value will be the Butcher table itself stored as a list of lists. The convergence order for each Runge Kutta method is also stored. We will construct the dictionary `Butcher_dict` one Butcher table at a time in the following sections."
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {
"execution": {
"iopub.execute_input": "2021-03-07T17:16:28.467913Z",
"iopub.status.busy": "2021-03-07T17:16:28.466544Z",
"iopub.status.idle": "2021-03-07T17:16:28.471177Z",
"shell.execute_reply": "2021-03-07T17:16:28.470508Z"
}
},
"outputs": [],
"source": [
"# Step 2a: Generating a Dictionary of Butcher Tables for Explicit Runge Kutta Techniques\n",
"\n",
"# Initialize the dictionary Butcher_dict\n",
"Butcher_dict = {}"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"### Step 2.a.i: Euler's Method [Back to [top](#toc)\\]\n",
"$$\\label{euler}$$\n",
"\n",
"[Forward Euler's method](https://en.wikipedia.org/w/index.php?title=Euler_method&oldid=896152463) is a first order Runge Kutta method. Euler's method obtains the solution $y(t)$ at time $t_{n+1}$ from $t_n$ via:\n",
"$$y_{n+1} = y_{n} + \\Delta tf(y_{n}, t_{n})$$\n",
"with the trivial corresponding Butcher table \n",
"$$\\begin{array}{c|c}\n",
"0 & \\\\ \\hline\n",
" & 1 \n",
"\\end{array}$$\n"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {
"execution": {
"iopub.execute_input": "2021-03-07T17:16:28.475988Z",
"iopub.status.busy": "2021-03-07T17:16:28.475345Z",
"iopub.status.idle": "2021-03-07T17:16:28.477769Z",
"shell.execute_reply": "2021-03-07T17:16:28.478357Z"
}
},
"outputs": [],
"source": [
"# Step 2.a.i: Euler's Method\n",
"\n",
"Butcher_dict['Euler'] = (\n",
"[[sp.sympify(0)],\n",
"[\"\", sp.sympify(1)]]\n",
", 1)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"### Step 2.a.ii: RK2 Heun's Method [Back to [top](#toc)\\]\n",
"$$\\label{rktwoheun}$$\n",
"\n",
"[Heun's method](https://en.wikipedia.org/w/index.php?title=Heun%27s_method&oldid=866896936) is a second-order RK method that obtains the solution $y(t)$ at time $t_{n+1}$ from $t_n$ via:\n",
"\\begin{align}\n",
"k_1 &= \\Delta tf(y_n, t_n), \\\\\n",
"k_2 &= \\Delta tf(y_n + k_1, t_n + \\Delta t), \\\\\n",
"y_{n+1} &= y_n + \\frac{1}{2}(k_1 + k_2) + \\mathcal{O}\\big((\\Delta t)^3\\big).\n",
"\\end{align}\n",
"with corresponding Butcher table\n",
"$$\\begin{array}{c|cc}\n",
" 0 & \\\\\n",
" 1 & 1 & \\\\ \\hline\n",
" & 1/2 & 1/2\n",
"\\end{array} $$\n"
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {
"execution": {
"iopub.execute_input": "2021-03-07T17:16:28.485084Z",
"iopub.status.busy": "2021-03-07T17:16:28.483946Z",
"iopub.status.idle": "2021-03-07T17:16:28.487509Z",
"shell.execute_reply": "2021-03-07T17:16:28.486801Z"
}
},
"outputs": [],
"source": [
"# Step 2.a.ii: RK2 Heun's Method\n",
"\n",
"Butcher_dict['RK2 Heun'] = (\n",
"[[sp.sympify(0)],\n",
"[sp.sympify(1), sp.sympify(1)],\n",
"[\"\", sp.Rational(1,2), sp.Rational(1,2)]]\n",
", 2)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"### Step 2.a.iii: RK2 Midpoint Method [Back to [top](#toc)\\]\n",
"$$\\label{rk2mp}$$\n",
"\n",
"[Midpoint method](https://en.wikipedia.org/w/index.php?title=Midpoint_method&oldid=886630580) is a second-order RK method that obtains the solution $y(t)$ at time $t_{n+1}$ from $t_n$ via:\n",
"\\begin{align}\n",
"k_1 &= \\Delta tf(y_n, t_n), \\\\\n",
"k_2 &= \\Delta tf(y_n + \\frac{2}{3}k_1, t_n + \\frac{2}{3}\\Delta t), \\\\\n",
"y_{n+1} &= y_n + \\frac{1}{2}k_2 + \\mathcal{O}\\big((\\Delta t)^3\\big).\n",
"\\end{align}\n",
"with corresponding Butcher table\n",
"$$\\begin{array}{c|cc}\n",
" 0 & \\\\\n",
" 1/2 & 1/2 & \\\\ \\hline\n",
" & 0 & 1\n",
"\\end{array} $$\n"
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {
"execution": {
"iopub.execute_input": "2021-03-07T17:16:28.494082Z",
"iopub.status.busy": "2021-03-07T17:16:28.493180Z",
"iopub.status.idle": "2021-03-07T17:16:28.495812Z",
"shell.execute_reply": "2021-03-07T17:16:28.496419Z"
}
},
"outputs": [],
"source": [
"# Step 2.a.iii: RK2 Midpoint (MP) Method\n",
"\n",
"Butcher_dict['RK2 MP'] = (\n",
"[[sp.sympify(0)],\n",
"[sp.Rational(1,2), sp.Rational(1,2)],\n",
"[\"\", sp.sympify(0), sp.sympify(1)]]\n",
", 2)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"### Step 2.a.iv: RK2 Ralston's Method [Back to [top](#toc)\\]\n",
"$$\\label{rk2ralston}$$\n",
"\n",
"Ralston's method (see [Ralston (1962)](https://www.ams.org/journals/mcom/1962-16-080/S0025-5718-1962-0150954-0/S0025-5718-1962-0150954-0.pdf), is a second-order RK method that obtains the solution $y(t)$ at time $t_{n+1}$ from $t_n$ via:\n",
"\\begin{align}\n",
"k_1 &= \\Delta tf(y_n, t_n), \\\\\n",
"k_2 &= \\Delta tf(y_n + \\frac{1}{2}k_1, t_n + \\frac{1}{2}\\Delta t), \\\\\n",
"y_{n+1} &= y_n + \\frac{1}{4}k_1 + \\frac{3}{4}k_2 + \\mathcal{O}\\big((\\Delta t)^3\\big).\n",
"\\end{align}\n",
"with corresponding Butcher table\n",
"$$\\begin{array}{c|cc}\n",
" 0 & \\\\\n",
" 2/3 & 2/3 & \\\\ \\hline\n",
" & 1/4 & 3/4\n",
"\\end{array} $$"
]
},
{
"cell_type": "code",
"execution_count": 6,
"metadata": {
"execution": {
"iopub.execute_input": "2021-03-07T17:16:28.503658Z",
"iopub.status.busy": "2021-03-07T17:16:28.502409Z",
"iopub.status.idle": "2021-03-07T17:16:28.506091Z",
"shell.execute_reply": "2021-03-07T17:16:28.505355Z"
}
},
"outputs": [],
"source": [
"# Step 2.a.iv: RK2 Ralston's Method\n",
"\n",
"Butcher_dict['RK2 Ralston'] = (\n",
"[[sp.sympify(0)],\n",
"[sp.Rational(2,3), sp.Rational(2,3)],\n",
"[\"\", sp.Rational(1,4), sp.Rational(3,4)]]\n",
", 2)\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"### Step 2.a.v: Kutta's Third-order Method [Back to [top](#toc)\\]\n",
"$$\\label{rk3}$$\n",
"\n",
"[Kutta's third-order method](https://en.wikipedia.org/w/index.php?title=List_of_Runge%E2%80%93Kutta_methods&oldid=896594269) obtains the solution $y(t)$ at time $t_{n+1}$ from $t_n$ via:\n",
"\\begin{align}\n",
"k_1 &= \\Delta tf(y_n, t_n), \\\\\n",
"k_2 &= \\Delta tf(y_n + \\frac{1}{2}k_1, t_n + \\frac{1}{2}\\Delta t), \\\\\n",
"k_3 &= \\Delta tf(y_n - k_1 + 2k_2, t_n + \\Delta t) \\\\\n",
"y_{n+1} &= y_n + \\frac{1}{6}k_1 + \\frac{2}{3}k_2 + \\frac{1}{6}k_3 + \\mathcal{O}\\big((\\Delta t)^4\\big).\n",
"\\end{align}\n",
"with corresponding Butcher table\n",
"\\begin{array}{c|ccc}\n",
" 0 & \\\\\n",
" 1/2 & 1/2 & \\\\\n",
" 1 & -1 & 2 & \\\\ \\hline\n",
" & 1/6 & 2/3 & 1/6\n",
"\\end{array}"
]
},
{
"cell_type": "code",
"execution_count": 7,
"metadata": {
"execution": {
"iopub.execute_input": "2021-03-07T17:16:28.514852Z",
"iopub.status.busy": "2021-03-07T17:16:28.513555Z",
"iopub.status.idle": "2021-03-07T17:16:28.517363Z",
"shell.execute_reply": "2021-03-07T17:16:28.516669Z"
}
},
"outputs": [],
"source": [
"# Step 2.a.v: Kutta's Third-order Method\n",
"\n",
"Butcher_dict['RK3'] = (\n",
"[[sp.sympify(0)],\n",
"[sp.Rational(1,2), sp.Rational(1,2)],\n",
"[sp.sympify(1), sp.sympify(-1), sp.sympify(2)],\n",
"[\"\", sp.Rational(1,6), sp.Rational(2,3), sp.Rational(1,6)]]\n",
", 3)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"### Step 2.a.vi: RK3 Heun's Method [Back to [top](#toc)\\]\n",
"$$\\label{rk3heun}$$\n",
"\n",
"[Heun's third-order method](https://en.wikipedia.org/w/index.php?title=List_of_Runge%E2%80%93Kutta_methods&oldid=896594269) obtains the solution $y(t)$ at time $t_{n+1}$ from $t_n$ via:\n",
"\n",
"\\begin{align}\n",
"k_1 &= \\Delta tf(y_n, t_n), \\\\\n",
"k_2 &= \\Delta tf(y_n + \\frac{1}{3}k_1, t_n + \\frac{1}{3}\\Delta t), \\\\\n",
"k_3 &= \\Delta tf(y_n + \\frac{2}{3}k_2, t_n + \\frac{2}{3}\\Delta t) \\\\\n",
"y_{n+1} &= y_n + \\frac{1}{4}k_1 + \\frac{3}{4}k_3 + \\mathcal{O}\\big((\\Delta t)^4\\big).\n",
"\\end{align}\n",
"\n",
"with corresponding Butcher table\n",
"\n",
"\\begin{array}{c|ccc}\n",
" 0 & \\\\\n",
" 1/3 & 1/3 & \\\\\n",
" 2/3 & 0 & 2/3 & \\\\ \\hline\n",
" & 1/4 & 0 & 3/4\n",
"\\end{array}\n"
]
},
{
"cell_type": "code",
"execution_count": 8,
"metadata": {
"execution": {
"iopub.execute_input": "2021-03-07T17:16:28.524600Z",
"iopub.status.busy": "2021-03-07T17:16:28.523786Z",
"iopub.status.idle": "2021-03-07T17:16:28.526890Z",
"shell.execute_reply": "2021-03-07T17:16:28.526238Z"
}
},
"outputs": [],
"source": [
"# Step 2.a.vi: RK3 Heun's Method\n",
"\n",
"Butcher_dict['RK3 Heun'] = (\n",
"[[sp.sympify(0)],\n",
"[sp.Rational(1,3), sp.Rational(1,3)],\n",
"[sp.Rational(2,3), sp.sympify(0), sp.Rational(2,3)],\n",
"[\"\", sp.Rational(1,4), sp.sympify(0), sp.Rational(3,4)]]\n",
", 3)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"### Step 2.a.vii: RK3 Ralton's Method [Back to [top](#toc)\\]\n",
"$$\\label{rk3ralston}$$\n",
"\n",
"Ralston's third-order method (see [Ralston (1962)](https://www.ams.org/journals/mcom/1962-16-080/S0025-5718-1962-0150954-0/S0025-5718-1962-0150954-0.pdf), obtains the solution $y(t)$ at time $t_{n+1}$ from $t_n$ via:\n",
"\n",
"\\begin{align}\n",
"k_1 &= \\Delta tf(y_n, t_n), \\\\\n",
"k_2 &= \\Delta tf(y_n + \\frac{1}{2}k_1, t_n + \\frac{1}{2}\\Delta t), \\\\\n",
"k_3 &= \\Delta tf(y_n + \\frac{3}{4}k_2, t_n + \\frac{3}{4}\\Delta t) \\\\\n",
"y_{n+1} &= y_n + \\frac{2}{9}k_1 + \\frac{1}{3}k_2 + \\frac{4}{9}k_3 + \\mathcal{O}\\big((\\Delta t)^4\\big).\n",
"\\end{align}\n",
"\n",
"with corresponding Butcher table\n",
"\n",
"\\begin{array}{c|ccc}\n",
" 0 & \\\\\n",
" 1/2 & 1/2 & \\\\\n",
" 3/4 & 0 & 3/4 & \\\\ \\hline\n",
" & 2/9 & 1/3 & 4/9\n",
"\\end{array}"
]
},
{
"cell_type": "code",
"execution_count": 9,
"metadata": {
"execution": {
"iopub.execute_input": "2021-03-07T17:16:28.535795Z",
"iopub.status.busy": "2021-03-07T17:16:28.534838Z",
"iopub.status.idle": "2021-03-07T17:16:28.537798Z",
"shell.execute_reply": "2021-03-07T17:16:28.538447Z"
}
},
"outputs": [],
"source": [
"# Step 2.a.vii: RK3 Ralton's Method\n",
"\n",
"Butcher_dict['RK3 Ralston'] = (\n",
"[[0],\n",
"[sp.Rational(1,2), sp.Rational(1,2)],\n",
"[sp.Rational(3,4), sp.sympify(0), sp.Rational(3,4)],\n",
"[\"\", sp.Rational(2,9), sp.Rational(1,3), sp.Rational(4,9)]]\n",
", 3)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"### Step 2.a.viii: Strong Stability Preserving Runge-Kutta (SSPRK3) Method [Back to [top](#toc)\\]\n",
"$\\label{ssprk3}$\n",
"\n",
"The [Strong Stability Preserving Runge-Kutta (SSPRK3)](https://en.wikipedia.org/wiki/List_of_Runge%E2%80%93Kutta_methods#Kutta's_third-order_method) method obtains the solution $y(t)$ at time $t_{n+1}$ from $t_n$ via:\n",
"\n",
"\\begin{align}\n",
"k_1 &= \\Delta tf(y_n, t_n), \\\\\n",
"k_2 &= \\Delta tf(y_n + k_1, t_n + \\Delta t), \\\\\n",
"k_3 &= \\Delta tf(y_n + \\frac{1}{4}k_1 + \\frac{1}{4}k_2, t_n + \\frac{1}{2}\\Delta t) \\\\\n",
"y_{n+1} &= y_n + \\frac{1}{6}k_1 + \\frac{1}{6}k_2 + \\frac{2}{3}k_3 + \\mathcal{O}\\big((\\Delta t)^4\\big).\n",
"\\end{align}\n",
"\n",
"with corresponding Butcher table\n",
"\n",
"\\begin{array}{c|ccc}\n",
" 0 & \\\\\n",
" 1 & 1 & \\\\\n",
" 1/2 & 1/4 & 1/4 & \\\\ \\hline\n",
" & 1/6 & 1/6 & 2/3\n",
"\\end{array}\n"
]
},
{
"cell_type": "code",
"execution_count": 10,
"metadata": {
"execution": {
"iopub.execute_input": "2021-03-07T17:16:28.544883Z",
"iopub.status.busy": "2021-03-07T17:16:28.543893Z",
"iopub.status.idle": "2021-03-07T17:16:28.546648Z",
"shell.execute_reply": "2021-03-07T17:16:28.546124Z"
}
},
"outputs": [],
"source": [
"# Step 2.a.viii: Strong Stability Preserving Runge-Kutta (SSPRK3) Method\n",
"\n",
"Butcher_dict['SSPRK3'] = (\n",
"[[0],\n",
"[sp.sympify(1), sp.sympify(1)],\n",
"[sp.Rational(1,2), sp.Rational(1,4), sp.Rational(1,4)],\n",
"[\"\", sp.Rational(1,6), sp.Rational(1,6), sp.Rational(2,3)]]\n",
", 3)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"### Step 2.a.ix: Classic RK4 Method [Back to [top](#toc)\\]\n",
"$$\\label{rkfour}$$\n",
"\n",
"The [classic RK4 method](https://en.wikipedia.org/w/index.php?title=Runge%E2%80%93Kutta_methods&oldid=894771467) obtains the solution $y(t)$ at time $t_{n+1}$ from $t_n$ via:\n",
"\n",
"\\begin{align}\n",
"k_1 &= \\Delta tf(y_n, t_n), \\\\\n",
"k_2 &= \\Delta tf(y_n + \\frac{1}{2}k_1, t_n + \\frac{\\Delta t}{2}), \\\\\n",
"k_3 &= \\Delta tf(y_n + \\frac{1}{2}k_2, t_n + \\frac{\\Delta t}{2}), \\\\\n",
"k_4 &= \\Delta tf(y_n + k_3, t_n + \\Delta t), \\\\\n",
"y_{n+1} &= y_n + \\frac{1}{6}(k_1 + 2k_2 + 2k_3 + k_4) + \\mathcal{O}\\big((\\Delta t)^5\\big).\n",
"\\end{align}\n",
"\n",
"with corresponding Butcher table\n",
"\n",
"$$\\begin{array}{c|cccc}\n",
" 0 & \\\\\n",
" 1/2 & 1/2 & \\\\ \n",
" 1/2 & 0 & 1/2 & \\\\\n",
" 1 & 0 & 0 & 1 & \\\\ \\hline\n",
" & 1/6 & 1/3 & 1/3 & 1/6\n",
"\\end{array} $$"
]
},
{
"cell_type": "code",
"execution_count": 11,
"metadata": {
"execution": {
"iopub.execute_input": "2021-03-07T17:16:28.554487Z",
"iopub.status.busy": "2021-03-07T17:16:28.553348Z",
"iopub.status.idle": "2021-03-07T17:16:28.555917Z",
"shell.execute_reply": "2021-03-07T17:16:28.556503Z"
}
},
"outputs": [],
"source": [
"# Step 2.a.vix: Classic RK4 Method\n",
"\n",
"Butcher_dict['RK4'] = (\n",
"[[sp.sympify(0)],\n",
"[sp.Rational(1,2), sp.Rational(1,2)],\n",
"[sp.Rational(1,2), sp.sympify(0), sp.Rational(1,2)],\n",
"[sp.sympify(1), sp.sympify(0), sp.sympify(0), sp.sympify(1)],\n",
"[\"\", sp.Rational(1,6), sp.Rational(1,3), sp.Rational(1,3), sp.Rational(1,6)]]\n",
", 4)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"### Step 2.a.x: RK5 Dormand-Prince Method [Back to [top](#toc)\\]\n",
"$$\\label{dp5}$$\n",
"\n",
"The fifth-order Dormand-Prince (DP) method from the RK5(4) family (see [Dormand, J. R.; Prince, P. J. (1980)](https://www.sciencedirect.com/science/article/pii/0771050X80900133?via%3Dihub)) Butcher table is:\n",
"\n",
"$$\\begin{array}{c|ccccccc}\n",
" 0 & \\\\\n",
" \\frac{1}{5} & \\frac{1}{5} & \\\\ \n",
" \\frac{3}{10} & \\frac{3}{40} & \\frac{9}{40} & \\\\\n",
" \\frac{4}{5} & \\frac{44}{45} & \\frac{-56}{15} & \\frac{32}{9} & \\\\ \n",
" \\frac{8}{9} & \\frac{19372}{6561} & \\frac{−25360}{2187} & \\frac{64448}{6561} & \\frac{−212}{729} & \\\\\n",
" 1 & \\frac{9017}{3168} & \\frac{−355}{33} & \\frac{46732}{5247} & \\frac{49}{176} & \\frac{−5103}{18656} & \\\\\n",
" 1 & \\frac{35}{384} & 0 & \\frac{500}{1113} & \\frac{125}{192} & \\frac{−2187}{6784} & \\frac{11}{84} & \\\\ \\hline\n",
" & \\frac{35}{384} & 0 & \\frac{500}{1113} & \\frac{125}{192} & \\frac{−2187}{6784} & \\frac{11}{84} & 0\n",
"\\end{array} $$"
]
},
{
"cell_type": "code",
"execution_count": 12,
"metadata": {
"execution": {
"iopub.execute_input": "2021-03-07T17:16:28.570783Z",
"iopub.status.busy": "2021-03-07T17:16:28.569438Z",
"iopub.status.idle": "2021-03-07T17:16:28.572740Z",
"shell.execute_reply": "2021-03-07T17:16:28.572151Z"
}
},
"outputs": [],
"source": [
"# Step 2.a.x: RK5 Dormand-Prince Method\n",
"\n",
"Butcher_dict['DP5'] = (\n",
"[[0],\n",
"[sp.Rational(1,5), sp.Rational(1,5)],\n",
"[sp.Rational(3,10),sp.Rational(3,40), sp.Rational(9,40)],\n",
"[sp.Rational(4,5), sp.Rational(44,45), sp.Rational(-56,15), sp.Rational(32,9)],\n",
"[sp.Rational(8,9), sp.Rational(19372,6561), sp.Rational(-25360,2187), sp.Rational(64448,6561), sp.Rational(-212,729)],\n",
"[sp.sympify(1), sp.Rational(9017,3168), sp.Rational(-355,33), sp.Rational(46732,5247), sp.Rational(49,176), sp.Rational(-5103,18656)],\n",
"[sp.sympify(1), sp.Rational(35,384), sp.sympify(0), sp.Rational(500,1113), sp.Rational(125,192), sp.Rational(-2187,6784), sp.Rational(11,84)],\n",
"[\"\", sp.Rational(35,384), sp.sympify(0), sp.Rational(500,1113), sp.Rational(125,192), sp.Rational(-2187,6784), sp.Rational(11,84), sp.sympify(0)]]\n",
", 5)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"### Step 2.a.xi: RK5 Dormand-Prince Method Alternative [Back to [top](#toc)\\]\n",
"$$\\label{dp5alt}$$\n",
"\n",
"The fifth-order Dormand-Prince (DP) method from the RK6(5) family (see [Dormand, J. R.; Prince, P. J. (1981)](https://www.sciencedirect.com/science/article/pii/0771050X81900103)) Butcher table is:\n",
"\n",
"$$\\begin{array}{c|ccccccc}\n",
" 0 & \\\\\n",
" \\frac{1}{10} & \\frac{1}{10} & \\\\\n",
" \\frac{2}{9} & \\frac{-2}{81} & \\frac{20}{81} & \\\\\n",
" \\frac{3}{7} & \\frac{615}{1372} & \\frac{-270}{343} & \\frac{1053}{1372} & \\\\\n",
" \\frac{3}{5} & \\frac{3243}{5500} & \\frac{-54}{55} & \\frac{50949}{71500} & \\frac{4998}{17875} & \\\\\n",
" \\frac{4}{5} & \\frac{-26492}{37125} & \\frac{72}{55} & \\frac{2808}{23375} & \\frac{-24206}{37125} & \\frac{338}{459} & \\\\\n",
" 1 & \\frac{5561}{2376} & \\frac{-35}{11} & \\frac{-24117}{31603} & \\frac{899983}{200772} & \\frac{-5225}{1836} & \\frac{3925}{4056} & \\\\ \\hline\n",
" & \\frac{821}{10800} & 0 & \\frac{19683}{71825} & \\frac{175273}{912600} & \\frac{395}{3672} & \\frac{785}{2704} & \\frac{3}{50}\n",
"\\end{array}$$"
]
},
{
"cell_type": "code",
"execution_count": 13,
"metadata": {
"execution": {
"iopub.execute_input": "2021-03-07T17:16:28.584591Z",
"iopub.status.busy": "2021-03-07T17:16:28.583355Z",
"iopub.status.idle": "2021-03-07T17:16:28.587317Z",
"shell.execute_reply": "2021-03-07T17:16:28.586555Z"
}
},
"outputs": [],
"source": [
"# Step 2.a.xi: RK5 Dormand-Prince Method Alternative\n",
"\n",
"Butcher_dict['DP5alt'] = (\n",
"[[0],\n",
"[sp.Rational(1,10), sp.Rational(1,10)],\n",
"[sp.Rational(2,9), sp.Rational(-2, 81), sp.Rational(20, 81)],\n",
"[sp.Rational(3,7), sp.Rational(615, 1372), sp.Rational(-270, 343), sp.Rational(1053, 1372)],\n",
"[sp.Rational(3,5), sp.Rational(3243, 5500), sp.Rational(-54, 55), sp.Rational(50949, 71500), sp.Rational(4998, 17875)],\n",
"[sp.Rational(4, 5), sp.Rational(-26492, 37125), sp.Rational(72, 55), sp.Rational(2808, 23375), sp.Rational(-24206, 37125), sp.Rational(338, 459)],\n",
"[sp.sympify(1), sp.Rational(5561, 2376), sp.Rational(-35, 11), sp.Rational(-24117, 31603), sp.Rational(899983, 200772), sp.Rational(-5225, 1836), sp.Rational(3925, 4056)],\n",
"[\"\", sp.Rational(821, 10800), sp.sympify(0), sp.Rational(19683, 71825), sp.Rational(175273, 912600), sp.Rational(395, 3672), sp.Rational(785, 2704), sp.Rational(3, 50)]]\n",
", 5)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"### Step 2.a.xii: RK5 Cash-Karp Method [Back to [top](#toc)\\]\n",
"$$\\label{ck5}$$\n",
"\n",
"The fifth-order Cash-Karp Method (see [J. R. Cash, A. H. Karp. (1980)](https://dl.acm.org/citation.cfm?doid=79505.79507)) Butcher table is:\n",
"\n",
"$$\\begin{array}{c|cccccc}\n",
" 0 & \\\\ \n",
"\t\\frac{1}{5} & \\frac{1}{5} & \\\\ \n",
"\t\\frac{3}{10} & \\frac{3}{40} & \\frac{9}{40} & \\\\ \n",
"\t\\frac{3}{5} & \\frac{3}{10} & \\frac{−9}{10} & \\frac{6}{5} & \\\\ \n",
"\t1 & \\frac{−11}{54} & \\frac{5}{2} & \\frac{−70}{27} & \\frac{35}{27} & \\\\ \n",
"\t\\frac{7}{8} & \\frac{1631}{55296} & \\frac{175}{512} & \\frac{575}{13824} & \\frac{44275}{110592} & \\frac{253}{4096} & \\\\ \\hline\n",
"\t & \\frac{37}{378} & 0 & \\frac{250}{621} & \\frac{125}{594} & 0 & \\frac{512}{1771} \n",
"\\end{array}$$\n",
"\n",
"\n"
]
},
{
"cell_type": "code",
"execution_count": 14,
"metadata": {
"execution": {
"iopub.execute_input": "2021-03-07T17:16:28.598595Z",
"iopub.status.busy": "2021-03-07T17:16:28.597930Z",
"iopub.status.idle": "2021-03-07T17:16:28.600613Z",
"shell.execute_reply": "2021-03-07T17:16:28.600084Z"
}
},
"outputs": [],
"source": [
"# Step 2.a.xii: RK5 Cash-Karp Method\n",
"\n",
"Butcher_dict['CK5'] = (\n",
"[[0],\n",
"[sp.Rational(1,5), sp.Rational(1,5)],\n",
"[sp.Rational(3,10),sp.Rational(3,40), sp.Rational(9,40)],\n",
"[sp.Rational(3,5), sp.Rational(3,10), sp.Rational(-9,10), sp.Rational(6,5)],\n",
"[sp.sympify(1), sp.Rational(-11,54), sp.Rational(5,2), sp.Rational(-70,27), sp.Rational(35,27)],\n",
"[sp.Rational(7,8), sp.Rational(1631,55296), sp.Rational(175,512), sp.Rational(575,13824), sp.Rational(44275,110592), sp.Rational(253,4096)],\n",
"[\"\",sp.Rational(37,378), sp.sympify(0), sp.Rational(250,621), sp.Rational(125,594), sp.sympify(0), sp.Rational(512,1771)]]\n",
", 5)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"### Step 2.a.xiii: RK6 Dormand-Prince Method [Back to [top](#toc)\\]\n",
"$$\\label{dp6}$$\n",
"\n",
"The sixth-order Dormand-Prince method (see [Dormand, J. R.; Prince, P. J. (1981)](https://www.sciencedirect.com/science/article/pii/0771050X81900103)) Butcher Table is\n",
"\n",
"\n",
"$$\\begin{array}{c|cccccccc}\n",
" 0 & \\\\\n",
" \\frac{1}{10} & \\frac{1}{10} & \\\\\n",
" \\frac{2}{9} & \\frac{-2}{81} & \\frac{20}{81} & \\\\\n",
" \\frac{3}{7} & \\frac{615}{1372} & \\frac{-270}{343} & \\frac{1053}{1372} & \\\\\n",
" \\frac{3}{5} & \\frac{3243}{5500} & \\frac{-54}{55} & \\frac{50949}{71500} & \\frac{4998}{17875} & \\\\\n",
" \\frac{4}{5} & \\frac{-26492}{37125} & \\frac{72}{55} & \\frac{2808}{23375} & \\frac{-24206}{37125} & \\frac{338}{459} & \\\\\n",
" 1 & \\frac{5561}{2376} & \\frac{-35}{11} & \\frac{-24117}{31603} & \\frac{899983}{200772} & \\frac{-5225}{1836} & \\frac{3925}{4056} & \\\\\n",
" 1 & \\frac{465467}{266112} & \\frac{-2945}{1232} & \\frac{-5610201}{14158144} & \\frac{10513573}{3212352} & \\frac{-424325}{205632} & \\frac{376225}{454272} & 0 & \\\\ \\hline\n",
" & \\frac{61}{864} & 0 & \\frac{98415}{321776} & \\frac{16807}{146016} & \\frac{1375}{7344} & \\frac{1375}{5408} & \\frac{-37}{1120} & \\frac{1}{10}\n",
"\\end{array}$$\n",
"\n"
]
},
{
"cell_type": "code",
"execution_count": 15,
"metadata": {
"execution": {
"iopub.execute_input": "2021-03-07T17:16:28.613220Z",
"iopub.status.busy": "2021-03-07T17:16:28.612528Z",
"iopub.status.idle": "2021-03-07T17:16:28.615243Z",
"shell.execute_reply": "2021-03-07T17:16:28.614681Z"
}
},
"outputs": [],
"source": [
"# Step 2.a.xiii: RK6 Dormand-Prince Method\n",
"\n",
"Butcher_dict['DP6'] = (\n",
"[[0],\n",
"[sp.Rational(1,10), sp.Rational(1,10)],\n",
"[sp.Rational(2,9), sp.Rational(-2, 81), sp.Rational(20, 81)],\n",
"[sp.Rational(3,7), sp.Rational(615, 1372), sp.Rational(-270, 343), sp.Rational(1053, 1372)],\n",
"[sp.Rational(3,5), sp.Rational(3243, 5500), sp.Rational(-54, 55), sp.Rational(50949, 71500), sp.Rational(4998, 17875)],\n",
"[sp.Rational(4, 5), sp.Rational(-26492, 37125), sp.Rational(72, 55), sp.Rational(2808, 23375), sp.Rational(-24206, 37125), sp.Rational(338, 459)],\n",
"[sp.sympify(1), sp.Rational(5561, 2376), sp.Rational(-35, 11), sp.Rational(-24117, 31603), sp.Rational(899983, 200772), sp.Rational(-5225, 1836), sp.Rational(3925, 4056)],\n",
"[sp.sympify(1), sp.Rational(465467, 266112), sp.Rational(-2945, 1232), sp.Rational(-5610201, 14158144), sp.Rational(10513573, 3212352), sp.Rational(-424325, 205632), sp.Rational(376225, 454272), sp.sympify(0)],\n",
"[\"\", sp.Rational(61, 864), sp.sympify(0), sp.Rational(98415, 321776), sp.Rational(16807, 146016), sp.Rational(1375, 7344), sp.Rational(1375, 5408), sp.Rational(-37, 1120), sp.Rational(1,10)]]\n",
", 6)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"### Step 2.a.xiv: RK6 Luther's Method [Back to [top](#toc)\\]\n",
"$$\\label{l6}$$\n",
"\n",
"Luther's sixth-order method (see [H. A. Luther (1968)](http://www.ams.org/journals/mcom/1968-22-102/S0025-5718-68-99876-1/S0025-5718-68-99876-1.pdf)) Butcher table is:\n",
"$$\\begin{array}{c|ccccccc}\n",
" 0 & \\\\\n",
" 1 & 1 & \\\\\n",
" \\frac{1}{2} & \\frac{3}{8} & \\frac{1}{8} & \\\\\n",
" \\frac{2}{3} & \\frac{8}{27} & \\frac{2}{27} & \\frac{8}{27} & \\\\\n",
" \\frac{(7-q)}{14} & \\frac{(-21 + 9q)}{392} & \\frac{(-56 + 8q)}{392} & \\frac{(336 - 48q)}{392} & \\frac{(-63 + 3q)}{392} & \\\\\n",
" \\frac{(7+q)}{14} & \\frac{(-1155 - 255q)}{1960} & \\frac{(-280 - 40q)}{1960} & \\frac{320q}{1960} & \\frac{(63 + 363q)}{1960} & \\frac{(2352 + 392q)}{1960} & \\\\\n",
" 1 & \\frac{(330 + 105q)}{180} & \\frac{2}{3} & \\frac{(-200 + 280q)}{180} & \\frac{(126 - 189q)}{180} & \\frac{(-686 - 126q)}{180} & \\frac{(490 - 70q)}{180} & \\\\ \\hline\n",
" & \\frac{1}{20} & 0 & \\frac{16}{45} & 0 & \\frac{49}{180} & \\frac{49}{180} & \\frac{1}{20}\n",
"\\end{array}$$\n",
"\n",
"where $q = \\sqrt{21}$."
]
},
{
"cell_type": "code",
"execution_count": 16,
"metadata": {
"execution": {
"iopub.execute_input": "2021-03-07T17:16:28.626801Z",
"iopub.status.busy": "2021-03-07T17:16:28.626106Z",
"iopub.status.idle": "2021-03-07T17:16:28.699010Z",
"shell.execute_reply": "2021-03-07T17:16:28.699497Z"
}
},
"outputs": [],
"source": [
"# Step 2.a.xiv: RK6 Luther's Method\n",
"\n",
"q = sp.sqrt(21)\n",
"Butcher_dict['L6'] = (\n",
"[[0],\n",
"[sp.sympify(1), sp.sympify(1)],\n",
"[sp.Rational(1,2), sp.Rational(3,8), sp.Rational(1,8)],\n",
"[sp.Rational(2,3), sp.Rational(8,27), sp.Rational(2,27), sp.Rational(8,27)],\n",
"[(7 - q)/14, (-21 + 9*q)/392, (-56 + 8*q)/392, (336 -48*q)/392, (-63 + 3*q)/392],\n",
"[(7 + q)/14, (-1155 - 255*q)/1960, (-280 - 40*q)/1960, (-320*q)/1960, (63 + 363*q)/1960, (2352 + 392*q)/1960],\n",
"[sp.sympify(1), ( 330 + 105*q)/180, sp.Rational(2,3), (-200 + 280*q)/180, (126 - 189*q)/180, (-686 - 126*q)/180, (490 - 70*q)/180],\n",
"[\"\", sp.Rational(1, 20), sp.sympify(0), sp.Rational(16, 45), sp.sympify(0), sp.Rational(49, 180), sp.Rational(49, 180), sp.Rational(1, 20)]]\n",
", 6)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"### Step 2.a.xv: RK8 Dormand-Prince Method [Back to [top](#toc)\\]\n",
"$$\\label{dp8}$$\n",
"\n",
"The eighth-order Dormand-Prince Method (see [Dormand, J. R.; Prince, P. J. (1981)](https://www.sciencedirect.com/science/article/pii/0771050X81900103)) Butcher table is:\n",
"\n",
"$$\\begin{array}{c|ccccccccc}\n",
" 0 & \\\\\n",
" \\frac{1}{18} & \\frac{1}{18} & \\\\\n",
" \\frac{1}{12} & \\frac{1}{48} & \\frac{1}{16} & \\\\\n",
" \\frac{1}{8} & \\frac{1}{32} & 0 & \\frac{3}{32} & \\\\\n",
" \\frac{5}{16} & \\frac{5}{16} & 0 & \\frac{-75}{64} & \\frac{75}{64} & \\\\\n",
" \\frac{3}{8} & \\frac{3}{80} & 0 & 0 & \\frac{3}{16} & \\frac{3}{20} & \\\\\n",
" \\frac{59}{400} & \\frac{29443841}{614563906} & 0 & 0 & \\frac{77736538}{692538347} & \\frac{-28693883}{1125000000} & \\frac{23124283}{1800000000} & \\\\\n",
" \\frac{93}{200} & \\frac{16016141}{946692911} & 0 & 0 & \\frac{61564180}{158732637} & \\frac{22789713}{633445777} & \\frac{545815736}{2771057229} & \\frac{-180193667}{1043307555} & \\\\\n",
" \\frac{5490023248}{9719169821} & \\frac{39632708}{573591083} & 0 & 0 & \\frac{-433636366}{683701615} & \\frac{-421739975}{2616292301} & \\frac{100302831}{723423059} & \\frac{790204164}{839813087} & \\frac{800635310}{3783071287} & \\\\\n",
" \\frac{13}{20} & \\frac{246121993}{1340847787} & 0 & 0 & \\frac{-37695042795}{15268766246} & \\frac{-309121744}{1061227803} & \\frac{-12992083}{490766935} & \\frac{6005943493}{2108947869} & \\frac{393006217}{1396673457} & \\frac{123872331}{1001029789} & \\\\\n",
" \\frac{1201146811}{1299019798} & \\frac{-1028468189}{846180014} & 0 & 0 & \\frac{8478235783}{508512852} & \\frac{1311729495}{1432422823} & \\frac{-10304129995}{1701304382} & \\frac{-48777925059}{3047939560} & \\frac{15336726248}{1032824649} & \\frac{-45442868181}{3398467696} & \\frac{3065993473}{597172653} & \\\\\n",
" 1 & \\frac{185892177}{718116043} & 0 & 0 & \\frac{-3185094517}{667107341} & \\frac{-477755414}{1098053517} & \\frac{-703635378}{230739211} & \\frac{5731566787}{1027545527} & \\frac{5232866602}{850066563} & \\frac{-4093664535}{808688257} & \\frac{3962137247}{1805957418} & \\frac{65686358}{487910083} & \\\\\n",
" 1 & \\frac{403863854}{491063109} & 0 & 0 & \\frac{-5068492393}{434740067} & \\frac{-411421997}{543043805} & \\frac{652783627}{914296604} & \\frac{11173962825}{925320556} & \\frac{-13158990841}{6184727034} & \\frac{3936647629}{1978049680} & \\frac{-160528059}{685178525} & \\frac{248638103}{1413531060} & 0 & \\\\\n",
" & \\frac{14005451}{335480064} & 0 & 0 & 0 & 0 & \\frac{-59238493}{1068277825} & \\frac{181606767}{758867731} & \\frac{561292985}{797845732} & \\frac{-1041891430}{1371343529} & \\frac{760417239}{1151165299} & \\frac{118820643}{751138087} & \\frac{-528747749}{2220607170} & \\frac{1}{4}\n",
"\\end{array}$$\n",
"\n"
]
},
{
"cell_type": "code",
"execution_count": 17,
"metadata": {
"execution": {
"iopub.execute_input": "2021-03-07T17:16:28.722301Z",
"iopub.status.busy": "2021-03-07T17:16:28.718749Z",
"iopub.status.idle": "2021-03-07T17:16:28.724997Z",
"shell.execute_reply": "2021-03-07T17:16:28.724423Z"
}
},
"outputs": [],
"source": [
"# Step 2.a.xv: RK8 Dormand-Prince Method\n",
"\n",
"Butcher_dict['DP8']=(\n",
"[[0],\n",
"[sp.Rational(1, 18), sp.Rational(1, 18)],\n",
"[sp.Rational(1, 12), sp.Rational(1, 48), sp.Rational(1, 16)],\n",
"[sp.Rational(1, 8), sp.Rational(1, 32), sp.sympify(0), sp.Rational(3, 32)],\n",
"[sp.Rational(5, 16), sp.Rational(5, 16), sp.sympify(0), sp.Rational(-75, 64), sp.Rational(75, 64)],\n",
"[sp.Rational(3, 8), sp.Rational(3, 80), sp.sympify(0), sp.sympify(0), sp.Rational(3, 16), sp.Rational(3, 20)],\n",
"[sp.Rational(59, 400), sp.Rational(29443841, 614563906), sp.sympify(0), sp.sympify(0), sp.Rational(77736538, 692538347), sp.Rational(-28693883, 1125000000), sp.Rational(23124283, 1800000000)],\n",
"[sp.Rational(93, 200), sp.Rational(16016141, 946692911), sp.sympify(0), sp.sympify(0), sp.Rational(61564180, 158732637), sp.Rational(22789713, 633445777), sp.Rational(545815736, 2771057229), sp.Rational(-180193667, 1043307555)],\n",
"[sp.Rational(5490023248, 9719169821), sp.Rational(39632708, 573591083), sp.sympify(0), sp.sympify(0), sp.Rational(-433636366, 683701615), sp.Rational(-421739975, 2616292301), sp.Rational(100302831, 723423059), sp.Rational(790204164, 839813087), sp.Rational(800635310, 3783071287)],\n",
"[sp.Rational(13, 20), sp.Rational(246121993, 1340847787), sp.sympify(0), sp.sympify(0), sp.Rational(-37695042795, 15268766246), sp.Rational(-309121744, 1061227803), sp.Rational(-12992083, 490766935), sp.Rational(6005943493, 2108947869), sp.Rational(393006217, 1396673457), sp.Rational(123872331, 1001029789)],\n",
"[sp.Rational(1201146811, 1299019798), sp.Rational(-1028468189, 846180014), sp.sympify(0), sp.sympify(0), sp.Rational(8478235783, 508512852), sp.Rational(1311729495, 1432422823), sp.Rational(-10304129995, 1701304382), sp.Rational(-48777925059, 3047939560), sp.Rational(15336726248, 1032824649), sp.Rational(-45442868181, 3398467696), sp.Rational(3065993473, 597172653)],\n",
"[sp.sympify(1), sp.Rational(185892177, 718116043), sp.sympify(0), sp.sympify(0), sp.Rational(-3185094517, 667107341), sp.Rational(-477755414, 1098053517), sp.Rational(-703635378, 230739211), sp.Rational(5731566787, 1027545527), sp.Rational(5232866602, 850066563), sp.Rational(-4093664535, 808688257), sp.Rational(3962137247, 1805957418), sp.Rational(65686358, 487910083)],\n",
"[sp.sympify(1), sp.Rational(403863854, 491063109), sp.sympify(0), sp.sympify(0), sp.Rational(-5068492393, 434740067), sp.Rational(-411421997, 543043805), sp.Rational(652783627, 914296604), sp.Rational(11173962825, 925320556), sp.Rational(-13158990841, 6184727034), sp.Rational(3936647629, 1978049680), sp.Rational(-160528059, 685178525), sp.Rational(248638103, 1413531060), sp.sympify(0)],\n",
"[\"\", sp.Rational(14005451, 335480064), sp.sympify(0), sp.sympify(0), sp.sympify(0), sp.sympify(0), sp.Rational(-59238493, 1068277825), sp.Rational(181606767, 758867731), sp.Rational(561292985, 797845732), sp.Rational(-1041891430, 1371343529), sp.Rational(760417239, 1151165299), sp.Rational(118820643, 751138087), sp.Rational(-528747749, 2220607170), sp.Rational(1, 4)]]\n",
", 8)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"# Step 3: Code validation against `MoLtimestepping.RK_Butcher_Table_Dictionary` NRPy+ module [Back to [top](#toc)\\]\n",
"$$\\label{code_validation}$$\n",
"\n",
"As a code validation check, we verify agreement in the dictionary of Butcher tables between\n",
"1. this tutorial and \n",
"2. the NRPy+ [MoLtimestepping.RK_Butcher_Table_Dictionary](../edit/MoLtimestepping/RK_Butcher_Table_Dictionary.py) module.\n",
"\n",
"We analyze all key/value entries in the dictionary for consistency."
]
},
{
"cell_type": "code",
"execution_count": 18,
"metadata": {
"execution": {
"iopub.execute_input": "2021-03-07T17:16:28.730871Z",
"iopub.status.busy": "2021-03-07T17:16:28.730202Z",
"iopub.status.idle": "2021-03-07T17:16:28.736101Z",
"shell.execute_reply": "2021-03-07T17:16:28.735368Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"The dictionaries match!\n"
]
}
],
"source": [
"# Step 3: Code validation against MoLtimestepping.RK_Butcher_Table_Dictionary NRPy+ module\n",
"import sys # Standard Python module for multiplatform OS-level functions\n",
"from MoLtimestepping.RK_Butcher_Table_Dictionary import Butcher_dict as B_dict\n",
"valid = True\n",
"for key, value in Butcher_dict.items():\n",
" if Butcher_dict[key] != B_dict[key]:\n",
" valid = False\n",
" print(key)\n",
"if valid == True and len(Butcher_dict.items()) == len(B_dict.items()):\n",
" print(\"The dictionaries match!\")\n",
"else:\n",
" print(\"ERROR: Dictionaries don't match!\")\n",
" sys.exit(1)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"# Step 4: 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_Dictionary.pdf](Tutorial-RK_Butcher_Table_Dictionary.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": 19,
"metadata": {
"execution": {
"iopub.execute_input": "2021-03-07T17:16:28.741473Z",
"iopub.status.busy": "2021-03-07T17:16:28.740462Z",
"iopub.status.idle": "2021-03-07T17:16:33.239118Z",
"shell.execute_reply": "2021-03-07T17:16:33.238544Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Created Tutorial-RK_Butcher_Table_Dictionary.tex, and compiled LaTeX file\n",
" to PDF file Tutorial-RK_Butcher_Table_Dictionary.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_Dictionary\")"
]
}
],
"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
}