Evaluación de fracciones continuas

Question

Estoy escribiendo un programa que calculará $\exp(z)$ .

Originalmente usé la serie Taylor, que funcionaba bien. Sin embargo, las fracciones continuas pueden converger más rápidamente que las series de potencias, así que decidí seguir ese camino.

Encontré esta fracción continua en múltiples lugares. Es de La aplicación de las fracciones continuas y sus generalizaciones a los problemas de la teoría de la aproximación de A. N. Khovanskii (1963), pg 114.

$${e^{z}=1+{\cfrac {2z}{2-z+{\cfrac {z^{2}}{6+{\cfrac {z^{2}}{10+{\cfrac {z^{2}}{14+\ddots }}}}}}}}}$$

Se puede representar como

$${e^{z}=1+\cfrac{2z}{2-z +}\cfrac{z^2/6}{1 +}\sum_{m=3}^{\infty}\left({\cfrac{{a_m}^{z^2}}{1}}\right)}$$

(Lo siento, es la primera vez que uso Mathjax).

Annie A.M. Cuyt, Vigdis Petersen, Brigitte Verdonk, Haakon Waadeland y William B. Jones. 2008. Manual de fracciones continuas para funciones especiales (1 ed.). Springer Publishing Company, Incorporated, pg 194.

Estoy calculando la fracción continua utilizando el algoritmo de Lentz modificado de Recetas numéricas en C: El Arte de la Computación Científica (ISBN 0-521-43105-5), pg 171.

La cuestión: sólo funciona para un pequeño conjunto de números. (Por lo que veo, [-30, 30]).

Así que mi pregunta es: ¿se espera esto? Soy relativamente nuevo en las fracciones continuas, así que mientras piense en Los capto, no estoy del todo seguro.

Dado un "generador", la biblioteca boost de C++ puede calcular fracciones continuas. Básicamente, cada llamada al "generador" devuelve el siguiente término de la FC.

Esto es lo que he utilizado (donde $z$ es la entrada y $m$ es el índice del término actual):

• Término 0:
  • A: $0$
  • B: $1$
• Término 1:
  • A: $2z$
  • B: $2-z$
• Término 2:
  • A: $z^2$
  • B: $1$
• Término 3..N:
  • A: $1 / (4 * (2m - 3) * (2m - 1))) * z^2$
  • B: $1$

Dado $z = 38.5$ El algoritmo de Lentz proporciona lo siguiente (cada línea es $f_j$ ):

-1.1095890410958904 1.3657233326736593 -1.8636200592602417 2.81700268519061 -4.711808282946167 8.70960008496205 -17.765408200530924 39.92118587023142 -98.65273729066934 267.59129323043885 -795.1270569517195 2582.9903072474267 -9154.51726459279 35324.19248620738 -148091.74385797494 673154.6572372171 -3.310843427214524e+06 1.758455218773504e+07 -1.0065603099002995e+08 6.197709770714103e+08 -4.097272719854685e+09 2.9029669729927063e+10 -2.2004739456295242e+11 1.7812408854416082e+12 -1.5389163725064727e+13 1.40311430366499e+14 -1.483409360122725e+15 8.317237533312957e+15 2.273093139324771e+16 1.964661705307002e+16 1.9877542316211204e+16 1.985862955594415e+16 1.9860080303623144e+16 1.985997463174032e+16 1.9859981941770148e+16 1.9859981460940616e+16 1.9859981491045844e+16 1.9859981489249772e+16 1.9859981489351984e+16 1.9859981489346428e+16 1.9859981489346716e+16 1.98599814893467e+16

Por lo tanto, converge en $1.98599814893467e+16$ cuando se supone que la respuesta real es: $\approx 52521552285925160$

Answer 1

Esto parece ser un problema computacional. En concreto, puede tratarse de un error de representación de redondeo, que suele producirse cuando intentamos representar un número con una cadena finita de dígitos. En este caso, los errores de redondeo, incluso en los dígitos decimales de los primeros números de la progresión mostrada en el PO, finalmente da lugar a un gran error de convergencia.

La fracción continua es correcta y funciona bien para cualquier $z $ (también hay que tener en cuenta que sólo es una variación de otras Fracciones continuas de Euler que dan $e^z $ ). Sin embargo, utilizando herramientas de cálculo estándar que trabajan con 16 dígitos de precisión, proporciona una muy buena aproximación sólo hasta aproximadamente $z=15$ . Para valores más altos de $z $ comienza a subestimar significativamente el verdadero valor de $e^z$ . En particular, este error no cambia si aumentamos la longitud de la fracción continua (es decir, si aumentamos $m $ ), lo que apoya la hipótesis de su naturaleza de redondeo.

Por ejemplo, calcular en un $16$ -de precisión, para $z=15\,$ tenemos una buena aproximación (error absoluto $0.0004\,\,\,$ ), pero para $z=20\,$ el error de subestimación es $\approx 4.9\,$ , elevando a $\approx 208585\,\,\,$ para $z=25\,$ y a $\approx 9218010394\,\,\,$ para $z=30\,$ .

Curiosamente, en estos primeros valores, la magnitud del error absoluto parece crecer exponencialmente: véase por ejemplo este ajuste lineal por WA para la relación entre $z $ y el logaritmo del error absoluto de subestimación, calculado para los valores $z=20\,$ , $z=25\,$ y $z=30\,$ (el análisis posterior del logaritmo del error para valores más altos de $z $ sin embargo, parece sugerir una tasa de crecimiento ligeramente inferior).