Lanczos Annäherung - Lanczos approximation

In der Mathematik ist die Lanczos-Näherung eine Methode zur numerischen Berechnung der Gammafunktion , die 1964 von Cornelius Lanczos veröffentlicht wurde. Sie ist eine praktische Alternative zur populäreren Stirling-Näherung zur Berechnung der Gammafunktion mit fester Genauigkeit.

Einführung

Die Lanczos-Näherung besteht aus der Formel

${\ displaystyle \ Gamma (z + 1) = {\ sqrt {2 \ pi}} {\ left (z + g + {\ tfrac {1} {2}} \ right)} ^ {z + 1/2} e ^ {- (z + g + 1/2)} A_ {g} (z)}$

für die Gammafunktion mit

${\ displaystyle A_ {g} (z) = {\ frac {1} {2}} p_ {0} (g) + p_ {1} (g) {\ frac {z} {z + 1}} + p_ {2} (g) {\ frac {z (z-1)} {(z + 1) (z + 2)}} + \ cdots.}$

Hier ist g eine Konstante , die willkürlich gewählt werden kann, vorbehaltlich der Einschränkung, dass Re ( z )>  1 /. 2 . Die von g abhängigen Koeffizienten p sind etwas schwieriger zu berechnen (siehe unten). Obwohl die hier angegebene Formel nur für Argumente in der rechten komplexen Halbebene gilt , kann sie durch die Reflexionsformel auf die gesamte komplexe Ebene erweitert werden .

${\ displaystyle \ Gamma (1-z) \; \ Gamma (z) = {\ pi \ over \ sin \ pi z}.}$

Die Reihe A ist konvergent und kann abgeschnitten werden, um eine Annäherung mit der gewünschten Genauigkeit zu erhalten. Durch die Wahl eines geeigneten g (typischerweise eine kleine ganze Zahl ist ), nur etwa 5-10 Glieder der Reihe nötig sind , um die Gamma - Funktion mit typischen zur Berechnung einzelner oder doppelter Fließkommagenauigkeit. Wenn ein festes g gewählt wird, können die Koeffizienten im Voraus berechnet werden und die Summe wird in die folgende Form umformuliert:

${\ displaystyle A_ {g} (z) = c_ {0} + \ sum _ {k = 1} ^ {N} {\ frac {c_ {k}} {z + k}}}$

Die Berechnung der Gammafunktion wird somit zu einer Frage der Auswertung nur einer kleinen Anzahl von Elementarfunktionen und der Multiplikation mit gespeicherten Konstanten. Die Lanczos-Näherung wurde von Numerical Recipes populär gemacht , wonach die Berechnung der Gammafunktion "nicht viel schwieriger wird als andere eingebaute Funktionen, die wir für selbstverständlich halten, wie sin  x oder e ^x ". Die Methode ist auch in der GNU Scientific Library , Boost , CPython und musl implementiert .

Koeffizienten

Die Koeffizienten sind gegeben durch

${\ displaystyle p_ {k} (g) = {\ frac {\ sqrt {2 \,}} {\ pi}} \ sum _ {\ ell = 0} ^ {k} C_ {2k + 1, \, 2 \ ell +1} \ left (\ ell - {\ tfrac {1} {2}} \ right)! {\ left (\ ell + g + {\ tfrac {1} {2}} \ right)} ^ {- (\ ell +1/2)} e ^ {\ ell + g + 1/2}}$

wo stellt das ( n , m ) -te Element der Matrix von Koeffizienten für die Chebyshev - Polynome , die berechnet werden können rekursiv aus diesen Identitäten: ${\ displaystyle C_ {n, m}}$

${\ displaystyle {\ begin {align} C_ {1, \, 1} & = 1 \\ [5px] C_ {2, \, 2} & = 1 \\ [5px] C_ {n + 1, \, 1 } & = - \, C_ {n-1, \, 1} & {\ text {for}} n & = 2,3,4 \, \ dots \\ [5px] C_ {n + 1, \, n + 1} & = 2 \, C_ {n, \, n} & {\ text {for}} n & = 2,3,4 \, \ dots \\ [5px] C_ {n + 1, \, m + 1 } & = 2 \, C_ {n, \, m} -C_ {n-1, \, m + 1} & {\ text {for}} n &> m = 1,2,3 \, \ dots \ end {ausgerichtet}}}$

