Si su función es conocida por ser dos veces diferenciable, utilice
f'(x) = (f(x + h) - f(x - h))/2h
que es de segundo orden preciso en h. Si solo es diferenciable una vez, use
f'(x) = (f(x + h) - f(x))/h (*)
que es de primer orden en h.
Esto es teoría. En la práctica, las cosas son bastante complicadas. Tomaré la segunda fórmula (primer orden) ya que el análisis es más simple. Haz el segundo orden uno como ejercicio.
La primera observación es que debe asegúrese de que (x + h) - x = h, de lo contrario obtendrá grandes errores. De hecho, f (x + h) yf (x) están cerca uno del otro (digamos 2.0456 y 2.0467), y cuando los restas, pierdes muchas cifras significativas (aquí es 0.0011, que tiene 3 cifras significativas menos que x). Entonces, cualquier error en h es probable que tenga un gran impacto en el resultado.
Así que, primer paso, arregle un candidato h (le mostraré en un minuto cómo elegirlo), y tome como h para su cálculo la cantidad h '= (x + h) - x. Si está utilizando un lenguaje como C, debe tener cuidado de definir h o x como volátil para que ese cálculo no se optimice.
A continuación, la elección de h. El error en (*) tiene dos partes: el error de truncamiento y el error de redondeo. El error de truncamiento es porque la fórmula no es exacta:
(f(x + h) - f(x))/h = f'(x) + e1(h)
donde e1(h) = h/2 * sup_{x in [0,h]} |f''(x)|
.
El error de redondeo proviene del hecho de que f (x + h) yf (x) están cerca uno del otro. Se puede estimar aproximadamente como
e2(h) ~ epsilon_f |f(x)/h|
donde epsilon_f
es la relativa precisión en el cálculo de f (x) (o f (x + h), que está cerca). Esto debe evaluarse a partir de su problema. Para funciones simples, epsilon_f
se puede tomar como la máquina épsilon. Para los más complicados, puede ser peor que eso en órdenes de magnitud.
Así que usted quiere h
que minimice e1(h) + e2(h)
. Conectar todo juntos y la optimización de los rendimientos en h
h ~ sqrt(2 * epsilon_f * f/f'')
que ha de estimarse a partir de su función. Puede tomar estimaciones aproximadas. En caso de duda, tome h ~ sqrt (epsilon) donde epsilon = precisión de la máquina. Para la elección óptima de h, la precisión relativa a la que se conoce la derivada es sqrt (epsilon_f), es decir. la mitad de las cifras significativas son correctas.
En resumen: demasiado pequeño a h => error de redondeo, demasiado grande a h => error de truncamiento.
Para la segunda fórmula orden, mismos rendimientos de cálculo
h ~ (6 * epsilon_f/f''')^(1/3)
y una precisión fraccionaria de (epsilon_f)^(2/3) para el derivado (que es típicamente una o dos cifras significativas mejor que el primero fórmula de orden, suponiendo doble precisión).
Si esto es demasiado impreciso, no dude en pedir más métodos, hay muchos trucos para obtener una mayor precisión. La extrapolación de Richardson es un buen comienzo para funciones suaves. Pero esos métodos generalmente se calculan bastantes veces, esto puede ser o no lo que usted desea si su función es compleja.
Si va a utilizar derivadas numéricas muchas veces en diferentes puntos, resulta interesante construir una aproximación de Chebyshev.
+1 por una buena respuesta :) – posdef
Estoy en lo cierto al pensar que cuanto más pequeño elijo dx, más precisa será mi respuesta. Si es así, ¿hay una constante MATLAB para un número real MUY pequeño? Algo como pi para 3.14 ... o i para sqrt (-1)? – lms
@codenoob: 'eps' te da un número muy pequeño. Sin embargo, para la mayoría de los propósitos prácticos, 0.0001 debería ser suficiente. Además, si trabajas con polinomios y usas 'polyder', no tienes que preocuparte por el tamaño de' dx'. – Jonas