Aproximación de Lanczos

Método numérico para calcular la función gamma.

En matemáticas , la aproximación de Lanczos es un método para calcular numéricamente la función gamma , publicado por Cornelius Lanczos en 1964. Es una alternativa práctica a la aproximación de Stirling, más popular , para calcular la función gamma con precisión fija.

Introducción

La aproximación de Lanczos consta de la fórmula

${\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)}$

para la función gamma, con

${\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 .}$

Aquí g es una constante real que puede elegirse arbitrariamente sujeta a la restricción de que Re( z + g + ⁠1/2⁠ ) ​​> 0. ^[1] Los coeficientes p , que dependen de g , son un poco más difíciles de calcular (ver más abajo). Aunque la fórmula indicada aquí solo es válida para argumentos en el semiplano complejo derecho , se puede extender a todo el plano complejo mediante la fórmula de reflexión ,

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

La serie A es convergente y puede truncarse para obtener una aproximación con la precisión deseada. Al elegir una g apropiada (normalmente un número entero pequeño), sólo se necesitan entre 5 y 10 términos de la serie para calcular la función gamma con la típica precisión de punto flotante simple o doble . Si se elige una g fija, los coeficientes se pueden calcular de antemano y, gracias a la descomposición en fracciones parciales , la suma se transforma en la siguiente forma:

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

Por tanto, calcular la función gamma se convierte en una cuestión de evaluar sólo un pequeño número de funciones elementales y multiplicarlas por constantes almacenadas. La aproximación de Lanczos fue popularizada por Numerical Recipes , según la cual calcular la función gamma se vuelve "no mucho más difícil que otras funciones integradas que damos por sentado, como sin  x o e ^x ". El método también se implementa en la Biblioteca Científica GNU , Boost , CPython y musl .

Coeficientes

Los coeficientes están dados por

${\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}}$

donde representa el ( n , m )ésimo elemento de la matriz de coeficientes de los polinomios de Chebyshev , que se puede calcular recursivamente a partir de estas identidades: ${\displaystyle C_{n,m}}$

${\displaystyle {\begin{aligned}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{aligned}}}$

Godfrey (2001) describe cómo obtener los coeficientes y también el valor de la serie truncada A como producto matricial . ^[2]

Derivación

Lanczos derivó la fórmula de la integral de Leonhard Euler

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

realizando una secuencia de manipulaciones básicas para obtener

${\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,}$

y derivar una serie para la integral.

Implementación sencilla

La siguiente implementación en el lenguaje de programación Python funciona para argumentos complejos y normalmente proporciona 13 decimales correctos. Tenga en cuenta que omitir los coeficientes más pequeños (en busca de velocidad, por ejemplo) da resultados totalmente inexactos; los coeficientes deben recalcularse desde cero para una expansión con menos términos.

desde  cmath  importar  pecado ,  sqrt ,  pi ,  exp""" Los coeficientes utilizados en el código son para cuando g = 7 y n = 9 Aquí hay algunos otros ejemplosg = 5 n = 5 p = [  1.0000018972739440364,  76.180082222642137322,  -86.505092037054859197,  24.012898581922685900,  -1.2296028490285820 771 ]g = 5 n = 7 p = [  1.0000000001900148240,  76.180091729471463483,  -86.505320329416767652,  24.014098240830910490,  -1.2317395724501553 875,  0.0012086509738661785061,  -5.3952393849531283785e-6 ]g = 8 n = 12 p = [  0.9999999999999999298,  1975.3739023578852322,  -4397.3823927922428918,  3462.6328459862717019,  -1156.9851431631167820,  154.53815050252775060,  -6.2536716123689161798,  0.034642762454736807441,  -7.4776171974442977377e-7,  6.3041253821852264261e-8,  -2.7405717035683877489e-8,  4.0486948817567609101e-9 ] """g  =  7 n  =  9 p  =  [  0.99999999999980993 ,  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 ):  si  abs ( z . imag )  <=  EPSILON :  z  =  z . retorno real  z def  gamma ( z ):  z  =  complejo ( z )  si  z . real  <  0.5 :  y  =  pi  /  ( sin ( pi  *  z )  *  gamma ( 1  -  z ))  # Fórmula de reflexión  else :  z  -=  1  x  =  p [ 0 ]  para  i  en el  rango ( 1 ,  len ( p )) :  x  +=  p [ i ]  /  ( z  +  i )  t  =  z  +  g  +  0.5  y  =  sqrt ( 2  *  pi )  *  t  **  ( z  +  0.5 )  *  exp ( - t )  *  x  return  drop_imag ( y ) """ El uso anterior de la reflexión (por lo tanto la estructura if-else) es necesario, aunque pueda parecer extraño, ya que permite extender la aproximación a valores de z donde Re(z) < 0,5, donde los Lanczos El método no es válido.imprimir ( gamma ( 1 )) imprimir ( gamma ( 5 ))  imprimir ( gamma ( 0.5 ))

Ver también

• La aproximación de Stirling
• La aproximación de Spouge

Referencias

1. Pugh, Glendon (2004). Un análisis de la aproximación de Lanczos Gamma (PDF) (Ph.D.).
2. Godfrey, Paul (2001). "Implementación de Lanczos de la función gamma". Numérica .

• Godofredo, Paul (2001). "Implementación de Lanczos de la función Gamma".
• Lanczos, Cornelio (1964). "Una aproximación de precisión de la función gamma". Revista de la Sociedad de Matemáticas Industriales y Aplicadas, Serie B: Análisis numérico . 1 : 86–96. Código Bib : 1964SJNA....1...86L. doi :10.1137/0701008. ISSN  0887-459X. JSTOR  2949767.
• Prensa, WH; Teukolsky, SA; Vetterling, WT; Flannery, BP (2007), "Sección 6.1. Función gamma", Recetas numéricas: el arte de la informática científica (3.ª ed.), Nueva York: Cambridge University Press, ISBN 978-0-521-88068-8
• Pugh, Glendon (2004). Un análisis de la aproximación de Lanczos Gamma (PDF) (tesis doctoral).
• Toth, Viktor (2005). "Calculadoras programables: la aproximación de Lanczos".
• Weisstein, Eric W. "Aproximación de Lanczos". MundoMatemático .