\documentclass[12pt,a4paper]{report}
\usepackage{latexsym,amssymb,german,a4wide,verbatim,xr}
\usepackage{fancyhdr} 
%\usepackage{pgffor}
\usepackage{enumitem}
\usepackage{amsmath}
\usepackage{verbatim}
\usepackage{mathtools}
\DeclarePairedDelimiter{\ceil}{\lceil}{\rceil}

% use labels from a different document
%\usepackage{xr}
% optional \externaldocument[blablbla]{...} to avoid multiple labels \ref{blablablaLABEL}
% is used to get the correct references from the external document
%\externaldocument{../vorles/vorles}

% include graphical environment
%\usepackage{graphics}

% include self defined matrix symbols
%\usepackage{matrix}

%\newcommand*{\pathToExercises}{../../Material_NumerikI_Matthias}
% use selfdefined commands
%\input{\pathToExercises/vorles/meinstil.def}


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Layout
%\addtolength{\topmargin}{-15ex}
%\addtolength{\textheight}{20ex}
\pagestyle{fancy}
\setlength{\headheight}{28pt}
\lhead{Technische Universität Berlin\\Numerische Mathematik I - WS 2019/2020}
\rhead{Prof.~V.~Mehrmann}
\cfoot{}

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

\newcommand{\CurrentWeekNb}{4}
\newcounter{taskNb}
\stepcounter{taskNb}

\begin{document}
%\thispagestyle{AbgabeDatum}
\begin{center}
%{\LARGE Numerische Mathematik}\\[1.5ex]
%%{\large \CurrentWeekNb. \"Ubungsblatt}\\ 
{\Large Lösungsskizzen für die \CurrentWeekNb. Tutorien}\\ 
{Woche 11.11-18.11.19}\\
%\\{Diese Aufgaben müßen spätestens am 24.10.19 in der Vorlesung abgegeben werden.}
%{Abgabe:} Spätestens am 24.10.19 in der Vorlesung\\[0.5cm]
\end{center}
%\hfill {Abgabe: Spätestens am 24.10.19 in der Vorlesung}\\[0.5cm]
%
%%\def\exerciseList{001, 002, 003}
%\begin{enumerate}
%
\begin{enumerate} 
\item 
\begin{enumerate}
\item[] 1a) \[
\begin{aligned}
((1\oplus x)\otimes(1\oplus x)) \ominus 1&= \big[[(1+x)(1+\epsilon_1)\times (1+x)(1+\epsilon_1)] (1+\epsilon_2)-1\big](1+\epsilon_3)\\
&=:f(\epsilon_1, \epsilon_2, \epsilon_3)
\end{aligned}
\]
mit $\epsilon_i\leq$\texttt{eps}, für $i=1,\ldots, 3$.

Wir wollen zuerst den absoluten Fehler $\Delta f(\epsilon_1, \epsilon_2, \epsilon_3)=|f(0,0,0)-f(\epsilon_1, \epsilon_2, \epsilon_3)|$ abschätzen, deswegen benutzen wir eine Taylor-Entwicklung 1. Ordnung:
\[
f(\epsilon_1, \epsilon_2, \epsilon_3)= f(0,0,0) + Df(0,0,0) [\epsilon_1, \epsilon_2, \epsilon_3]^T  + \mathcal{O}(\texttt{eps}^2)
\]
wobei $Df$ der Gradient
\[
\begin{aligned}
Df(e_1, e_2, e_3)=[&2(1+x)^2(1+e_1)(1+e_2)(1+e_3), \big((1+x)(1+e_1)\big)^2(1+e_3),\\
&\big((1+x)(1+e_1)\big)^2(1+e_2)-1
]
\end{aligned}
\]
von $f$ ist. Deswegen finden wir $Df(0,0,0)=[2(1+x)^2, (1+x)^2, (1+x)^2-1
]$ und 
\[
\Delta f(\epsilon_1, \epsilon_2, \epsilon_3) \approx |Df(0,0,0)[\epsilon_1, \epsilon_2, \epsilon_3]^T|\leq || Df(0,0,0) ||_1 \cdot \texttt{eps}
\]
Wir berechnen
\[
||Df(0,0,0)||_1 = 3(x+1)^2 + |(1+x)^2-1|.
\]
Dann ist der relative Fehler
\[
\frac{\Delta f(\epsilon_1, \epsilon_2, \epsilon_3)}{f(0,0,0)}=\left[\frac{3(x+1)^2}{ |(1+x)^2-1|}+1\right] \texttt{eps}.
\]
Man sieht, dass der erste Term in der Summe groß wird wenn $x\approx$ \texttt{eps}. Hier tritt also Auslöschung auf.
%uns sehen dass $||Df(0,0,0)||_1$ beschränkt bleibt, wenn $x\to 0$.
\item[] 1b) Die gleiche Methode liefert hier:
\[
\begin{aligned}
(x \otimes x)\oplus(2\otimes x) &= [(x \times x)(1+\epsilon_1)\oplus(2\times x)(1+\epsilon_2)](1+\epsilon_3)\\
&=:f(\epsilon_1, \epsilon_2, \epsilon_3)
\end{aligned}
\]
$\Rightarrow Df(0,0,0)=[x^2, 2x, x^2+2x]$\\
$\Rightarrow ||Df(0,0,0)||_1 = x^2 + 2|x| +|x^2+2x|$\\
%\leq 2x^2 +4|x|\longrightarrow 0$, wenn $x\longrightarrow 0$.

