diff --git a/Lectures_my/NumMet/2016/Lecture9/lecture9.tex b/Lectures_my/NumMet/2016/Lecture9/lecture9.tex index 943c763..ed85bb2 100644 --- a/Lectures_my/NumMet/2016/Lecture9/lecture9.tex +++ b/Lectures_my/NumMet/2016/Lecture9/lecture9.tex @@ -92,7 +92,7 @@ \alert{linear} ordinary differential equations: \begin{overlayarea}{\textwidth}{.7\textheight} \begin{description} - \item<1-2>[differential equation] Any equation of the form + \item<1->[differential equation] Any equation of the form \begin{equation*} F(x; y, y', y^{(2)}, \dots, y^{(n)}) = 0 \end{equation*} @@ -104,9 +104,9 @@ \begin{equation*} F(x; y, \dots, y^{(n)}) = A(x) \cdot \vec{y}(x) + r(x) \end{equation*} - with $A(x) \in \mathbb{R}^{m\times n}$, $\vec{y}(x) = (y(x), \dots y^{(n)}(x))^T$ + with $A(x) \in \mathbb{R}^{m\times (n+1)}$, $\vec{y}(x) = (y(x), \dots y^{(n)}(x))^T$ and $r(x) \in \mathbb{R}^m$, then - $\Phi$ is a \alert{linear} system of ODEs. + $F$ describes a \alert{linear} system of ODEs. } \only<2>{ \item[ordinary] The quality \alert{ordinary} hinges on the fact that $y$ only @@ -125,18 +125,25 @@ \end{frame} \begin{frame}{Explicit linear ODEs of order $n=1$/Initial value problem} + Let's use a well-known linear ODE of 2nd order as an example: \begin{equation*} - y' = f(x; y) + y' = f(x; y) = -y(x) \end{equation*} - Initial value: $y_0 \equiv y(x_0)$ provided externally. + Solving it, regardless of method, poses an \alert{initial value problem}: + \begin{equation} + \frac{y(x)}{y(x_0)} = \exp(-(x - x_0)) + \end{equation} + We need to know the initial value $y_0\equiv y(x_0)$ to solve this problem. + It must be provided \alert{externally}. \end{frame} \begin{frame}{Single Step Methods} - We want to increment by one \alert{step}, i.e. from point $n$ to point $n + 1$. + We want to increment by one \alert{step}, i.e. from the point $(x_n, y_n)$ to the + point $(x_n + h, y_{n + 1})$. The increment is defined by the incremental function $\Phi$: \begin{equation*} - y_{n+1} = y_n + \Phi(f; x_0, h; y_n, y_{n+1}); + y_{n+1} = y_n + \Phi(F; x_n, h; y_n, y_{n+1}); \end{equation*} with a step size $h$. This is an implicit method. If $\Phi$ does not depend on $y_{n+1}$, we have an explicit method. @@ -153,18 +160,18 @@ \vspace{-2\bigskipamount} \begin{align*} y' & = f(x, y) & - y_{n+1} & = y_n + \Phi(f; x_0, h; y_n) + y_{n+1} & = y_n + \Phi(f; x_n, h; y_n) \end{align*} \end{block} Standard procedure for an unknown function: Taylor expand! \begin{align*} - y(x_n + h) & = y(x_n) + h y'(x_n) + \frac{h^2}{2} y''(x_n) + O{h^3} + y(x_n + h) & = y(x_n) + h y'(x_n) + \frac{h^2}{2} y''(x_n) + O(h^3) \end{align*} Keeping only terms of order $h$, and substituting $y'(x_n) = f(x_n = x_0 + n\cdot h, y_n)$ we obtain \begin{block}{Euler method} \begin{align*} - y_{n+1} & = y_n + \underbrace{h f(x_n, y_n)}_{=\Phi_\text{Euler}(f; x_0, h; y_n)} + y_{n+1} & = y_n + \underbrace{h\,f(x_n, y_n)}_{=\Phi_\text{Euler}(f; x_n, h; y_n)} \end{align*} \end{block} \end{frame} @@ -213,6 +220,7 @@ \end{axis} \end{tikzpicture} } + \begin{center}$h$\end{center} \end{column} \end{columns} \end{overlayarea} @@ -226,7 +234,7 @@ y_{n+1} = y_n + \int_{x_n}^{x_n + h} \mathrm{d}x\, f(x, y(x))\,. \end{equation*} \vspace{-\medskipamount} - The \alert{$s$}-stage Runge-Kutta method then reads + The \alert{$s$}-stage Runge-Kutta method then estimates this integral \begin{equation*} y_{n+1} = y_n + \sum_{i=1}^{\alert{s}} \hl{b_i} k_i \end{equation*} @@ -265,36 +273,255 @@ $f$ \item the LTE is of order $p$, i.e.: $\text{LTE} = O(h^{p + 1})$ \item in general $s > p$ - \item up to $s = 4$ there are are methods known with $s \geq p$ + \item up to $s = 4$ there are methods known with $s \geq p$ \end{itemize} \end{frame} \begin{frame}{Explicit Single Step Methods: Runge-Kutta (cont'd)} + \begin{block}{Euler's method} \begin{columns} \begin{column}{.5\textwidth} - \begin{block}{Euler's method} \begin{tabular}{c|c} $0$ & \\ \hline & $b_1$ \\ \end{tabular} - \end{block} \end{column} \begin{column}{.5\textwidth} - \begin{tabular}{c|ccccc} - $0$ & \\ - $c_2$ & $a_{2,1}$ \\ - $c_3$ & $a_{3,1}$ & $a_{3,2}$ \\ - $\vdots$ & $\vdots$ & & $\ddots$ \\ - $c_s$ & $a_{s,1}$ & $a_{s,2}$ & $\dots$ & $a_{s,s-1}$ \\ - \hline - & $b_1$ & $b_2$ & $\dots$ & $b_{s-1}$ & $b_s$ \\ - \end{tabular} + \begin{itemize} + \item method of stage 1 and order 1 + \end{itemize} \end{column} \end{columns} + \end{block} + + \begin{block}{Ralston's method} + \begin{columns} + \begin{column}{.5\textwidth} + \begin{tabular}{c|cc} + $0$ & \\ + $2/3$ & $2/3$ \\ + \hline + & $1/4$ & $3/4$ \\ + \end{tabular} + \end{column} + \begin{column}{.5\textwidth} + \begin{itemize} + \item method of stage 2 and order 2 + \end{itemize} + \end{column} + \end{columns} + \end{block} + + \begin{block}{RK4 (\emph{The} Runge-Kutta method)} + \begin{columns} + \begin{column}{.5\textwidth} + \begin{tabular}{c|ccccc} + $0$ & \\ + $1/2$ & $1/2$ \\ + $1/2$ & $0$ & $1/2$ \\ + $1$ & $0$ & $0$ & $1$ & $0$ \\ + \hline + & $1/6$ & $1/3$ & $1/3$ & $1/6$ \\ + \end{tabular} + \end{column} + \begin{column}{.5\textwidth} + \begin{itemize} + \item method of stage 4 and order 4 + \end{itemize} + \end{column} + \end{columns} + \end{block} +\end{frame} + +\begin{frame}{Global Truncation Error/Accumulated Error} + \begin{columns} + \begin{column}{.7\textwidth} + \resizebox{\textwidth}{!}{ + \begin{tikzpicture} + \begin{axis}[% + axis x line=center, + axis y line=left, + ymin=-1.2040,ymax=-0.59784, + xmin=0.0,xmax=+0.5, + samples=160 + ] + \addplot[ultra thick,gray,domain=0:+1] { + ln(0.5) - x + }; + % Euler, h = 0.125 + \addplot[thick,blue,mark=x,only marks,mark size=3pt] coordinates { + (0.125, -0.82668) % ln(0.43750) + (0.250, -0.96022) % ln(0.38281) + (0.375, -1.09370) % ln(0.33496) + }; + % Runge-Kutta, h = 0.125 + \addplot[thick,red,mark=+,only marks,mark size=3pt] coordinates { + (0.125, -0.81814) % ln(0.44125) + (0.250, -0.94314) % ln(0.38940) + (0.375, -1.06810) % ln(0.34364) + }; + \end{axis} + \end{tikzpicture} + } + \end{column} + \begin{column}{.3\textwidth} + \begin{itemize} + \item[{\color{blue} x}] Euler method\\ + $h = 0.125$\\ + order $p = 1$ + \item[{\color{red} +}] Runge-Kutta 4\\ + $h = 0.125$\\ + order $p = 4$ + \end{itemize} + \end{column} + \end{columns} + The differences in accumulated trunction errors is clearly visible +\end{frame} + +\begin{frame}{Wastefulness of Single-Step Methods} + In order to achieve a single step, the previously shown methods + of stage $s$ evaluate the function $f(x; y)$ in at least $s$ points.\\ + \vfill + However, all the information obtained from these evaluations is + thrown away before the next step.\\ + \vfill + It would be efficient to reuse this information, which gives + rise to the \alert{multiple step methods}. \end{frame} \begin{frame}{Multiple Step Methods} + An explicit \alert{multiple step method} of \alert{length $m$} follows from: + \begin{equation*} + y_{n+1} = -\sum_{j=0}^{m - 1} \hl{a_{j}} \, y_{n - j} + h\,\sum_{j=0}^{m - 1} \hl{b_{j}} \, f(x_{n - j}, y_{n - j}) + \end{equation*} + + Again, one strives to choose the coefficients $\lbrace \hl{a_i}\rbrace$ and $\lbrace \hl{b_i}\rbrace$ + such that the LTE is of high order in $h$. +\end{frame} + +\begin{frame}{Adams-Bashforth methods} + In Adams-Bashforth methods one chooses $\hl{a_0} = -1$, and $\hl{a_i} = 0$ for $i\geq 1$. + The coefficients $\lbrace \hl{b_i} \rbrace$ are chosen such that $y(x)$ is interpolated + in the last $m$ steps: + \begin{equation*} + \hl{b_{m - j - 1}} = \frac{(-1)^j}{j!\, (m - j - 1)!} \int_0^1 \prod_{i=0,i\neq j}^{s - 1} \mathrm{d}u\, (u + i)\,. + \end{equation*} + for $j = 0, \dots, m - 1$.\\[\medskipamount] + + The order $p$ of any Adams-Bashforth method is $p=m$. +\end{frame} + +\begin{frame}{Stiff ODEs} + \begin{block}{Example ODE} + \begin{equation*} + y' = -15\, y\qquad \text{with} \quad y(x_0) = y_0 = 1.0\,. + \end{equation*} + \end{block} + \begin{overlayarea}{\textwidth}{.7\textheight} + \only<1>{ + \vfill + \includegraphics[width=\textwidth]{StiffEquationNumericalSolvers}\\ + ~\hfill{\tiny shamelessly taken from Wikipedia} + } + \only<2>{ + \begin{itemize} + \item a stiff ODE forces you to go to \alert{unreasonably small} + step size $h$ in order to have reasonable convergence of the method + \vfill + \item in general, \alert{explicit} methods demand smaller step sizes + than computationally reasonable + \vfill + \item implicit methods solve this problem at the expense of additional + evaluations of the generation function $F$ when solving for the + increment + \end{itemize} + } + \end{overlayarea} +\end{frame} + +\begin{frame}{Implicit Euler Method} + (Also known as the backward Euler method) + \begin{equation*} + y_{n + 1} = y_n + h\cdot f(x_{n+1}, y_{n+1}) + \end{equation*} + ~\\ + \vspace{-2\bigskipamount} + \begin{block}{Example with stiff ODE} + \begin{equation*} + y' = -y,\qquad \text{with} \quad y_0 = y(x_0) = 0.5 + \end{equation*} + ~\\ + \begin{overlayarea}{\textwidth}{.5\textheight} + \vspace{-2\bigskipamount} + \only<1,3>{ + \begin{columns} + \begin{column}{.6\textwidth} + \resizebox{\textwidth}{!}{ + \begin{tikzpicture} + \begin{axis}[% + axis x line=center, + axis y line=left, + ymin=0.30,ymax=0.51, + xmin=0.00,xmax=0.50, + samples=160 + ] + \addplot[ultra thick,black,domain=0:+1] { + 0.5 * exp(-x) + }; + % Euler, h = 0.125 + \addplot[thick,blue,mark=x,only marks,mark size=3pt] coordinates { + (0.125, 0.43750) + (0.250, 0.38281) + (0.375, 0.33496) + }; + \only<3>{ + % Implict Euler, h = 0.125 + \addplot[thick,red,mark=+,only marks,mark size=3pt] coordinates { + (0.125, 0.444) + (0.250, 0.395) + (0.375, 0.351) + }; + } + \end{axis} + \end{tikzpicture} + } + \end{column} + \begin{column}{.4\textwidth} + \begin{itemize} + \item[{\color{blue} x}] regular Euler\\ + $h = 0.125$ + \item<3>[{\color{red} +}] implicit Euler\\ + $h = 0.125$ + \end{itemize} + \end{column} + \end{columns} + } + \only<2>{ + \begin{align*} + y_{n + 1} & = y_n + h\, f(x_{n+1}, y_{n +1})\\ + & = y_n - h\, y_{\hl{n + 1}}\\ + \Rightarrow + y_{n + 1} & = \frac{y_n}{1 + h} + \end{align*} + This yields the analytic result in the limit $h\to 0$ ($h = h_0 / n$, $n \to \infty$), since: + \begin{equation*} + \lim_{n\to \infty} \frac{y_0}{(1 + \frac{x}{n})^n} = y_0 \exp(-x) + \end{equation*} + } + \end{overlayarea} + \end{block} +\end{frame} + +\begin{frame}{Next time on Ordinary Differential Equations} + \begin{itemize} + \item \alert{implicit} Runge-Kutta methods + \vfill + \item how to turn an $n$-order ODE into an $n$-size \alert{system} of + first-order ODEs + \vfill + \item a glimpse at \alert{boundary problems} + \end{itemize} \end{frame} \backupbegin