Para bajar por una escalera de mano sin saltar, no puede faltar ningún peldaño. Si falta uno, tienes que dar un rodeo. La eliminación de Gauss es esa escalera: cada pivote es un peldaño, y la LU existe si puedes bajar hasta el final sin intercambiar filas. Hoy veremos cómo saber de antemano si falta algún peldaño (hay un criterio con determinantes), qué rodeo se da cuando falta (\(PA=LU\)) y, sobre todo, por qué todo este esfuerzo compensa: contaremos las operaciones.
desliza para ver todo →
Con una matriz \(1000\times1000\), resolver un sistema cuesta lo mismo con Gauss que con LU. Pero si hay que resolver 100 sistemas con la misma matriz, repetir Gauss cada vez cuesta unas 77 veces más multiplicaciones que factorizar una vez y reutilizar. Con 1000 sistemas, unas 251 veces más. Las cuentas están más abajo.
En el post anterior dejamos dos preguntas abiertas. La primera: ¿qué matrices tienen LU? Vimos que \(\begin{pmatrix}0 & 1\\ 1 & 1\end{pmatrix}\) no la tiene, aunque es invertible. La segunda: ¿cuánto se ahorra de verdad? Dijimos «del orden de \(n^3\) frente a \(n^2\)», y hoy lo precisamos.
Los menores principales deciden
Los menores principales de una matriz \(A\) de orden \(n\) son los determinantes de sus esquinas superiores izquierdas:
\[\begin{gathered}\Delta_1=a_{11},\quad\Delta_2=\begin{vmatrix}a_{11} & a_{12}\\ a_{21} & a_{22}\end{vmatrix},\\ \dots,\quad\Delta_n=\det A.\end{gathered}\]Sea \(A\) una matriz cuadrada invertible. Son equivalentes:
- \(A\) tiene factorización LU.
- \(A\) se puede escalonar con Gauss sin intercambiar filas.
- \(\Delta_k\neq0\) para todo \(k=1,\dots,n\).
Y en ese caso, los pivotes son \(u_{11}=\Delta_1\) y \(u_{kk}=\dfrac{\Delta_k}{\Delta_{k-1}}\) para \(k\ge2\).
Por qué. \((2)\Rightarrow(1)\) es la receta del post 11. Para las otras dos implicaciones hay una observación clave: la esquina \(k\times k\) de \(LU\) es el producto de las esquinas \(k\times k\) de \(L\) y de \(U\). Es porque el elemento \((i,j)\) de \(LU\) es \(\sum_m\ell_{im}u_{mj}\), y \(\ell_{im}=0\) si \(m>i\): para \(i,j\le k\) solo intervienen índices \(m\le k\).
- \((1)\Rightarrow(3)\): si \(A=LU\), entonces \(\Delta_k=\det L_k\det U_k=1\cdot u_{11}\cdots u_{kk}\). Como \(\det A=u_{11}\cdots u_{nn}\neq0\), ningún \(u_{ii}\) es cero, así que ningún \(\Delta_k\) lo es.
- \((3)\Rightarrow(2)\): supón que ya has hecho \(k-1\) pasos de Gauss sin intercambios. Las operaciones \(f_i-\ell f_j\) con \(j<i\) no cambian los menores principales (dentro de la esquina \(k\times k\), solo suman a una fila un múltiplo de otra de más arriba). Y tras esos pasos, la esquina \(k\times k\) es triangular, con diagonal \(u_{11},\dots,u_{kk}\). Así que \(\Delta_k=u_{11}\cdots u_{kk}\), es decir, \(u_{kk}=\Delta_k/\Delta_{k-1}\neq0\). El siguiente pivote nunca es cero. \(\blacksquare\)
¿Tienen LU \(A=\begin{pmatrix}2 & 1 & 1\\ 4 & 3 & 3\\ -2 & 1 & 4\end{pmatrix}\) y \(C=\begin{pmatrix}1 & 2 & 3\\ 2 & 4 & 5\\ 1 & 3 & 4\end{pmatrix}\)?
¿Y si \(A\) no es invertible? El teorema habla de matrices invertibles. Una matriz singular puede tener LU o no. Por ejemplo, la del ejercicio 11.7 con \(a=2\) la tiene, con un cero en el último pivote de \(U\). Para resolver sistemas, lo que interesa es el caso invertible.
Cuando falta un peldaño: \(PA=LU\)
Una matriz de permutación \(P\) es la identidad con las filas reordenadas. Multiplicar \(PA\) reordena las filas de \(A\) de la misma forma. Por ejemplo, \(P=\begin{pmatrix}0 & 1 & 0\\ 1 & 0 & 0\\ 0 & 0 & 1\end{pmatrix}\) intercambia las filas 1 y 2.
Para toda matriz cuadrada invertible \(A\) existe una matriz de permutación \(P\) tal que \(PA\) tiene factorización LU:
\[PA=LU.\]La idea es simple: Gauss siempre puede seguir intercambiando filas (si \(A\) es invertible, siempre hay un elemento no nulo por debajo del pivote que falta). Si sabes de antemano qué intercambios vas a necesitar, los haces todos al principio, y la matriz resultante \(PA\) ya se escalona sin intercambios.
Cómo se usa: \(AX=B\iff PAX=PB\iff LUX=PB\). Es decir, se reordenan las componentes de \(B\) igual que las filas, y después se resuelve como siempre. Y para el determinante: cada intercambio cambia el signo, así que \(\det A=(-1)^s\,u_{11}\cdots u_{nn}\), donde \(s\) es el número de intercambios.
Truco práctico. No hace falta adivinar los intercambios. Escalona normalmente apuntando multiplicadores y, si en algún momento intercambias dos filas, intercambia también los multiplicadores que ya hubieras apuntado en esas dos filas. Al final, \(P\) recoge todos los intercambios hechos.
Calcula \(PC=LU\) para \(C=\begin{pmatrix}1 & 2 & 3\\ 2 & 4 & 5\\ 1 & 3 & 4\end{pmatrix}\).
Resuelve \(AX=(3,2,9)^t\) con \(A=\begin{pmatrix}0 & 2 & 1\\ 1 & 1 & 0\\ 2 & 4 & 3\end{pmatrix}\).
Solución: \((1,1,1)\). Comprobación en el sistema original: \(0+2+1=3\) ✔, \(1+1+0=2\) ✔, \(2+4+3=9\) ✔. Y \(\det A=-1\cdot2\cdot2=-4\).
Intercambiar también cuando no hace falta
Hay una razón más para permutar, y es la que usan todos los programas: los ordenadores redondean. Un pivote muy pequeño, aunque no sea cero, obliga a usar multiplicadores enormes, y los errores de redondeo se disparan. La solución habitual es el pivoteo parcial: en cada columna se sube a la posición del pivote el elemento de mayor valor absoluto.
Resuelve \(\begin{cases}0{,}0001\,x+y=1\\ x+y=2\end{cases}\) redondeando cada resultado a 3 cifras significativas.
La solución exacta es \(x=\tfrac{10000}{9999}\approx1{,}0001\) e \(y\approx0{,}9999\). Redondeadas a 3 cifras, las dos valen \(1{,}00\).
El problema no era que \(A\) fuese mala (\(\det A\approx-1\), está lejos de ser singular), sino el orden de las operaciones. Con el pivote pequeño, el error al redondear \(y\) se multiplica por \(10000\) al despejar \(x\).
Contando multiplicaciones
Contamos multiplicaciones y divisiones, que es lo tradicional. Las sumas y restas son aproximadamente tantas como multiplicaciones, así que, si las cuentas también, todo se multiplica por 2 y las conclusiones no cambian.
Factorizar: \(\frac{n^3-n}{3}\)
En el paso \(k\) de Gauss hay \(n-k\) filas por debajo del pivote. Para cada una hay que calcular su multiplicador (1 división) y actualizar los \(n-k\) elementos que quedan a la derecha (\(n-k\) multiplicaciones). En total, \((n-k)(n-k+1)\) operaciones en el paso \(k\). Sumando y llamando \(j=n-k\):
\[\sum_{j=1}^{n-1}j(j+1)=\frac{(n-1)\,n\,(n+1)}{3}=\frac{n^3-n}{3}\approx\frac{n^3}{3}.\](La fórmula de la suma la demuestras en el ejercicio 12.6.)
Cada \(B\): \(n^2\)
- \(LY=B\): en la fila \(i\) hay \(i-1\) multiplicaciones (no hay división, porque la diagonal de \(L\) son unos). Total: \(\frac{n(n-1)}{2}\).
- \(UX=Y\): en la fila \(i\) hay \(n-i\) multiplicaciones y 1 división. Total: \(\frac{n(n+1)}{2}\).
Sumando: \(\frac{n(n-1)}{2}+\frac{n(n+1)}{2}=n^2\). Y Gauss completo para un sistema cuesta lo mismo que factorizar más una resolución: \(\frac{n^3-n}{3}+n^2\).
Con \(k\) sistemas de la misma matriz
\[\begin{aligned}\text{Gauss cada vez:}&\quad k\left(\tfrac{n^3-n}{3}+n^2\right)\\ \text{LU una vez:}&\quad \tfrac{n^3-n}{3}+k\,n^2\end{aligned}\]Con \(n=1000\), factorizar cuesta unos \(3{,}3\cdot10^8\) multiplicaciones y cada \(B\) solo \(10^6\). De ahí salen los cocientes del dibujo inicial: con \(k=100\), \(3{,}34\cdot10^{10}\) frente a \(4{,}33\cdot10^8\), unas 77 veces menos.
FACTORIZAR LU
\(\approx n^3/3\)
Una sola vez por matriz. Después, cada \(B\) cuesta \(n^2\).
CALCULAR \(A^{-1}\)
\(\approx n^3\)
Unas tres veces más que la LU. Y después, \(A^{-1}B\) cuesta también \(n^2\): no gana nada, y es más sensible al redondeo.
CRAMER
\(\approx n^4/3\)
\(n+1\) determinantes, cada uno calculado con Gauss. Con cofactores, sería del orden de \(n!\): impracticable ya con \(n=20\) (\(20!\approx2{,}4\cdot10^{18}\)).
Para hacerse una idea, imagina un ordenador que hiciera mil millones (\(10^9\)) de multiplicaciones por segundo, y \(n=1000\). La LU tardaría unas tres décimas de segundo, y cada \(B\) nuevo, una milésima. Cramer con determinantes por Gauss, unos \(\tfrac{10^{12}/3}{10^9}\approx333\) segundos, más de cinco minutos. Y Cramer con cofactores necesitaría del orden de \(1000!\) operaciones, un número de más de 2500 cifras: ni con todos los ordenadores del mundo durante la edad del universo.
Cuando \(n\) es de verdad grande
- SimulaciónMillones de incógnitasCuando se discretiza un problema físico (una estructura, un fluido, un campo eléctrico), salen sistemas con muchísimas incógnitas. Ahí la diferencia entre \(n^3\) y \(n^2\) es la diferencia entre posible e imposible.
- Matrices bandaAprovechar los cerosEn esos problemas cada incógnita solo se relaciona con sus vecinas, y la matriz tiene casi todo ceros. Si los elementos no nulos están a distancia como mucho \(b\) de la diagonal, Gauss sin intercambios no rellena fuera de esa banda, y la LU cuesta del orden de \(n\,b^2\) en lugar de \(n^3\) (ejercicio 12.8).
- Cálculo numéricoPivoteoLos programas de cálculo científico resuelven sistemas generales con \(PA=LU\) y pivoteo parcial, precisamente por el problema del redondeo del ejemplo 4.
Ocho ejercicios sobre existencia y coste
Intenta cada uno sin abrir nada. Si te atascas, abre la pista. La solución, solo al final.
¿Tiene LU?
Usa los menores principales para decidir si cada matriz tiene LU:
Ver pista
Calcula \(\Delta_1\), \(\Delta_2\) y \(\Delta_3\). En c), ¿es invertible? ¿Se aplica el teorema?
Ver solución
a) \(\Delta_1=2\), \(\Delta_2=2+4=6\), \(\Delta_3=68\). Todos no nulos: sí. Pivotes: \(2\), \(3\) y \(\tfrac{68}{6}=\tfrac{34}{3}\).
b) \(\Delta_1=1\), \(\Delta_2=0\), \(\Delta_3=-1\). Es invertible y \(\Delta_2=0\): no.
c) \(\det=0\): no es invertible y el teorema no se aplica. De hecho, sí tiene LU: \(f_2-2f_1\) da \(U=\begin{pmatrix}3 & 1\\ 0 & 0\end{pmatrix}\) y \(L=\begin{pmatrix}1 & 0\\ 2 & 1\end{pmatrix}\). El cero aparece en el último pivote, cuando ya no hay nada que anular.
Pivotes sin escalonar
De una matriz \(3\times3\) sabes que sus menores principales son \(\Delta_1=2\), \(\Delta_2=6\) y \(\Delta_3=-12\). Sin conocer la matriz, ¿tiene LU? ¿Cuál es la diagonal de \(U\)? ¿Y \(\det A\)?
Ver pista
\(u_{kk}=\Delta_k/\Delta_{k-1}\).
Ver solución
Todos no nulos (y \(\det A=\Delta_3\neq0\)): tiene LU. Pivotes: \(u_{11}=2\), \(u_{22}=\tfrac62=3\) y \(u_{33}=\tfrac{-12}{6}=-2\). \(\det A=\Delta_3=-12=2\cdot3\cdot(-2)\) ✔.
LU con parámetro
Sea \(A=\begin{pmatrix}1 & 2 & 0\\ 2 & a & 1\\ 0 & 1 & 1\end{pmatrix}\). a) ¿Para qué valores de \(a\) es invertible y tiene LU? b) Para \(a=4\), ¿es invertible? Calcula una factorización \(PA=LU\).
Ver pista
a) \(\Delta_2=a-4\) y \(\Delta_3=\det A=a-5\). b) Con \(a=4\), el segundo pivote es cero.
Ver solución
a) \(\Delta_1=1\), \(\Delta_2=a-4\) y \(\Delta_3=1\cdot(a-1)-2\cdot(2-0)=a-5\). Invertible con LU si y solo si \(a\neq4\) y \(a\neq5\).
b) \(\det A=4-5=-1\neq0\): invertible, pero sin LU. Con \(f_2-2f_1\), la fila 2 queda \((0,0,1)\): pivote cero. Intercambiamos las filas 2 y 3 (y su multiplicador, que pasa a la fila 3):
\[PA=\begin{pmatrix}1 & 2 & 0\\ 0 & 1 & 1\\ 2 & 4 & 1\end{pmatrix},\]\[L=\begin{pmatrix}1 & 0 & 0\\ 0 & 1 & 0\\ 2 & 0 & 1\end{pmatrix},\]\[U=\begin{pmatrix}1 & 2 & 0\\ 0 & 1 & 1\\ 0 & 0 & 1\end{pmatrix}.\]Comprobación del determinante: \(\det A=(-1)^1\cdot1\cdot1\cdot1=-1\) ✔.
Resolver con permutación
Calcula \(PA=LU\) para \(A=\begin{pmatrix}0 & 1 & 1\\ 2 & 1 & 0\\ 1 & 1 & 1\end{pmatrix}\) y resuelve \(AX=(2,3,3)^t\).
Ver pista
\(\Delta_1=0\). Intercambia las filas 1 y 2, y no olvides reordenar también \(B\).
Ver solución
Con \(P\) intercambiando las filas 1 y 2, \(PA=\begin{pmatrix}2 & 1 & 0\\ 0 & 1 & 1\\ 1 & 1 & 1\end{pmatrix}\). Con \(f_3-\tfrac12f_1\) queda \((0,\tfrac12,1)\), y con \(f_3-\tfrac12f_2\), \((0,0,\tfrac12)\):
\[L=\begin{pmatrix}1 & 0 & 0\\ 0 & 1 & 0\\ \tfrac12 & \tfrac12 & 1\end{pmatrix},\]\[U=\begin{pmatrix}2 & 1 & 0\\ 0 & 1 & 1\\ 0 & 0 & \tfrac12\end{pmatrix}.\]\(PB=(3,2,3)\). \(LY=PB\): \(y_1=3\), \(y_2=2\), \(y_3=3-\tfrac32-1=\tfrac12\). \(UX=Y\): \(x_3=1\), \(x_2=2-1=1\), \(2x_1=3-1\), así que \(x_1=1\). \(X=(1,1,1)\).
Comprobación: \(0+1+1=2\) ✔, \(2+1+0=3\) ✔, \(1+1+1=3\) ✔.
Haciendo números
Hay que resolver 50 sistemas con la misma matriz \(100\times100\). Cuenta las multiplicaciones y divisiones: a) repitiendo Gauss cada vez; b) con una LU. c) ¿Cuántas veces más rápida es la LU?
Ver pista
Con \(n=100\): \(\frac{n^3-n}{3}=333\,300\) y \(n^2=10\,000\).
Ver solución
a) \(50\cdot(333\,300+10\,000)=17\,165\,000\).
b) \(333\,300+50\cdot10\,000=833\,300\).
c) \(17\,165\,000/833\,300\approx20{,}6\) veces.
La suma del coste
Demuestra por inducción que \(\displaystyle\sum_{j=1}^{m}j(j+1)=\frac{m(m+1)(m+2)}{3}\) para todo \(m\ge1\). Deduce que factorizar una matriz \(n\times n\) cuesta \(\frac{n^3-n}{3}\) multiplicaciones y divisiones, y compruébalo contando a mano para \(n=3\).
Ver pista
Para el paso de inducción, suma \((m+1)(m+2)\) a los dos lados y saca \((m+1)(m+2)\) factor común. Después toma \(m=n-1\).
Ver solución
Caso \(m=1\): \(1\cdot2=2=\frac{1\cdot2\cdot3}{3}\) ✔.
Paso: si vale para \(m\), entonces \(\sum_{j=1}^{m+1}j(j+1)=\frac{m(m+1)(m+2)}{3}+(m+1)(m+2)=(m+1)(m+2)\left(\frac m3+1\right)=\frac{(m+1)(m+2)(m+3)}{3}\), que es la fórmula para \(m+1\). \(\blacksquare\)
Con \(m=n-1\): \(\frac{(n-1)n(n+1)}{3}=\frac{n^3-n}{3}\).
Para \(n=3\): en el paso 1, dos filas con 1 división y 2 multiplicaciones cada una (6 operaciones). En el paso 2, una fila con 1 división y 1 multiplicación (2 operaciones). Total: 8, y \(\frac{27-3}{3}=8\) ✔.
El pivote pequeño
Resuelve \(\begin{cases}0{,}0002\,x+y=1\\ x+2y=3\end{cases}\) redondeando cada resultado intermedio a 3 cifras significativas: a) sin intercambiar filas; b) con pivoteo parcial. c) Compara con la solución exacta, \(x\approx1{,}0004\), \(y\approx0{,}9998\).
Ver pista
a) El multiplicador es \(1/0{,}0002=5000\). Al redondear \(2-5000\) y \(3-5000\) a 3 cifras, ¿qué queda?
Ver solución
a) La segunda ecuación queda \((2-5000)\,y=3-5000\), es decir, \(-4998\,y=-4997\). Con 3 cifras, los dos se redondean a \(-5000\), así que \(y=1{,}00\). Entonces \(x=\dfrac{1-1{,}00}{0{,}0002}=0\). Totalmente equivocado.
b) Intercambiamos (el 1 es mayor que \(0{,}0002\)). El multiplicador es \(0{,}0002\), y la nueva segunda ecuación es \((1-0{,}0004)\,y=1-0{,}0006\). Con 3 cifras, \(1{,}00\,y=0{,}999\), así que \(y=0{,}999\) y \(x=3-2\cdot0{,}999=3-1{,}998=1{,}002\), que con 3 cifras es \(1{,}00\).
c) Con pivoteo: \(x=1{,}00\) e \(y=0{,}999\). Errores de \(0{,}0004\) y \(0{,}0008\), del tamaño del redondeo. Sin pivoteo, el error en \(x\) es del 100 %.
Una barra con un millón de puntos
Para calcular con precisión la temperatura de una barra (como en el ejercicio 11.8), la dividimos en tantos tramos que hay \(n=10^6\) puntos interiores. La matriz es tridiagonal, con \(2\) en la diagonal y \(-1\) justo encima y justo debajo. a) Si la tratáramos como una matriz cualquiera, ¿cuántas multiplicaciones costaría la LU? b) Demuestra que, aprovechando que es tridiagonal, la LU cuesta \(2(n-1)\) operaciones, y cada sistema \(3n-2\). c) Con un ordenador imaginario de \(10^9\) multiplicaciones por segundo, ¿cuánto tardaría cada opción?
Ver pista
b) En cada paso de Gauss solo hay una fila que tiene algo debajo del pivote, y en esa fila solo cambia un elemento (el de la diagonal). Y \(L\) y \(U\) son bidiagonales.
Ver solución
a) \(\frac{n^3-n}{3}\approx\frac{10^{18}}{3}\approx3{,}3\cdot10^{17}\).
b) En el paso \(k\) solo la fila \(k+1\) tiene un elemento no nulo bajo el pivote. Se calcula su multiplicador (1 división) y se actualiza solo su elemento diagonal (1 multiplicación): el de su derecha no cambia, porque encima de él, en la fila \(k\), hay un cero. Son 2 operaciones por paso y \(n-1\) pasos: \(2(n-1)\). Para cada \(B\), \(LY=B\) tiene 1 multiplicación por fila salvo la primera (\(n-1\) en total), y \(UX=Y\) tiene 1 división por fila y 1 multiplicación por fila salvo la última (\(2n-1\)). Total: \(3n-2\). \(\blacksquare\)
c) Matriz general: \(3{,}3\cdot10^{17}/10^9\approx3{,}3\cdot10^8\) segundos, unos diez años y medio. Tridiagonal: unos \(2\cdot10^6+3\cdot10^6=5\cdot10^6\) operaciones, unas cinco milésimas de segundo. Aprovechar la estructura de la matriz física lo cambia todo.
Las trampas de este post
- Mirar solo \(\det A\). Para la LU hacen falta todos los menores principales no nulos. \(\det A\neq0\) solo garantiza \(PA=LU\).
- Usar cualquier menor. Son los de las esquinas superiores izquierdas, no otros menores de orden \(k\).
- Intercambiar filas y olvidar los multiplicadores. Si intercambias dos filas a mitad de camino, los multiplicadores ya apuntados en ellas también se intercambian.
- Olvidar reordenar \(B\). Con \(PA=LU\), el sistema es \(LUX=PB\), no \(LUX=B\).
- Olvidar el signo del determinante. \(\det A=(-1)^s\,u_{11}\cdots u_{nn}\), con \(s\) el número de intercambios.
- Creer que calcular \(A^{-1}\) ahorra trabajo. Cuesta unas tres veces más que la LU, y después cada \(B\) cuesta lo mismo que con la LU.
Lo imprescindible del post
Existencia
- \(A\) invertible tiene LU \(\iff\Delta_1,\dots,\Delta_n\neq0\).
- Pivotes: \(u_{kk}=\Delta_k/\Delta_{k-1}\).
- Siempre: \(PA=LU\), y \(AX=B\iff LUX=PB\).
- \(\det A=(-1)^s\,u_{11}\cdots u_{nn}\).
- En la práctica: pivoteo parcial (el mayor en valor absoluto).
Coste (multiplicaciones y divisiones)
- Factorizar: \(\frac{n^3-n}{3}\approx\frac{n^3}{3}\).
- Cada \(B\): \(n^2\).
- \(k\) sistemas: \(\frac{n^3}{3}+kn^2\) frente a \(k\left(\frac{n^3}{3}+n^2\right)\).
- Inversa: \(\approx n^3\). Cramer: \(\approx\frac{n^4}{3}\).
- Tridiagonal: del orden de \(n\).
Por qué los ordenadores usan LU
La LU existe exactamente cuando ningún menor principal se anula, y cuando alguno se anula basta reordenar las filas: \(PA=LU\). Contar operaciones explica todo lo demás: factorizar cuesta del orden de \(n^3/3\), cada sistema nuevo solo \(n^2\), y ni la inversa ni Cramer pueden competir. Con eso, el Tema 2 queda completo.
En el último post del tema lo juntamos todo: una batería de ejercicios resueltos, una tabla para elegir la herramienta adecuada en cada caso y un simulacro de examen.
Si no falta ningún peldaño, baja directo. Si falta, reordena y baja igual.
Siguiente: ejercicios resueltos del Tema 2 →Como la escalera de mano: lo importante no es lo larga que sea, sino que no le falte ningún peldaño.
Todas las Etiquetas:
Comentarios
Todavía no hay comentarios. ¡Sé el primero en comentar!