\[
\Rightarrow \frac{\Delta f(\epsilon_1, \epsilon_2, \epsilon_3)}{f(0,0,0)}= \left[ \frac{x^2 + 2|x|}{|x^2+2x|}+1 \right]\texttt{eps} \approx 2\texttt{eps}
\]
wenn $x\approx$ \texttt{eps}.

Das Algorithmus b) ist also hier besser (keine Auslöschung).

\item[] 2a)  Analog:
\[
\begin{aligned}
(\text{pseudo-exp}(x)\ominus 1) \oslash x &= (\exp(x)(1+\epsilon_1)-1) (1+\epsilon_2) / x \cdot (1+\epsilon_3) \\
&=:f(\epsilon_1,\epsilon_2,\epsilon_3)
\end{aligned}
\]
$\Rightarrow Df(0,0,0)=[\exp(x)/x, (\exp(x)-1)/x,(\exp(x)-1)/x ]$\\
$\Rightarrow ||Df(0,0,0)||_1=2|\exp(x)-1|/|x| +\exp(x)/|x|$% \longrightarrow \infty \text{, wenn } x\longrightarrow 0$.
\[
\Rightarrow \frac{\Delta f(\epsilon_1, \epsilon_2, \epsilon_3)}{f(0,0,0)}=\left[ 2 +\frac{\exp(x)}{\exp(x)-1} \right] \texttt{eps}
\]

Hier tritt also Auslöschung auf, wenn $x$ klein ist, da der 2. Term der Summe groß wird.

\item[] 2b) Mit der gleichen Methode kommt hier:\\
$\Rightarrow ||Df(0,0,0)||_1\leq 2+\frac{3}{2}|x|+\frac{3}{6}|x|^2$.
\\
...\\
Keine Auslöschung.

\end{enumerate}

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%\item  
%$x_0=f(a,b)=-\frac{b}{a}$\\
%Berechnet man $x_0$ anhand von $f$ mit leicht gestörten Daten $\tilde{a}=a+\Delta a$, $\tilde{b}=b + \Delta b$, findet man $\tilde{x}_0 = - \frac{b+\Delta b}{a+\Delta a}$. Der absolute Fehler ist also
%\[
%\begin{aligned}
%|\tilde{x}_0 - x_0|\approx |\frac{\partial f}{\partial a} \Delta a + \frac{\partial f}{\partial b} \Delta b | &= |\frac{b}{a^2} \Delta a - \frac{1}{a} \Delta b| =  |x_0||\frac{1}{a} \Delta a - \frac{1}{b} \Delta b| \\ &\leq \max\{\frac{1}{a}, \frac{1}{b}\} (|\Delta a|+| \Delta b|) |x_0|
%% \leq  |\frac{b}{a}| |\frac{\Delta a}{a}| + |\frac{b}{a}| |\frac{\Delta b}{b}|.
%\end{aligned}
%\]
%Man sieht, dass die Konstante $C:=\max\{\frac{1}{a}, \frac{1}{b}\}$ groß sein kann, wenn $a$ oder $b$ klein sind. 
%%
\item  Man sieht, dass die 4. Stelle trotzdem falsch sein kann, sogar wenn man das Cosinus mit doppelter Genauigkeit berechnet hat.

