Laboratorio Métodos Numéricos Diferenciación Numérica Objetivo 1. Analizar los errores incurridos al realizar diferenciaciones numéricas simples: diferencia diferencia hacia atras, adelante y central. 2. Entender Entender los efectos efectos que el “ruido” “ruido” produce sobre las derivadas numéricas. numéricas.
Procedimiento Escriba el siguiente guión: // Este Este scri script pt sirv sirve e inves investi tiga gar r el uso uso de las las deri derivad vadas as numer numeric icas as sobr sobre e datos datos // """"muestreados """"muestreados"""" """" // Empe Empeza zamo mos s por por gene genera rar r punt puntos os para para la func funcio ion n y = x *cos(x) h = 0.1; // paso paso de mues muestr treo eo x = 0:h: 0:h:4; 4; sdev sdev = 0.1; 0.1; y = x .*cos(x) cos(x) + rand(x, rand(x,"no "norma rmal") l") *sdev; // La deri deriva vada da exac exacta ta de esta esta func funcio ion n es faci facil l de hall hallar ar: : // y’’ = cos(x) - x *sin(x) dyexa dyexact cta a = cos( cos(x) x)-x -x . *sin(x); // Ahora Ahora calc calcul ulamo amos s la deriv derivad ada a usan usando do la dife difere renc ncia ia finit finita a haci hacia a atras atras // Busc Buscam amos os el larg largo o de x prim primer ero: o: nl = max(siz max(size(x) e(x)); ); // Aplic Aplicam amos os la dife diferen renci cia a haci hacia a atra atras. s. Notes Notese e que que no es posib posible le (sin (sin exten extende der r los los dato dato // deriv derivar ar el prim primer er punto punto. . Saca Sacamo mos s la divis divisio ion n del del ciclo ciclo para para hace hacerl rlo o mas mas efici eficien ente. te. for i = 2:n 2:nl dyatras dyatras(1, (1,i) i) = y(i)-y( y(i)-y(i-1 i-1); ); end; dyatras dyatras = dyatra dyatras/h; s/h;
// Aplicamos la diferencia hacia adelante ahora. Ahora no es posible encontrar la // derivada en el ultimo punto. Teniendo el resultado anterior, existe una manera // mas eficiente de hacer lo que sigue? for i = 1:nl-1 dyadelante(1,i) = y(i+1)-y(i); end; dyadelante = dyadelante/h; // Ahora encontramos la diferencia central. Ni el primer punto ni el final pueden // ser hallados. for i = 2:nl-1 dycentral(1,i) = y(i+1)-y(i-1); end; dycentral = dycentral/(2 *h); // Como las derivadas no fueron calculadas en todos los puntos y queremos graficarlas // todas para compararlas a las reales, entonces vamos a """"rellenar"""" los puntos fa // con los datos verdaderos. dyatras(1) = dyexacta(1); dyadelante(nl) = dyexacta(nl); dycentral(1) = dyexacta(1); dycentral(nl) = dyexacta(nl); // Ahora podemos graficar y ver las diferencias con respecto al resultado // ideal. scf(1) clf(1) plot(x,dyexacta,"y",x,dyatras,"r",x,dyadelante,"g",x,dycentral,"b"); xtitle("Comparacion de las derivadas","x","dy/dx"); legend("Exacta","Atras","Adelante","Central"); // Ahora podemos graficar los errores. scf(2) clf(2) plot(x,dyatras-dyexacta,"r",x,dyadelante-dyexacta,"g",x,dycentral-dyexacta,"b"); xtitle("Errores de las derivadas","x","Error dy/dx"); legend("Atras","Adelante","Central"); // Grafique la funcion real y la que tienen la bulla scf(3) clf(3) plot(x,x . *cos(x),"g",x,y,"r"); xtitle("Efecto del ruido gaussiano") legend("Sin ruido","Con ruido");
Ruido y su efecto Ahora añadiremos “ruido” a la “señal” y . Asumiremos que el ruido es gaussiano (tiene una distribución normal) con valor esperado de 0 (promedio cero) y desviación estándard 0,1. Para hacerlo usamos la instrucción rand(N,’normal’) que genera un vector de dimensión N aleatorio con distribución normal unitaria. Multiplicamos el mismo por 0.1 para cambiar su desviación estándar. Modificamos la linea que calcula el vector y , añadiendo el vector aleatorio: y = x .*cos(x) + rand(x,’normal’) *sdev;
y añadimos al final del programa instrucciones para visualizar los resultados.
Preguntas 1. ¿Qué método de diferenciación resulta mas “exacto” para este caso? 2. Compare los errores de la derivada en la presencia de ruido con los errores en y . ¿Qué efecto tuvo el ruido sobre la derivada?
