reorganize chapter 7

14 files changed, 2589 insertions, 0 deletions

diff --git a/buch/chapters/070-orthogonalitaet/Makefile.inc b/buch/chapters/070-orthogonalitaet/Makefile.inc
new file mode 100644
index 0000000..80bb54b
--- /dev/null
+++ b/buch/chapters/070-orthogonalitaet/Makefile.inc
@@ -0,0 +1,14 @@
+#
+# Makefile.inc -- Makefile dependencies for chapter 7
+#
+# (c) 2022 Prof Dr Andreas Müller, OST Ostschweizer Fachhochschule
+#
+
+CHAPTERFILES = $(CHAPTERFILES)						\
+	chapters/070-orthogonalitaet/orthogonal.tex			\
+	chapters/070-orthogonalitaet/legendredgl.tex			\
+	chapters/070-orthogonalitaet/bessel.tex				\
+	chapters/070-orthogonalitaet/jacobi.tex				\
+	chapters/070-orthogonalitaet/sturm.tex				\
+	chapters/070-orthogonalitaet/gaussquadratur.tex			\
+	chapters/070-orthogonalitaet/chapter.tex
diff --git a/buch/chapters/070-orthogonalitaet/bessel.tex b/buch/chapters/070-orthogonalitaet/bessel.tex
new file mode 100644
index 0000000..3e9412a
--- /dev/null
+++ b/buch/chapters/070-orthogonalitaet/bessel.tex
@@ -0,0 +1,91 @@
+%
+% Besselfunktionen also orthogonale Funktionenfamilie
+%
+\section{Bessel-Funktionen als orthogonale Funktionenfamilie}
+\rhead{Bessel-Funktionen}
+Auch die Besselfunktionen sind eine orthogonale Funktionenfamilie.
+Sie sind Funktionen differenzierbaren Funktionen $f(r)$ für $r>0$
+mit $f'(r)=0$ und für $r\to\infty$ nimmt $f(r)$ so schnell ab, dass
+auch $rf(r)$ noch gegen $0$ strebt.
+Das Skalarprodukt ist
+\[
+\langle f,g\rangle
+=
+\int_0^\infty r f(r) g(r)\,dr,
+\]
+als Operator verwenden wir
+\[
+A = \frac{d^2}{dr^2} + \frac{1}{r}\frac{d}{dr} + s(r),
+\]
+wobei $s(r)$ eine beliebige integrierbare Funktion sein kann.
+Zunächst überprüfen wir, ob dieser Operator wirklich selbstadjungiert ist.
+Dazu rechnen wir
+\begin{align}
+\langle Af,g\rangle
+&=
+\int_0^\infty
+r\,\biggl(f''(r)+\frac1rf'(r)+s(r)f(r)\biggr) g(r)
+\,dr
+\notag
+\\
+&=
+\int_0^\infty rf''(r)g(r)\,dr
++
+\int_0^\infty f'(r)g(r)\,dr
++
+\int_0^\infty s(r)f(r)g(r)\,dr.
+\notag
+\intertext{Der letzte Term ist symmetrisch in $f$ und $g$, daher
+ändern wir daran weiter nichts.
+Auf das erste Integral kann man partielle Integration anwenden und erhält}
+&=
+\biggl[rf'(r)g(r)\biggr]_0^\infty
+-
+\int_0^\infty f'(r)g(r) + rf'(r)g'(r)\,dr
++
+\int_0^\infty f'(r)g(r)\,dr
++
+\int_0^\infty s(r)f(r)g(r)\,dr.
+\notag
+\intertext{Der erste Term verschwindet wegen der Bedingungen an die
+Funktionen $f$ und $g$.
+Der erste Term im zweiten Integral hebt sich gegen das
+zweite Integral weg.
+Der letzte Term ist das Skalarprodukt von $f'$ und $g'$.
+Somit ergibt sich
+}
+&=
+-\langle f',g'\rangle
++
+\int_0^\infty s(r) f(r)g(r)\,dr.
+\label{buch:integrale:orthogonal:besselsa}
+\end{align}
+Vertauscht man die Rollen von $f$ und $g$, erhält man das Gleiche, da im
+letzten Ausdruck~\eqref{buch:integrale:orthogonal:besselsa} die Funktionen
+$f$ und $g$ symmetrische auftreten.
+Damit ist gezeigt, dass der Operator $A$ selbstadjungiert ist.
+Es folgt nun, dass Eigenvektoren des Operators $A$ automatisch
+orthogonal sind.
+
+Eigenfunktionen von $A$ sind aber Lösungen der Differentialgleichung
+\[
+\begin{aligned}
+&&
+Af&=\lambda f
+\\
+&\Rightarrow\qquad&
+f''(r) +\frac1rf'(r) + s(r)f(r) &= \lambda f(r)
+\\
+&\Rightarrow\qquad&
+r^2f''(r) +rf'(r)+ (-\lambda r^2+s(r)r^2)f(r) &= 0
+\end{aligned}
+\]
+sind.
+
+Durch die Wahl $s(r)=1$ wird der Operator $A$ zum Bessel-Operator
+$B$ definiert in
+\eqref{buch:differentialgleichungen:bessel-operator}.
+Die Lösungen der Besselschen Differentialgleichung zu verschiedenen Werten
+des Parameters müssen also orthogonal sein, insbesondere sind die
+Besselfunktion $J_\nu(r)$ und $J_\mu(r)$ orthogonal wenn $\mu\ne\nu$ ist.
+
diff --git a/buch/chapters/070-orthogonalitaet/chapter.tex b/buch/chapters/070-orthogonalitaet/chapter.tex
new file mode 100644
index 0000000..4c6019f
--- /dev/null
+++ b/buch/chapters/070-orthogonalitaet/chapter.tex
@@ -0,0 +1,31 @@
+%
+% chapter.tex -- Spezielle Funktionen definiert durch Integrale
+%
+% (c) 2021 Prof Dr Andreas Müller, Hochschule Rapperswil
+%
+% !TeX spellcheck = de_CH
+\chapter{Orthogonalität
+\label{buch:chapter:orthogonalitaet}}
+\lhead{Orthogonalität}
+\rhead{}
+\input{chapters/070-orthogonalitaet/orthogonal.tex}
+\input{chapters/070-orthogonalitaet/jacobi.tex}
+\input{chapters/070-orthogonalitaet/legendredgl.tex}
+\input{chapters/070-orthogonalitaet/bessel.tex}
+\input{chapters/070-orthogonalitaet/sturm.tex}
+\input{chapters/070-orthogonalitaet/gaussquadratur.tex}
+
+\section{TODO}
+\begin{itemize}
+\item Jacobi-Polynome
+\item Tschebyscheff-Polynome
+\end{itemize}
+
+\section*{Übungsaufgaben}
+\rhead{Übungsaufgaben}
+\aufgabetoplevel{chapters/070-orthogonalitaet/uebungsaufgaben}
+\begin{uebungsaufgaben}
+%\uebungsaufgabe{0}
+%\uebungsaufgabe{1}
+\end{uebungsaufgaben}
+
diff --git a/buch/chapters/070-orthogonalitaet/gaussquadratur.tex b/buch/chapters/070-orthogonalitaet/gaussquadratur.tex
new file mode 100644
index 0000000..870c8a8
--- /dev/null
+++ b/buch/chapters/070-orthogonalitaet/gaussquadratur.tex
@@ -0,0 +1,485 @@
+%
+% Anwendung: Gauss-Quadratur
+%
+\section{Anwendung: Gauss-Quadratur}
+\rhead{Gauss-Quadratur}
+Orthogonale Polynome haben eine etwas unerwartet Anwendung in einem
+von Gauss erdachten numerischen Integrationsverfahren.
+Es basiert auf der Beobachtung, dass viele Funktionen sich sehr
+gut durch Polynome approximieren lassen.
+Wenn man also sicherstellt, dass ein Verfahren für Polynome
+sehr gut funktioniert, darf man auch davon ausgehen, dass es für
+andere Funktionen nicht allzu schlecht sein wird.
+
+\subsection{Interpolationspolynome}
+Zu einer stetigen Funktion $f(x)$ auf dem Intervall $[-1,1]$ 
+ist ein Polynome vom Grad $n$ gesucht, welches in den Punkten
+$x_0<x_1<\dots<x_n$ die Funktionswerte $f(x_i)$ annimmt.
+Ein solches Polynom $p(x)$ hat $n+1$ Koeffizienten, die aus dem
+linearen Gleichungssystem der $n+1$ Gleichungen $p(x_i)=f(x_i)$ 
+ermittelt werden können.
+
+Das Interpolationspolynom $p(x)$ lässt sich abera uch direkt 
+angeben.
+Dazu konstruiert man zuerst die Polynome
+\[
+l_i(x)
+=
+\frac{
+(x-x_0)(x-x_1)\cdots\widehat{(x-x_i)}\cdots (x-x_n)
+}{
+(x_i-x_0)(x_i-x_1)\cdots\widehat{(x_i-x_i)}\cdots (x_i-x_n)
+}
+\]
+vom Grad $n$, wobei der Hut bedeutet, dass diese Faktoren
+im Produkt wegzulassen sind.
+Die Polynome $l_i(x)$ haben die Eigenschaft
+\[
+l_i(x_j) = \delta_{ij}
+=
+\begin{cases}
+1&\qquad i=j\\
+0&\qquad\text{sonst}.
+\end{cases}
+\]
+Die Linearkombination
+\[
+p(x) = \sum_{i=0}^n f(x_i)l_i(x)
+\]
+ist dann ein Polynom vom Grad $n$, welches am den Stellen $x_j$
+die Werte
+\[
+p(x_j) 
+=
+\sum_{i=0}^n f(x_i)l_i(x_j)
+=
+\sum_{i=0}^n f(x_i)\delta_{ij}
+=
+f(x_j)
+\]
+hat, das Polynome $p(x)$ ist also das gesuchte Interpolationspolynom.
+
+\subsection{Integrationsverfahren auf der Basis von Interpolation}
+Das Integral einer stetigen Funktion $f(x)$ auf dem Intervall $[-1,1]$
+kann mit Hilfe des Interpolationspolynoms approximiert werden.
+Wenn $|f(x)-p(x)|<\varepsilon$ ist im Intervall $[-1,1]$, dann gilt
+für die Integrale
+\[
+\biggl|\int_{-1}^1 f(x)\,dx -\int_{-1}^1p(x)\,dx\biggr|
+\le
+\int_{-1}^1 |f(x)-p(x)|\,dx
+\le
+2\varepsilon.
+\]
+Ein Interpolationspolynom mit kleinem Fehler liefert also auch
+eine gute Approximation für das Integral.
+
+Da das Interpolationspolynome durch die Funktionswerte $f(x_i)$
+bestimmt ist, muss auch das Integral allein aus diesen Funktionswerten
+berechnet werden können.
+Tatsächlich ist
+\begin{equation}
+\int_{-1}^1 p(x)\,dx
+=
+\int_{-1}^1 \sum_{i=0}^n f(x_i)l_i(x)\,dx
+=
+\sum_{i=0}^n f(x_i)
+\underbrace{\int_{-1}^1
+l_i(x)\,dx}_{\displaystyle = A_i}.
+\label{buch:integral:gaussquadratur:eqn:Aidef}
+\end{equation}
+Das Integral von $f(x)$ wird also durch eine mit den Zahlen $A_i$
+gewichtete Summe
+\[
+\int_{-1}^1 f(x)\,dx
+\approx
+\sum_{i=1}^n f(x_i)A_i
+\]
+approximiert.
+
+\subsection{Integrationsverfahren, die für Polynome exakt sind}
+Ein Polynom vom Grad $2n$ hat $2n+1$ Koeffizienten.
+Um das Polynom durch ein Interpolationspolynom exakt wiederzugeben,
+braucht man $2n+1$ Stützstellen.
+Andererseits gilt
+\[
+\int_{-1}^1 a_{2n}x^{2n} + a_{2n-1}x^{2n-1} + \dots + a_2x^2 + a_1x + a_0\,dx
+=
+\int_{-1}^1 a_{2n}x^{2n} + a_{2n-2}x^{2n-2}+\dots +a_2x^2 +a_0\,dx,
+\]
+das Integral ist also bereits durch die $n+1$ Koeffizienten mit geradem
+Index bestimmt.
+Es sollte daher möglich sein, aus $n+1$ Funktionswerten eines beliebigen
+Polynoms vom Grad höchstens $2n$ an geeignet gewählten Stützstellen das
+Integral exakt zu bestimmen.
+
+\begin{beispiel}
+Wir versuchen dies für quadratische Polynome durchzuführen, also 
+für $n=1$.
+Gesucht sind also zwei Werte $x_i$, $i=0,1$ und Gewichte $A_i$, $i=0,1$
+derart, dass für jedes quadratische Polynome $p(x)=a_2x^2+a_1x+a_0$ 
+das Integral durch
+\[
+\int_{-1}^1 p(x)\,dx
+=
+A_0 p(x_0) + A_1 p(x_1)
+\]
+gebeben ist.
+Indem wir für $p(x)$ die Polynome $1$, $x$, $x^2$ und $x^3$ einsetzen,
+erhalten wir vier Gleichungen
+\[
+\begin{aligned}
+p(x)&=\rlap{$1$}\phantom{x^2}\colon& 2       &= A_0\phantom{x_0}+ A_1     \\
+p(x)&=x^{\phantom{2}}\colon& 0       &= A_0x_0   + A_1x_1  \\
+p(x)&=x^2\colon& \frac23 &= A_0x_0^2 + A_1x_1^2\\
+p(x)&=x^3\colon& 0       &= A_0x_0^3 + A_1x_1^3.
+\end{aligned}
+\]
+Dividiert man die zweite und vierte Gleichung in der Form
+\[
+\left.
+\begin{aligned}
+A_0x_0 &= -A_1x_1\\
+A_0x_0^2 &= -A_1x_1^2
+\end{aligned}
+\quad
+\right\}
+\quad
+\Rightarrow
+\quad
+x_0^2=x_1^2
+\quad
+\Rightarrow
+\quad
+x_1=-x_0.
+\]
+Indem wir dies in die zweite Gleichung einsetzen, finden wir 
+\[
+0 = A_0x_0 + A_1x_1 = A_0x_1 -A_1x_0 = (A_0-A_1)x_0
+\quad\Rightarrow\quad
+A_0=A_1.
+\]
+Aus der ersten Gleichung folgt jetzt
+\[
+2= A_0+A_1 = 2A_0 \quad\Rightarrow\quad A_0 = 1.
+\]
+Damit bleiben nur noch die Werte von $x_i$ zu bestimmen, was 
+mit Hilfe der zweiten Gleichung geschehen kann:
+\[
+\frac23 = A_0x_0^2 + A_1x_1^2 = 2x_0^2
+\quad\Rightarrow\quad
+x_0 = \frac{1}{\sqrt{3}}, x_1 = -\frac{1}{\sqrt{3}}
+\]
+Damit ist das Problem gelöst: das Integral eines Polynoms vom Grad 3
+im Interval $[-1,1]$ ist exakt gegeben durch
+\[
+\int_{-1}^1 p(x)\,dx
+=
+p\biggl(-\frac{1}{\sqrt{3}}\biggr)
++
+p\biggl(\frac{1}{\sqrt{3}}\biggr).
+\]
+Das Integral kann also durch nur zwei Auswertungen des Polynoms
+exakt bestimmt werden.
+
+Im Laufe der Lösung des Gleichungssystems wurden die Gewichte $A_i$
+mit bestimmt.
+Es ist aber auch möglich, die Gewichte zu bestimmen, wenn man die
+Stützstellen kennt.
+Nach \eqref{buch:integral:gaussquadratur:eqn:Aidef}
+sind sie die $A_i$ gegeben als Integrale der Polynome
+$l_i(x)$, die im vorliegenden Fall linear sind:
+\begin{align*}
+l_0(x)
+&=
+\frac{x-x_1}{x_0-x_1}
+=
+\frac{x-\frac1{\sqrt{3}}}{-\frac{2}{\sqrt{3}}}
+=
+\frac12(1-\sqrt{3}x)
+\\
+l_1(x)
+&=
+\frac{x-x_0}{x_1-x_0}
+=
+\frac{x+\frac1{\sqrt{3}}}{\frac{2}{\sqrt{3}}}
+=
+\frac12(1+\sqrt{3}x)
+\end{align*}
+Diese haben die Integrale
+\[
+\int_{-1}^1\frac12(1\pm\sqrt{3}x)\,dx
+=
+\int_{-1}^1 \frac12\,dx
+=
+1,
+\]
+da das Polynom $x$ verschwindendes Integral hat.
+Dies stimmt mit $A_0=A_1=1$ überein.
+\label{buch:integral:beispiel:gaussquadraturn1}
+\end{beispiel}
+
+Das eben vorgestellt Verfahren kann natürlich auf beliebiges $n$
+verallgemeinert werden.
+Allerdings ist die Rechnung zur Bestimmung der Stützstellen und
+Gewichte sehr mühsam.
+
+\subsection{Stützstellen und Orthogonalpolynome}
+Sei $R_n=\{p(X)\in\mathbb{R}[X] \mid \deg p\le n\}$ der Vektorraum
+der Polynome vom Grad $n$.
+
+\begin{satz}
+\label{buch:integral:satz:gaussquadratur}
+Sei $p$ ein Polynom vom Grad $n$, welches auf allen Polynomen in $R_{n-1}$
+orthogonal sind.
+Seien ausserdem $x_0<x_1<\dots<x_n$ Stützstellen im Intervall $[-1,1]$ 
+und $A_i\in\mathbb{R}$ Gewichte derart dass
+\[
+\int_{-1}^1 f(x)\,dx =
+\sum_{i=0}^n A_if(x_i)
+\]
+für jedes Polynom $f$ vom Grad höchstens $2n-1$, dann sind die Zahlen
+$x_i$ die Nullstellen des Polynoms $p$.
+\end{satz}
+
+\begin{proof}[Beweis]
+Sei $f(x)$ ein beliebiges Polynom vom Grad $2n-1$.
+Nach dem Polynomdivisionsalgorithmus gibt es
+Polynome $q,r\in R_{n-1}$ derart, dass $f=qp+r$.
+Dann ist das Integral von $f$ gegeben durch
+\[
+\int_{-1}^1 f(x)\,dx
+=
+\int_{-1}^1q(x) p(x)\,dx + \int_{-1}^1 r(x)\,dx
+=
+\langle q,p\rangle + \int_{-1}^1 r(x)\,dx.
+\]
+Da $p\perp R_{n-1}$ folgt insbesondere, dass $\langle q,p\rangle=0$.
+
+Da die Integrale auch aus den Werten in den Stützstellen berechnet
+werden können, muss auch
+\[
+0
+=
+\int_{-1}^1 q(x)p(x)\,dx
+=
+\sum_{i=0}^n q(x_i)p(x_i)
+\]
+für jedes beliebige Polynom $q\in R_{n-1}$ gelten.
+Da man für $q$ die Interpolationspolynome $l_j(x)$ verwenden
+kann, den Grad $n-1$ haben, folgt
+\[
+0
+=
+\sum_{i=0}^n
+l_j(x_i)p(x_i)
+=
+\sum_{i=0}^n \delta_{ij}p(x_i),
+\]
+die Stützstellen $x_i$ müssen also die Nullstellen des Polynoms
+$p(x)$ sein.
+\end{proof}
+
+Der Satz~\ref{buch:integral:satz:gaussquadratur} begründet das
+{\em Gausssche Quadraturverfahren}.
+Die in Abschnitt~\ref{buch:integral:section:orthogonale-polynome}
+bestimmten Legendre-Polynome $P_n$ haben die im Satz
+verlangte Eigenschaft,
+dass sie auf allen Polynomen geringeren Grades orthogonal sind.
+Wählt man die $n$ Nullstellen von $P_n$ als Stützstellen, erhält man 
+automatisch ein Integrationsverfahren, welches für Polynome vom Grad
+$2n-1$ exakt ist.
+
+\begin{beispiel}
+Das Legendre-Polynom $P_2(x) = \frac12(3x^2-1)$ hat die
+Nullstellen $x=\pm1/\sqrt{3}$, dies sind genau die im Beispiel
+auf Seite~\pageref{buch:integral:beispiel:gaussquadraturn1} befundenen
+Sützstellen.
+\end{beispiel}
+
+\subsection{Fehler der Gauss-Quadratur}
+Das Gausssche Quadraturverfahren mit $n$ Stützstellen berechnet
+Integrale von Polynomen bis zum Grad $2n-1$ exakt.
+Für eine beliebige Funktion kann man die folgende Fehlerabschätzung
+angeben \cite[theorem 7.3.4, p.~497]{buch:numal}.
+
+\begin{satz}
+Seien $x_i$ die Stützstellen und $A_i$ die Gewichte einer
+Gaussschen Quadraturformel mit $n+1$ Stützstellen und sei $f$
+eine auf dem Interval $[-1,1]$ $2n+2$-mal stetig differenzierbare
+Funktion, dann ist der $E$ Fehler des Integrals
+\[
+\int_{-1}^1 f(x)\,dx = \sum_{i=0}^n A_i f(x_i) + E
+\]
+gegeben durch
+\begin{equation}
+E = \frac{f^{(2n+2)}(\xi)}{(2n+2)!}\int_{-1}^1 l(x)^2\,dx,
+\label{buch:integral:gaussquadratur:eqn:fehlerformel}
+\end{equation}
+wobei $l(x)=(x-x_0)(x-x_1)\dots(x-x_n)$  und $\xi$ ein geeigneter
+Wert im Intervall $[-1,1]$ ist.
+\end{satz}
+
+Dank dem Faktor $(2n+2)!$ im Nenner von
+\eqref{buch:integral:gaussquadratur:eqn:fehlerformel}
+geht der Fehler für grosses $n$ sehr schnell gegen $0$.
+Man kann auch zeigen, dass die mit Gauss-Quadratur mit $n+1$
+Stützstellen berechneten Näherungswerte eines Integrals einer
+stetigen Funktion $f(x)$ für $n\to\infty$ immer gegen den wahren
+Wert des Integrals konvergieren.
+
+\begin{table}
+\def\u#1{\underline{#1}}
+\centering
+\begin{tabular}{|>{$}c<{$}|>{$}r<{$}|>{$}r<{$}|}
+\hline
+           n & \text{Gauss-Quadratur} & \text{Trapezregel} \\
+\hline
+\phantom{0}2 & 0.\u{95}74271077563381 & 0.\u{95}63709682242596 \\
+\phantom{0}4 & 0.\u{95661}28333449730 & 0.\u{956}5513401768598 \\
+\phantom{0}6 & 0.\u{9566114}812034364 & 0.\u{956}5847489712136 \\
+\phantom{0}8 & 0.\u{956611477}5028123 & 0.\u{956}5964425360520 \\
+          10 & 0.\u{9566114774905}637 & 0.\u{9566}018550715587 \\
+          12 & 0.\u{956611477490518}7 & 0.\u{9566}047952369826 \\
+          14 & 0.\u{95661147749051}72 & 0.\u{9566}065680717177 \\
+          16 & 0.\u{956611477490518}7 & 0.\u{9566}077187127541 \\
+          18 & 0.\u{956611477490518}3 & 0.\u{9566}085075898731 \\
+          20 & 0.\u{956611477490518}4 & 0.\u{9566}090718697414 \\
+\hline
+      \infty & 0.9566114774905183 & 0.9566114774905183 \\
+\hline
+\end{tabular}
+\caption{Integral von $\sqrt{1-x^2}$ zwischen $-\frac12$ und $\frac12$ 
+berechnet mit Gauss-Quadratur und der Trapezregel, aber mit zehnmal
+so vielen Stützstellen.
+Bereits mit 12 Stützstellen erreicht die Gauss-Quadratur
+Maschinengenauigkeit, die Trapezregel liefert auch mit 200 Stützstellen
+nicht mehr als 4 korrekte Nachkommastellen.
+\label{buch:integral:gaussquadratur:table0.5}}
+\end{table}
+
+%\begin{table}
+%\def\u#1{\underline{#1}}
+%\centering
+%\begin{tabular}{|>{$}c<{$}|>{$}r<{$}|>{$}r<{$}|}
+%\hline
+%           n & \text{Gauss-Quadratur} & \text{Trapezregel} \\
+%\hline
+%\phantom{0}2 & 1.\u{5}379206741571556 & 1.\u{5}093105464758343 \\
+%\phantom{0}4 & 1.\u{51}32373472933831 & 1.\u{51}13754509594814 \\
+%\phantom{0}6 & 1.\u{512}1624557410367 & 1.\u{51}17610879524799 \\
+%\phantom{0}8 & 1.\u{51207}93479994321 & 1.\u{51}18963282632112 \\
+%          10 & 1.\u{51207}13859966004 & 1.\u{51}19589735776959 \\
+%          12 & 1.\u{512070}5317779943 & 1.\u{51}19930161260693 \\
+%          14 & 1.\u{5120704}334802813 & 1.\u{5120}135471596636 \\
+%          16 & 1.\u{5120704}216176006 & 1.\u{5120}268743889558 \\
+%          18 & 1.\u{5120704}201359081 & 1.\u{5120}360123137213 \\
+%          20 & 1.\u{5120704199}459651 & 1.\u{5120}425490275837 \\
+%\hline
+%      \infty & 1.5120704199172947 & 1.5120704199172947 \\
+%\hline
+%\end{tabular}
+%\end{table}
+
+%\begin{table}
+%\def\u#1{\underline{#1}}
+%\centering
+%\begin{tabular}{|>{$}c<{$}|>{$}r<{$}|>{$}r<{$}|}
+%\hline
+%           n & \text{Gauss-Quadratur} & \text{Trapezregel} \\
+%\hline
+%\phantom{0}2 & 1.\u{}6246862220133462 & 1.\u{5}597986803933712 \\
+%\phantom{0}4 & 1.\u{5}759105515463101 & 1.\u{56}63563456168101 \\
+%\phantom{0}6 & 1.\u{5}706630058381434 & 1.\u{56}77252866190838 \\
+%\phantom{0}8 & 1.\u{56}94851106536780 & 1.\u{568}2298707696152 \\
+%          10 & 1.\u{56}91283195332679 & 1.\u{568}4701957758742 \\
+%          12 & 1.\u{56}90013806299465 & 1.\u{568}6030805941198 \\
+%          14 & 1.\u{5689}515434853885 & 1.\u{568}6841603070025 \\
+%          16 & 1.\u{5689}306507843050 & 1.\u{568}7372230731711 \\
+%          18 & 1.\u{5689}214761291217 & 1.\u{568}7738235496322 \\
+%          20 & 1.\u{56891}73062385982 & 1.\u{568}8001228530786 \\
+%\hline
+%      \infty & 1.5689135396691616 & 1.5689135396691616 \\
+%\hline
+%\end{tabular}
+%\end{table}
+
+\begin{table}
+\def\u#1{\underline{#1}}
+\centering
+\begin{tabular}{|>{$}c<{$}|>{$}r<{$}|>{$}r<{$}|}
+\hline
+           n & \text{Gauss-Quadratur} & \text{Trapezregel} \\
+\hline
+\phantom{0}2 & 1.\u{}6321752373234928 & 1.\u{5}561048774629949 \\
+\phantom{0}4 & 1.\u{57}98691557134743 & 1.\u{5}660124134617943 \\
+\phantom{0}6 & 1.\u{57}35853681692993 & 1.\u{5}683353001877542 \\
+\phantom{0}8 & 1.\u{57}19413565928206 & 1.\u{5}692627503425400 \\
+          10 & 1.\u{57}13388119633434 & 1.\u{5}697323578543481 \\
+          12 & 1.\u{57}10710489948883 & 1.\u{570}0051217458713 \\
+          14 & 1.\u{570}9362135398341 & 1.\u{570}1784766276063 \\
+          16 & 1.\u{570}8621102742815 & 1.\u{570}2959121005231 \\
+          18 & 1.\u{570}8186779483588 & 1.\u{570}3793521168343 \\
+          20 & 1.\u{5707}919411931615 & 1.\u{570}4408749735932 \\
+\hline
+      \infty & 1.5707367072605671 & 1.5707367072605671 \\
+\hline
+\end{tabular}
+\caption{Integral von $\sqrt{1-x^2}$ zwischen $-0.999$ und $0.999$ 
+berechnet mit Gauss-Quadratur und der Trapezregel, aber mit zehnmal
+so vielen Stützstellen.
+Wegen der divergierenden Steigung des Integranden bei $\pm 1$ tun
+sich beide Verfahren sehr schwer. 
+Trotzdem erreich die Gauss-Quadrator 4 korrekte Nachkommastellen
+mit 20 Stütztstellen, während die Trapezregel auch mit 200 Stützstellen
+nur 3 korrekte Nachkommastellen findet.
+\label{buch:integral:gaussquadratur:table0.999}}
+\end{table}
+
+\begin{figure}
+\centering
+\includegraphics{chapters/060-integral/gq/gq.pdf}
+\caption{Approximationsfehler des
+Integrals~\eqref{buch:integral:gaussquadratur:bspintegral}
+in Abhängigkeit von $a$.
+Die Divergenz der Ableitung des Integranden an den Intervallenden
+$\pm 1$ führt zu schlechter Konvergenz des Verfahrens, wenn $a$
+nahe an $1$ ist.
+\label{buch:integral:gaussquadratur:fehler}}
+\end{figure}
+
+Zur Illustration der Genauigkeit der Gauss-Quadratur berechnen wir
+das Integral
+\begin{equation}
+\int_{-a}^a \sqrt{1-x^2}\,dx
+=
+\arcsin a + a \sqrt{1-a^2}
+\label{buch:integral:gaussquadratur:bspintegral}
+\end{equation}
+mit Gauss-Quadratur einerseits und dem Trapezverfahren
+andererseits.
+Da Gauss-Quadratur mit sehr viel weniger Sützstellen auskommt,
+berechnen wir die Trapeznäherung mit zehnmal so vielen Stützstelln.
+In den Tabellen~\ref{buch:integral:gaussquadratur:table0.5}
+und
+\ref{buch:integral:gaussquadratur:table0.999}
+sind die Resultate zusammengestellt.
+Für $a =\frac12$ zeigt
+Tabelle~\ref{buch:integral:gaussquadratur:table0.5}
+die sehr schnelle Konvergenz der Gauss-Quadratur, schon mit
+12 Stützstellen wird Maschinengenauigkeit erreicht.
+Das Trapezverfahren dagegen erreicht auch mit 200 Stützstellen nur
+4 korrekte Nachkommastellen.
+
+An den Stellen $x=\pm 1$ divergiert die Ableitung des Integranden
+des Integrals \eqref{buch:integral:gaussquadratur:bspintegral}.
+Da grösste und kleinste Stützstelle der Gauss-Quadratur immer
+deutlich vom Rand des Intervalls entfernt ist, kann das Verfahren
+diese ``schwierigen'' Stellen nicht erkennen.
+Tabelle~\ref{buch:integral:gaussquadratur:table0.999} zeigt, wie
+die Konvergenz des Verfahrens in diesem Fall sehr viel schlechter ist.
+Dies zeigt auch der Graph in
+Abbildung~\ref{buch:integral:gaussquadratur:fehler}.
+
+\subsection{Skalarprodukte mit Gewichtsfunktion}
+
diff --git a/buch/chapters/070-orthogonalitaet/images/Makefile b/buch/chapters/070-orthogonalitaet/images/Makefile
new file mode 100644
index 0000000..e3a988a
--- /dev/null
+++ b/buch/chapters/070-orthogonalitaet/images/Makefile
@@ -0,0 +1,14 @@
+#
+# Makefile -- Bilder zum Kapitel über durch Integrale definierte spezielle
+#             Funktionen
+#
+# (c) 2021 Prof Dr Andreas Müller, OST Ostschweizer Fachhochschule
+#
+all:	legendre.pdf orthogonal.pdf
+
+legendrepaths.tex:	legendre.m
+	octave legendre.m
+legendre.pdf:	legendre.tex legendrepaths.tex
+	pdflatex legendre.tex
+orthogonal.pdf:	orthogonal.tex legendrepaths.tex
+	pdflatex orthogonal.tex
diff --git a/buch/chapters/070-orthogonalitaet/images/legendre.m b/buch/chapters/070-orthogonalitaet/images/legendre.m
new file mode 100644
index 0000000..8e8317d
--- /dev/null
+++ b/buch/chapters/070-orthogonalitaet/images/legendre.m
@@ -0,0 +1,64 @@
+#
+# legendre.m
+#
+# (c) 2021 Prof Dr Andreas Müller, OST Ostschweizer Fachhochschule
+#
+pkg load miscellaneous
+global N;
+N = 30;
+
+function retval = legendrepath(fn, n, name)
+	global N;
+	m = n * N;
+	c = legendrepoly(n)
+	x = (-m:m)/m;
+	v = polyval(c, x);
+	fprintf(fn, "\\def\\legendre%s{\n", name)
+	fprintf(fn, "\t   ({%.5f*\\dx},{%.5f*\\dy})", x(1), v(1));
+	for i = (2:length(v))
+		fprintf(fn, "\n\t-- ({%.5f*\\dx},{%.5f*\\dy})", x(i), v(i));
+		
+	endfor
+	fprintf(fn, "\n}\n");
+	ci = polyint(conv(c, c))
+polyval(ci, 1)
+	normalization = sqrt(polyval(ci, 1) - polyval(ci, -1))
+	fprintf(fn, "\\def\\normalization%s{%.5f}\n", name, normalization);
+endfunction
+
+function retval = legendreprodukt(fn, a, b, name)
+	global N;
+	n = max(a, b);
+	m = n * N;
+	pa = legendrepoly(a);
+	pb = legendrepoly(b);
+	p = conv(pa, pb);
+	x = (-m:m)/m;
+	v = polyval(p, x);
+	fprintf(fn, "\\def\\produkt%s{\n", name)
+	fprintf(fn, "\t   ({%.5f*\\dx},{%.5f*\\dy})", x(1), v(1));
+	for i = (2:length(v))
+		fprintf(fn, "\n\t-- ({%.5f*\\dx},{%.5f*\\dy})", x(i), v(i));
+	endfor
+	fprintf(fn, "\n}\n");
+endfunction
+
+fn = fopen("legendrepaths.tex", "w");
+legendrepath(fn,  1, "one");
+legendrepath(fn,  2, "two");
+legendrepath(fn,  3, "three");
+legendrepath(fn,  4, "four");
+legendrepath(fn,  5, "five");
+legendrepath(fn,  6, "six");
+legendrepath(fn,  7, "seven");
+legendrepath(fn,  8, "eight");
+legendrepath(fn,  9, "nine");
+legendrepath(fn, 10, "ten");
+
+legendreprodukt(fn, 4, 7, "ortho");
+legendreprodukt(fn, 4, 4, "vier");
+legendreprodukt(fn, 7, 7, "sieben");
+
+fclose(fn);
+
+
diff --git a/buch/chapters/070-orthogonalitaet/images/legendre.pdf b/buch/chapters/070-orthogonalitaet/images/legendre.pdf
new file mode 100644
index 0000000..a893c26
--- /dev/null
+++ b/buch/chapters/070-orthogonalitaet/images/legendre.pdf
diff --git a/buch/chapters/070-orthogonalitaet/images/legendre.tex b/buch/chapters/070-orthogonalitaet/images/legendre.tex
new file mode 100644
index 0000000..8409da0
--- /dev/null
+++ b/buch/chapters/070-orthogonalitaet/images/legendre.tex
@@ -0,0 +1,99 @@
+%
+% legendre.tex -- plots of legendre polynomials
+%
+% (c) 2021 Prof Dr Andreas Müller, OST Ostschweizer Fachhochschule
+%
+\documentclass[tikz]{standalone}
+\usepackage{amsmath}
+\usepackage{times}
+\usepackage{txfonts}
+\usepackage{pgfplots}
+\usepackage{csvsimple}
+\usetikzlibrary{arrows,intersections,math}
+\begin{document}
+\def\skala{6.5}
+\input{legendrepaths.tex}
+\begin{tikzpicture}[>=latex,thick,scale=\skala]
+
+\definecolor{fone}{rgb}{1,0,0}
+\definecolor{ftwo}{rgb}{1.0,0,0.8}
+\definecolor{fthree}{rgb}{0.8,0,1}
+\definecolor{ffour}{rgb}{0,0,1}
+\definecolor{ffive}{rgb}{0,0.8,1}
+\definecolor{fsix}{rgb}{0,1,1}
+\definecolor{fseven}{rgb}{0,0.6,0}
+\definecolor{feight}{rgb}{0.2,1,0.6}
+\definecolor{fnine}{rgb}{0.6,0.8,0.2}
+\definecolor{ften}{rgb}{1,0.4,0}
+
+\def\dx{1}
+\def\dy{0.25}
+
+\def\achsen{
+	\draw[->] ({-1-(0.1/\skala)},0) -- ({1+(0.3/\skala)},0)
+		coordinate[label={$x$}];
+	\draw[->] (0,{-(\dy)-(0.1/\skala)}) -- (0,{(\dy)+(0.3/\skala)})
+		coordinate[label={right:$y$}];
+	\foreach \x in {-1,-0.9,...,1.001}{
+		\draw ({\dx*\x},{-0.05/\skala}) -- ({\dx*\x},{0.05/\skala});
+	}
+	\foreach \x in {-1,-0.5,0.5,1}{
+		\draw ({-0.05/\skala},{\dy*\x}) -- ({0.05/\skala},{\dy*\x});
+		\node at ({\dx*\x},{-0.05/\skala}) [below] {$\mathstrut\x$};
+		\node at ({-0.05/\skala},{\dy*\x}) [left] {$\mathstrut\x$};
+	}
+}
+
+\begin{scope}[yshift=0cm]
+	\node[color=fone] at (-0.70,{-0.9*\dy}) [right] {$n=1\mathstrut$};
+	\node[color=ftwo] at (-0.90,{0.9*\dy}) [right] {$n=2\mathstrut$};
+	\draw[line width=1.4pt,color=fone] \legendreone;
+	\draw[line width=1.4pt,color=ftwo] \legendretwo;
+	\achsen
+\end{scope}
+
+\begin{scope}[yshift=-0.6cm]
+	\node[color=fthree] at (-0.96,{-0.9*\dy}) [right] {$n=3\mathstrut$};
+	\node[color=ffour] at (-0.96,{0.9*\dy}) [right] {$n=4\mathstrut$};
+	\draw[line width=1.4pt,color=fthree] \legendrethree;
+	\draw[line width=1.4pt,color=ffour] \legendrefour;
+	\achsen
+\end{scope}
+
+\begin{scope}[yshift=-1.2cm]
+	\node[color=ffive] at (-0.96,{-0.9*\dy}) [right] {$n=5\mathstrut$};
+	\node[color=fsix] at (-0.96,{0.9*\dy}) [right] {$n=6\mathstrut$};
+	\draw[line width=1.4pt,color=ffive] \legendrefive;
+	\draw[line width=1.4pt,color=fsix] \legendresix;
+	\achsen
+\end{scope}
+
+\begin{scope}[yshift=-1.8cm]
+	\node[color=fseven] at (-0.96,{-0.9*\dy}) [right] {$n=7\mathstrut$};
+	\node[color=feight] at (-0.96,{0.9*\dy}) [right] {$n=8\mathstrut$};
+	\draw[line width=1.4pt,color=fseven] \legendreseven;
+	\draw[line width=1.4pt,color=feight] \legendreeight;
+	\achsen
+\end{scope}
+
+\begin{scope}[yshift=-2.4cm]
+	\node[color=fnine] at (-0.96,{-0.9*\dy}) [right] {$n=9\mathstrut$};
+	\node[color=ften] at (-0.96,{0.9*\dy}) [right] {$n=10\mathstrut$};
+	\draw[line width=1.4pt,color=fnine] \legendrenine;
+	\draw[line width=1.4pt,color=ften] \legendreten;
+	\achsen
+\end{scope}
+
+%\draw[line width=1.4pt,color=ftwo] \legendretwo;
+%\draw[line width=1.4pt,color=fthree] \legendrethree;
+%\draw[line width=1.4pt,color=ffour] \legendrefour;
+%\draw[line width=1.4pt,color=ffive] \legendrefive;
+%\draw[line width=1.4pt,color=fsix] \legendresix;
+%\draw[line width=1.4pt,color=fseven] \legendreseven;
+%\draw[line width=1.4pt,color=feight] \legendreeight;
+%\draw[line width=1.4pt,color=fnine] \legendrenine;
+%\draw[line width=1.4pt,color=ften] \legendreten;
+
+\end{tikzpicture}
+\end{document}
+
diff --git a/buch/chapters/070-orthogonalitaet/images/orthogonal.pdf b/buch/chapters/070-orthogonalitaet/images/orthogonal.pdf
new file mode 100644
index 0000000..960c4ff
--- /dev/null
+++ b/buch/chapters/070-orthogonalitaet/images/orthogonal.pdf
diff --git a/buch/chapters/070-orthogonalitaet/images/orthogonal.tex b/buch/chapters/070-orthogonalitaet/images/orthogonal.tex
new file mode 100644
index 0000000..8600281
--- /dev/null
+++ b/buch/chapters/070-orthogonalitaet/images/orthogonal.tex
@@ -0,0 +1,79 @@
+%
+% orthogonal.tex -- plots of legendre polynomials
+%
+% (c) 2021 Prof Dr Andreas Müller, OST Ostschweizer Fachhochschule
+%
+\documentclass[tikz]{standalone}
+\usepackage{amsmath}
+\usepackage{times}
+\usepackage{txfonts}
+\usepackage{pgfplots}
+\usepackage{csvsimple}
+\usetikzlibrary{arrows,intersections,math}
+\begin{document}
+\def\skala{6}
+\input{legendrepaths.tex}
+\begin{tikzpicture}[>=latex,thick,scale=\skala]
+
+\definecolor{fone}{rgb}{1,0,0}
+\definecolor{ftwo}{rgb}{1.0,0,0.8}
+\definecolor{fthree}{rgb}{0.8,0,1}
+\definecolor{ffour}{rgb}{0,0,1}
+\definecolor{ffive}{rgb}{0,0.8,1}
+\definecolor{fsix}{rgb}{0,1,1}
+\definecolor{fseven}{rgb}{0,0.6,0}
+\definecolor{feight}{rgb}{0.2,1,0.6}
+\definecolor{fnine}{rgb}{0.6,0.8,0.2}
+\definecolor{ften}{rgb}{1,0.4,0}
+
+\def\dx{1}
+\def\Dy{3}
+\def\dy{3}
+
+\begin{scope}
+\clip (-1,-0.6) rectangle (1,1);
+
+%\pgfmathparse{\Dy/(\normalizationfour*\normalizationfour)}
+%\xdef\dy{\pgfmathresult}
+\fill[color=ffour!70,opacity=0.5] (-1,0) -- \produktvier -- (1,0) -- cycle;
+
+%\pgfmathparse{\Dy/(\normalizationeight*\normalizationeight)}
+%\xdef\dy{\pgfmathresult}
+\fill[color=fseven!70,opacity=0.5] (-1,0) -- \produktsieben -- (1,0) -- cycle;
+
+%\pgfmathparse{\Dy/(\normalizationfour*\normalizationeight)}
+%\xdef\dy{\pgfmathresult}
+\fill[color=red!50,opacity=0.5] (-1,0) -- \produktortho -- (1,0) -- cycle;
+
+%\pgfmathparse{\Dy/\normalizationfour}
+%\xdef\dy{\pgfmathresult}
+%\draw[line width=1.4pt,color=ffour] \legendrefour;
+%
+%\pgfmathparse{\Dy/\normalizationeight}
+%\xdef\dy{\pgfmathresult}
+%\draw[line width=1.4pt,color=feight] \legendreeight;
+
+%\pgfmathparse{\Dy/(\normalizationfour*\normalizationeight)}
+%\xdef\dy{\pgfmathresult}
+\draw[line width=1.4pt,color=red] \produktortho;
+
+\end{scope}
+
+\draw[->] ({-1-(0.1/\skala)},0) -- ({1+(0.3/\skala)},0)
+	coordinate[label={$x$}];
+\draw[->] (0,{-{0.2*\Dy}-(0.1/\skala)}) -- (0,{1+(0.3/\skala)})
+	coordinate[label={right:$y$}];
+\foreach \x in {-1,-0.9,...,1.001}{
+	\draw ({\dx*\x},{-0.1/\skala}) -- ({\dx*\x},{0.1/\skala});
+}
+\foreach \y in {-0.2,-0.1,0.1,0.2,0.3}{
+	\draw ({-0.1/\skala},{\Dy*\y}) -- ({0.1/\skala},{\Dy*\y});
+	\node at ({-0.1/\skala},{\Dy*\y}) [left] {$\mathstrut\y$};
+}
+\foreach \x in {-1,-0.5,0.5,1}{
+	\node at ({\dx*\x},{-0.1/\skala}) [below] {$\mathstrut\x$};
+}
+
+\end{tikzpicture}
+\end{document}
+
diff --git a/buch/chapters/070-orthogonalitaet/jacobi.tex b/buch/chapters/070-orthogonalitaet/jacobi.tex
new file mode 100644
index 0000000..fb7d5ff
--- /dev/null
+++ b/buch/chapters/070-orthogonalitaet/jacobi.tex
@@ -0,0 +1,9 @@
+%
+% jacobi.tex
+%
+% (c) 2022 Prof Dr Andreas Müller, OST Ostsdchweizer Fachhochschule
+%
+\section{Jacobi-Polynome
+\label{buch:integrale:subsection:jacobi-polynome}}
+\rhead{Jacobi-Polynome}
+
diff --git a/buch/chapters/070-orthogonalitaet/legendredgl.tex b/buch/chapters/070-orthogonalitaet/legendredgl.tex
new file mode 100644
index 0000000..12555b8
--- /dev/null
+++ b/buch/chapters/070-orthogonalitaet/legendredgl.tex
@@ -0,0 +1,369 @@
+%
+% legendredgl.tex
+%
+% (c) 2021 Prof Dr Andreas Müller, OST Ostschweizer Fachhochschule
+%
+\section{Orthogonale Polynome und Differentialgleichungen}
+\rhead{Differentialgleichungen orthogonaler Polynome}
+Legendre hat einen ganz anderen Zugang zu den nach ihm benannten
+Polynomen gefunden.
+Er hat sie gefunden als die Lösungen einer speziellen Differentialgleichungen.
+In diesem Abschnitt sollen diese Funktionen mit der Potenzreihen-Methode
+wiedergefunden werden.
+Dabei stellt sich heraus, dass diese Polynome auch Eigenfunktionen eines
+selbstadjungierten Differentialgoperator sind.
+Die Orthogonalität wird dann aus einer Verallgemeinerung der bekannten
+Eingeschaft folgen, dass Eigenvektoren einer symmetrischen Matrix zu 
+verschiedenen Eigenwerten orthogonal sind.
+
+\subsection{Legendre-Differentialgleichung}
+Die {\em Legendre-Differentialgleichung} ist die Differentialgleichung
+\begin{equation}
+(1-x^2) y'' - 2x y' + n(n+1) y = 0
+\label{buch:integral:eqn:legendre-differentialgleichung}
+\end{equation}
+für eine Funktion $y(x)$ auf dem Intervall $[-1,1]$.
+
+Sei $y(x)$ eine Lösung der Differentialgleichung
+\eqref{buch:integral:eqn:legendre-differentialgleichung}.
+Setzt man $y_s(x)=y(-x)$ in die Differentialgleichung ein, erhält
+man
+\[
+(1-x^2)y_s''(x) - 2x y'_s(x) + n(n+1)y_s(x)
+=
+(1-x^2)y''(-x) +2x y(-x) +n(n+1)y(-x).
+\]
+Ersetzt man $t=-x$, dann wird daraus
+\[
+(1-x^2)y''(t) -2t y(t) + n(n+1) y(t) = 0
+\]
+aus der Differentialgleichung
+\eqref{buch:integral:eqn:legendre-differentialgleichung}.
+Insbesondere ist die gespiegelte Funktion $y_s(x)$ ebenfalls
+eine Lösung der Differentialgleichung.
+
+Ist $y(x)$ eine Lösung der Differentialgleichung ist, dann lässt
+sie sich in die Summe einer geraden und einer ungeraden Funktion
+\[
+\left.
+\begin{aligned}
+y_g(x) &= \frac{y(x)+y(-x)}{2}\\
+y_u(x) &= \frac{y(x)-y(-x)}{2}
+\end{aligned}
+\quad
+\right\}
+\quad
+\Rightarrow
+\quad
+y(x) = y_g(x) + y_u(x)
+\]
+zerlegen, die als Linearkombinationen der beiden Lösungen
+$y(x)$ und $y_s(x)$ ebenfalls Lösungen der Differentialgleichung
+sind.
+
+\subsection{Potenzreihenlösung}
+Wir suchen eine Lösung in Form einer Potenzreihe um $x=0$ und 
+verwenden dazu den Ansatz
+\[
+y(x) = a_0+a_1x+a_2x^2+ \dots = \sum_{k=0}^\infty a_kx^k.
+\]
+\begin{align*}
+(1-x^2) \sum_{k=2}^\infty k(k-1)a_kx^{k-2}
+-2x\sum_{k=0}^\infty ka_kx^{k-1}
++
+n(n+1)\sum_{k=0}^\infty  a_kx^k
+&=
+0
+\\
+\sum_{k=0}^\infty (k+2)(k+1)a_{k+2}x^k
+-
+\sum_{k=2}^\infty k(k-1)a_kx^k
+-
+2\sum_{k=1}^\infty ka_kx^k
++
+n(n+1)\sum_{k=0}^\infty  a_kx^k
+&=
+0
+\end{align*}
+Die Koeffizienten zur Potenz $k$ sind daher
+\begin{align}
+k&=0:
+&
+0&=
+2a_2+n(n+1)a_0
+\notag
+\\
+&&
+a_2&=-\frac{n(n+1)}{2}a_0
+\notag
+\\
+k&=1:
+&
+0&=
+6a_3-2a_1+n(n+1)a_1
+\notag
+\\
+&&
+a_3&= \frac{2-n(n+1)}{6}a_1
+\notag
+\\
+k&>1:
+&
+0&=
+(k+2)(k+1)a_{k+2} -k(k-1)a_k -2ka_k +n(n+1) a_k
+\notag
+\\
+&&
+a_{k+2}
+&=
+\frac{ k(k+1)-n(n+1) }{(k+2)(k+1)}
+a_k
+\label{buch:integral:legendre-dgl:eqn:akrek}
+\end{align}
+Wenn $a_1=0$ und $a_0\ne 1$ ist, dann ist die Funktion $y(x)$ gerade,
+alle ungeraden Koeffizienten verschwinden.
+Ebenso verschwinden alle geraden Koeffizienten, wenn $a_0=0$ und $a_1\ne 0$.
+Für jede Lösung $y(x)$ der Differentialgleichung ist
+$y_g(x)$ ein Lösung mit $a_1=0$ und $y_u(x)$ eine Lösung mit $a_0=0$.
+Wir können die Diskussion der Lösungen daher auf gerade oder ungerade
+Lösungen einschränken.
+
+Gesucht ist jetzt eine Lösung in Form eines Polynoms.
+In diesem Fall müssen die Koeffizienten $a_k$ ab einem
+gewissen Index verschwinden.
+Dies tritt nach \eqref{buch:integral:legendre-dgl:eqn:akrek} genau
+dann auf, wenn der Zähler für ein $k$ verschwindet.
+Folglich gibt es genau dann Polynomlösungen der Differentialgleichungen,
+wenn $n$ eine natürlich Zahl ist.
+Ausserdem ist die Lösung ein Polynom $\bar{P}_n(x)$ vom Grad $n$.
+Das Polynom soll wieder so normiert sein, dass $\bar{P}_n(1)=1$ ist.
+
+Die Lösungen der Differentialgleichungen können jetzt explizit
+berechnet werden.
+Zunächst ist $\bar{P}_0(x)=1$ und $\bar{P}_1(x)=x$.
+Für $n=2$ setzen wir zunächst $a_0=1$ und $a_1=0$ und erhalten
+\[
+y(x)
+=
+1 + \frac{0(0+1) - 2(2+1)}{(0+2)(0+1)}a_0 x^2
+=
+1
+-3x^2
+\qquad\text{oder}\qquad
+\bar{P}_3(x) = \frac12(3x^2-1).
+\]
+Für $n=3$ starten wir von $a_1=1$ und $a_0=0$, was zunächst $a_2=0$
+impliziert.
+Für $a_3$ finden wir
+\[
+a_3=\frac{1(1+1)-3(3+1)}{(1+2)(1+1)} = -\frac53
+\qquad\Rightarrow\qquad
+y(x) = x-\frac53x^3
+\qquad\Rightarrow\qquad
+\bar{P}_3(x) = \frac12(5x^3-3x).
+\]
+Dies stimmt überein mit den früher gefundenen Ausdrücken für
+die Legendre-Polynome.
+
+Die Potenzreihenlösung zeigt zwar, dass es für jedes $n\in\mathbb{N}$
+eine Polynomlösung $\bar{P}_n(x)$ vom Grad $n$ gibt.
+Dies kann aber nicht erklären, warum die so gefundenen Polynome
+orthogonal sind.
+
+\subsection{Eigenfunktionen}
+Die Differentialgleichung
+\eqref{buch:integral:eqn:legendre-differentialgleichung}
+Kann mit dem Differentialoperator
+\[
+D = \frac{d}{dx}(1-x^2)\frac{d}{dx}
+\]
+als
+\[
+Dy + n(n+1)y = 0
+\]
+geschrieben werden.
+Tatsächlich ist
+\[
+Dy
+=
+\frac{d}{dx} (1-x^2) \frac{d}{dy}
+=
+\frac{d}{dx} (1-x^2)y'
+=
+(1-x^2)y'' -2x y'.
+\]
+Dies bedeutet, dass die Lösungen $\bar{P}_n(x)$ Eigenfunktionen
+des Operators $D$ zum Eigenwert $n(n+1)$ sind:
+\[
+D\bar{P}_n = -n(n+1) \bar{P}_n.
+\]
+
+\subsection{Orthogonalität von $\bar{P}_n$ als Eigenfunktionen}
+Ein Operator $A$ auf Funktionen heisst {\em selbstadjungiert}, wenn
+für zwei beliebige Funktionen $f$ und $g$ gilt
+\[
+\langle Af,g\rangle = \langle f,Ag\rangle
+\]
+gilt.
+Im vorliegenden Zusammenhang möchten wir die Eigenschaft nutzen,
+dass Eigenfunktionen eines selbstadjungierten Operatores zu verschiedenen
+Eigenwerten orthogonal sind.
+Dazu seien $Df = \lambda f$ und $Dg=\mu g$ und wir rechnen
+\begin{equation}
+\renewcommand{\arraycolsep}{2pt}
+\begin{array}{rcccrl}
+\langle Df,g\rangle &=& \langle \lambda f,g\rangle &=& \lambda\phantom{)}\langle f,g\rangle
+&\multirow{2}{*}{\hspace{3pt}$\biggl\}\mathstrut-\mathstrut$}\\
+=\langle f,Dg\rangle &=& \langle f,\mu g\rangle &=& \mu\phantom{)}\langle f,g\rangle&
+\\[2pt]
+\hline
+         0           & &                        &=& (\lambda-\mu)\langle f,g\rangle&
+\end{array}
+\label{buch:integrale:eqn:eigenwertesenkrecht}
+\end{equation}
+Da $\lambda-\mu\ne 0$ ist, muss $\langle f,g\rangle=0$ sein.
+
+Der Operator $D$ ist selbstadjungiert, d.~h.
+für zwei beliebige zweimal stetig differenzierbare Funktion $f$ und $g$
+auf dem Intervall $[-1,1]$ gilt
+\begin{align*}
+\langle Df,g\rangle
+&=
+\int_{-1}^1 (Df)(x) g(x) \,dx
+\\
+&=
+\int_{-1}^1
+\biggl(\frac{d}{dx} (1-x^2)\frac{d}{dx}f(x)\biggr) g(x)
+\,dx
+\\
+&=
+\underbrace{
+\biggl[
+\biggl((1-x^2)\frac{d}{dx}f(x)\biggr) g(x)
+\biggr]_{-1}^1
+}_{\displaystyle = 0}
+-
+\int_{-1}^1
+\biggl((1-x^2)\frac{d}{dx}f(x)\biggr) \frac{d}{dx}g(x)
+\,dx
+\\
+&=
+-
+\int_{-1}^1
+\biggl(\frac{d}{dx}f(x)\biggr) \biggl((1-x^2)\frac{d}{dx}g(x)\biggr)
+\,dx
+\\
+&=
+-
+\underbrace{
+\biggl[
+f(x) \biggl((1-x^2)\frac{d}{dx}g(x)\biggr)
+\biggr]_{-1}^1}_{\displaystyle = 0}
++
+\int_{-1}^1
+f(x) \biggl(\frac{d}{dx}(1-x^2)\frac{d}{dx}g(x)\biggr)
+\,dx
+\\
+&=
+\langle f,Dg\rangle.
+\end{align*}
+Dies beweist, dass $D$ selbstadjungiert ist.
+Da $\bar{P}_n$ Eigenwerte des selbstadjungierten Operators $D$ zu
+den verschiedenen Eigenwerten $-n(n+1)$ sind, folgt auch, dass
+die $\bar{P}_n$ orthogonale Polynome vom Grad $n$ sind, die die 
+gleiche Standardierdisierungsbedingung wie die Legendre-Polyonome
+erfüllen, also ist $\bar{P}_n(x)=P_n(x)$.
+
+\subsection{Legendre-Funktionen zweiter Art}
+%Siehe Wikipedia-Artikel \url{https://de.wikipedia.org/wiki/Legendre-Polynom}
+%
+Die Potenzreihenmethode liefert natürlich auch Lösungen der
+Legendreschen Differentialgleichung, die sich nicht als Polynome
+darstellen lassen.
+Ist $n$ gerade, dann liefern die Anfangswerte $a_0=0$ und $a_1=1$ 
+eine ungerade Funktion, die Folge der Koeffizienten bricht
+aber nicht ab, vielmehr ist
+\begin{align*}
+a_{k+2}
+&=
+\frac{k(k+1)}{(k+1)(k+2)}a_k
+=
+\frac{k}{k+2}a_k.
+\end{align*}
+Durch wiederholte Anwendung dieser Rekursionsformel findet man
+\[
+a_{k}
+=
+\frac{k-2}{k}a_{k-2}
+=
+\frac{k-2}{k}\frac{k-4}{k-2}a_{k-4}
+=
+\frac{k-2}{k}\frac{k-4}{k-2}\frac{k-6}{k-4}a_{k-6}
+=
+\dots
+=
+\frac{1}{k}a_1.
+\]
+Die Lösung hat daher die Reihenentwicklung
+\[
+Q_0(x) = x+\frac13x^3 + \frac15x^5 + \frac17x^7+\dots
+=
+\frac12\log \frac{1+x}{1-x}
+=
+\operatorname{artanh}x.
+\]
+Die Funktion $Q_0(x)$ heisst {\em Legendre-Funktion zweiter Art}.
+
+Für $n=1$ wird die Reihenentwicklung $a_0=1$ und $a_1=0$ etwas
+interessanter.
+Die Rekursionsformel für die Koeffizienten ist
+\[
+a_{k+2}
+=
+\frac{k(k+1)-2}{(k+1)(k+2)} a_k.
+\qquad\text{oder}\qquad
+a_k
+=
+\frac{(k-1)(k-2)-2}{k(k-1)}
+a_{k-2}
+\]
+Man erhält der Reihe nach
+\begin{align*}
+a_2 &= \frac{-2}{2\cdot 1} a_0 = -1
+\\
+a_3 &= 0
+\\
+a_4 &= \frac{3\cdot 2-2}{4\cdot 3} a_2 = \frac{4}{4\cdot 3}a_2 = \frac13a_2 = -\frac13
+\\
+a_5 &= 0
+\\
+a_6 &= \frac{5\cdot 4-2}{6\cdot 5}a_4 = \frac{18}{6\cdot 5}a_4 = -\frac15
+\\
+a_7 &= 0
+\\
+a_8 &= \frac{7\cdot 6-2}{8\cdot 7}a_6 = \frac{40}{8\cdot 7} = -\frac17
+\\
+a_9 &= 0
+\\
+a_{10} &= \frac{9\cdot 8-2}{10\cdot 9}a_8 = \frac{70}{10\cdot 9} = -\frac19,
+\end{align*}
+woraus sich die Reihenentwicklung
+\begin{align*}
+y(x)
+&=
+-x^2 -\frac13x^4 -\frac15x^6 - \frac17x^8 -\frac19x^{10}-\dots
+\\
+&=
+-x\biggl(x+\frac13x^3 + \frac15x^5 + \frac17x^7 + \frac19x^9+\dots\biggr)
+=
+-x\operatorname{artanh}x.
+\end{align*}
+Die {\em Legendre-Funktionen zweiter Art} $Q_n(x)$  werden allerdings
+so definiert, dass gewisse Rekursionsformeln für die Legendre-Polynome,
+die wir hier nicht hergeleitet haben, auch für die $Q_n(x)$ gelten.
+In dieser Normierung muss statt des eben berechneten $y(x)$ die Funktion
+\[
+Q_1(x) = x \operatorname{artanh}x-1
+\]
+verwendet werden.
+
diff --git a/buch/chapters/070-orthogonalitaet/orthogonal.tex b/buch/chapters/070-orthogonalitaet/orthogonal.tex
new file mode 100644
index 0000000..2b7bf41
--- /dev/null
+++ b/buch/chapters/070-orthogonalitaet/orthogonal.tex
@@ -0,0 +1,725 @@
+%
+% orthogonal.tex
+%
+% (c) 2021 Prof Dr Andreas Müller, OST Ostschweizer Fachhochschule
+%
+\section{Orthogonale Funktionenfamilien
+\label{buch:orthogonalitaet:section:orthogonale-funktionen}}
+\rhead{Orthogonale Funktionenfamilien}
+Die Fourier-Theorie basiert auf der Idee, Funktionen durch 
+Funktionenreihen mit Summanden zu bilden, die im Sinne eines
+Skalarproduktes orthogonal sind, welches mit Hilfe eines Integrals
+definiert sind.
+Solche Funktionenfamilien treten jedoch auch als Lösungen von
+Differentialgleichungen.
+Besonders interessant wird die Situation, wenn die Funktionen 
+Polynome sind.
+
+%
+% Skalarprodukt
+%
+\subsection{Skalarprodukt}
+Der reelle Vektorraum $\mathbb{R}^n$ trägt das Skalarprodukt
+\[
+\langle\;\,,\;\rangle
+\colon
+\mathbb{R}^n \times \mathbb{R}^n \to \mathbb{R}
+:
+(x,y)\mapsto \langle x, y\rangle = \sum_{k=1}^n x_iy_k,
+\]
+welches viele interessante Anwendungen ermöglicht.
+Eine orthonormierte Basis macht es zum Beispiel besonders leicht,
+eine Zerlegung eines Vektors in dieser Basis zu finden.
+In diesem Abschnitt soll zunächst an die Eigenschaften erinnert
+werden, die zu einem nützlichen 
+
+\subsubsection{Eigenschaften eines Skalarproduktes}
+Das Skalarprodukt erlaubt auch, die Länge eines Vektors $v$
+als $|v| = \sqrt{\langle v,v\rangle}$ zu definieren.
+Dies funktioniert natürlich nur, wenn die Wurzel auch immer
+definiert ist, d.~h.~das Skalarprodukt eines Vektors mit sich
+selbst darf nicht negativ sein.
+Dazu dient die folgende Definition.
+
+\begin{definition}
+Sei $V$ ein reeller Vektorraum.
+Eine bilineare Abbildung
+\[
+\langle\;\,,\;\rangle
+\colon
+V\times V
+\to
+\mathbb{R}
+:
+(u,v) \mapsto \langle u,v\rangle.
+\]
+heisst {\em positiv definit}, wenn für alle Vektoren $v \in V$ mit
+$v\ne 0 \Rightarrow \langle v,v\rangle > 0$ 
+Die {\em Norm} eines Vektors $v$ ist
+$|v|=\sqrt{\langle v,v\rangle}$.
+\end{definition}
+
+Damit man mit dem Skalarprodukt sinnvoll rechnen kann, ist ausserdem
+erforderlich, dass es eine einfache Beziehung zwischen 
+$\langle x,y\rangle$ und $\langle y,x\rangle$ gibt.
+
+\begin{definition}
+Ein {\em Skalarprodukt} auf einem reellen Vektorraum $V$ ist eine
+positiv definite, symmetrische bilineare Abbildung
+\[
+\langle\;\,,\;\rangle
+\colon
+V\times V
+\to
+\mathbb{R}
+:
+(u,v) \mapsto \langle u,v\rangle.
+\]
+\end{definition}
+
+Das Skalarprodukt $\langle u,v\rangle=u^tv$ auf dem Vektorraum 
+$\mathbb{R}^n$ erfüllt die Definition ganz offensichtlich,
+sie führt auf die Komponentendarstellung
+\[
+\langle u,v\rangle = u^tv = \sum_{k=1}^n u_iv_i.
+\]
+Weitere Skalarprodukte ergeben ergeben sich mit jeder symmetrischen,
+positiv definiten Matrix $G$ und der Definition
+$\langle u,v\rangle_G=u^tGv$.
+Ein einfacher Spezialfall tritt auf, wenn $G$ eine Diagonalmatrix
+$\operatorname{diag}(w_1,\dots,w_n)$
+mit positiven Einträgen $w_i>0$ auf der Diagonalen ist.
+In diesem Fall schreiben wir
+\[
+\langle u,v\rangle_w
+=
+u^t\operatorname{diag}(w_1,\dots,w_n)v
+=
+\sum_{k=1}^n u_iv_i\,w_i
+\]
+und nennen $\langle \;\,,\;\rangle_w$ das {\em gewichtete Skalarprodukt}
+mit {\em Gewichten $w_i$}.
+
+\subsubsection{Skalarprodukte auf Funktionenräumen}
+Das Integral ermöglicht jetzt, ein Skalarprodukt auf dem reellen
+Vektorraum der stetigen Funktionen auf einem Intervall zu definieren.
+
+\begin{definition}
+Sei $V$ der reelle Vektorraum $C([a,b])$ der reellwertigen, stetigen
+Funktion auf dem Intervall $[a,b]$.
+Dann ist 
+\[
+\langle\;\,,\;\rangle
+\colon
+C([a,b]) \times C([a,b]) \to \mathbb{R}
+:
+(f,g) \mapsto \langle f,g\rangle = \int_a^b f(x)g(x)\,dx.
+\]
+ein Skalarprodukt.
+\end{definition}
+
+Die Definition ist offensichtlich symmetrisch in $f$ und $g$ und
+aus den Eigenschaften des Integrals ist klar, dass das Produkt
+bilinear ist:
+\begin{align*}
+\langle \lambda_1 f_1+\lambda_2f_2,g\rangle
+&=
+\int_a^b (\lambda_1f_(x) +\lambda_2f_2(x))g(x)\,dx
+=
+\lambda_1\int_a^b f_1(x) g(x)\,dx
++
+\lambda_2\int_a^b f_2(x) g(x)\,dx
+\\
+&=
+\lambda_1\langle f_1,g\rangle
++
+\lambda_2\langle f_2,g\rangle.
+\end{align*}
+Ausserdem ist es positiv definit, denn wenn $f(x_0) \ne 0$ ist,
+dann gibt es wegen der Stetigkeit von $f$ eine Umgebung
+$U=[x_0-\varepsilon,x_0+\varepsilon]$, derart, dass $|f(x)| > \frac12|f(x_0)|$
+ist für alle $x\in U$.
+Somit ist das Integral
+\[
+\langle f,f\rangle
+=
+\int_a^b |f(x)|^2\,dx
+\ge
+\int_{x_0-\varepsilon}^{x_0+\varepsilon} |f(x)|^2\,dx
+\ge
+\int_{x_0-\varepsilon}^{x_0+\varepsilon} \frac14|f(x_0)|^2\,dx
+=
+\frac{1}{4}|f(x_0)|^2\cdot 2\varepsilon
+=
+\frac{|f(x_0)|^2\varepsilon}{2}
+>0,
+\]
+was beweist, dass $\langle\;,\;\rangle$ positiv definit und damit
+ein Skalarprodukt ist.
+
+Die Definition kann noch etwas verallgemeinert werden, indem 
+die Funktionswerte nicht überall auf dem Definitionsbereich 
+gleich gewichtet werden. 
+
+\begin{definition}
+Sei $w\colon [a,b]\to \mathbb{R}^+$ eine positive, stetige Funktion,
+dann ist
+\[
+\langle\;\,,\;\rangle_w
+\colon
+C([a,b]) \times C([a,b]) \to \mathbb{R}
+:
+(f,g) \mapsto \langle f,g\rangle_w = \int_a^b f(x)g(x)\,w(x)\,dx.
+\]
+das {\em gewichtete Skalarprodukt} mit {\em Gewichtsfunktion $w(x)$}.
+\end{definition}
+
+\subsubsection{Gram-Schmidt-Orthonormalisierung}
+In einem reellen Vektorraum $V$ mit Skalarprodukt $\langle\;\,,\;\rangle$
+kann aus einer beleibigen Basis $b_1,\dots,b_n$ mit Hilfe des 
+Gram-Schmidtschen Orthogonalisierungsverfahrens immer eine
+orthonormierte Basis $\tilde{b}_1,\dots,\tilde{b}_n$ Basis
+gewonnen werden.
+Es stellt sicher, dass für alle $k\le n$ gilt
+\[
+\langle b_1,\dots,b_k\rangle
+=
+\langle \tilde{b}_1,\dots,\tilde{b}_k\rangle.
+\]
+Zur Vereinfachung der Formeln schreiben wir $v^0=v/|v|$ für einen zu
+$v$ parallelen Einheitsvektor.
+Die Vektoren $\tilde{b}_i$ können mit Hilfe der Formeln
+\begin{align*}
+\tilde{b}_1
+&=
+(b_1)^0
+\\
+\tilde{b}_2
+&=
+\bigl(
+b_2
+-
+\langle \tilde{b}_1,b_2\rangle \tilde{b}_1
+\bigr)^0
+\\
+\tilde{b}_3
+&=
+\bigl(
+b_3
+-
+\langle \tilde{b}_1,b_3\rangle \tilde{b}_1
+-
+\langle \tilde{b}_2,b_3\rangle \tilde{b}_2
+\bigr)^0
+\\
+&\;\vdots
+\\
+\tilde{b}_n
+&=
+\bigl(
+b_n
+-
+\langle \tilde{b}_1,b_n\rangle \tilde{b}_1
+-
+\langle \tilde{b}_2,b_n\rangle \tilde{b}_2
+-\dots
+-
+\langle \tilde{b}_{n-1},b_n\rangle \tilde{b}_{n-1}
+\bigr)^0
+\end{align*}
+iterativ berechnet werden.
+Dieses Verfahren lässt sich auch auf Funktionenräume anwenden.
+
+Die Normierung ist nicht unbedingt nötig und manchmal unangenehm,
+da die Norm unschöne Quadratwurzeln einführt.
+Falls es genügt, eine orthogonale Basis zu finden, kann darauf
+verzichtet werden, bei der Orthogonalisierung muss aber berücksichtigt
+werden, dass die Vektoren $\tilde{b}_i$ jetzt nicht mehr Einheitslänge
+haben.
+Die Formeln
+\begin{align*}
+\tilde{b}_0
+&=
+b_0
+\\
+\tilde{b}_1
+&=
+b_1
+-
+\frac{\langle b_1,\tilde{b}_0\rangle}{\langle \tilde{b}_0,\tilde{b}_0\rangle}\tilde{b}_0
+\\
+\tilde{b}_2
+&=
+b_2
+-
+\frac{\langle b_2,\tilde{b}_0\rangle}{\langle \tilde{b}_0,\tilde{b}_0\rangle}\tilde{b}_0
+-
+\frac{\langle b_2,\tilde{b}_1\rangle}{\langle \tilde{b}_1,\tilde{b}_1\rangle}\tilde{b}_1
+\\
+&\;\vdots
+\\
+\tilde{b}_n
+&=
+b_n
+-
+\frac{\langle b_n,\tilde{b}_0\rangle}{\langle \tilde{b}_0,\tilde{b}_0\rangle}\tilde{b}_0
+-
+\frac{\langle b_n,\tilde{b}_1\rangle}{\langle \tilde{b}_1,\tilde{b}_1\rangle}\tilde{b}_1
+-
+\dots
+-
+\frac{\langle b_n,\tilde{b}_{n-1}\rangle}{\langle \tilde{b}_{n-1},\tilde{b}_{n-1}\rangle}\tilde{b}_{n-1}.
+\end{align*}
+berücksichtigen dies.
+
+\subsubsection{Selbstadjungierte Operatoren und Eigenvektoren}
+Symmetrische Matrizen spielen eine spezielle Rolle in der
+endlichdimensionalen linearen Algebra, weil sie sich immer
+mit einer orthonormierten Basis diagonalisieren lassen.
+In der vorliegenden Situation undendlichdimensionaler Vektorräume
+brauchen wir eine angepasste Definition.
+
+\begin{definition}
+Eine lineare Selbstabbildung $A\colon V\to V$
+eines Vektorrraums mit Skalarprodukt
+heisst {\em selbstadjungiert}, wenn für alle Vektoren $u,v\in V$
+heisst $\langle Au,v\rangle = \langle u,Av\rangle$.
+\end{definition}
+
+Es ist wohlbekannt, dass Eigenvektoren einer symmetrischen Matrix
+zu verschiedenen Eigenwerten orthogonal sind.
+Der Beweis ist direkt übertragbar, wir halten das Resultat hier
+für spätere Verwendung fest.
+
+\begin{satz}
+Sind $f$ und $g$ Eigenvektoren eines selbstadjungierten Operators $A$
+zu verschiedenen Eigenwerten $\lambda$ und $\mu$, dann sind $f$ und $g$
+orthogonal.
+\end{satz}
+
+\begin{proof}[Beweis]
+Im vorliegenden Zusammenhang möchten wir die Eigenschaft nutzen,
+dass Eigenfunktionen eines selbstadjungierten Operatores zu verschiedenen
+Eigenwerten orthogonal sind.
+Dazu seien $Df = \lambda f$ und $Dg=\mu g$ und wir rechnen
+\begin{equation*}
+\renewcommand{\arraycolsep}{2pt}
+\begin{array}{rcccrl}
+\langle Df,g\rangle &=& \langle \lambda f,g\rangle &=& \lambda\phantom{)}\langle f,g\rangle
+&\multirow{2}{*}{\hspace{3pt}$\biggl\}\mathstrut-\mathstrut$}\\
+=\langle f,Dg\rangle &=& \langle f,\mu g\rangle &=& \mu\phantom{)}\langle f,g\rangle&
+\\[2pt]
+\hline
+         0           & &                        &=& (\lambda-\mu)\langle f,g\rangle&
+\end{array}
+\end{equation*}
+Da $\lambda-\mu\ne 0$ ist, muss $\langle f,g\rangle=0$ sein.
+\end{proof}
+
+\begin{beispiel}
+Sei $C^1([0,2\pi], \mathbb{C})=C^1(S^1,\mathbb{C})$
+der Vektorraum der $2\pi$-periodischen differenzierbaren Funktionen mit
+dem Skalarprodukt 
+\[
+\langle f,g\rangle
+=
+\frac{1}{2\pi}\int_0^{2\pi} \overline{f(t)}g(t)\,dt
+\]
+enthält die Funktionen $e_n(t) = e^{int}$.
+Der Operator
+\[
+D=i\frac{d}{dt}
+\]
+ist selbstadjungiert, denn mit Hilfe von partieller Integration erhält man
+\[
+\langle Df,g\rangle
+=
+\frac{1}{2\pi}
+\int_0^{2\pi}
+\underbrace{
+\overline{i\frac{df(t)}{dt}}
+}_{\uparrow}
+\underbrace{g(t)}_{\downarrow}
+\,dt
+=
+\underbrace{
+\frac{-i}{2\pi}
+\biggl[
+\overline{f(t)}g(t)
+\biggr]_0^{2\pi}
+}_{\displaystyle=0}
++
+\frac{1}{2\pi}
+\int_0^{2\pi}
+\overline{f(t)}i\frac{dg(t)}{dt}
+\,dt
+=
+\langle f,Dg\rangle
+\]
+unter Ausnützung der $2\pi$-Periodizität der Funktionen.
+
+Die Funktionen $e_n(t)$ sind Eigenfunktionen des Operators $D$, denn
+\[
+De_n(t) = i\frac{d}{dt}e^{int} = -n e^{int} = -n e_n(t).
+\]
+Nach obigem Satz sind die Eigenfunktionen von $D$ orthogonal.
+\end{beispiel}
+
+Das Beispiel illustriert, dass orthogonale Funktionenfamilien
+ein automatisches Nebenprodukt selbstadjungierter Operatoren sind.
+
+%%
+%% Besselfunktionen also orthogonale Funktionenfamilie
+%%
+%\subsection{Bessel-Funktionen als orthogonale Funktionenfamilie}
+%Auch die Besselfunktionen sind eine orthogonale Funktionenfamilie.
+%Sie sind Funktionen differenzierbaren Funktionen $f(r)$ für $r>0$
+%mit $f'(r)=0$ und für $r\to\infty$ nimmt $f(r)$ so schnell ab, dass
+%auch $rf(r)$ noch gegen $0$ strebt.
+%Das Skalarprodukt ist
+%\[
+%\langle f,g\rangle
+%=
+%\int_0^\infty r f(r) g(r)\,dr,
+%\]
+%als Operator verwenden wir
+%\[
+%A = \frac{d^2}{dr^2} + \frac{1}{r}\frac{d}{dr} + s(r),
+%\]
+%wobei $s(r)$ eine beliebige integrierbare Funktion sein kann.
+%Zunächst überprüfen wir, ob dieser Operator wirklich selbstadjungiert ist.
+%Dazu rechnen wir
+%\begin{align}
+%\langle Af,g\rangle
+%&=
+%\int_0^\infty
+%r\,\biggl(f''(r)+\frac1rf'(r)+s(r)f(r)\biggr) g(r)
+%\,dr
+%\notag
+%\\
+%&=
+%\int_0^\infty rf''(r)g(r)\,dr
+%+
+%\int_0^\infty f'(r)g(r)\,dr
+%+
+%\int_0^\infty s(r)f(r)g(r)\,dr.
+%\notag
+%\intertext{Der letzte Term ist symmetrisch in $f$ und $g$, daher
+%ändern wir daran weiter nichts.
+%Auf das erste Integral kann man partielle Integration anwenden und erhält}
+%&=
+%\biggl[rf'(r)g(r)\biggr]_0^\infty
+%-
+%\int_0^\infty f'(r)g(r) + rf'(r)g'(r)\,dr
+%+
+%\int_0^\infty f'(r)g(r)\,dr
+%+
+%\int_0^\infty s(r)f(r)g(r)\,dr.
+%\notag
+%\intertext{Der erste Term verschwindet wegen der Bedingungen an die
+%Funktionen $f$ und $g$.
+%Der erste Term im zweiten Integral hebt sich gegen das
+%zweite Integral weg.
+%Der letzte Term ist das Skalarprodukt von $f'$ und $g'$.
+%Somit ergibt sich
+%}
+%&=
+%-\langle f',g'\rangle
+%+
+%\int_0^\infty s(r) f(r)g(r)\,dr.
+%\label{buch:integrale:orthogonal:besselsa}
+%\end{align}
+%Vertauscht man die Rollen von $f$ und $g$, erhält man das Gleiche, da im
+%letzten Ausdruck~\eqref{buch:integrale:orthogonal:besselsa} die Funktionen
+%$f$ und $g$ symmetrische auftreten.
+%Damit ist gezeigt, dass der Operator $A$ selbstadjungiert ist.
+%Es folgt nun, dass Eigenvektoren des Operators $A$ automatisch
+%orthogonal sind.
+%
+%Eigenfunktionen von $A$ sind aber Lösungen der Differentialgleichung
+%\[
+%\begin{aligned}
+%&&
+%Af&=\lambda f
+%\\
+%&\Rightarrow\qquad&
+%f''(r) +\frac1rf'(r) + s(r)f(r) &= \lambda f(r)
+%\\
+%&\Rightarrow\qquad&
+%r^2f''(r) +rf'(r)+ (-\lambda r^2+s(r)r^2)f(r) &= 0
+%\end{aligned}
+%\]
+%sind.
+%
+%Durch die Wahl $s(r)=1$ wird der Operator $A$ zum Bessel-Operator
+%$B$ definiert in
+%\eqref{buch:differentialgleichungen:bessel-operator}.
+%Die Lösungen der Besselschen Differentialgleichung zu verschiedenen Werten
+%des Parameters müssen also orthogonal sein, insbesondere sind die
+%Besselfunktion $J_\nu(r)$ und $J_\mu(r)$ orthogonal wenn $\mu\ne\nu$ ist.
+%
+%
+% Orthogonale Polynome
+%
+\subsection{Orthogonale Polynome
+\label{buch:integral:subsection:orthogonale-polynome}}
+Die Polynome $1,x,x^2,\dots,x^n$ bilden eine Basis des Vektorraums
+der Polynome vom Grad $\le n$.
+Bezüglich des Skalarproduktes
+\[
+\langle p,q\rangle
+=
+\int_{-1}^1 p(x)q(x)\,dx
+\]
+sind sie jedoch nicht orthogonal, denn es ist
+\[
+\langle x^i,x^j\rangle
+=
+\int_{-1}^1 x^{i+j}\,dx
+=
+\biggl[\frac{x^{i+j+1}}{i+j+1}\biggr]_{-1}^1
+=
+\begin{cases}
+\displaystyle
+\frac{2}{i+j+1}&\qquad\text{$i+j$ gerade}\\
+              0&\qquad\text{$i+j$ ungerade}.
+\end{cases}
+\]
+Wir können daher das Gram-Schmidtsche Orthonormalisierungsverfahren
+anwenden, um eine orthogonale Basis von Polynomen zu finden, was
+wir im Folgenden tun wollen.
+
+% XXX Orthogonalisierungsproblem so formulieren, dass klar wird,
+% XXX dass man ein "Normierungskriterium braucht.
+
+Da wir auf die Normierung verzichten, brauchen wir ein anderes
+Kriterium, welches die Polynome eindeutig festlegen kann.
+Wir bezeichnen das Polynom vom Grad $n$, das bei diesem Prozess
+entsteht, mit $P_n(x)$ und legen willkürlich aber traditionskonform
+fest, dass $P_n(1)=1$ sein soll.
+
+Das Skalarprodukt berechnet ein Integral eines Produktes von zwei
+Polynomen über das symmetrische Interval $[-1,1]$.
+Ist die eine gerade und die andere ungerade, dann ist das
+Produkt eine ungerade Funktion und das Skalarprodukt verschwindet.
+Sind beide Funktionen gerade oder ungerade, dann ist das Produkt
+gerade und das Skalarprodukt ist im Allgmeinen von $0$ verschieden.
+Dies zeigt, dass es tatsächlich etwas zu Orthogonalisieren gibt.
+
+Die ersten beiden Funktionen sind das konstante Polynom $1$ und
+das Polynome $x$.
+Nach obiger Beobachtung ist das Skalarprodukt $\langle 1,x\rangle=0$,
+also ist $P_1(x)=x$.
+Die Graphen der entstehenden Polynome sind in
+Abbildung~\ref{buch:integral:orthogonal:legendregraphen}
+dargestellt.
+\begin{figure}
+\centering
+\includegraphics{chapters/070-orthogonalitaet/images/legendre.pdf}
+\caption{Graphen der Legendre-Polynome $P_n(x)$ für $n=1,\dots,10$.
+\label{buch:integral:orthogonal:legendregraphen}}
+\end{figure}
+
+\begin{lemma}
+Die Polynome $P_{2n}(x)$ sind gerade, die Polynome $P_{2n+1}(x)$ sind
+ungerade Funktionen von $x$.
+\end{lemma}
+
+\begin{proof}[Beweis]
+Wir verwenden vollständige Induktion nach $n$.
+Wir wissen bereits, dass $P_0(x)=1$ und $P_1(x)=x$ die verlangten
+Symmetrieeigenschaften haben.
+Im Sinne der Induktionsannahme nehmen wir daher an, dass die
+Symmetrieeigenschaften für $P_k(x)$, $k<n$, bereits bewiesen sind.
+$P_n(x)$ entsteht jetzt durch Orthogonalisierung nach der Formel
+\[
+P_n(x)
+=
+x^n
+-
+\langle P_{n-1},x^n\rangle P_{n-1}(x)
+-
+\langle P_{n-2},x^n\rangle P_{n-2}(x)
+-\dots-
+\langle P_1,x^n\rangle P_1(x)
+-
+\langle P_0,x^n\rangle P_0(x).
+\]
+Die Skalarprodukte
+$\langle P_{n-1},x^n\rangle$,
+$\langle P_{n-3},x^n\rangle$, $\dots$ verschwinden alle, so dass
+$P_n(x)$ eine Linearkombination der Funktionen $x^n$, $P_{n-2}(x)$,
+$P_{n-4}(x)$ ist, die die gleiche Parität wie $x^n$ haben.
+Also hat auch $P_n(x)$ die gleiche Parität, was das Lemma beweist.
+\end{proof}
+
+Die Ortogonalisierung von $x^2$ liefert daher
+\[
+p(x) = x^2
+-
+\frac{\langle x^2,P_0\rangle}{\langle P_0,P_0\rangle} P_0(x)
+=
+x^2 - \frac{\int_{-1}^1x^2\,dx}{\int_{-1}^11\,dx}
+=
+x^2 - \frac{\frac{2}{3}}{2}=x^2-\frac13
+\]
+Dieses Polynom erfüllt die Standardisierungsbedingung noch 
+nicht den $p(1)=\frac23$.
+Daraus leiten wir ab, dass
+\[
+P_2(x) = \frac12(3x^2-1)
+\]
+ist.
+
+Für $P_3(x)$ brauchen wir nur die Skalaprodukte
+\[
+\left.
+\begin{aligned}
+\langle x^3,P_1\rangle
+&=
+\int_{-1}^1  x^3\cdot x\,dx
+=
+\biggl[\frac15x^5\biggr]_{-1}^1
+=
+\frac25
+\qquad
+\\
+\langle P_1,P_1\rangle
+&=
+\int_{-1}^1 x^2\,dx
+=
+\frac23
+\end{aligned}
+\right\}
+\qquad
+\Rightarrow
+\qquad
+p(x) = x^3 - \frac{\frac25}{\frac23}x=x^3-\frac{3}{5}x
+\]
+Die richtige Standardisierung ergibt sich,
+indem man durch $p(1)=\frac25$ dividiert, also
+\[
+P_2(x) = \frac12(5x^3-3x).
+\]
+
+Die Berechnung weiterer Polynome verlangt, dass Skalarprodukte
+$\langle x^n,P_k\rangle$ berechnet werden müssen, was wegen
+der zunehmend komplizierten Form von $P_k$ etwas mühsam ist.
+Wir berechnen den Fall $P_4$.
+Dazu muss das Polynom $x^4$ um eine Linearkombination von
+$P_2$ und $P_0(x)=1$ korrigiert werden.
+Die Skalarprodukte sind
+\begin{align*}
+\langle x^4, P_0\rangle
+&=
+\int_{-1}^1 x^4\,dx = \frac25
+\\
+\langle P_0,P_0\rangle
+&=
+\int_{-1}^1 \,dx = 2
+\\
+\langle x^4,P_2\rangle
+&=
+\int_{-1}^1 \frac32x^6-\frac12 x^4\,dx
+=
+\biggl[\frac{3}{14}x^7-\frac{1}{10}x^5\biggr]_{-1}^1
+=
+\frac6{14}-\frac15
+=
+\frac8{35}
+\\
+\langle P_2,P_2\rangle
+&=
+\int_{-1}^1 \frac14(3x^2-1)^2\,dx
+=
+\int_{-1}^1 \frac14(9x^4-6x^2+1)\,dx
+=
+\frac14(\frac{18}{5}-4+2)
+=\frac25.
+\end{align*}
+Daraus folgt für $p(x)$
+\begin{align*}
+p(x)
+&=
+x^4
+-
+\frac{\langle x^4,P_2\rangle}{\langle P_2,P_2\rangle}P_2(x)
+-
+\frac{\langle x^4,P_0\rangle}{\langle P_0,P_0\rangle}P_0(x)
+\\
+&=
+x^4
+-\frac47 P_2(x) - \frac15 P_0(x)
+\\
+&=
+x^4 - \frac{6}{7}x^2 + \frac{3}{35}
+\end{align*}
+mit $p(1)=\frac{8}{35}$, so dass man
+\[
+P_4(x) =
+\frac18(35x^4-30x^2+3)
+\]
+setzen muss.
+
+\begin{figure}
+\centering
+\includegraphics{chapters/070-orthogonalitaet/images/orthogonal.pdf}
+\caption{Orthogonalität der Legendre-Polynome $P_4(x)$ ({\color{blue}blau})
+und $P_7(x)$ ({\color{darkgreen}grün}).
+Die blaue Fläche ist die Fläche unter dem Graphen 
+von $P_4(x)^2$, $P_4(x)$ muss durch die Wurzel aus diesem Flächeninhalt
+geteilt werden, um ein Polynome mit Norm $1$ zu erhalten.
+Für die grüne Fläche ist es $P_7(x)$.
+Die rote Kurve ist der Graph der Funktion $P_4(x)\cdot P_7(x)$,
+die rote Fläche ist deren Integral, sie ist $0$, d.~h.~die beiden
+Funktionen sind orthogonal.
+\label{buch:integral:orthogonal:legendreortho}}
+\end{figure}
+
+\begin{table}
+\centering
+\renewcommand{\arraystretch}{1.2}
+\begin{tabular}{|>{$}c<{$}|>{$}l<{$}|}
+\hline
+n&P_n(x)\\
+\hline
+ 0&1
+\\
+ 1&x
+\\
+ 2&\frac12(3x^2-1)
+\\
+ 3&\frac12(5x^3-3x)
+\\
+ 4&\frac18(35x^4-30x^2+3)
+\\
+ 5&\frac18(63x^5-70x^3+15x)
+\\
+ 6&\frac1{16}(231x^6-315x^4+105x^2-5)
+\\
+ 7&\frac1{16}(429x^7-693x^5+315x^3-35x)
+\\
+ 8&\frac1{128}(6435x^8-12012x^6+6930x^4-1260x^2+35)
+\\
+ 9&\frac1{128}(12155x^9-25740x^7+18018x^5-4620x^3+315x)
+\\
+10&\frac1{256}(46189x^{10}-109395x^8+90090x^6-30030x^4+3465x^2-63)
+\\[2pt]
+\hline
+\end{tabular}
+\caption{Die Legendre-Polynome $P_n(x)$ für $n=0,1,\dots,10$ sind
+orthogonale Polynome vom Grad $n$, die den Wert $P_n(1)=1$ haben.
+\label{buch:integral:table:legendre-polynome}}
+\end{table}
+
+
+
+Die so konstruierten Polynome heissen die {\em Legendre-Polynome}.
+Durch weitere Durchführung des Verfahrens liefert die Polynome in
+Tabelle~\ref{buch:integral:table:legendre-polynome}.
+Die Graphen sind in Abbildung~\ref{buch:integral:orthogonal:legendregraphen}
+dargestellt.
+Abbildung~\ref{buch:integral:orthogonal:legendreortho} illustriert, 
+dass die die beiden Polynome $P_4(x)$ und $P_7(x)$ orthogonal sind.
+Das Produkt $P_4(x)\cdot P_7(x)$ hat Integral $=0$.
+
diff --git a/buch/chapters/070-orthogonalitaet/sturm.tex b/buch/chapters/070-orthogonalitaet/sturm.tex
new file mode 100644
index 0000000..c8ee11a
--- /dev/null
+++ b/buch/chapters/070-orthogonalitaet/sturm.tex
@@ -0,0 +1,609 @@
+%
+% sturm.tex
+%
+% (c) 2022 Prof Dr Andreas Müller, OST Ostschweizer Fachhochschule
+%
+\section{Das Sturm-Liouville-Problem
+\label{buch:integrale:subsection:sturm-liouville-problem}}
+\rhead{Das Sturm-Liouville-Problem}
+Sowohl bei den Bessel-Funktionen wie bei den Legendre-Polynomen
+konnte die Orthogonalität der Funktionen dadurch gezeigt werden,
+dass sie als Eigenfunktionen eines bezüglich eines geeigneten
+Skalarproduktes selbstadjungierten Operators erkannt wurden.
+
+\subsection{Differentialgleichung}
+Das klassische Sturm-Liouville-Problem ist das folgende Eigenwertproblem.
+Gesucht sind Lösungen der Differentialgleichung
+\begin{equation}
+((p(x)y'(x))' + q(x)y(x) = \lambda w(x) y(x)
+\label{buch:integrale:eqn:sturm-liouville}
+\end{equation}
+auf dem Intervall $(a,b)$, die zusätzlich die Randbedingungen
+\begin{equation}
+\begin{aligned}
+k_a y(a) + h_a p(a) y'(a) &= 0 \\
+k_b y(b) + h_b p(b) y'(b) &= 0
+\end{aligned}
+\label{buch:integrale:sturm:randbedingung}
+\end{equation}
+erfüllen, wobei $|k_i|^2 + |h_i|^2\ne 0$ mit $i=a,b$.
+Weitere Bedingungen an die Funktionen $p(x)$, $q(x)$, $w(x)$  sowie die
+Lösungsfunktionen $y(x)$ sollen später geklärt werden.
+
+\subsection{Das verallgemeinerte Eigenwertproblem für symmetrische Matrizen}
+Ein zu \eqref{buch:integrale:eqn:sturm-liouville} analoges Eigenwertproblem
+für Matrizen ist das folgende verallgemeinerte Eigenwertproblem.
+Das gewohnte Eigenwertproblem verwendet die Matrix $B=E$.
+
+\begin{definition}
+Seien $A$ und $B$ $n\times n$-Matrizen.
+$v$ heisst verallgemeinerter Eigenvektor zum Eigenwert $\lambda$,
+wenn
+\[
+Av = \lambda Bv.
+\]
+\end{definition}
+
+Für symmetrische Matrizen lässt sich dieses Problem auf ein 
+Optimierungsproblem reduzieren.
+
+\begin{satz}
+Seien $A$ und $B$ symmetrische $n\times n$-Matrizen und sei ausserdem
+$B$ positiv definit.
+Ist $v$ ein Vektor, der die Grösse
+\[
+f(v)=\frac{v^tAv}{v^tBv}
+\]
+maximiert, ist ein verallgemeinerter Eigenvektor für die Matrizen $A$
+und $B$.
+\end{satz}
+
+\begin{proof}[Beweis]
+Sei $\lambda = f(v)$ der maximale Wert und $u\ne 0$ ein beliebiger Vektor. 
+Da $v$ die Grösse $f(v)$ maximiert, muss die Ableitung
+von $f(u+tv)$ für $t=0$ verschwinden.
+Um diese Ableitung zu berechnen, bestimmen wir zunächst die Ableitung
+von $(v+tu)^tM(v+tu)$  an der Stelle $t=0$ für eine beliebige
+symmetrische Matrix:
+\begin{align*}
+\frac{d}{dt}
+(v+tu)^tM(v+tu)
+&=
+\frac{d}{dt}\bigl(
+v^tv + t(v^tMu+u^tMv) + t^2 u^tMu
+\bigr)
+=
+v^tMu+u^tMv + 2tv^tMv
+\\
+\frac{d}{dt}
+(v^t+tu^t)M(v+tu)
+\bigg|_{t=0}
+&=
+v^tMu+u^tMv
+=
+2v^tMu
+\end{align*}
+Dies wenden wir jetzt auf den Quotenten $\lambda(v+tu)$ an.
+\begin{align*}
+\frac{d}{dt}f(v+tu)\bigg|_{t=0}
+&=
+\frac{d}{dt}
+\frac{(v+tu)^tA(v+tu)}{(v+tu)^tB(v+tu)}\bigg|_{t=0}
+\\
+&=
+\frac{2u^tAv(v^tBv) - 2u^tBv(v^tAv)}{(v^tBv)^2}
+=
+\frac{2}{v^tBv}
+u^t
+\biggl(
+Av - \frac{v^tAv}{v^tBv} Bv
+\biggr)
+\\
+&=
+2
+\frac{
+u^t(
+Av - \lambda Bv
+)
+}{v^tBv}
+\end{align*}
+Da $v$ ein Maximum von $\lambda(v)$ ist, verschwindet diese Ableitung
+für alle Vektoren $u$, somit gilt
+\[
+u^t(Av-\lambda Bv)=0
+\]
+für alle $u$, also auch $Av=\lambda Bv$.
+Dies beweist, dass $v$ ein verallgemeinerter Eigenvektor zum
+Eigenwert $\lambda$ ist.
+\end{proof}
+
+\begin{satz}
+Verallgemeinerte Eigenvektoren $u$ und $v$ von $A$ und $B$
+zu verschiedenen Eigenwerten erfüllen $u^tBv=0$.
+\end{satz}
+
+\begin{proof}
+Seien $\lambda$ und $\mu$ die Eigenwerte, also $Au=\lambda Bu$
+und $Av=\mu Bv$.
+Wie in \eqref{buch:integrale:eqn:eigenwertesenkrecht}
+berechnen wir das Skalarprodukt auf zwei Arten
+\[
+\renewcommand{\arraycolsep}{2pt}
+\begin{array}{rcccrl}
+ u^tAv &=&u^t\lambda Bv &=& \lambda\phantom{\mathstrut-\mu)} u^tBv
+	&\multirow{2}{*}{\hspace{3pt}$\bigg\}\mathstrut-\mathstrut$}\\
+=v^tAu &=&v^t\mu Bu     &=&  \mu\phantom{)}u^tBv         &\\
+\hline
+     0 & &              &=& (\lambda - \mu)u^tBv.        &
+\end{array}
+\]
+Da die Eigenwerte verschieden sind, ist $\lambda-\mu\ne 0$, es folgt, 
+dass $u^tBv=0$ sein muss.
+\end{proof}
+
+Verallgemeinerte Eigenwerte und Eigenvektoren verhalten sich also
+ganz analog zu den gewöhnlichen Eigenwerten und Eigenvektoren.
+Da $B$ positiv definit ist, ist $B$ auch invertierbar.
+Zudem kann $B$ zur Definition des verallgemeinerten Skalarproduktes
+\[
+\langle u,v\rangle_B = u^tBv
+\]
+verwendet werden.
+Die Matrix 
+\[
+\tilde{A} = B^{-1}A
+\]
+ist bezüglich dieses Skalarproduktes selbstadjungiert, denn es gilt
+\[
+\langle\tilde{A}u,v\rangle_B
+=
+(B^{-1}Au)^t Bv
+=
+u^tA^t(B^{-1})^tBv
+=
+u^tAv
+=
+u^tBB^{-1}Av
+=
+\langle u,\tilde{A}v\rangle.
+\]
+Das verallgemeinerte Eigenwertproblem für symmetrische Matrizen
+ist damit ein gewöhnliches Eigenwertproblem für selbstadjungierte
+Matrizen des Operators $\tilde{A}$ bezüglich des verallgemeinerten
+Skalarproduktes $\langle\,\;,\;\rangle_B$.
+
+\subsection{Der Operator $L_0$ und die Randbedingung}
+Die Differentialgleichung kann auch in Operatorform geschrieben werden.
+Dazu schreiben wir
+\[
+L_0 
+=
+-\frac{d}{dx}p(x)\frac{d}{dx}.
+\]
+Bezüglich des gewöhnlichen Skalarproduktes
+\[
+\langle f,g\rangle
+=
+\int_a^b f(x)g(x)\,dx
+\]
+für Funktionen auf dem Intervall $[a,b]$ ist der Operator $L_0$
+tatsächlich selbstadjungiert.
+Mit partieller Integration rechnet man nach:
+\begin{align}
+\langle f,L_0g\rangle
+&=
+\int_a^b f(x) \biggl(-\frac{d}{dx}p(x)\frac{d}{dx}\biggr)g(x)\,dx
+\notag
+\\
+&=
+-\int_a^b f(x) \frac{d}{dx}\bigl( p(x) g'(x) \bigr)\,dx
+\notag
+\\
+&=
+-\biggl[f(x) p(x)g'(x)\biggr]_a^b
++
+\int_a^b f'(x) p(x) g'(x) \,dx
+\notag
+\\
+\langle L_0f,g\rangle
+&=
+-\biggl[f'(x)p(x)g(x)\biggr]_a^b
++
+\int_a^b f'(x) p(x) g'(x) \,dx.
+\notag
+\intertext{Die beiden Skalarprodukte führen also auf das gleiche
+Integral, sie unterscheiden sich nur um die Randterme}
+\langle f,L_0g\rangle
+-
+\langle L_0f,g\rangle
+&=
+-f(b)p(b)g'(b) + f(a)p(a)g'(a)
++f'(b)p(b)g(b) - f'(a)p(a)g(a)
+\label{buch:integrale:sturm:sabedingung}
+\\
+&=
+-
+p(b)
+\left|\begin{matrix}
+f(b) &g(b)\\
+f'(b)&g'(b)
+\end{matrix}\right|
++
+p(a)
+\left|\begin{matrix}
+f(a) &g(a)\\
+f'(a)&g'(a)
+\end{matrix}\right|
+\notag
+\\
+&=
+-
+\left|\begin{matrix}
+f(b) &g(b)\\
+p(b)f'(b)&p(b)g'(b)
+\end{matrix}\right|
++
+\left|\begin{matrix}
+f(a) &g(a)\\
+p(a)f'(a)&p(a)g'(a)
+\end{matrix}\right|.
+\notag
+\end{align}
+Um zu erreichen, dass der Operator selbstadjunigert wird, muss 
+sichergestellt werden, dass entweder $p$ oder die Determinanten
+an den Intervallenden verschwinden.
+Dies passiert genau dann, wenn die Vektoren 
+\[
+\begin{pmatrix}
+f(a)\\
+p(a)f'(a)
+\end{pmatrix}
+\text{\;und\;}
+\begin{pmatrix}
+g(a)\\
+p(a)g'(a)
+\end{pmatrix}
+\]
+linear abhängig sind.
+In zwei Dimensionen bedeutet lineare Abhängigkeit, dass es
+eine nichttriviale Linearform mit Koeffizienten $h_a, k_a$ gibt,
+die auf beiden Vektoren verschwindet.
+Ausgeschrieben bedeutet dies, dass die Randbedingung
+\eqref{buch:integrale:sturm:randbedingung}
+erfüllt sein muss.
+
+\subsection{Skalarprodukt}
+Das Ziel der folgenden Abschnitte ist, das Sturm-Liouville-Problem als
+Eigenwertproblem für einen selbstadjungierten Operator in einem 
+Funktionenraum mit einem geeigneten Skalarprodukt zu finden.
+
+Wir haben bereits gezeigt, dass die Randbedingung
+\eqref{buch:integrale:sturm:randbedingung} sicherstellt, dass der
+Operator $L_0$ für das Standardskalarprodukt selbstadjungiert ist.
+Dies entspricht der Symmetrie der Matrix $A$.
+
+Die Komponente $q(x)$ stellt keine besonderen Probleme, denn
+\[
+\langle f,qg\rangle
+=
+\int_a^b f(x)q(x)g(x)\,dx
+=
+\langle qf,g\rangle.
+\]
+Der Operator $f(x) \mapsto q(x)f(x)$, der eine Funktion mit 
+der Funktion $q(x)$ multipliziert, ist also ebenfalls symmetrisch.
+Dasselbe gilt für einen Operator, der mit $w(x)$ multipliziert.
+Da $w(x)$ eine positive Funktion ist, ist der Operator $f(x)\mapsto w(x)f(x)$
+sogar positiv definit.
+Dies entspricht der Matrix $B$.
+Nach den Erkenntnissen des vorangegangenen Abschnittes ist das
+verallgemeinerte Eigenwertproblem daher ein Eigenwertproblem
+für einen modifizierten Operator bezüglich eines alternativen
+Skalarproduktes.
+
+Als Skalarprodukt muss also das Integral
+\[
+\langle f,g\rangle_w
+=
+\int_a^b f(x)g(x)w(x)\,dx
+\]
+mit der Gewichtsfunktion $w(x)$ verwendet werden.
+Damit dies ein vernünftiges Skalarprodukt ist, muss $w(x)>0$ im
+Innerend es Intervalls sein.
+
+\subsection{Der Vektorraum $H$}
+Damit können wir jetzt die Eigenschaften der in Frage kommenden
+Funktionen zusammenstellen.
+Zunächst müssen sie auf dem Intervall $[a,b]$ definiert sein und
+das Integral
+\[
+\int_a^b |f(x)|^2 w(x)\,dx < \infty
+\]
+muss existieren.
+Wir bezeichnen den Vektorraum der Funktionen, deren Quadrat mit
+der Gewichtsfunktion $w(x)$ integrierbar sind, mit
+$L^2([a,b],w)$.
+
+Damit auch $\langle qf,f\rangle_w$ und $\langle L_0f,0f\rangle_w$
+wohldefiniert sind, müssen zusätzlich die Integrale
+\[
+\int_a^b |f(x)|^2 q(x) w(x)\,dx
+\qquad\text{und}\qquad
+\int_a^b |f'(x)|^2 p(x) w(x)\,dx
+\]
+existieren.
+Wir setzen daher
+\[
+H
+=
+\biggl\{
+f\in L^2([a,b],w)\;\bigg|\;
+\int_a^b |f'(x)|^2p(x)w(x)\,dx<\infty,
+\int_a^b |f(x)|^2q(x)w(x)\,dx<\infty
+\biggr\}.
+\]
+
+\subsection{Differentialoperator}
+Das verallgemeinerte Eigenwertproblem für $A$ und $B$ ist ein
+gewöhnliches Eigenwertproblem für die Operator $\tilde{A}=B^{-1}A$
+bezüglich des modifizierten Skalarproduktes.
+Das Sturm-Liouville-Problem ist also ein Eigenwertproblem im
+Vektorraum $H$ mit dem Skalarprodukt $\langle\,\;,\;\rangle_w$.
+Der Operator
+\[
+L = \frac{1}{w(x)} \biggl(-\frac{d}{dx} p(x)\frac{d}{dx} + q(x)\biggr)
+\]
+heisst der {\em Sturm-Liouville-Operator}.
+Eine Lösung des Sturm-Liouville-Problems ist eine Funktion $y(x)$ derart,
+dass 
+\[
+Ly = \lambda y,
+\]
+$\lambda$ ist der zu $y(x)$ gehörige Eigenwert.
+Der Operator ist definiert auf Funktionen des im vorangegangenen Abschnitt
+definierten Vektorraumes $H$.
+
+\subsection{Beispiele}
+Die meisten der früher vorgestellten Funktionenfamilien stellen sich
+als Lösungen eines geeigneten Sturm-Liouville-Problems heraus.
+Alle Eigenschaften aus der Sturm-Liouville-Theorie gelten daher
+automatisch für diese Funktionenfamilien.
+
+\subsubsection{Trigonometrische Funktionen}
+Die trigonometrischen Funktionen sind Eigenfunktionen des Operators
+$d^2/dx^2$, also eines Sturm-Liouville-Operators mit $p(x)=1$, $q(x)=0$
+und $w(x)=0$.
+Auf dem Intervall $(-\pi,\pi)$ können wir die Randbedingungen
+\bgroup
+\renewcommand{\arraycolsep}{2pt}
+\[
+\begin{aligned}
+&
+\begin{array}{lclclcl}
+k_{-\pi}          &=&1,&\qquad&h_{-\pi}          &=&0\\
+k_{\phantom{-}\pi}&=&1,&\qquad&h_{\phantom{-}\pi}&=&0
+\end{array}
+\;\bigg\}
+&&\Rightarrow&
+\begin{array}{lcl}
+y(-\pi)          &=&0\\
+y(\phantom{-}\pi)&=&0\\
+\end{array}
+\;\bigg\}
+&\quad\Rightarrow&
+y(x) &= B\sin nx
+\\
+&
+\begin{array}{lclclcl}
+k_{-\pi}          &=&0,&\qquad&h_{-\pi}          &=&1\\
+k_{\phantom{-}\pi}&=&0,&\qquad&h_{\phantom{-}\pi}&=&1
+\end{array}
+\;\bigg\}
+&&\Rightarrow&
+\begin{array}{lcl}
+y'(-\pi)          &=&0\\
+y'(\phantom{-}\pi)&=&0\\
+\end{array}
+\; \bigg\}
+&\quad\Rightarrow&
+y(x) &= A\cos nx
+\end{aligned}
+\]
+\egroup
+verwenden.
+Die Orthogonalität der Sinus- und Kosinus-Funktionen folgt jetzt
+ganz ohne weitere Rechnung.
+
+An dieser Lösung ist nicht ganz befriedigend, dass die trigonometrischen
+Funktionen nicht mit einer einzigen Randbedingung gefunden werden können.
+Der Ausweg ist, periodische Randbedingungen zu verlangen, also
+$y(-\pi)=y(\pi)$ und $y'(-\pi)=y'(\pi)$.
+Dann ist wegen
+\begin{align*}
+\langle f,L_0g\rangle - \langle L_0f,g\rangle
+&=
+-f(\pi)g'(\pi)+f(-\pi)g'(-\pi)+f'(\pi)g(\pi)-f'(-\pi)g(-\pi)
+\\
+&=
+-f(\pi)g'(\pi)+f(\pi)g'(\pi)+f'(\pi)g(\pi)-f'(\pi)g(\pi)
+=0
+\end{align*}
+die Bedingung~\eqref{buch:integrale:sturm:sabedingung}
+ebenfalls erfüllt, $L_0$ ist in diesem Raum selbstadjungiert.
+
+\subsubsection{Bessel-Funktionen}
+Der Bessel-Operator \eqref{buch:differentialgleichungen:bessel-operator}
+hat die Form eines Sturm-Liouville-Operators
+\[
+x^2\frac{d^2}{dx^2} + x\frac{d}{dx} + x^2
+=
+\frac{d}{dx} x^2 \frac{d}{dx} + x^2
+\]
+mit $p(x)=x^2$, $q(x)=x^2$.
+
+XXX TODO: Faktor 2 fehlt.
+
+\subsubsection{Tschebyscheff-Polynome}
+Die Tschebyscheff-Polynome sind Lösungen der
+Tschebyscheff-Differentialgleichung
+\[
+(1-x^2)y'' -xy' = n^2y
+\]
+auf dem Intervall $(-1,1)$.
+Auch dieses Problem kann als Sturm-Liouville-Problem formuliert
+werden mit
+\begin{align*}
+w(x) &= \frac{1}{\sqrt{1-x^2}} \\
+p(x) &= \sqrt{1-x^2} \\
+q(x) &= 0
+\end{align*}
+Das zugehörige Sturm-Liouville-Eigenwertproblem ist
+\[
+\frac{d}{dx}\sqrt{1-x^2}\frac{d}{dx} y(x)
+=
+\lambda \frac{1}{\sqrt{1-x^2}} y(x).
+\]
+Führt man die Ableitungen auf der linken Seite aus, entsteht die
+Gleichung
+\begin{align*}
+\sqrt{1-x^2} y''(x) -\frac{x}{\sqrt{1-x^2}}y'(x)
+&=  \lambda \frac{1}{\sqrt{1-x^2}} y(x)
+\intertext{Multiplikation mit $\sqrt{1-x^2}$ ergibt}
+(1-x^2)
+y''(x) 
+-
+xy'(x)
+&=
+\lambda y(x).
+\end{align*}
+Es folgt, dass die Tschebyscheff-Polynome orthogonal sind 
+bezüglich des Skalarproduktes
+\[
+\langle f,g\rangle = \int_{-1}^1 f(x)g(x)\frac{dx}{\sqrt{1-x^2}}.
+\]
+
+\subsubsection{Jacobi-Polynome}
+TODO
+
+\subsubsection{Hypergeometrische Differentialgleichungen}
+%\url{https://encyclopediaofmath.org/wiki/Hypergeometric_equation}
+Auch die Eulersche hypergeometrische Differentialgleichung
+lässt sich in die Form eines Sturm-Liouville-Operators
+bringen.
+Dazu setzt man
+\begin{align*}
+p(z)
+&=
+z^c(z-1)^{a+b+1-c}
+\\
+q(z)
+&=
+-abz^{c-1}(z-1)^{a+b-c}
+\\
+w(z)
+&=
+z^{c-1}(z-1)^{a+b-c}.
+\end{align*}
+Setzt man dies in den Sturm-Liouville-Operator ein, erhält man
+\begin{equation}
+L
+=
+-\frac{d}{dz}p(z)\frac{d}{dz} + q(z)
+=
+-p(z)\frac{d^2}{dz^2}
+-p'(z)\frac{d}{dz}
++q(z)
+\label{buch:orthgonalitaet:eqn:hypersturm}
+\end{equation}
+Wir brauchen also
+\begin{align*}
+p'(z)
+&=
+cz^{c-1}(z-1)^{a+b+1-c}
++
+(a+b+1-c)
+z^c
+(z-1)^{a+b-c}
+\\
+&=
+\bigl(
+c(z-1)+
+(a+b+1-c)z
+\bigr)
+\cdot
+z^{c-1}(z-1)^{a+b-c}
+\\
+&=
+-
+\bigl(
+c-(a+b+1)z
+\bigr)
+\cdot
+z^{c-1}(z-1)^{a+b-c}.
+\end{align*}
+Einsetzen in~\eqref{buch:orthgonalitaet:eqn:hypersturm} liefert
+\begin{align*}
+L
+%=
+%-\frac{d}{dz}p(z)\frac{d}{dz}+q(z)
+&=
+-z^c(z-1)^{a+b+1-c} \frac{d^2}{dz^2}
++
+w(z)
+(c-(a+b+1)z)
+\frac{d}{dz}
+-
+abw(z)
+\\
+&=
+w(z)
+\biggl(
+-
+z(z-1)
+\frac{d^2}{dz^2}
++
+(c-(a+b+1)z)
+\frac{d}{dz}
+-ab
+\biggr)
+\\
+&=
+w(z)
+\biggl(
+z(1-z)
+\frac{d^2}{dz^2}
++
+(c-(a+b+1)z)
+\frac{d}{dz}
+-ab
+\biggr).
+\end{align*}
+Die Klammer auf der rechten Seite ist tatsächlich die linke Seite der
+eulerschen hypergeometrischen Differentialgleichung.
+
+Die hypergeometrische Funktion $\mathstrut_2F_1(a,b;c;z)$ ist ein
+Eigenvektor des Operators $L$ zum Eigenwert $\lambda$.
+Sei jetzt $w(z)$ eine Eigenfunktion zum Eigenwert $\lambda\ne 0$,
+also
+\[
+z(1-z)w''(z) + (c-(a+b+1)z)w'(z) - ab w(z) = \lambda w(z).
+\]
+Kann man $a$ und $b$ so in $a_1$ und $b_1$ ändern, dass $a+b=a_1+b_1$
+gleich bleiben aber das Produkt den Wert $a_1b_1=ab-\lambda$?
+$a_1$ und $b_1$ sind die Lösungen der quadratischen Gleichung
+\[
+x^2 - (a+b)x + ab-\lambda = 0.
+\]
+Alle Eigenfunktionen des Operators $L$ sind also hypergeometrische
+Funktion $\mathstrut_2F_1$.
+
+Da die Gewichtsfunktion $w(z)$ bei der Ersetzung $a\to a_1$ und $b\to b_1$
+sich nicht ändert ($w(z)$ hängt nur von der Summe $a+b$ ab, welche sich
+nicht ändert), sind die beide beiden Eigenfunktionen bezüglich
+des Skalarproduktes mit der Gewichtsfunktion $w(z)$ orthogonal.
+
+
+
+
+
+