Godfrey (2001) beschreibt, wie man die Koeffizienten und auch den Wert der verkürzten Reihe A als Matrixprodukt erhält .

Ableitung

Lanczos abgeleitet die Formel aus Leonhard Euler ist integral

${\ displaystyle \ Gamma (z + 1) = \ int _ {0} ^ {\ infty} t ^ {z} \, e ^ {- t} \, dt,}$

Durchführen einer Folge grundlegender Manipulationen, um zu erhalten

${\ displaystyle \ Gamma (z + 1) = (z + g + 1) ^ {z + 1} e ^ {- (z + g + 1)} \ int _ {0} ^ {e} {\ Big ( } v (1- \ log v) {\ Big)} ^ {z} v ^ {g} \, dv,}$

und Ableiten einer Reihe für das Integral.

Einfache Implementierung

Die folgende Implementierung in der Programmiersprache Python funktioniert für komplexe Argumente und gibt normalerweise 15 korrekte Dezimalstellen an. Beachten Sie, dass das Weglassen der kleinsten Koeffizienten nicht zu einer schnelleren, aber etwas weniger genauen Implementierung führt. Die Koeffizienten müssen für eine Erweiterung mit weniger Termen von Grund auf neu berechnet werden.

from cmath import sin, sqrt, pi, exp

p = [676.5203681218851
    ,-1259.1392167224028
    ,771.32342877765313
    ,-176.61502916214059
    ,12.507343278686905
    ,-0.13857109526572012
    ,9.9843695780195716e-6
    ,1.5056327351493116e-7
    ]

EPSILON = 1e-07
def drop_imag(z):
    if abs(z.imag) <= EPSILON:
        z = z.real
    return z

def gamma(z):
    z = complex(z)
    if z.real < 0.5:
        y = pi / (sin(pi * z) * gamma(1 - z))  # Reflection formula
    else:
        z -= 1
        x = 0.99999999999980993
        for (i, pval) in enumerate(p):
            x += pval / (z + i + 1)
        t = z + len(p) - 0.5
        y = sqrt(2 * pi) * t ** (z + 0.5) * exp(-t) * x
    return drop_imag(y)
"""
The above use of the reflection (thus the if-else structure) is necessary, even though 
it may look strange, as it allows to extend the approximation to values of z where 
Re(z) < 0.5, where the Lanczos method is not valid.
"""

print(gamma(1))
print(gamma(5))   
print(gamma(0.5))

Siehe auch

• Stirlings Annäherung
• Spouges Annäherung

Verweise

• Godfrey, Paul (2001). "Lanczos Implementierung der Gammafunktion" .
• Lanczos, Cornelius (1964). "Eine präzise Annäherung der Gammafunktion". Zeitschrift der Gesellschaft für industrielle und angewandte Mathematik, Reihe B: Numerische Analyse . 1 : 86–96. Bibcode : 1964SJNA .... 1 ... 86L . doi : 10.1137 / 0701008 . ISSN   0887-459X . JSTOR   2949767 .
• Drücken Sie, WH; Teukolsky, SA; Vetterling, WT; Flannery, BP (2007), "Abschnitt 6.1. Gammafunktion" , Numerische Rezepte: Die Kunst des wissenschaftlichen Rechnens (3. Aufl.), New York: Cambridge University Press, ISBN   978-0-521-88068-8
• Pugh, Glendon (2004). Eine Analyse der Lanczos-Gamma-Näherung (PDF) (Dissertation).
• Toth, Viktor (2005). "Programmierbare Taschenrechner: Die Lanczos-Näherung" .
• Weisstein, Eric W. "Lanczos Approximation" . MathWorld .