commit bbe3d07a2a6baf4df5b567b715cf494add1047c1
parent 644a7400c4df8a6bc08dff2944b05603406abd83
Author: Orangerot <purple@orangerot.dev>
Date: Tue, 15 Apr 2025 00:28:16 +0200
feat: numerik summary
Diffstat:
| A | numerik/summary.typ | | | 653 | +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ |
1 file changed, 653 insertions(+), 0 deletions(-)
diff --git a/numerik/summary.typ b/numerik/summary.typ
@@ -0,0 +1,653 @@
+#set page("a4", margin: 1cm, numbering: "1/1", header: rect(width: 100%, stroke:
+(bottom: 1pt), [Numerik Zusammenfassung #h(1fr) Gero Beckmann]))
+
+#let TODO(..rest, content) = text(stroke: red, ..rest, [TODO: #content])
+
+== Kondition und Normen
+
+Problem: Wie wirken sich Störungen der Eingabegrößen (hier A und b) auf die
+Lösung unabhängig vom gewählten Algorithmus aus?
+
+=== Zugehörige Matrixnorm
+
+$
+norm(A) &:= sup_(x != 0) norm(A x)/norm(x) = sup_(x !=0) norm(A x/norm(x)) =
+sup_(norm(x)=1) norm(A x) = max_(norm(x) = 1) norm(A x) \
+
+norm(A^(-1)) &= sup_(x != 0) norm(A^(-1) x)/ norm(x) =^(A^(-1)x = 1) sup_(z !=
+0) norm(z)/norm(A z) = max_(norm(z) = 1) 1/norm(A z) = 1 / (min_(norm(z) = 1)
+ norm(A z))
+$
+
+ #table(columns: 4, align: center + horizon,
+[Spaltensummennorm], $ norm(A)_1 = max_(m=1,...,N) sum_(n=1)^N abs(a_"nm") $,
+$ norm(x)_1 = sum_(n=1)^N abs(x_n) $, [1-Norm],
+[Zeilensummennorm], $ norm(A)_1 = max_(n=1,...,N) sum_(m=1)^N abs(a_"nm") $,
+$ max_(n=1,...,N) abs(x_n) $, [Maximumsnorm],
+[Spektralnorm], $ norm(A)_2 = sqrt("größter EW von " A^T A) $, $ norm(x)_2 = sqrt(sum_(n=1)^N x_n^2) $, [euklidische Norm]
+ )
+
+ Eigenschaften: positive Definitheit $norm(x) >= 0$, Dreiecksungleichung $norm(x
+ + y ) <= norm(x) + norm(y)$, Homogenität $norm(lambda x) = abs(lambda) norm(x)$
+
+$norm(A x) <= norm(A) norm(x)$, $norm(I_N) = 1$, Submiltiplikativität $norm(A B)
+<= norm(A) norm(B)$, $norm(x) = sqrt(angle.l x"," x angle.r)$
+
+=== Störungen
+
+#let cond = "cond";
+
+$
+x - tilde(x) &= A^(-1) b - A^(-1) tilde(b) = A^(-1)(b - tilde(b)) \
+<=> norm(x - tilde(x)) &= norm(A^(-1) (b - tilde(b))) <= norm(A^(-1)) norm(b - tilde(b))
+&& "(absolute Abweichung)" \
+<=> norm(x - tilde(x)) / norm(x) &= norm(A^(-1) ( b - tilde(b))) / norm(x) <=
+(norm(b) norm(A^(-1))) / norm(x) norm(b - tilde(b)) / norm(b) && "(relative Abweichung)" \
+<=> norm(x - tilde(x)) / norm(x) &<=^(norm(b) = norm(A x) <= norm(A) norm(x))
+underbrace(norm(A) norm(A^(-1)), =: cond(A)) norm(b - tilde(b)) / norm(b) $
+
+=== Konditionszahl einer Matrix
+
+$
+#rect($cond(A) := norm(A) norm(A^(-1))$) >=^"Submiltiplikativit" norm(A A^(-1))
+= norm(I_N) = 1
+$
+
+$
+cond_2(A) = (max{abs(lambda) : lambda "EW von " A})/(max{abs(lambda) : lambda
+"EW von " A}), cond(A) = (max_(norm(y) = 1) norm(A y))/(min_(norm(z) = 1) norm(
+ A z)), cond(A) = cond(alpha A), alpha in bb(R) "\\" {0}
+$
+
+Seien relative Abweichung beschränkt
+$
+norm(A - tilde(A)) / norm(A) <= epsilon_A, " " && norm(b- tilde(b)) / norm(b) <= epsilon_b
+$
+Dann gilt im Fall $epsilon_A cond(A) < 1$ die Abschätzung
+$
+norm(x - tilde(x))/norm(x) <= cond(A) / (1 - epsilon_A cond(A)) (epsilon_A + epsilon_b)
+$
+
+$cond(alpha I_N) = 1$, für orthogonale Matrizen Q: $ norm(Q x)_2^2 = x^T Q^T T
+x = x^T x = norm(x)_2^2, norm(Q)_2 = 1 => cond_2(Q) = 10 = norm(Q) norm(Q^(-1)) =
+norm(Q)^2$
+
+#pagebreak()
+
+== QR-Zerlegung
+
+Geg. $A in bb(R)^(M times N), M >= N, "voller Rang" (=N)$\
+Ges. orthogonales $Q in bb(R)^(M times N)$ und $R = mat(dash(R);0) in bb(R)^(M
+times N)$, rechte obere Dreiecksmatrix $dash(R) in bb(R)^(N times N)$
+
+#rect()[
+=== Householder-Transformation
+
+Householdermatrix Q ist symmetrisch, orthogonal und eine Spiegelung eines
+Vektors v auf ein Vielfaches des Vektors u $(=e^1)$ mit Spiegelachse $perp w$ (Householdervektor).
+
+$
+Q = I_N 2 w w^T => Q v = sigma e^1 = -"sign"(v_1) norm(v)_2 e^1 = v - 2 w^T v w
+#h(1cm) (2M " Ops")
+$
+
+Q wird eigl. nie ausgerechnet, sondern nur w gespeichert und QA wird berechnet
+aus $(Q a^1, Q a^2, ...)$ Spaltenvektoren von A. (2M(N-1)) Ops
+
+$
+"(Householdervektor) " w = (v - u)/norm(v - u)_2 = (v - sigma e^1)/norm(v - sigma e^1)_2
+#h(1cm) sigma = -"sign"(v_1)norm(v)_2 = cases(-&norm(v)_2 ", falls " v_1 > 0,
+&norm(v)_2 ", falls " v_1 <= 0)
+$
+]
+
+Berechne Householder-Transformation des ersten Spaltenvektors $a^1$ von $A^((n))$
+
+$
+w = (a^1 - sigma e^1) / norm(a^1 - sigma e^1)
+#h(1cm) ,
+Q a^1 = r^1 = sigma e^1 = -"sign"(a^1) norm(a^1)_2 e^1
+#h(1cm) ,
+Q^((n)) A^((n))
+= mat(
+r_(11), r_(12)^T;
+, A^((n+1))
+)
+= mat(
+sigma, r_(12)^T;
+0, A^((n+1))
+)
+$
+
+Durch die Householder-Transformation bekommen wir die Matrix für die nächste
+Iteration.
+
+Householdervektor $tilde(w)^((n+1))$ ist eine Dimension kleiner. Um diesen an
+das vorherige Ergebnis zu multiplizieren $w^((n+1)) = vec(0, tilde(w)^((n+1)))$
+sodass $Q^((n+1)) = mat(1, 0; 0, tilde(Q)^((n+1)))$ (alternativ: es wird auf
+$e^n$-Vektor gespiegelt)
+
+Insgesamt: in $M N^2-1/3 N^3$ Ops #h(2cm) $Q_N ... Q_2 Q_1 A =: R
+#h(1cm)
+Q_1 ... Q_N =: Q$
+
+#grid(columns: 2*(1fr,), align: left,
+[
+=== Lineares Ausgleichsproblem $norm(A x - b)_2 = min!$
+Berechne $A = Q R$ \
+Berechne $Q^T b = vec(c, d)$ \
+Löse $dash(R) x = c$ (Rückwärtssubstitution)
+],
+[
+=== Lösen von LGS $A x = b$ (falls $M=N$)
+$
+A x &= b \
+<=> Q R x &= b \
+<=> R x &= Q^T b \
+<=> R x &= c " (Rückwärtssubstitution)"
+
+$
+],
+)
+
+#v(-1cm)
+=== QR-Algorithmus (Eigenwertberechnung)
+
+1. setze $A_0 = A$ und $k = 0$
+2. zerlege $A_k = Q_k R_k$
+3. berchne $A_(k+1) = R_k Q_k$, erhöhe k um 1 und gehe zu 2 #h(2cm)
+ $O(N^3)$ Ops pro Iteration
+
+$k -> oo: A_k -> R = mat(lambda_1, *, *, *; , lambda_2, *, *; ,, dots.down, *;
+0,,,lambda_N), Q_k -> I_N$ #h(2cm) Konvergenz: $abs(lambda_N)/abs(lambda_(N-1))$
+
+Benutze *Hessenbergform* für weniger Ops pro Iteration $O(N^2)$. Kann durch N-2
+Householder-Transformation berechnet werden $Q^T A Q = H$
+
+Householder-Transformation $tilde(Q)_2 = I_(N-1) - 2 tilde(w)^2 (tilde(w)^2)^T
+in bb(R)^((N-1) times (N-1))$ mit $tilde(Q)_2 vec(a_(21), a_(31), dots.v, a_(N
+1)) = vec(*, 0, dots.v. 0) in bb(R)^(N-1), Q_2 = I_N - 2 w w^T =
+ mat(1,;,tilde(Q)_2), w^2 = vec(0, tilde(w)^2)$
+
+#v(-.5cm)
+$
+A^((1)) = Q_2 A Q_2 \
+Q_(N-1) A^((N-3)) Q_(N-1) = Q_(N-1) ... Q_2 A Q_2 ... Q_(N-1) = Q^T A Q = H
+#h(1cm)
+, Q = Q_2 ... Q_(N-1)
+$
+
+Wenn wir eine Hessenbergmatrix QR-Zerlegung $H=Q R$, dann ist die Matrix
+$dash(H) = R Q$ ebenfalls Hessenbergmatrix, wegen $dash(H) = R Q = Q^T Q R Q =
+Q^T H Q$
+
+*QR-Algorithmus mit Shift* $mu_k = h_(NN)^((k))$ für schnellere Konvergenz
+$abs(tilde(h)_(n+1,n)^((k))) = O(abs((lambda_(n+1) - mu)/(lambda_n - mu)))$
+
+1. setze $H_0 = H, k = 0$ und wähle Toleranz $epsilon$
+2. zerlege $H_k - h_(NN)^((k')) I_N = Q_k R_K$
+3. berechne $H_(k+1) = R_k Q_k + h_(NN)^((k)) I_N$
+4. ist $abs(h_(N,N-1r^((k+1)))) <= epsilon (abs(h_(NN)^((k+1)) +
+ abs(h_(N-1,N-1)^((k+1))))$ so akzeptiere $h_(NN)^((k+1))$ als EW. Sonst
+ erhöhe k um 1 und gehe zu 2.
+
+
+#pagebreak()
+
+== Newton-Verfahren
+
+Problem: Für eine vorgegebene (nichtlineare) Funktion $f: D subset bb(R)^N ->
+bb(R)^N$ finde $x^* in D$ mit $f(x^*) = 0_N in bb(R)^N$
+
+Sei $x^0$ in der Nähe der Lösung $x^*$. Eine Taylor-Entwicklung liefert
+#v(-.8cm)
+$
+ && 0 = f(x^*) = & f(x^0) + f'(x^0)(x^* - x^0) + overbrace(O(norm(x^0 - x^*)^2),"vernachlässigen" approx 0) \
+~> &&0 = & f(x^0) + f'(x^0)(x^*-x^0) \
+ &&" "= & f(x^0) + f'(x^0)x^* - f'(x^0)x^0 \
+<=>&& f'(x^0)x^0 - f(x^0)= & f'(x^0) x^* \
+<=>&& f'(x^0)^(-1) (f'(x^0) x^0 - f(x^0)) = & x^* \
+<=>&& x^* approx & x^0 - f'(x^0)^(-1) f(x^0) \
+ =>&& x^(k+1) = & x^k - f'(x^k)^(-1) f(x^k) \
+ & ""#place()[Berechnung der inversen Jacobi-Matrix vermeiden] \
+<=>&& x^(k+1) = & -f'(x^k)^(-1) f(x^k) \
+<=>&& f'(x^k) underbrace((x^(k+1) - x^k), =d^k) = & - f(x^k) \
+$
+
+1. Wähle Startwert $x^0 in D$ und Toleranz $epsilon$. Setze $k=0$
+2. Löse $f'(x^k)d^k = -f(x^k)$ (LGS durch LR-Zerlegung) \
+ $x^(k+1) = x^k + d^k$
+3. Falls ($norm(d^k) < epsilon$): STOP \
+ Sonst erhöhe k um 1 und gehe zu 2.
+
+Bem.: Man wählt nicht $norm(f(x^k)) < epsilon$ als STOP-Kriterium, weil das
+Problem $A f(x) = 0$ (A regulär) dann eine andere Lösung hätte.
+$A f(x) = 0$ ändert Iterierten nicht, deshalb ist $norm(d^k) < epsilon$ besser,
+weil so bei $A f(x) = 0$ und $f(x) = 0$ die Lösung und Iterierten gleich bleiben.
+
+==== lokale quadratische Konvergenz
+
+Es gibt eine Kugel $kappa := kappa_p (x^*) = { x in bb(R)^N: norm(x-x^*) <= p }
+subset D$ sodass $x^*$ die einzige Nullstelle von f in $kappa$0ist. Zudem liegen
+die Folgeglieder $x^(k+1) = x^k - f'(x^;)^(-1) f(x^k)$ für jeden Startwert $x^0
+in kappa$ ebenfalls in $kappa$ und es gilt $lim_(k->oo) x^k = x^*$.
+ Weiter existiert eine Konstante $c > 0$ mit $norm(x^* - x^k) <= c norm(x^* -
+ x^(k-1))^2$ für $k in bb(N)$
+
+==== Homotopie-Verfahren: guten Startwerd suchen
+$g: bb(R)^N -> bb(R)^N$ eine (einfache) Funktion mit bekannter Nullstelle
+$x^0$. Definiere $f_s (x) = (1-s) g(x) + s f (x)$, wählen $0 < s_1 < s_2 < ...
+< s_M = 1$ und lösen $f_(s_n) (x) = 0$ mit Newton mit Startvektor $x^(m-1)$ mit
+$f_(s_(m-1)) (x^(m-1) = 0$ für $m=1,...,M$
+
+==== Dämpfungsstrategie: Verbesserung der Konvergenz
+
+Bestimme $t_k in (0,1]$ sodass für $x^(k+1) = x^k - t_k f'(x^k)^(-1) f(x^k)$ das
+Residuum kleiner wird $norm(f(x^(k+1))) < norm(f(x^k))$
+
+==== Differenzenqoutient
+Wenn $f'(x^k)$ aufwendig zu berechnen, wähle $A_k approx f'(x^k)$ durch
+$(A_k)_(n m) = 1/h (f_n (x^k + h ^m) - f_n (x^k))$
+
+==== Vereinfachtes Newton-Verfahren
+
+pro Iteration $f'$ auswerten und LR-Zerlegen ist teuer. Stattdessen $A approx
+f'(x^0)$ alse einmal $f'$ auswerten und LR-Zerlegen. $F(x) = x - A^(-1) f (x)$
+Nurnoch lineare Konvergenz.
+
+==== Globale Minimierung
+
+Wir suchen $x^* in bb(R)^N$ mit $f(x^*) <= f(x)$ für alle $x in bb(R)^N$. dh
+insbesodere grad $f(x^*) = 0_N$. Wir berechnen Nullstelle von $gradient f(x) =
+("grad" f(x))^T$ mit Newton $x^(k+1) = x^k - H_f (x^k)^(-1) gradient f(x^k)$.
+$H_f$ Hessematrix.
+
+==== Nichlineare ausgleichrrechnung $norm(f(x))_2 = min!$
+
+Wähle Startwert $x^0 in bb(R)^N$ für $k=0,1,...$\
+1. löse lineares Ausglechsproblem $norm(f'(x^k) d^k + f(x^k))_2 = min$
+2. setze $x^(k+1) = x^k + d^k$
+
+#pagebreak()
+
+== Interpolation
+
+Problem: Suche zu den Stützstellen $(x_0,f_0), ..., (x_N, f_N)$ eine Funktion p
+mit $p(x_n) = f_n "für alle " n=0,...,N$
+
+=== Polynom-Interpolation
+
+
+
+Problem: Suche zu N+1 Stützstellen $(x_0, f_0), ..., (x_0, f_0)$ ein Polynom p
+vom Grad $<= N (p in bb(P)_N)$, welches das Interpolationsproblem (s.o.) lösst:
+$p(x_n) = f_n " für alle " n = 0,...,N$,
+
+==== Eindeutigkeit des Interolationspolynoms zu einem Interpolationsproblem:
+
+Nehmen wir an es gäbe zwei Polynome $p,q in bb(P)_N$ vom Grad $<= N$ welche das
+gleiche Interpolationsproblem lösen. $p(x_n) = q(x_n) = f_n "für " n = 0,...,N$.
+Dann wäre $p-q in bb(P)_N$ ein Polynom vom Grad $<= N$, welcehs $N+1$
+Nullstellen an den Stützstellen $x_n$ hat: $p(x_n) - q(x_n) = f_n - f_n = 0 "für alle " n=0,...,N$.
+
+Nach dem Fundamentalsatz der Algebra kann ein Polynom vom Grad $=N$ nur maximal
+N Nullstellen haben außer es ist das konstante 0-Polynom $f=0$.
+
+Damit das Polynom $p-q$ vom Grad $<= N $ also $N+1$ Nullstellen haben kann, muss
+nach dem Fundamentalsatz der Algebra $p - q = 0 <=> p = q$ das konstante
+0-Polynom sein, was bedeutet, das p und q identisch sein müssen.
+
+==== Lagrange-Basis
+
+Das n-te Lagrnge-Polynom $L_n$ zu den Stztzstellen $x_0, ..., x_N$ ist definiert
+als:
+
+$
+p(x) = sum_(n=0)^N f_n L_n (x) #h(1cm)
+, L_n (x) = sum_(m=0,m != n)^N (x-x_m)/(x_n-x_m) = sigma_(n m) = cases(1 ", falls "
+x_m = x_n, 0 ", falls " x_m != x_n)
+
+$
+
+Lagrnge-Polynome $L_n$ bilden aufgrund von Skalarprodukt $angle.l p,q angle.r =
+sum_(n=0)^N p(x_n) q(x_n)$ eine orthonormalbasis und linearkombination
+($bb(P)_N$).
+
+===== Kondition
+
+Frage: Wie wirken sich Störungen der Eingabegrößen (Stützwerte $f_n$) auf die
+Lösung der Interpolation aus?
+
+$
+p(x) = sum_(n=0)^N f_n L_n (x) #h(2cm)
+tilde(p)(x) = sum_(n=0)^N tilde(f_n) L_n (x) \
+$
+$
+&abs(p(x) - tilde(p)(x)) = abs(sum_(n=0)^N (f_n - tilde(f_n)) L_n (x))
+<= sum_(n=0)^N abs(f_n - tilde(f_n)) abs(L_n (x))
+<= max_(n=0,...,N) abs(f_n - tilde(f_n)) sum_(n=0)^N abs(L_n (x)) \
+<=> &max_(x in [a,b]) abs(p(x) - tilde(p)(x)) <= underbrace( max_(x in [a,b]) sum_(n=0)^N
+abs(L_n (x)), =: Lambda_N >= 1) max_(n=0,...,N) abs(f_n - tilde(f_n))
+$
+
+Lebesgue-Konstante $Lambda_N$ ist invariant under Affinen Transformationen daher
+nur von relativer Lage der Stützstellen $x_n$ zueinander abhängig.
+
+#let annotate(..args) = {
+ box(place(bottom + center, dy: -10pt, stack(math.script(..args), line(angle: 90deg, length: 10pt))))
+ sym.wj
+ h(0pt, weak: true)
+}
+
+$
+#[==== Newton-Darstellung] #h(2cm)
+p(x) &= a_0 + a_1 w_1(x) + ... + #annotate([Leitkoeffizient])a_N w_N (x) #h(1cm)
+&& w_n (x) = product_(m=0)^(n-1) (x - x_m)
+\
+"Aus Darstellung folgt "
+p(x) &= p_(0,N-1)(x) + a_N w_N (x) && p_(n,n)(x) = f_n = f_(n,n) \
+&= p_(0,N)(x) = f_(0,0) + f_(0,1) w_1(x) + ... + f_(0,N) w_N (x) #h(1cm) &&"mit " a_n = f_(0,n)
+$
+
+#grid(columns: 2, align: center + horizon,
+$
+#[*Lemma von Aitken*] #h(1cm)
+p_(n,k)(x) &= ((x_n - x) p_(n+1,k)(x) - (x_k-x)p_(n,k-1))/(x_n - x_k) \
+=> "Für Koeffizienten" f_(n,k) &= (f_(,k-1) - f_(n+1,k)) / (x_n - x_k),
+#h(.5cm) "(Dividierte Differenzen)"
+$,
+$
+#h(.5cm) mat(delim: #none,
+f_(n,k-1);
+,arrow.br;
+f_(n+1,k), arrow.r, f_(n,k);
+)
+mat(delim: #none,
+"2(N-1) Additionen";
+"N-1 Divisionen"
+)
+$
+)
+#grid(columns: 2,align: left+horizon, gutter: .5cm,
+[ *Interpolationsfehler* ],
+$
+f(x) - p(x) = w_(N+1)(x) (f^(N+1)(xi))/(N+1)! #h(1cm),
+xi = xi(x) in (a,b) "Zwischenstelle" \
+$
+)
+#v(-.3cm)
+#grid(columns: 2,align: left+horizon, gutter: .5cm,
+[*Maximaler Interpolationsfehler/Approximationsgüte*],
+$
+max_(x in [a,b]) abs(f(x) - p(x)) <= max_(x in [a,b]) abs(w_(N+1)(x)) max_(x in [a,b]) abs(f^(N+1)(x))/(N+1)!
+$
+)
+
+==== Tschebyscheff-Interpolation
+
+Tschebyscheff-Interpolationsformel
+
+$
+p(x) = 1/2 c_0 + c_1 T_1(x) + ... + c_N T_N (x)
+$
+
+mit Koeffizienten $c_m = 2/(N+1) sum_(n=0)^N f_n cos(m (2n +1)/(2N + 2) pi)$ mit
+$(N+1)^2$ Multiplikationen
+
+Tschebyscheff-Polynom vom Grad N ist auf $[-1,1]$ gegeben durch
+
+#grid(
+columns: 2, align: center + horizon,
+$
+T_0(x) &= 1 \
+T_1(x) &= x \
+T_(N+1)(x) &= 2 dot T_N (x) - T_(N-1)(x) $,
+$
+"oder "
+T_N (x) = cos(N arccos(x))
+$
+)
+
+*Tschebyscheff-Stützstellen*: $x_n = cos((2n+1)/(2N + 2) pi) $ für alle $n=0,...,N$
+
+*Min-Max-Eigenschaft*: $max_(x in [-1,1]) abs(w_(N+1)(x))$ wird genau für
+Tschebyscheff-Stützstellen minimal mit $2^(-N)$
+
+
+mit FFT lässt sich Koeffizient $c_m$ in $O(N log N)$ Multiplikationen berechnen,
+wenn $N+1 = 2^k$ eine 2-er Potenz ist.
+
+*Clenshaw-Algorithmus* zur Auswertung $p(x)$ mit N+2 Multiplikationen und 2N Additionen
+$
+d_(N+2) &= d_(N+1) = 0 \
+d_n &= c_n + 2x dot d_(n+1) - d_(n+2) " für " n=N,N-1,...,0 \
+"Dann gilt "& p(x) = 1/2 (d_0 - d_2)
+$
+
+#table(columns: 3,
+table.header([Zusammenfassung], [Bestimmung], [Auswertung]),
+[Monom],
+[LGS mit $1/3 N^3$ Ops],
+[Horner-Schema mit N Ops],
+[Lagrange],
+[direkt gegeben],
+[zu aufwendig],
+[Newton],
+[Dividierte Differenzen mit $N^2$ Ops],
+[Horner-Schema mit $2N$ Ops],
+[Tschebyscheff],
+[FFT mit $O(N log N)$ Ops],
+[Clenshaw-Algorithmus mit $2N$ Ops],
+)
+
+nur in Newton-Darstellung können zusätzliche Stützstellen hinzugefüngt werden,
+wobei nicht alle Koeffizienten nei berechnet werden müssen.
+
+=== Spline-Interpolation
+
+Problem: Suche zeimal stetig differenzierbare Funktion $s in C^2([a,b], bb(R))$
+mit mit $s(x_n) = y_n$ für $n=0,...,N$, soduss folgendes Integral minimal ist:
+
+$
+integral_a^b abs(s''(x))^2 d x <= integral_a^b abs(s''(x) + epsilon h''(x))^2 d x
+\
+<=> integral_a^b s''(x)h''(x) d x
+= sum_(n=1)^N integral_(x_(n-1))^(x_n) s''(x)h''(x) d x
+= s''(b)h'(b) - s''(a)h'(a) =^! 0
+$
+für alle $epsilon in bb(R)$ und $C^2$-Funktionen $h: [a,b] -> bb(R)$ mit $h(x_n)
+= 0$ für alle $n=0,...,N$
+
+Aus dieser Gleichung ergeben sich derei Lösungen und damit drei Typen von Splines:
+
+#align(center, table(stroke: none, columns: 2, align: left+bottom,
+$s'(a) = v_0 " und " s'(b) = v_N$, [eingespannt / hermitisch],
+$ s''(a) = 0 = s''(b)$, [natürlich],
+$s'(a) = s'(b) " und " s''(a) = s''(b)$, [periodisch (sinnvoll für $overbrace(s(a) = s(b), y_0 = y_N)$)]
+))
+
+not-a-knot-Bedingung falls Ableidung an Randpunkten unbekannt: $s'''_1(x_1) =
+s'''_2(x_1)$ und $s'''_(N-1) = s'''_N (x_(N-1))$
+
+#pagebreak()
+==== kubischer interpolierender Spline
+
+Eine $C^2$-Funktion $s: [a,b] -> bb(R)$ heißt kubischer interpolierender Spline
+zu den Stützpunkten $(x_0, y_0), ..., (x_N, y_N)$ mit der Eigenschaft $a = x_0 <
+... < x_N = b$, falls
+
+1. $s|_([x_(n-1),x_n]) in bb(P)_3$ für alle $n=1,...,N$
+2. $s(x_n) = y_n$ für alle $n=0,...,N$
+
+$
+s|_([x_(n-1),x_n]) = s_n (x) =
+underbracket(y_(n-1) + (x - x_(n-1)) overbrace(y_(n-1,n), = (y_(n-1) - y_n)/(x_(n-1) - x_n)), "interpolierender Anteil") +
+underbracket((x-x_(n-1))(x-x_n)[alpha_n (x - x_(n-1)) + beta_n (x-x_n)],
+"glättender Anteil")
+$
+
+Daraus erbeben sich 2N Unbekannte $alpha_n, beta_n$. Diese erhält man durch die
+2N - 2 Glattheitsbedingungen + 2 Randbedingungen (s.o.). Also erhält man
+$alpha_n, beta_n$ aus LGS mit dim 2N.
+
+$
+cases(reverse: #true,
+s'_n (x_n) = s'_(n+1)(x_n),
+s''_n (x_n) = s''_(n+1)(x_n)
+) " für " n=1,...,N-1 => 2N-2 " Bedingungen"
+$
+
+Die Berechnung der Unbekannten können wir vereinfachen, indem wir die $alpha_n,
+beta_n$ durch N+1 Unbekannte Momente $gamma_n$ darstellen:
+
+$
+cases(reverse: #true,
+s''_n (x_(n-1)) &= gamma_(n-1),
+s''_n (x_n) &= gamma_n
+) " damit " s''_n (x_n) = s''_(n-1) (x_n) = gamma_n
+$
+
+Dies sind 2N Gleichungen für die N+1 Unbekannten $gamma_0, ..., gamma_N$. Durch
+den Ansatz sind N-1 Glattheitsbedingungen an die zweite Ableitung automatisch
+erfüllt.
+
+Mit den Gitterweiten $h_n := x_n - x_(n-1)$ erhalten wir
+$#h(1cm) alpha_n = 1/(6 h_n) (gamma_(n-1) + 2 gamma_n) #h(1cm)
+beta_n = -1/(6 h_n) (2 gamma_(n-1) + gamma_n)$
+
+Eingesetzt in $s'$ in der Glattheitsbedingung $s'_n (x_n) = s'_(n+1) (x_n)$ folgt
+
+$
+h_n / 6 (gamma_(n-1) + 2 gamma_n) + h_(n+1) / 6 ( 2 gamma_n + gamma_(n+1)) =
+gamma_(n,n+1) - gamma_(n-1,n) =: d
+$
+
+mit Einschränkung auf eingespannte kubische Splines $s'_1(x_0) = 0_0, s'_N(x_N)
+= v_N$ erhalten wir außerdem
+
+$
+h_1/6 ( 2 gamma_0 + gamma_1) = y_(0,1) - v_0 =: d_0
+#h(2cm)
+h_N / 6 (gamma_(N-1) + 2 gamma_N) = v_N - y_(N-1,N) =: d_N
+\
+"Insgesamt ein LGS: "
+1/6 underbrace( mat(
+2 h_1, h_1;
+h_1, 2(h_1 + h_2), h_2;
+,dots.down, dots.down, dots.down;
+,,h_(N-1), 2(h_(N-1) + h_N), h_N;
+,,,h_N, 2h_N;
+), #[
+strikt diagonaldominante symmetrische Tridiagonalmatrix (spd, regulär) \
+$=>$ mit Cholesky in $O(N)$ Ops lösbar.
+])
+vec(gamma_0, gamma_1, dots.v, gamma_(N-1), gamma_N)
+=
+vec(d_0, d_1, dots.v, d_(N-1), d_N)
+$
+
+===== Kondition
+
+$tilde(s)(x)$ der Spline zu gestörten Daten $(x_n, tilde(y)_n)$, Lagrange-Spline
+$l_n (x_m) = cases(1 ", falls " n=m, 0 ", falls" n != m), #h(1cm) l'_n (a) = 0, l'_n (b) =0$
+
+$
+abs(s(x) - tilde(s)(x)) <= sum_(n=0)^N abs(y_n - tilde(y)_n) norm(l_n (x))
+<= sum_(n=0)^N abs(l_n (x)) max_(m=0,...,N) abs(y_m - tilde(y)_m)
+<= Lambda_N max_(n=0,...,N) abs(x_m - tilde(y)_m)
+$
+
+$
+"mit " Lambda_N := max_(x in [a,b]) sum_(n=0)^N abs(l_n (x)) " also "
+max_(x in [a,b]) abs(s(x) - tilde(s)(x)) <= Lambda_N max_(m=0,...,N) abs(y_m - tilde(y)_n)
+$
+
+Bei äquidisanten Unterteilung gilt $Lambda_N <= 2$ für alle $N in bb(N)$
+
+#pagebreak()
+
+
+== Numerische Integration
+
+Problem: Berechne für eine gegebene Funktion $f: [a,b] -> bb(R)$ das
+Riemann-Integral.
+
+=== Eigenschaften von Integralen
+
+#table(columns: 2*(auto, 1fr,), stroke: none,
+[Additiv:],[ $ integral_a^b f(x) d
+x = integral_a^c f(x) d x + integral_c^b f(x) d x$],
+[Linear:],[ $I(lambda f + mu g) = lambda I(f) + mu I(g)$],
+)
+Monoton: falls $f >= g$ auf $[a,b]$ dann $integral_a^b f(x) d x >= integral_a^b g(x) d x$
+ \ $abs(I(f)) <= I(abs(f)) <=> cond_1(f) = I(abs(f)) / abs(I(f)) >= 1$
+
+Norm: $norm(f)_1 = I(abs(f))$
+
+=== Quadraturformeln
+
+$
+integral_a^b f(x) d x approx integral_a^b p(x) d x = integral_a^b sum_(n=0)^N f(x_n)
+L_n (x) = sum_(n=0)^N integral_a^b L_n (x) d x f(x) = (b-a) sum_(k=1)^s b_k f( a +
+c_k (b-a))
+$
+
+#table(columns: 3,
+[Rechteckregel],
+$s = 1, b_1 = 1, c_1 = 0, p = 1$,
+$I(f) approx (b - a) f(a)$,
+[Mittelpunktregel],
+$s = 1, b_1 = 1, c_1 = 1/2, && p = 2$,
+$I(f) approx (b-a) f((a+b)/2)$,
+[Trapezregel],
+$s = 2, b_1 = b_2 = 1/2, c_1 = 0, c_1 = 1, p = 2$,
+$I(f) approx (b-a) (f(a) + f(b))/2$,
+[Simpsonregel],
+$s = 3, b_1 = b_3 = 1/6, b_2 = 4/6, c_1 = 0, c_2 = 1/2, c_3 = 1, p = 4$,
+$I(f) approx (b-a) 1/6 (f(a) + 4 f((a+b)/2) + f(b))$,
+)
+
+Eine QF $(b_k, c_k)_(k=1,...,s)$ besitzt die Ordnung p, falls sie für alle
+Polynome vom Grad $<= p - 1$ das Integral exakt berechnet, wobei p maximal ist.
+
+#TODO[BEWEIS ZU ORDNUNGBEDINGUNGEN!!!!!]\
+
+*Ordnungsbedingungen*: Eine QF besitzt genau dann die Ordnung p, falls
+
+$
+1 / q = sum_(k=1)^s b_k c_k^(q-1) "für alle " q = 1, ..., p " aber nicht mehr für " q = p + 1 " gilt"
+$
+
+Für eine QF mit vorgegebenen Knoten $c_1 < ... c_s$ können die Gewichte genau
+dann mi $b_k = integral_0^1 underbrace(L_k (x), #place(center, $L_k (x) = sum_(j=1,j != k)^s
+(x-c_j)/(c_k - c_j)$)) d x$ eindeutig bestimmt werden wenn $p >= s$.
+
+*Symmetrische QF*:
+
+$
+b_k = b_(s+1-k) && c_k = 1 - c_(s+1-k) "also symmetrisch um "1/2" verteilt"
+$
+
+Die Ordnung p einer symmetrischen QF ist gerade. #TODO[Beweis]
+
+Für symmetrische Knoten $c_k$ und $p >= s$ sind die eindeutigen Gewichte $b_k$
+ebenfalls symmetrisch.
+
+=== QF mit erhöhter Ordnung
+
+Die maximale Ordnung einer QF ist 2s (insb. bei Gauß-QF)
+
+Bew.: $angle.l M,M angle.r = integral_0^1 M(x)^2 d x > 0, deg(M) = s$
+
+*Gauß-QF*: Ordnung 2s gegeben durch $c_k = 1/2(1+gamma_k)$ mit
+$gamma_1,...,gamma_s$ die Nullstellen des Legendre-Polynoms vom Grad s.
+Die Gewichte $b_k$ einer Gauß-QF sind positiv. Gauß-QF sind symmetrisch, weil
+Legendre-Nullstellen $gamma_1, ..., gamma_s$ symmetrisch zum Punkt 0 im Interval
+$(-1, 1)$ und so Knoten symmetrisch zu Punkt $1/2$ im Interval $(0,1)$.
+
+*Summierte QF*: Um Fehler zu verkleinern zerlege $[a,b]$ in Teilintervalle
+$[x_(n-1), x_n], n=1,...,N$ mit Intervallängen $x_n - x_(n-1) = h_n$
+
+$
+integral_a^b f(x) d x approx sum_(n=1)^N h_n sum_(k=1)^s b_k (x_(n-1) + c_k h_n)
+$
