MAESTRÍA EN INVESTIGACIÓN E INNOVACIÓN DE TECNOLOGÍAS DE LA INFORMACIÓN Y COMUNICACIÓN LABORATORY 1: NUMERICAL AND DATA-INTENSIVE COMPUTING (COURSE 2012/2013) SYSTEM OF LINEAR EQUATIONS PROFESOR: DR. CARLOS SANTACRUZ José Luis Carrillo Medina Cesar Naranjo Hidalgo
1. Ejercicios 1) En matemáticas la propiedad propiedad asociativa asociativa de la suma hace que el orden en la cual las operaciones se realizan no importa siempre y cuando la secuencia de operaciones no cambie. Sin embargo esto no es verdadero cuando se implementa en un computador, esto se debe a su representación en punto flotante. Un ejemplo donde el redondeo produce el error y no cumple la propiedad asociativa de la es: Consideremos la siguiente serie infinita
Diferentes resultados se obtienen, si la serie es sumada en orden ascendente y descendente: PROGRAMA
function y=seriepi(n) function y=seriepi(n) format long long; ; j = 1; while (j while (j <= n) auxAsc = 0; auxDes = 0; i = 1; while (i while (i <= j) auxAsc = auxAsc + 1/i^4; auxDes = auxDes + 1/(j-i+1)^4; 1/(j-i+1)^4; i = i + 1; end Res(j) = auxAsc - auxDes; j = j + 1; end j=1:1:n; plot (j,Res)
Cuando n = 10000, obtenemos
Preguntas: a) Porque la propiedad asociativa no cumple?
En una computadora el almacenamiento de un número solo puede hacerse con una cantidad finita de bits. Esta cantidad está determinada por la máquina en la cual se hará la representación. El número de bits generalmente se llama palabra y estas van desde ocho bits hasta 64. Generalmente una palabra almacena un número, sin embargo a veces es necesaria mas de una palabra para almacenar ciertos números. Por ejemplo, si se quiere almacenar números enteros en una palabra de 16 bits, el primero de estos representa el signo del número (un cero es signo más y un uno un signo menos). Los 15 bits restantes pueden usarse para guardar números binarios de 0 a 32767 (2 15 -1). Es decir una palabra de 16 bits puede contener un número cualquiera entre -32768 a +32767. Para número reales se emplea la representación binaria llamada de punto flotante: 0. d 1d2d3d4d5d6d7d8x2d’ 1d’ 2d’ 3d’ donde di con i = 1…8 representan la mantisa, y d’ j con j= 1...7 la característica. El primer bit del exponente representa el signo de este. La propiedad asociativa no se cumple por lo anteriormente indicado ya que los computadores representan los números puntos flotantes (con muchos decimales), por medio de aproximaciones, a través del redondeo y/o truncamiento que realiza, donde redondear un número
quiere decir reducir el número de cifras manteniendo un valor parecido. Los tipos de datos float y double tienen defectos cuando representan resultados de ciertas operaciones aritméticas entre números reales con total precisión, en algunos casos solo son capaces de representar aproximaciones a ellos ya que las operaciones se realizan en base dos, al igual que nosotros que no podemos representar 1/3 con números en base diez con total exactitud, por lo que usarlos para hacer ciertas operaciones aritméticas, así como acumularlas, puede dar lugar a errores de redondeo y precisión en el resultado final, por lo que usar este tipo de datos no son los más adecuado para trabajar en una aplicación que calcule precios. El error de redondeo es el que resulta de reemplazar un número por su forma de punto flotante, es decir, por su representación en una máquina concreta. b) Que orden es más seguro? Porque? Es más seguro hacer los cálculos en orden ascendente debido a que los resultados iniciales parciales son los que aportan más significativamente al resultado y los errores de aproximación no se ven reflejados tan profundamente como los que se llevaría acabo si se procesa el resultado en orden descendente, es decir desde un inicio se puede empezar con redondeos y/o truncamientos, por lo que requieren de una mejor representación.
c) Escriba una función llamada sumSeriesA.m con el número de términos a ser tomados en cuenta en el argumento. La salida de la función debería ser la primera de n términos sumados en orden ascend**ente. (Prototipo function n= sumSeriesA(n)). function n=sunSeriesA(nIt er) format long; j = nIter; while (i <= j) auxAsc = auxAsc + 1/i^4; i = i + 1; end
4
%Calculo de la serie 1/i %Orden Ascendente %Definir variables tipo entero %largo
%Acumulador de sumas parciales
n = auxAsc;
d) Escriba la función llamada sumSeriesD.m pero en este caso sume los términos en orden descendente. function n = sunSeriesD(nIte r) format long; i = 1; j = nIter; while (i <= j) auxDes = auxDes + 1/j^4; j = j - 1;
%Calculo de la serie 1/i 4 %Orden Ascendente %Definir variables tipo entero %largo
%Acumulador de sumas parciales
end sumaD = auxDes:
e)
Escriba un programa en matlab usando funciones para comparar los resultados. function y=sumSeriesPI(n )
%Cálculo de suma de la serie PI
format long; j = 1; auxAsc = 0; auxDes = 0; while (j <= n) auxAsc = sumSeriesA(j); auxDes = sumSeriesD(j); Res(j) = auxAsc - auxDes; j = j + 1; end j=1:1:n; plot (j,Res)
%Cálculo en orden ascendente de la %suma de la serie PI %Cálculo en orden descendente de la %suma de la serie PI %Diferencias entre cálculos en %orden ascendente y descendente
f) Encuentre el valor de n que se obtiene de diferentes resultados. Expliqué porque se obtiene estos resultados. No. Iteraciones 100
Figura
Observación La función:
Al comparar la serie calculada en forma ascendente y descendente se puede indicar que en el rango de 0 a 12, no habiendo errores de redondeo, debido a que la capacidad de representación de puntos flotantes de la memoria contiene los valores generados. A partir de 13, se produce un error, debido a que la capacidad de representación de punto flotante del computador no lo puede representar y asigna valores próximos, lo que produce errores de redondeo o truncamiento, inesperados.
500
Existen ruidos propios del sistema, como que se trunca en el orden ascendente y en descendente.
1000
Lo mismo que en el anterior caso
5000
Empieza la inestabilidad sistema en 3500
del
10000
Se puede indicar que el sistema es aceptable hasta los 5000, en donde las variabilidades son grandes
15000
A partir de 5000 se puede indicar que el sistema no es estable, las variabilidades son muy grandes
De lo anterior se puede indicar que de acuerdo a las graficas obtenidas se representa una señal oscilatoria y las señales oscilatorias no cumplen con la propiedad asociativa.
Consejos: a) Imprima los resultados con 20 dígitos significativos b) Grafique las diferencias entre ambos valores obtenidos c) Grafique los diferencias entre cada valor y el real Para graficar las diferencias utilizamos la siguiente función: function y=compararPIconsumSeriePI(n) format long;
j = 1; auxAsc = 0; auxDes = 0; while (j <= n) ResAsc(j) = (pi^4/90)-sumSeriesA(j); ResDes(j) = (pi^4/90)-sumSeriesD(j); j = j + 1; end ResAsc(j-1) ResDes(j-1) j=1:1:n; subplot(1,2,1); plot (j,ResAsc); subplot(1,2,2); plot (j,ResDes);
(pi^4/90)-sumSeriesA(j); %Calculo del error en Orden Ascendente %Error en 20 iteraciones 4.496830281475184e-007 (pi^4/90)-sumSeriesA(j); %Calculo del error en Orden Descendente %Error en 20 iteraciones 4.496830279254738e-007
De estos gráficos se puede indicar que el valor real con respecto a la serie calculada tiende a converger hacia el infinito, es decir mientras más iteraciones converge la serie, tanto en orden ascendente como descendente.
2) Analíticamente las expresiones pueden no ser equivalentes cuando son evaluadas en un computador. Como ejemplo de esta situación. Consideremos la ecuación cuadrática:
Donde las raíces de esta ecuación están dadas por la siguiente expresión equivalente: (1)
(2)
En estas ecuaciones la suma se hace delicada y se prueba el error de redondeo cuando b>0 y sobre la condición |ac| << b2. Si cualquiera de los dos a o c son pequeños, entonces una de las raíces implicara la resta de b a partir de una cantidad igual (El discriminante). El correcto camino para calcular las raíces es:
(3)
Las raíces son:
Preguntas: a) Muestre la equivalencia, entre los tres expresiones previas Dada la solución de un sistema de ecuación cuadrática (1):
Racionalizamos multiplicando el numerador y dividiendo para el mismo: Multiplico por el conjugado
Aplico la ley distributiva entre términos Hago la resta de terminos
De (1) obtenemos (2):
Para demostrar la siguiente equivalencia tenemos que la función signo de un número real es sgn(x), (En MatLab es sign(x)). Donde su dominio de definición es R y su conjunto imagen es: {-1;0;1}.
Remplazó:
por
b) Escriba la función cuadratic1.m, quadratic2.m y quadratic3.m con los parámetros a, b y c como argumentos y las raíces usando diferentes expresiones de salida.
function [x1,x2]=cuadratic1(a,b,c) d = sqrt(b*b - 4*a*c); x1 = (-b + r)/(2 * a); x2 = (-b - r)/(2 * a); end
function [x1,x2]=cuadratic2(a,b,c) d = sqrt(b*b - 4*a*c); x1 = (2 * c)/(-b - r); x2 = (2 * c)/(-b + r); end
function [x1,x2]=cuadratic3(a,b,c) d = sqrt(b*b - 4*a*c); q = -1/2 * (b + sgn(b) * r); x1 = c/q; x2 = q/a; end
c) Encuentre un ejemplo donde los resultados sean diferentes. Explique los resultados. Se realizo la prueba con datos arbitrarios para obtener las raíces de la ecuación: Ax2 + Bx + C Probando los resultados con A = 1; B = 100; C = 10 >> probar(1,100,10) EC:(1):
x1 = -0.100100200501402
EC:(2):
y1 = -0.100100200501404
EC:(3):
z1 = -0.100100200501404
x2 = -99.899899799498598 y2 = -99.899899799501242 z2 = -99.899899799498598
De los resultados obtenidos y pruebas realizadas se puede indicar que vasta que el coeficiente b sea mucho mas grande que la multiplicación de a * b se obtienen raíces desiguales.
d) Cuál es el resultado correcto? y Porque?. El resultado correcto es de la función cuadratic1, en la cual no existen divisiones de cantidades muy pequeñas, en las ecuaciones cuadratic2 y cuadratic3 existen divisiones pequeñas, además si b 2 es mucho mayor que 4ac, b2 >> 4ac, por lo que la resta entre -b y la raíz seria de números casi iguales lo que haría que el error tienda a cero.