\item 
Wir betrachten Zahlen $a, b, c$ die nicht unbedingt Maschinenzahlen sind. D.h $gl(a)=a(1+\epsilon_1)$, $gl(b)=b(1+\epsilon_2)$, $gl(c)=c(1+\epsilon_3)$ wobei $|\epsilon_1|, |\epsilon_2|, |\epsilon_3|\leq \texttt{eps}$.

Rückwärtsanalyse unter Vernachlässigung der Terme höhere Ordnung in $\texttt{eps}$:
\[(gl(a) \oplus gl(b)) \otimes gl(c) = (a(1+\epsilon_1)+b(1+\epsilon_2))(1+\epsilon_4) \times c(1+\epsilon_3) (1+\epsilon_5)\approx (\tilde{a}+\tilde{b}) \times \tilde{c}\]
mit $\tilde{a}=a(1+\epsilon_1)$, $\tilde{b}=b(1+\epsilon_2)$, $\tilde{c}=c(1+\epsilon_3+\epsilon_4+\epsilon_5)$, wobei $|\epsilon_4|, |\epsilon_5|\leq \texttt{eps}$.

Deswegen können die absolute Rückwärtsfehler wir folgt abgeschätzt werden:
\[
\begin{aligned}
&|\tilde{a}-a| =|a-a(1+\epsilon_1) | \leq |a| \texttt{eps} \\
&\text{analog}:\\
&|\tilde{b}-b| =|b||\epsilon_1|\leq |b| \texttt{eps} \\
&|\tilde{c}-c| = |c| |\epsilon_3+\epsilon_4+\epsilon_5| \leq 3|c| \texttt{eps} 
\end{aligned}
\]
Die relative Rückwärtsfehler sind dementsprechend
\[
\frac{|\tilde{a}-a|}{|a|} \leq \texttt{eps}, \qquad \frac{|\tilde{b}-b|}{|b|} \leq  \texttt{eps}, \qquad \frac{|\tilde{c}-c|}{|c|} \leq 3 \texttt{eps} .
\]
Man kann also sagen, dass die Berechnung rückwärtsstabil ist, da die relativen Rückwärtsfehler die Größenordnung $\texttt{eps}$ haben.



Vorwärtsanalyse unter Vernachlässigung der Terme höhere Ordnung in $\texttt{eps}$:
\[
(gl(a) \oplus gl(b)) \otimes gl(c) \approx (a+b)c + (a\epsilon_1+b\epsilon_2)c +(a+b)c (\epsilon_3+\epsilon_4+\epsilon_5)
\]
Der relative Vorwärtsfehler ist dann
\[
\begin{aligned}
\frac{|(gl(a) \oplus gl(b)) \otimes gl(c) - (a + b)c|}{|(a + b)c|} &= \frac{| (a\epsilon_1+b\epsilon_2) +(a+b) (\epsilon_3+\epsilon_4+\epsilon_5)|}{|(a + b)|}\\
&\leq \left(\frac{|(a\epsilon_1+b\epsilon_2)|}{|a+b|} \right) +3\texttt{eps}\\
&\leq \frac{|a|+|b|}{|a+b|} \texttt{eps}+3\texttt{eps}.
\end{aligned}
\]
Also die Berechnung ist nicht vorwärststabil, wenn $|a+b|\ll |a\epsilon_1+b\epsilon_2|$.

