¿Cuáles son otros métodos más rápidos para calcular $e^x$ hasta cualquier número de decimales aparte de usar la fórmula de la serie de Taylor?
Respuestas
¿Demasiados anuncios?
jwarzech
Puntos
2769
Si quisiera implementar una aproximación de precisión arbitraria a $e^x$ desde "cero" (construyendo sobre una biblioteca numérica de precisión arbitraria como GMP ) con una convergencia más rápida que las aproximaciones de series/polinomios de Taylor, diseñaría alrededor de aproximaciones de fracción continua/función racional .
Así como la precisión dada por la conocida serie de Taylor puede ser extendida según se requiera usando términos adicionales, así también la expansión continua de la fracción dará cualquier precisión requerida usando términos adicionales convergentes .
En un primer momento puede parecer que la evaluación de la aproximación de fracciones continuas conlleva una desventaja significativa en comparación con las series de potencia. Se pueden añadir términos adicionales a las sumas parciales de las series de potencia, para obtener aproximaciones secuenciales de izquierda a derecha. Hay una forma de lograr el mismo efecto con las fracciones continuas, a saber, evaluar los numeradores y denominadores parciales en los convergentes a través de una relación de recurrencia .
Hay un par de ideas importantes que vale la pena considerar incluso si se utilizara otro método para evaluar la función exponencial en un intervalo limitado, como por ejemplo $[0,1]$ como por ejemplo series de poder o incluso tablas de búsqueda. Uno de ellos es Simetría de la función exponencial, y el otro es un método de reducción de rango .
Aunque la función exponencial no es par (o impar), satisface una relación:
$$e^{-x} = 1/e^x$$
que implica que es una simple transformación racional de una función uniforme:
$$e^x = (1 + x*f(x^2))/(1 - x*f(x^2))$$
que puede resolverse para una serie de potencia o una aproximación de fracción continua de f(z). Esta es la idea de Simetría como se aplica a la función exponencial. No sólo reduce la evaluación de la función exponencial en toda la línea real a sólo la mitad positiva (o negativa), sino que también "contrae" a la mitad el número de términos necesarios para una exactitud dada (conservando sólo $x^2$ términos en la expansión).
Aunque una fracción continua o una aproximación de series de Taylor a la función exponencial puede converger para valores reales arbitrarios (o complejos)), la tasa de convergencia no es uniforme. La serie de Taylor converge más rápido cuanto más se acerca el argumento a cero. Así pues, una reducción de rango es importante para expresar la función exponencial a un valor real arbitrario en términos de un valor en algún intervalo como $[0,1]$ . Por ejemplo, si la aritmética de punto flotante es binaria (radix 2), puede ser especialmente conveniente utilizar la conocida ley de los exponentes con múltiplos de $ \ln 2$ :
$$e^x = 2^k * e^r \ \text {where} \ x = r + k* \ln 2$$
que permite $e^x$ para ser evaluada, hasta un cambio en el exponente binario, usando una aproximación que converge rápidamente sobre $[0, \ln 2]$ .
Combinando estas dos ideas (simetría, reducción del rango) la velocidad de convergencia puede limitarse al intervalo $[0, \ln 2/2]$ . Limitar el intervalo de evaluación puede permitirle determinar de antemano cuántos términos de la fracción continuada o de la expansión de la serie de Taylor deben conservarse para obtener una precisión deseada. Cuando esto se conoce, la evaluación puede hacerse más eficientemente (recursión hacia atrás para las fracciones continuas o el método de Horner para las series/polinomios de Taylor truncados) que si nos viéramos obligados a introducir continuamente más términos hasta alcanzar la precisión deseada.
Añadido:
La fracción continua "más rápida" que tenía en mente es la fórmula (11.1.2) aquí:
Manual de fracciones continuas para funciones especiales (Cuyt et al, 2008)]
http://books.google.com/books?id=DQtpJaEs4NIC&dq=exponential+función+continuado+fracción+convergencia
Su referencia es a este clásico moderno:
La aplicación de las fracciones continuas y sus generalizaciones a los problemas de la teoría de la aproximación
A.N. Khovanskii, 1963 (P.Noordhoff)
Un ordenado sitio basado en la webMathematic de Andreas Lauschke presenta algunas ideas para acelerar la convergencia de las fracciones continuas por "contracciones". La dirección IP cambia de vez en cuando y no puede ser usada como un enlace dentro de StackExchange, pero puedes buscarla en Google:
[Aceleración de la convergencia con contracciones canónicas de fracciones continuas: Función exponencial -- Andreas Lauschke Consulting]
Tengo algunas notas sobre sus fórmulas (derivadas por contracción de la señalada anteriormente) si eso sería útil. Se ha aportado algún material conexo a Demostraciones de Wolfram .
Computar la constante $ \ln 2$
Generaciones de estudiantes de matemáticas han sido introducidos al concepto de convergencia condicional versus absoluta por el ejemplo de la series armónicas alternas :
$$ \ln 2 = 1 - 1/2 + 1/3 - 1/4 + ...$$
Por supuesto que esta serie, derivada de una expansión de la serie de potencia de $ \ln x$ alrededor de $x = 1$ tiene una convergencia tan lenta que, incluso si combinamos pares de términos consecutivos:
$$ \ln 2 = 1/2 + 1/12 + 1/30 + ...$$
la serie absolutamente convergente resultante es inútil para obtener valores de precisión arbitrarios de $ \ln 2$ . Por comodidad, las primeras siete sumas parciales de esto son:
0.50000000...
0.58333333...
0.61666666...
0.63452381...
0.64563492...
0.65870518... Desde $ \ln 2$ es 0,69314718..., tenemos un error de alrededor de un tercio de unidad en el primer decimal. En otras palabras, no hay mucha más convergencia que un decimal correcto.
Por lo tanto, hace un sorprendente contraste con la agradable convergencia de una continua expansión de la fracción del mismo valor:
$$ \ln 2 = \cfrac {1}{1 + \cfrac {1}{2 + \cfrac {1}{3 + \cfrac {4}{4 + \cfrac {4}{5 + \cfrac {9}{6 + \cfrac {9}{7 + ...}}}}}}}$$
Los primeros siete convergentes son:
0.66666666...
0.70000000...
0.69230769...
0.69333333...
0.69312169...
0.69315245... El error aquí es aproximadamente media unidad en el quinto decimal.
Johan
Puntos
1007
La forma teóricamente más rápida parece ser utilizar la iteración de Newton para reducir el problema a computar la función del logaritmo y luego utilizar un algoritmo basado en la media aritmética-geométrica para computar el logaritmo. Vea esta página de wikipedia. http://en.wikipedia.org/wiki/Computational_complexity_of_mathematical_operations
En la práctica, la serie de Taylor debería funcionar bien, dada una buena implementación. La siguiente página web http://numbers.computation.free.fr/Constants/constants.html tiene un ejemplo de implementación del uso de la serie Taylor para calcular la e. Afirman que se tardó 0,02 segundos en calcular la e hasta mil decimales en un procesador PentiumII de 450 MHz.
Vincent
Puntos
5027
Alex Andronov
Puntos
178
Una forma muy razonable de hacerlo es usar
$e^x = 1 + x/1! + x^2/2! + x^3/3!...$
ya que converge muy rápido.
Un "truco" que leí en alguna biblioteca para calcular $e^1$ es calcular las series para $x=1/2^k$ donde $k$ es la raíz de la cantidad de bits requeridos, entonces tienes que cuadrar el resultado k veces para obtener tu número final, esto evita que los números grandes/pequeños hagan el resultado demasiado inexacto.
Editar: También este es exactamente el mismo tema como Algoritmo de refinamiento iterativo para computar exp(x) con precisión arbitraria
Tal vez debería ser fusionado / cerrado?
Peter Taylor
Puntos
5221
No sé si esto es necesariamente el camino más rápido, pero es divertido. Si quieres saber cómo funciona, Google.
package com.akshor.pjt33.math;
import com.akshor.pjt33.math.symbalg.Rational;
import java.math.BigInteger;
public abstract class Spigot
{
/\*\* The relevant components of the prefix matrix. \*/
Rational h00, h01=Rational.ZERO, h11=Rational.ONE;
BigInteger i=BigInteger.ZERO;
boolean verify=true;
final int base;
private StringBuilder sb=new StringBuilder();
public Spigot(BigInteger multiplier) {
this(multiplier, 10);
}
public Spigot(BigInteger multiplier, int base) {
h00=new Rational(multiplier, BigInteger.ONE);
this.base=base;
}
/\*\*
\* Advances this spigot algorithm, producing as many digits as each consumption permits.
\* @param beta An upper bound (positive) on the magnitude of all subsequent term ratios.
\* @return The number of digits output.
\*/
public int advance(Rational beta) {
if (verify && beta.signum()0) throw new IllegalArgumentException(Et\_t+" exceeds bound "+beta);
h00=h00.mul(Et\_t);
i=i.add(BigInteger.ONE);
// Produce (if we can: without a converging bound we can't).
int rv=0;
if (beta.compareTo(Rational.ONE)0; i++) numdigits-=advance(beta);
return i;
}
protected abstract Rational termRatio();
public String toString() {
return sb.toString();
}
/\*\*
\* A spigot implementation which evaluates any hypergeometric sum for which you can
\* provide a term bound.
\* E.g. PI = 3 F(1/2, 1, 1, 8/5 ; 3/5, 4/3, 5/3 | 2/27) with term limit 2/27.
\*/
public static class Hypergeometric extends Spigot
{
private final Rational\[\] a, b;
protected final Rational z;
public Hypergeometric(Rational\[\] a, Rational\[\] b, Rational z) {
this(a,b,z,BigInteger.ONE);
}
public Hypergeometric(Rational\[\] a, Rational\[\] b, Rational z, BigInteger mul) {
super(mul);
this.a=a.clone();
this.b=new Rational\[b.length+1\];
this.b\[0\]=new Rational(1,1);
System.arraycopy(b,0,this.b,1,b.length);
this.z=z;
}
protected Rational termRatio() {
Rational rv=z;
for (Rational num : a) rv=rv.mul(num.add(i));
for (Rational num : b) rv=rv.div(num.add(i));
return rv;
}
}
public static class Exp extends Hypergeometric
{
public Exp(Rational z) {
super(new Rational\[0\], new Rational\[0\], z);
}
// We rely on the fact that the term ratio for Exp always decreases.
public int output(int numdigits) {
int i=0;
for (; numdigits>0; i++) numdigits-=advance(termRatio().abs());
return i;
}
}
public static void main(String\[\] args) {
Exp e=new Exp(Rational.ONE);
int count=e.output(100);
System.out.println(e);
System.out.println("Iterations required: "+(count+1));
Exp e2=new Exp(new Rational(2,1));
count=e2.output(100);
System.out.println(e2);
System.out.println("Iterations required: "+(count+1));
Exp em2=new Exp(new Rational(-2,1));
count = em2.output(100);
System.out.println(em2);
System.out.println("Iterations required: "+(count+1));
}
}La implementación de la clase Racional se deja como un ejercicio. Y sí, este código es un argumento para la sobrecarga limitada del operador.
