el método de cuasi minimización de residuos.
• bicgstab: aplica el método de gradiente biconjugado estabilizado.
Todas estas rutinas permiten reiniciar las iteraciones, incorporar criterios de tolerancia, número máximo de iteraciones, etcétera. Para mayores detalles de cada método, consultar por ejemplo Axelsson [4] y la ayuda en Matlab®.
1.3 Análisis del error
Consideremos la posibilidad de que haya errores en los datos de entrada del problema y, por lo tanto, necesitamos estimar la magnitud del error en la solución de las ecuaciones. Para estos fines, la matriz inversa A-1 es fundamental en la determinación del error de la solución.
Cuando los pequeños errores en la matriz de coeficientes A producen grandes cambios en la solución, entonces se dice que el sistema de ecuaciones está mal condicionado. Esto normalmente se debe a que las filas de la matriz A son casi linealmente dependientes entre sí.
Consideremos el caso cuando haya error en el lado derecho de las ecuaciones, es decir:
y cuya solución es x+δx. Luego se tiene que
y además se cumple
Por lo tanto, de las ecuaciones anteriores se deduce que una cota superior para el error relativo viene dada por (Moler, [1])
En donde κ(A) es el número de condición2 de la matriz A. Si κ >> 1, se dice que la matriz A está mal condicionada. Por ejemplo, si κ(A) es 108, un error relativo de 10-8 en el lado derecho, ¡podría producir un error tan grande como 100% en el error relativo de la solución!
Un punto importante que conviene destacar aquí es que la característica de mal condicionamiento de un matriz a menudo proviene del mismo usuario que la genera en forma inadvertida. Una manera fácil de producir matrices mal condicionadas consiste en tener ecuaciones con coeficientes con órdenes de magnitud distintos entre una ecuación y otra. No existen procedimientos numéricos que nos garanticen que una matriz no presente este problema.
Una solución que se puede plantear cuando el sistema es mal condicionado, es aplicar alguna de las rutinas iterativas que se han mencionado en la sección 1.2, ya que así no se utiliza el procedimiento de inversión de la matriz (que es el que genera normalmente los problemas de mal condicionamiento y ambigüedad de los resultados), sino un método iterativo que no trabaja con la inversa de la matriz original.
1.4 Problemas propuestos
1) Considere las matrices:
Ingrese las matrices al ambiente Matlab® e indique los resultados al aplicar las siguientes operaciones:
i) x = A(3,:)
ii) y = B(:,1)
iii) C = A(2:4;1:2)
iv) D = diag(A)
v) E = B’
vi) C = diag(diag(A))
vii) A(3,:)*B(:,2)
viii) A(2,3)*B(4,1)
ix) eig(A)
x) [U, R] = eig(A)
xi) A+4*ones(4,4)-eye(4)
2) Considere los objetos:
Encuentre lo que se obtendría al efectuar las siguientes operaciones en Matlab®:
i) A*B
ii) C = diag(diag(A))
iii) A(1:2,3:4)*B(3:4)
iv) A(3,:)*B
v) A(:,2)’*B
vi) A(2,3)*B(1)
vii) B’* diag(A)
viii) [U, R] = eig(A*A’)
ix) A-5*ones(4,4)+3*eye(4)
3) Para las matrices A de los dos problemas anteriores, descompóngalas de manera que se puedan expresar de la siguiente forma: A =D + L + U, donde D es la diagonal principal de A, L una matriz triangular inferior y U una matriz triangular superior. Utilice los comandos diag, tril y triu de Matlab®.
4) ¿Cuál de las siguientes matrices es singular?
Si la matriz no es singular, encuentre su inversa.
5) Encuentre la pendiente e intercepto de la recta de mínimos cuadrados que pasa por los puntos: (1, 1), (2, 1.5) y (3, 2.2). Grafique sus resultados usando plot.
6) ¿Tienen solución los siguientes sistemas de ecuaciones? Encuéntrela si ella existe.
7) Resuelva las siguientes ecuaciones lineales usando eliminación de Gauss-Jordan con diferentes opciones de pivoteo.
8) Emplee el método de Jacobi y el de Gauss-Seidel para resolver el sistema de ecuaciones del problema anterior, caso b); en ambos casos use como punto de partida x(0) = [1,1,0.4]T.
9) Considere el sistema lineal de ecuaciones A·x = b con
a) Resuelva el sistema anterior en Matlab® usando: x = A\b, operador que realiza la eliminación gaussiana. ¿Se obtiene el mismo resultado que al aplicar: x=inv(A)*b?
b) Defina un nuevo vector del lado derecho como: b1 = [30.00002; 42.00006; 22.00004] y calcule la nueva solución x1 utilizando ‘\’ y la función inv. ¿Hay diferencias en los resultados? ¿Esperaba usted que hubiera diferencias en la solución numérica x?
c) Compare los errores relativos al modificar el lado derecho3: norm(b1-b)/norm(b)
