Exercicio práctico de derivación numérica con Sage

En Sage non existe un comando para calcular a derivada numérica. Porén é moi doado construir unha función que o faga por nós:

#script numerical_derivative.sage
def dn2(x,y): #derivada numerica de segundo orden
    d = vector(RDF, len(y))
    d[0] = (-y[2]+4*y[1]-3*y[0])/(x[2]-x[0])
    d[1:len(y)-1] = [(y[i+1]-y[i-1])/(x[i+1]-x[i-1]) for i in [1..len(y)-2]]
    end = len(y)-1
    d[end] = (3*y[end]-4*y[end-1]+y[end-2])/(x[end]-x[end-2])
    return d

Obsérvese que a fórmula centrada non se pode usar no primeiro e derradeiro punto. Para aproximar a derivada neles usamos as fórmulas cara adiante e cara atrás. Ademais, a orde de erro é O(h²) en todas elas, para obter unha aproximación da mesma precisión en tódolos puntos.
Agora, podemos aproximar o valor da derivada dunha función, por exemplo, do seno no intervalo [0, π/4]:

sage: h = 0.01
sage: x = [0,h..pi/4]
sage: d = dn2(x, map(sin, x))

Se representamos a derivada exacta do seno xunto coa aproximación, veremos a súa similitude:

sage: plot(cos, 0, pi/4, legend_label='cos(x)') + \
....: line(zip(x,d), legend_label='dn2(x,sin(x))', \
....: color='green', marker='+', linestyle='')

dn2-sage

Se calculamos o erro máximo, vermos que é pequeno:

sage: max([abs(cos(i)-j) for i, j in zip(x,d)])
3.333216667877892e-05

E como as fórmulas son da orde O(h²), ao reducir o paso por 10, o erro máximo redúcese por cen:

sage: x = [0,0.001..pi/4]
sage: d = dn2(x, map(sin, x));
sage: max([abs(cos(i)-j) for i, j in zip(x,d)])
3.33333216806508e-7
Advertisements
Esta entrada foi publicada en Cálculo numérico, Derivación numérica, Matlab/Octave. Ligazón permanente.

2 Responses to Exercicio práctico de derivación numérica con Sage

  1. Pingback: Derivación numérica | fundmat

  2. Pingback: Exercicios das clases prácticas | fundmat

Deixar unha resposta

introduce os teu datos ou preme nunha das iconas:

Logotipo de WordPress.com

Estás a comentar desde a túa conta de WordPress.com. Sair / Cambiar )

Twitter picture

Estás a comentar desde a túa conta de Twitter. Sair / Cambiar )

Facebook photo

Estás a comentar desde a túa conta de Facebook. Sair / Cambiar )

Google+ photo

Estás a comentar desde a túa conta de Google+. Sair / Cambiar )

Conectando a %s