Das kann man einfacher sehen anhand der Faustregel: \begin{center}
(vorwärtsfehler) $\tilde\leq$ (Kondition)$\times$ (Rückwärtsfehler),\end{center} wobei die Kondition wie folgt gerechnet wird: $f(a,b,c)=(a+b)c$,
\[
\partial_a f(a,b,c) = c = C_a f(a,b,c),  \qquad \partial_b f = c= C_b f(a,b,c), \qquad \partial_c f = a+b = C_c f(a,b,c)
\]
mit $C_a=C_b=\frac{1}{a+b}$, $C_c=\frac{1}{c}$. Also der relative Vorwärtsfehler kann wie folgt abgeschätzt werden:
\[
\frac{|f(a,b,c)-f(\tilde{a},\tilde{b},\tilde{c})|}{|f(a,b,c)|} \tilde\leq \frac{1}{|a+b|} |a|\texttt{eps}+\frac{1}{|a+b|}|b|\texttt{eps}+\frac{1}{|c|}3|c|\texttt{eps}.
\]
Man sieht, dass es groß wird, wenn $a\approx -b$ (in diesem Fall ist die Kondition schlecht!).









\item Optionale extra Aufgabe:\\
Betrachten Sie das Anfangswertproblem $y'(t)=y(t), y(0)=y_0$.
\begin{enumerate}
\item Untersuchen Sie das explizite Euler-Verfahrens mit konstanter Schrittweite $h=1/n$ auf einem Rechner mit Maschinengenauigkeit \texttt{eps}. Geben Sie eine explizite Formel für $\tilde{u}_n:=\tilde{u}_h(t_n)$ in  Abhängigkeit von 
$y_0$ an, wobei wir davon ausgehen, dass jede Addition/Multiplikation einen relativen Fehler in der Größenordnung 
\texttt{eps} liefern. Terme der Größenordnung \texttt{eps}$^2$, \texttt{eps}$^3$, ... dürfen vernachlässigt werden. 
\item Gilt immer noch für kleine $h$, $1\gg h \approx \texttt{eps}$, dass $\tilde{u}_h(1)\approx y(1) = y_0 e^1$? 
\end{enumerate}
Hinweis: Es gilt die Formel $(1 + \frac{x}{n})^n =e^x(1+\mathcal{O}(h)).$ \\[0.8cm]
Lösungsskizzen zur optionalen Aufgabe:
\begin{enumerate}
\item $1\oplus h=(1+h)(1+\epsilon_1)$ mit $|\epsilon_1|\leq \texttt{eps}$. Explizit Euler gibt hier $u_j=u_{j-1}+h u_{j-1}=(1+h)u_{j-1}$. Wegen Rundungsfehler haben wir:
\[
\begin{aligned}
\tilde{u}_1&=(1\oplus h) \otimes y_0 =(1+h)(1+\epsilon_1) y_0 (1+\epsilon_2) \approx y_0 (1+h) (1+\epsilon_1+\epsilon_2)\\
\tilde{u}_2&=(1\oplus h) \otimes \tilde{u}_1 =(1+h)(1+\epsilon_1) \tilde{u}_1 (1+\epsilon_3) \approx y_0 (1+h)^2(1+2 \epsilon_1+\epsilon_2+\epsilon_3)\\
\tilde{u}_3&=(1\oplus h) \otimes \tilde{u}_2 =(1+h)(1+\epsilon_1) \tilde{u}_2 (1+\epsilon_4) \approx y_0 (1+h)^3(1+3\epsilon_1+\epsilon_2+\epsilon_3+\epsilon_4)\\
\vdots \\
\tilde{u}_n& \approx y_0 (1+h)^n (1+n\epsilon_1+\sum\limits_{j=2}^{n+1} \epsilon_j)=y_0 e^1 (1+\mathcal{O}(h)) (1+n\epsilon_1+\sum\limits_{j=2}^n \epsilon_j)
\end{aligned}
\]
\item Benutzt man $n=\frac{1}{h}$, findet man $\tilde{u}_n \approx y_0 e^1 (1+\mathcal{O}(h)) (1+\frac{\epsilon_1}{h}+\sum\limits_{j=2}^n \epsilon_j)$.
Wenn $h$ so klein wird, dass $h \approx \texttt{eps}$, dann ist $\frac{\epsilon_1}{h}$ nicht mehr vernachlässigbar und deswegen kann der Faktor $(1+\frac{\epsilon_1}{h}+\sum\limits_{j=2}^n \epsilon_j)$ nicht mehr mit 1 annähert werden. Die Annäherung $\tilde{u}_n\approx y_0 e^1$ gilt also nicht mehr.
\end{enumerate}
\end{enumerate}
%...
\end{document}


