aboutsummaryrefslogtreecommitdiffstats
diff options
context:
space:
mode:
authorPatrik Müller <patrik.mueller@ost.ch>2022-07-15 11:40:55 +0200
committerPatrik Müller <patrik.mueller@ost.ch>2022-07-15 11:40:55 +0200
commit3cb2fa354f814fa98474610dac744281285dafc6 (patch)
treeeecf8ed0e258f58670213f97dc469509e128ae10
parentCorrect Makefile, add text to gamma.tex, separate python-scripts for each image (diff)
downloadSeminarSpezielleFunktionen-3cb2fa354f814fa98474610dac744281285dafc6.tar.gz
SeminarSpezielleFunktionen-3cb2fa354f814fa98474610dac744281285dafc6.zip
First version of section 'Gauss Quadratur', fix to gamma_approx.py when z=0
-rw-r--r--buch/papers/laguerre/quadratur.tex148
-rw-r--r--buch/papers/laguerre/references.bib11
-rw-r--r--buch/papers/laguerre/scripts/gamma_approx.py18
-rw-r--r--buch/papers/laguerre/scripts/rel_error_range.py2
4 files changed, 151 insertions, 28 deletions
diff --git a/buch/papers/laguerre/quadratur.tex b/buch/papers/laguerre/quadratur.tex
index 851fe8a..7cbae48 100644
--- a/buch/papers/laguerre/quadratur.tex
+++ b/buch/papers/laguerre/quadratur.tex
@@ -5,25 +5,57 @@
%
\section{Gauss-Quadratur
\label{laguerre:section:quadratur}}
- {\large \color{red} TODO: Einleitung und kurze Beschreibung Gauss-Quadratur}
-
-Siehe Abschnitt~\ref{buch:orthogonalitaet:section:gauss-quadratur}
+Die Gauss-Quadratur ist ein numerisches Integrationsverfahren,
+welches die Eigenschaften von orthogonalen Polynomen ausnützt.
+Herleitungen und Analysen der Gauss-Quadratur können im
+Abschnitt~\ref{buch:orthogonalitaet:section:gauss-quadratur} gefunden werden.
+Als grundlegende Idee wird die Beobachtung,
+dass viele Funktionen sich gut mit Polynomen approximieren lassen,
+verwendet.
+Stellt man also sicher,
+dass ein Verfahren gut für Polynome gut funktioniert,
+sollte es auch für andere Funktionen nicht schlecht funktionieren.
+Es wird ein Polynom verwendet,
+welches an den Punkten $x_0 < x_1 < \ldots < x_n$
+die Funktionwerte~$f(x_i)$ annimmt.
+Als Resultat kann das Integral via eine gewichtete Summe der Form
\begin{align}
\int_a^b f(x) w(x) \, dx
\approx
\sum_{i=1}^n f(x_i) A_i
\label{laguerre:gaussquadratur}
\end{align}
+berechnet werden.
+Die Gauss-Quadratur ist exakt für Polynome mit Grad $2n -1$,
+wenn ein Interpolationspolynom von Grad $n$ gewählt wurde.
\subsection{Gauss-Laguerre-Quadratur
\label{laguerre:subsection:gausslag-quadratur}}
-Die Gauss-Quadratur kann auch auf Skalarprodukte mit Gewichtsfunktionen
-ausgeweitet werden.
-In unserem Falle möchten wir die Gauss Quadratur auf die Laguerre-Polynome
-$L_n$ ausweiten.
-Diese sind orthogonal im Intervall $(0, \infty)$ bezüglich
-der Gewichtsfunktion $e^{-x}$.
-Gleichung~\eqref{laguerre:laguerrequadratur} lässt sich wie folgt umformulieren:
+Wir möchten nun die Gauss-Quadratur auf die Berechnung
+von uneigentlichen Integralen erweitern,
+spezifisch auf das Interval $(0, \infty)$.
+Mit dem vorher beschriebenen Verfahren ist dies nicht direkt möglich.
+Mit einer Transformation die das unendliche Intervall $(a, \infty)$ mit
+\begin{align*}
+x
+=
+a + \frac{1 - t}{t}
+\end{align*}
+auf das Intervall $[0, 1]$ transformiert.
+Für unser Fall gilt $a = 0$.
+Das Integral eines Polynomes in diesem Intervall ist immer divergent,
+darum müssen wir sie mit einer Funktion multiplizieren,
+die schneller als jedes Polynom gegen $0$ geht,
+damit das Integral immer noch konvergiert.
+Die Laguerre-Polynome $L_n$ bieten hier Abhilfe,
+da ihre Gewichtsfunktion $e^{-x}$ schneller
+gegen $0$ konvergiert als jedes Polynom.
+% In unserem Falle möchten wir die Gauss Quadratur auf die Laguerre-Polynome
+% $L_n$ ausweiten.
+% Diese sind orthogonal im Intervall $(0, \infty)$ bezüglich
+% der Gewichtsfunktion $e^{-x}$.
+Gleichung~\eqref{laguerre:gaussquadratur} lässt sich wie folgt
+umformulieren:
\begin{align}
\int_{0}^{\infty} f(x) e^{-x} dx
\approx
@@ -43,20 +75,93 @@ l_i(x_j)
\delta_{ij}
=
\begin{cases}
-1 & i=j \\
-0 & \text{sonst.}
+1 & i=j \\
+0 & \text{.}
\end{cases}
+% .
\end{align*}
-Laut \cite{abramowitz+stegun} sind die Gewichte
-\begin{align}
+die Lagrangschen Interpolationspolynome.
+Laut \cite{hildebrand2013introduction} können die Gewicht mit
+\begin{align*}
A_i
+ & =
+-\frac{C_{n+1} \gamma_n}{C_n \phi'_n(x_i) \phi_{n+1} (x_i)}
+\end{align*}
+berechnet werden.
+$C_i$ entspricht dabei dem Koeffizienten von $x^i$
+des orthogonalen Polynoms $\phi_n(x)$, $\forall i =0,\ldots,n$ und
+\begin{align*}
+\gamma_n
+=
+\int_0^\infty w(x) \phi_n^2(x)\,dx
+\end{align*}
+dem Normalisierungsfaktor.
+Wir setzen nun $\phi_n(x) = L_n(x)$ und
+nutzen den Vorzeichenwechsel der Laguerrekoeffizienten aus,
+damit erhalten wir
+\begin{align*}
+A_i
+ & =
+-\frac{C_{n+1} \gamma_n}{C_n L'_n(x_i) L_{n+1} (x_i)}
+\\
+ & = \frac{C_n}{C_{n-1}} \frac{\gamma_{n-1}}{L_{n-1}(x_i) L'_n(x_i)}
+.
+\end{align*}
+Für Laguerre-Polynome gilt
+\begin{align*}
+\frac{C_n}{C_{n-1}}
+=
+-\frac{1}{n}
+\quad \text{und} \quad
+\gamma_n
=
-\frac{x_i}{(n + 1)^2 \left[ L_{n + 1}(x_i)\right]^2}
+1
+.
+\end{align*}
+Daraus folgt
+\begin{align}
+A_i
+&=
+- \frac{1}{n L_{n-1}(x_i) L'_n(x_i)}
+.
+\label{laguerre:gewichte_lag_temp}
+\end{align}
+Nun kann die Rekursionseigenschaft der Laguerre-Polynome
+\begin{align*}
+x L'_n(x)
+&=
+n L_n(x) - n L_{n-1}(x)
+\\
+&= (x - n - 1) L_n(x) + (n + 1) L_{n+1}(x)
+\end{align*}
+umgeformt werden und da $x_i$ die Nullstellen von $L_n(x)$ sind,
+folgt
+\begin{align*}
+x_i L'_n(x_i)
+&=
+- n L_{n-1}(x_i)
+\\
+&=
+ (n + 1) L_{n+1}(x_i)
+.
+\end{align*}
+Setzen wir das nun in \eqref{laguerre:gewichte_lag_temp} ein ergibt sicht
+\begin{align}
+\nonumber
+A_i
+&=
+\frac{1}{x_i \left[ L'_n(x_i) \right]^2}
+\\
+&=
+\frac{x_i}{(n+1)^2 \left[ L_{n+1}(x_i) \right]^2}
.
\label{laguerre:quadratur_gewichte}
\end{align}
\subsubsection{Fehlerterm}
+Die Gauss-Laguerre-Quadratur mit $n$ Stützstellen berechnet Integrale
+von Polynomen bis zum Grad $2n - 1$ exakt.
+Für beliebige Funktionen kann eine Fehlerabschätzung angegeben werden.
Der Fehlerterm $R_n$ folgt direkt aus der Approximation
\begin{align*}
\int_0^{\infty} f(x) e^{-x} \, dx
@@ -66,16 +171,15 @@ Der Fehlerterm $R_n$ folgt direkt aus der Approximation
und \cite{abramowitz+stegun} gibt ihn als
\begin{align}
R_n
-=
+ & =
+\frac{f^{(2n)}(\xi)}{(2n)!} \int_0^\infty l(x)^2 e^{-x}\,dx
+\\
+ & =
\frac{(n!)^2}{(2n)!} f^{(2n)}(\xi)
,\quad
0 < \xi < \infty
\label{laguerre:lag_error}
\end{align}
an.
-
-{
-\large \color{red}
-TODO:
-Noch mehr Text / bessere Beschreibungen in allen Abschnitten
-}
+Der Fehler ist also abhängig von der $2n$-ten Ableitung
+der zu integrierenden Funktion.
diff --git a/buch/papers/laguerre/references.bib b/buch/papers/laguerre/references.bib
index e12e218..2371922 100644
--- a/buch/papers/laguerre/references.bib
+++ b/buch/papers/laguerre/references.bib
@@ -4,6 +4,17 @@
% (c) 2020 Autor, Hochschule Rapperswil
%
+@book{hildebrand2013introduction,
+ title={Introduction to Numerical Analysis: Second Edition},
+ author={Hildebrand, F.B.},
+ isbn={9780486318554},
+ series={Dover Books on Mathematics},
+ url={https://books.google.ch/books?id=ic2jAQAAQBAJ},
+ year={2013},
+ publisher={Dover Publications},
+ pages = {389}
+}
+
@book{abramowitz+stegun,
added-at = {2008-06-25T06:25:58.000+0200},
address = {New York},
diff --git a/buch/papers/laguerre/scripts/gamma_approx.py b/buch/papers/laguerre/scripts/gamma_approx.py
index 208f770..9f9dae7 100644
--- a/buch/papers/laguerre/scripts/gamma_approx.py
+++ b/buch/papers/laguerre/scripts/gamma_approx.py
@@ -1,7 +1,5 @@
from pathlib import Path
-import matplotlib as mpl
-import matplotlib.pyplot as plt
import numpy as np
import scipy.special
@@ -58,6 +56,8 @@ def laguerre_gamma_shifted(z, x=None, w=None, n=8, target=11):
def laguerre_gamma_opt_shifted(z, x=None, w=None, n=8):
+ if z == 0.0:
+ return np.infty
x, w = _prep_zeros_and_weights(x, w, n)
n = len(x)
@@ -73,6 +73,8 @@ def laguerre_gamma_opt_shifted(z, x=None, w=None, n=8):
def laguerre_gamma_simple(z, x=None, w=None, n=8):
+ if z == 0.0:
+ return np.infty
x, w = _prep_zeros_and_weights(x, w, n)
z += 0j
res = np.sum(x ** (z - 1) * w)
@@ -81,6 +83,8 @@ def laguerre_gamma_simple(z, x=None, w=None, n=8):
def laguerre_gamma_mirror(z, x=None, w=None, n=8):
+ if z == 0.0:
+ return np.infty
x, w = _prep_zeros_and_weights(x, w, n)
z += 0j
if z.real < 1e-3:
@@ -90,8 +94,8 @@ def laguerre_gamma_mirror(z, x=None, w=None, n=8):
return laguerre_gamma_simple(z, x, w)
-def eval_laguerre_gamma(z, x=None, w=None, n=8, func="simple", **kwargs):
- x, w = _prep_zeros_and_weights(x, w, n)
+def eval_laguerre_gamma(z, x=None, w=None, n=8, func="simple", **kwargs):
+ x, w = _prep_zeros_and_weights(x, w, n)
if func == "simple":
f = laguerre_gamma_simple
elif func == "mirror":
@@ -104,4 +108,8 @@ def eval_laguerre_gamma(z, x=None, w=None, n=8, func="simple", **kwargs):
def calc_rel_error(x, y):
- return (y - x) / x
+ mask = np.abs(x) != np.infty
+ rel_error = np.zeros_like(y)
+ rel_error[mask] = (y[mask] - x[mask]) / x[mask]
+ rel_error[~mask] = 0.0
+ return rel_error
diff --git a/buch/papers/laguerre/scripts/rel_error_range.py b/buch/papers/laguerre/scripts/rel_error_range.py
index 7d017a7..7c74d76 100644
--- a/buch/papers/laguerre/scripts/rel_error_range.py
+++ b/buch/papers/laguerre/scripts/rel_error_range.py
@@ -5,7 +5,7 @@ if __name__ == "__main__":
import gamma_approx as ga
- N = 1000
+ N = 1001
xmin = -5
xmax = 5
ns = np.arange(2, 12, 2)