This nimib example document shows how to:
- apply (or not) a custom style (latex.css); you are looking at default style, for custom style click here
- using latex rendering of equations (thanks to katex)
- having as output an html table (using standard markdown converter to table)
The document itself shows how to use numericalnim to integrate an ODE.
Using NumericalNim
Example of usage of numericalnim.
import math, numericalnim
ODE
Example 3 from Paul's Online Notes
We want to solve the IVP (Initial Value Problem) for this linear first order differential equation:
$$y' - y = - \frac{1}{2} e^{\frac{t}{2}} \sin(5t) + 5e^{\frac{t}{2}}\cos(5t) \qquad y(0) = 0$$
This ODE has an analytical solution:
$$y(t) = e^{\frac{t}{2}}\sin(5t)$$
We want to find the approximation to the solution and compare it with the analytical solution at $t=5$.
We will use two fixed timestep methods to find the solution:
- Heun's 2nd order method (
Heun2
) - classic 4th order Runge-Kutta (
RK4
)
As timestamps we will use $h=0.1, 0.05, 0.01$.
First we translate the ODE as $y'=f(y,t)$ with $f$:
proc f(t, y: float, ctx: NumContext[float, float]): float =
y - 0.5*exp(0.5*t)*sin(5*t) + 5*exp(0.5*t)*cos(5*t)
proc y(t: float): float =
## analytical solution
exp(0.5*t)*sin(5*t)
let y0 = 0.0
## we will not be using the NumContext object
var ctx = newNumContext[float, float]()
## first derivative that will be used
echo "y'(0): ", f(0, y0, ctx) # "$y'(0)$" will not be converted to latex (katex has a protection not to look into code)
## expected solution
let y5 = y(5)
echo "y(5): ", y5
y'(0): 5.0 y(5): -1.612374396254655
We want the solution to be available for every point in $[0, 5]$ with timestep $h=0.05$
let tspan = arange(0.0, 5.0, 0.05, includeEnd=true)
echo tspan[0 .. 2], " . . . ", tspan[^3 .. ^1]
@[0.0, 0.05, 0.1] . . . @[4.9, 4.95, 5.0]
We compute the solution according to our selected 2 methods and 3 timesteps:
let
options1 = newODEoptions(dt = 0.1)
options2 = newODEoptions(dt = 0.05)
options3 = newODEoptions(dt = 0.01)
let
(t1hn, y1hn) = solveODE(f, y0, tspan, options = options1, integrator = "Heun2")
(t2hn, y2hn) = solveODE(f, y0, tspan, options = options2, integrator = "Heun2")
(t3hn, y3hn) = solveODE(f, y0, tspan, options = options3, integrator = "Heun2")
(t1rk, y1rk) = solveODE(f, y0, tspan, options = options1, integrator = "RK4")
(t2rk, y2rk) = solveODE(f, y0, tspan, options = options2, integrator = "RK4")
(t3rk, y3rk) = solveODE(f, y0, tspan, options = options3, integrator = "RK4")
## all returned time are the same as timespan (is this always true for fixed timestep methods?)
assert t1hn == tspan
assert t2hn == tspan
assert t3hn == tspan
assert t1rk == tspan
assert t2rk == tspan
assert t3rk == tspan
## solutions
echo "Heun2:"
echo fmt"{y1hn[^1]:.7f}"
echo fmt"{y2hn[^1]:.7f}"
echo fmt"{y3hn[^1]:.7f}"
echo "RK4:"
echo fmt"{y1rk[^1]:.7f}"
echo fmt"{y2rk[^1]:.7f}"
echo fmt"{y3rk[^1]:.7f}"
echo "Analytical:"
echo fmt"{y5:.7f}"
Heun2: -1.1273196 -1.4821550 -1.6068528 RK4: -1.6131927 -1.6124233 -1.6123745 Analytical: -1.6123744
As expected Heun is not very accurate even with smaller timesteps, while Runge-Kutta is very reliable even with bigger timesteps.
To compute percentage error of each method at $y(5)$ I will use:
proc pe(yApprox: float): string = fmt"{(100.0*abs(yApprox - y5) / abs(y5)):.3f}%"
echo pe y1hn[^1] ## used like this
30.083%
The following table is built as a Markdown table in a nbText
block.
Table 1. Percentage errors
Timestep | Heun2 | RK4 |
---|---|---|
$h=0.1$ | 30.083% | 0.051% |
$h=0.05$ | 8.076% | 0.003% |
$h=0.01$ | 0.342% | 0.000% |
and here is the code for this Markdown block:
fmt"""
The following table is built as a Markdown table in a `nbText` block.
**Table 1.** Percentage errors
Timestep | Heun2 | RK4
---------|--------------:|-------------:
$h=0.1$ | {pe y1hn[^1]} | {pe y1rk[^1]}
$h=0.05$ | {pe y2hn[^1]} | {pe y2rk[^1]}
$h=0.01$ | {pe y3hn[^1]} | {pe y3rk[^1]}
and here is the code for this Markdown block:
""" # right alignment of numbers does not seem to work. is this an issue of latex.css?