Práctica 3a - Resolución de sistemas lineales. Métodos directos
Contents
1. A Ejercicio 6 Calcular las normas l l y l2 del vec tor x 1 1 2 de R y las normas l y l de la matriz 1 2 1 A 0 3 1 5 1 1 Ejercicio 7 Escribir sentencias Fortran para calcu lar las normas vectoriales lo ly y l2 y las normas ma triciales lo y l Ayuda considerar el uso de las funcio nes intr nsecas SUM ABS MAXVAL y DOT_PRODUCT Estimaci n a posteriori del error N mero de condici n En la pr ctica una soluci n num rica X de un sistema de ecuaciones lineales Ax b tendr cierto error x X el cual por supuesto no puede ser calculado Ahora bien siempre podemos calcular el vector residuo r b AX Puesto que r 0 implica que X es la soluci n exacta del sistema resulta natural sugerir que si r es peque a entonces X es una soluci n aproximada precisa Sin embargo esto no es siempre el caso La exactitud de la soluci n calculada depende tambi n del n mero de condici n k A AINATH de la matriz A del sistema de acuerdo a la siguiente desigualdad Ix Ixl Esta desigualdad nos dice que el error relativo en una soluci n num rica puede ser tan grande como k A veces su residuo relativo De este modo solo si A 1 el error relativo y el residuo relativo son del mismo tama o y el residuo podr ser utilizado con seguridad como una estimaci n del error En esta situaci n se dice que la matriz A del sistema es una matriz bien condicionada Si por el c2. An lisis Num rico I PRACTICA 3a Resoluci n de sistemas lineales M todos directos Se abre el tel n y se ven dos sistemas lineales incompatibles C mo se llama la obra Kramer vs Kramer Eliminaci n gaussiana Factorizaci n LU Sea Ax b un sistema de n ecuaciones lineales con n inc gnitas donde A es una matriz cuadrada de orden n de elementos reales y no singular por lo cual el sistema admite una soluci n y sta es nica Los m todos num ricos directos para resolver dicho sistema son aquellos que en ausencia de errores de redondeo obtienen la soluci n exacta en un n mero finito de pasos El m todo fundamental es el procedimiento de eliminaci n gaussiana con la estrategia de pivo teo parcial la cual permite el intercambio de filas Dicho m todo muestra que la matriz A admite una factorizaci n de la forma A PLU Aqu L es una matriz triangular inferior con elemen tos diagonales iguales a la unidad y cuyos elementos bajo la diagonal son los multiplicadores utilizados durante la eliminaci n U es una matriz triangular su perior cuyos elementos son los coeficientes del sistema final equivalente siendo los elementos de la diagonal no nulos Si los intercambios de fila durante el pro ceso de eliminaci n son registrados en un vector p cuyo valor inicial es 1 2 n entonces la matriz de permutaci n P se obtiene a partir de la matriz identidad Z de orden n permutando sucesivamente la co
3. LAPACK y su interfaz LAPACK95 como sigue lPara otros sistemas el comando puede variar gfortran mall o programa programa f90 A I usr local include lapack95 llepaek llapackos5 Ejercicio 3 Utilice las rutina LA_GESV para resol ver los sistemas A X B para las matrices de coefi cientes A dadas por 0 1 2 2 1 0 a 1 2 3 b 0 3 2 2 3 2 1 2 4 1 2 1 c 2 0 2 1 2 3 y t rminos independientes dados por 2 4 1 B 4 12 0 5 17 1 Explicitar la factorizaci n PLU de las matrices A Matrices especiales Cuando la matriz A del sis tema tiene ciertas propiedades especiales es posible re ducir ya sea el costo computacional del procedimiento de eliminaci n gaussiana el costo de almacenamiento de la matriz o ambos En primer lugar para algunas clases de matrices la eliminaci n gaussiana es num ricamente estable s n ninguna estrategia de pivoteo a saber cuando A es diagonalmente dominante por columnas o filas esto es z laj l gt X la para cada 1 2 n 14 Jail gt X la para cada 1 2 n IA A es sim trica A A y definida positiva esto es xt Ax gt 0 para todo x 0 Por otra parte para otras matrices conocidas como matrices de banda la presencia de un gran n mero de elementos nulos en un patr n regular permite sim plificar el algoritmo de factorizaci n Entre ellas se encuentran en particular las matrices tridiagonales Estas variantes ser n
4. discutidas a continuaci n Matrices sim tricas definidas positivas M to do de Cholesky Si la matriz de coeficientes A es sim trica y definida positiva existe una factorizaci n de la forma A L Lt donde L es una matriz trian gular inferior con elementos positivos en la diagonal Dicha factorizaci n determinada por el m todo de Cholesky no requiere intercambio de filas y es estable frente a los errores de redondeo Esta factorizaci n es utilizada para resolver el sistema e implementa da en LAPACK95 por la subrutina LA_POSV donde A es almacenada en un arreglo bidimensional como es usual aunque por la simetr a solo se accede al tri ngulo superior o inferior Pr ctica 3a An lisis Num rico I Ejercicio 4 Resolver los sistemas de ecuaciones AX B donde 3 48 6 7 7 A 4 6 52 1 6 6 7 1 56 10 5 7 6 10 74 56 112 168 71 142 213 B 69 138 207 80 160 240 102 204 306 Explicitar la factorizaci n de A Matrices tridiagonales Las matrices tridiagona les tienen todos sus elementos nulos por fuera de la diagonal principal y la subdiagonal y superdiagonal adyacentes En tal caso los elementos relevantes de la matriz pueden ser entonces almacenados en tres arreglos unidimensionales independientes uno por cada diagonal Esto permite que s lo se requiera un espacio de almacenamiento O 3n en vez del arreglo n completo y adem s el costo computacional de la factorizaci n resulta ser O 9n lo cual implica q
5. ente Por ejemplo si como haremos trabajamos en doble precisi n la sentencia apropiada es USE la _precision ONLY WP gt DP F95_LAPACK quien define las interfaces expl citas de las subrutinas a utilizar Por ejemplo si vamos a utilizar la subrutina LA_GESV la sentencia apropiada es USE F95_LAPACK ONLY la_gesv La invocaci n de la subrutina LA_GESV procede en tonces como sigue donde indicamos con corchetes los argumentos opcionales CALL la_gesv A B IPIV p INFO info0 Como dato de entrada A es un arreglo bidimensional que almacena la matriz de coeficientes del sistema y B un arreglo unidimensional o bidimensional que contiene los t rminos independientes La subrutina devuelve en A los factores L y U de la factorizaci n y en B las soluciones de los sistemas Opcionalmente devuelve tambi n el vector de permutaciones p con el cual puede construirse la matriz de permutaci n P que completa la factorizaci n y una clave entera de error info cuyo valor igual a cero implica que la resoluci n del sistema se efectu con xito y otros valores alg n tipo de error 3 Una discusi n completa de los argumentos necesarios para invocar cada subrutina de LAPACK95 est dada en el manual de usuario el cual puede ser consultado on line en http www netlib org lapack95 1ug95 Finalmente para compilar el programa fuente de bemos informar explicitamente al compilador que utilizaremos la biblioteca de rutinas
6. inaci n gaussiana en vez de de programar nuestras propias subrutinas utilizare mos la biblioteca de rutinas LAPACK Linear Algebra PACKAage la cual es una colecci n de subrutinas para resolver num ricamente problemas matem ticos que se enmarcan en el campo del lgebra lineal junto con su interfaz LAPACK95 la cual explota las carac ter sticas propias de Fortran 95 asignaci n din mica de memoria arreglos de tama o ajustable en la lista Pr ctica 3a An lisis Num rico I de argumentos y argumentos opcionales Para la re soluci n de sistemas lineales LAPACK95 proporciona un amplio conjunto de subrutinas dependiendo de la naturaleza particular de la matriz de coeficientes En particular la subrutina LA_GESV permite resolver un sistema de ecuaciones lineales para una matriz A general almacenada en un arreglo bidimensional y uno o varios t rminos independientes almacenados en un arreglo bidimensional B seg n el esquema discuti do La cuesti n de utilizar las rutinas de LAPACK95 consta de dos partes c mo llamar a las subrutinas en nuestro programa y como compilar el mismo A continuaci n discutimos stos dos puntos El uso de las rutinas de LAPACK95 requiere de la invocaci n de dos m dulos a trav s de la sentencia USE LA PRECISION quien define la precisi n de los tipos de datos reales ya sea de simple o doble preci si n via el par metro WP cuya asignaci n ser SP DP respectivam
7. l sistema 0 4 1 m 9 1 illa sl 6 2 2 A E 1 Ejercicio 2 El siguiente ejercicio muestra que la es trategia de pivoteo parcial es crucial para mejorar la estabilidad num rica del procedimiento de eliminaci n gaussiana El sistema 1 566 z 1 569 2 436 z2 11 018 tiene la soluci n exacta x 10 1 Utilice aritm tica decimal de cuatro d gitos para resolverlo sin y con pivoteo parcial Indique el porque de la diferencia en las soluciones num ricas obtenidas 0 0003 0 3454 IS La estabilidad num rica del m todo de eliminaci n gaussiana depende del tama o del llamado factor de cre cimiento definido como la raz n de las magnitudes del mayor coeficiente de U al mayor coeficiente de A Sin pi voteo este factor puede ser arbitrariamente grande y por lo tanto la eliminaci n gaussiana sin pivoteo es num rica mente inestable Sin embargo a n con pivoteo parcial este factor puede ser tan grande como 2 el cual es enorme a n para modestos valores de n y por lo tanto el m todo estrictamente hablando es num ricamente inestable Sin embargo la experiencia ha mostrado que ste compor tamiento del factor de crecimiento es extremadamente raro En la pr ctica este factor se mantiene moderado En consecuencia el m todo de eliminaci n gaussiana con pivoteo parcial puede considerarse num ricamente estable en la pr ctica Para la implementaci n en una computadora del procedimiento de elim
8. lumna por la columna j p i parai 1 2 n Por ejemplo si para n 4 p 3 4 3 4 entonces SO 0 0 SO 0 1 0 0 ooo Conocida la factorizaci n el sistema de ecuaciones lineales es equivalente a PLUx b el cual conduce a dos sistemas triangulares Ly Ptb Ux y As para obtener el vector soluci n x se resuelve en primer lugar el primer sistema por sustituci n hacia adelante para obtener y y luego se resuelve el segundo sistema por sustituci n hacia atr s IS La factorizaci n LU de una matriz de orden n es un procedimiento finito de n 1 pasos cuyo costo compu tacional medido por el n mero de operaciones aritm ticas necesarias llevarla a cabo es O n En tanto que el costo computacional de la resoluci n de los dos sistemas trian gulares es O n el cual es un orden de magnitud menor que el costo de la factorizaci n Notar tambi n que en la factorizaci n el costo computacional de la determinaci n de los pivotes requiere de O n comparaciones el cual es entonces despreciable frente al costo de las operaciones aritm ticas IS Si tenemos varios sistemas de ecuaciones con la misma matriz de coeficientes A AX bi AX2 ba AXm bm los mismos pueden ser tratados simult neamente como el sistema matricial AX B donde X x1 x2 xm y B bi1 b2 bm son matrices rectangulares n x m Ejercicio 1 Resolver a mano por eliminaci n gaus siana con pivoteo parcial e
9. ontrario k A gt 1 el intervalo en que se puede situar el error relativo es muy amplio y por lo tanto no se puede asegurar que la soluci n num rica es precisa a n cuando el residuo sea peque o En esta situaci n se dice que la matriz A es una matriz mal condicionada Ejercicio 8 Mostrar que para toda matriz no singu lar k 4 gt 1 Pr ctica 3a An lisis Num rico I Ejercicio 9 El ejemplo can nico de matrices mal condicionadas son las matrices de Hilbert La matriz de Hilbert H de orden n est definida por my 1 n A ee P y es una matriz sim trica definida positiva El siguien te fragmento de c digo Fortran permite calcular su inversa T H t 1 1 REAL n n WP i 1 DO 3 2 n ela E I e o BRABIIAS Sl ae darle po REAL 1 43 1 3 1 3 1 WP END DO DO i 2 n DO j 1 n AE aeS e E E BEARD WEE REAL 1 3 1 x 1 1 1 1 WP END DO END DO Utilice dicho c digo para calcular el n mero de con dici n K de las diez primeras matrices de Hilbert Ejercicio 10 Sea Ax b el sistema lineal dado por 0 89 0 53 0 36 ES os a EA o i Dada la soluci n aproximada 11 5 20 calcule el vector residual y estime el error relativo Sabiendo que la soluci n exacta es x 1 1 calcule el error exacto Compare los resultados Estimaci n a priori del error An lisis pertur bativo En la secci n anterior hemos supuesto que A y b est n libres de e
10. orrespondiente est acotada aproximadamente en la norma l por x Ixlloo lt 2uk 4 donde u es la unidad de redondeo de la m quina y suponemos que A est bien condicionada Pr ctica 3a
11. rror Sin embargo si el sistema lineal proviene de la modelizaci n de un problema f sico debemos esperar que los coeficientes del sistema est n sujetos a error bien porque provienen de otros c lculos o de mediciones Adem s a n por la sola causa de la representaci n de los n meros en pun to flotante en la computadora A y b se encuentran afectados de errores El n mero de condici n juega tambi n aqu un papel importante en el an lisis de tales errores Ejercicio 11 Dado el sistema Ax b suponga que el t rmino b se perturba por una cantidad db Mostrar que el efecto de esta perturbaci n sobre la soluci n ser x x con x acotada por 16 lt y l xl A Ton Si ahora es la matriz A la que es perturbada por una cantidad 04 mostrar que el efecto de tal perturbaci n sobre la soluci n ser x x con x acotada por 19x 19 4 ETA Ejercicio 12 Sea 6 13 17 A 13 29 38 17 38 50 Si se desea resolver el sistema Ax b donde b 1 9358 2 3412 3 5718 tiene sus elementos redondeados a 5 d gitos qu puede decirse del error relativo en la soluci n En general si A y b est n perturbados de modo que A 0A x x b db entonces 9x A Cea E lx 1 CASAN AD IAL bll 7 siempre que 194 14 lt 1 Ejercicio 13 Mostrar que si A y b son conocidos exactamente al ser representados como n mero de puntos flotante en la computadora la perturbaci n x c
12. ue el procedimiento es muy eficiente respecto al caso ge neral LAPACK95 proporciona la subrutina LA_GTSV para resolver sistemas de ecuaciones lineales con ma trices de coeficientes tridiagonales almacenados de este modo Ejercicio 5 Resolver los sistemas de ecuaciones AX B donde 220000 154000 048400 A ag g e 000643 000049 y 4 8 12 10 20 30 16 32 48 B gt 40 60 13 26 39 13 26 39 Normas vectoriales y matriciales Toda solu ci n num rica de un sistema lineal Ax b debe considerarse como una soluci n aproximada X debido a los errores de redondeo introducidos en el c mputo de la misma Para medir cuan buena es la aproxima ci n X a la soluci n exacta x se recurre al concepto de norma vectorial La distancia entre los vectores X y x es entonces definida como la norma de la di ferencia de los mismos lx X En la pr ctica las normas vectoriales en R de uso frecuente son siendo X 21 2 Zn un vector de R a norma l xl m xi lt i lt n1 25 norma ly x Xalil norma la x 2 Pret ae E Asimismo se define una norma matricial natural su bordinada a la respectiva norma vectorial siendo entonces para A aij una matriz de R norma loo Alloo m x1 lt i lt n 110451 norma h A m xi lt jcn 110551 norma la All 4 A donde el radio espectral de una matriz A de n x n se define como p A m x A A autovalor de
Download Pdf Manuals
Related Search