From b4ded112e21249774920681e00456b2a5c2b39de Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Andreas=20M=C3=BCller?= Date: Mon, 7 Jun 2021 09:03:26 +0200 Subject: add numerical computation of lambert W function --- buch/chapters/020-exponential/lambertw.tex | 42 ++++++++++++++++++++++++++++++ 1 file changed, 42 insertions(+) (limited to 'buch/chapters/020-exponential/lambertw.tex') diff --git a/buch/chapters/020-exponential/lambertw.tex b/buch/chapters/020-exponential/lambertw.tex index a7a882c..410d145 100644 --- a/buch/chapters/020-exponential/lambertw.tex +++ b/buch/chapters/020-exponential/lambertw.tex @@ -319,6 +319,48 @@ W(-cbe^{ac}) Die Gleichung hat eine Lösung wenn $-cbe^{ac} > -1/e$ ist. \end{proof} +\subsection{Numerische Berechnung +\label{buch:subsection:lambertberechnung}} +Die $W$-Funktionen sind nur dann nützlich, wenn man sie effizient +berechnen kann. +Leider ist sie nicht Teil der C- oder C++-Standardbibliothek, +man muss sich also mit einer spezialisierten Bibliothek oder einer +eigenen Implementation behelfen. + +\subsubsection{Berechnung mit dem Newton-Algorithmus} +Für $x>-1$ ist die Funktion $W(x)$ ist die Umkehrfunktion der +streng monoton wachsenden und konvexen Funktion $f(x)=xe^x$. +In dieser Situation konvergiert der Newton-Algorithmus zur Bestimmung +der Nullstelle $x=W_0(y)$ von $f(x)-y$ für alle Werte von $y>-1/e$. +Für $W_{-1}(y)$ ist die Situation etwas komplizierter, da für +$x<-1$ die Funktion $f(x)$ nicht konvex ist. + +Ausgehend vom Startwert $x_0$ ist die Iterationsfolge definiert +durch +\[ +x_{n+1} += +x_n - \frac{f(x_n) - y}{f'(x_n)} += +x_n - \frac{x_ne^{x_n}-y}{(1+x_n)e^{x_n}}. +\] +Die Theorie verspricht, dass die Folge quadratisch konvergiert, wenn +der Startwert $x_0$ genügend genau ist. +Für $W_0(y)$ scheint $x_0=\log(1+y)$ ein guter Startwert zu sein, für +$W_{-1}(y)$ funktioniert $x_0=\log(-y)$. + +Die Steigung des Graphen der Funktion $f(x)$ ist für grosse positive +Werte von $x$ sehr gross und für grosse negative Werte von $x$ sehr +klein, was die Konvergenz stark beeinträchtigen kann. +An der Stelle $x=-1$ mit dem Wert $f(-1)=-1/e$ ist die Steigung $0$ +und die Konvergenz des Newton-Algorithmus ist nur noch linear. +Für mittelgrosse Werte von $y$ weg von $-1/e$ kann $W_0(y)$ oder $W_{-1}(y)$ +mit einer Genauigkeit $10^{-15}$ innert weniger als 10 Iterationen +bestimmt werden. + +\subsubsection{GNU scientific library} +Die Lambert $W$-Funktionen $W_0(x)$ und $W_{-1}(x)$ sind auch in der +GNU scientific library \cite{buch:library:gsl} implementiert. % % Verfolgungskurven -- cgit v1.2.1