1. Educación

Factorización LU de matrices cuadradas
Antes de comenzar a ver este tema convendría hacer un repaso de las nociones de
matrices:
 Operaciones con matrices (Suma de matrices, producto por un escalar).
 Matrices cuadradas, matriz inversa de una matriz.
 Matriz fila, matriz columna.
 Matriz diagonal, matriz triangular.
En concreto, sobre este último aspecto tenemos:
a) Matriz diagonal:
 d11 0

0 d 22
D 
 0
0

0
 0
0
0
d33
0
0 

0 
0 

d 44 
b) Matriz triangular superior (Upper triangular):
 u11 u12

0 u22
U 
 0
0

0
 0
u13
u23
u33
0
u14 

u24 
u34 

u44 
c) Matriz triangular inferior (Lower triangular):
 l11 0

l21 l22
L
 l31 l32

 l41 l42
0
0
l33
l43
0

0
0

l44 
En los tres casos citados el determinante de la matriz es el producto de los elementos
de la diagonal.

Matriz estrictamente dominante diagonalmente.
Se dice que una matriz es estrictamente dominante diagonalmente si cumple:
n
aii   aij
i 1
i j
Las matrices estrictamente dominantes diagonalmente son no singulares (su
determinante es no nulo), además pueden resolverse por eliminación gaussiana sin
intercambio de filas y por tanto pueden siempre factorizarse en la forma A = L.U.
Ejemplo de matriz estrictamente dominante diagonalmente:
6
2
A
1

 1
2 1 1
4 1 0 
1 4 1

0 1 3 
Para una matriz cuadrada A es muy útil descomponerla como producto de una matriz
triangular inferior por una triangular superior:
A=L.U
La factorización L.U de una matriz cuadrada tiene la siguiente ventaja:
Al operar con un sistema en la forma A . x = b, si expresamos A = L . U tendremos:
L.U.x=b
Y ahora podemos poner:
U . x = L-1. b
Y teniendo en cuenta que U es triangular y L-1 también es triangular, la solución del
sistema se limita a:
 a ) L1. b  y
U . x  L1. b 
 b) U . x  y
Es decir, se resuelven dos sistemas triangulares (más simple y exacto que resolver uno
no triangular): en a) se obtiene y, finalmente en b) se obtiene x.
La factorización LU de una matriz cuadrada A se podría hacer por el método de Gauss
(consultar las pag. 173-174 de “Curso de Métodos Numéricos” de Virginia Muto), sin
embargo hay otros métodos como el de Doolittle, el de Crout, y el de Cholesky.
 l11 0

l
l
L.U   21 22
 l31 l32

 l41 l42
0 0   u11

0 0   0

l33 0   0
 
l43 l44   0
 a11 a12
a
a
  21 22
 a31 a32

 a41 a42
u12
u22
0
0
a13
a23
a33
a43
u13 u14 
u23 u24 

u33 u34 

0 u44 
a14 
a24 
a34 

a44 
Observad arriba que conocida A, tendríamos 16 ecuaciones con 20 incógnitas: lij,
uij, sin embargo, nosotros podemos dar cualquier valor a 4 de ellas. Así hay dos métodos
que se suelen utilizar: a) Método de Doolittle: se toman los cuatro elementos de la
diagonal de L como ‘1’, lii=1 . b) Método de Crout: se toman los cuatro elementos de la
diagonal de Ucomo ‘1’, uii=1.
Nosotros aquí seguiremos el método de Doolittle, con un ejemplo. Sea dada A,
entonces es sencillo hallar L y U, sin más que resolver un sistema:
 u11 u12
 0 u
22
U 
 0
0

0
 0
1 0
l
1
L   21
 l31 l32

 l41 l42
0 0  6
0 0   2

1 0  1
 
l43 1   1
2
4
1
0
u13
u23
u33
0
1
1
4
1
u14 
u24 
u34 

u44 
 1
0 
A
 1

3 
La forma de operar “manualmente” es:
- 1) Considerar el desarrollo de los elementos de la primera fila de A.
- 2) Considerar el desarrollo de los elementos de la primera columna de A.
- 3) Considerar el desarrollo de los elementos de la segunda fila de A.
- Etc. etc.
Es decir, u11  6, u12  2, u13  1, u14  1  uij 
Ahora continuamos con:
l21u11  2, l31u11  1, l41u11  1  l j1 
a j1
u11
1
1
1
Así obtenemos: l21  , l31  , l41   .
3
6
6
a1 j
l11
Hasta ahora tenemos:
2
 6
u u
U   11 22
 0
0
 0
0

0
 1
 1/ 3 1
L
 1/ 6 l32

 1/ 6 l42
0
0
1
l43
0  6
0   2

0  1
 
1   1
1 1 
u23 u24 
u33 u34 

0 u44 
2
1  1
4
1
0 
A
1
4
 1

0 1
3 
Seguiríamos con la siguiente fila de U:
l21u12  l22u22  a22  u22 
1
(a22  l21u12 )
l22
l21u13  l22u23  a23  u23 
1
(a23  l21u13 )
l22

 
En realidad para obtener la segunda fila de U sirve la siguiente fórmula:
u2 j 
i 1
1 

a

l2 k ukj  (j=1, …, 4)

2j

l22 
k 1

Para una fila i-ésima de U serviría:
i 1
1

uij   aij   lik ukj 
lii 
k 1

De una manera análoga, la columna i-ésima de L se calcula por la fórmula:
l ji 
i 1
1

a

l jk uki 

ji

uii 
k 1

ALGORITMO DE FACTORIZACIÓN DE DOOLITTLE.
Entrada: Elementos de A. Tomaremos inicialmente L= eye(4), U = zeros(4).
Salida: Matrices L y U, producto L*U.
Paso 1: Para j = 1, …, n tomar:
a
uij  1 j (primera fila de U).
l11
a
l j1  j1 (primera columna de L).
u11
Paso 2: Para i = 2, 3,…, n-1 seguir paso 3.
Paso 3: Para j = i+1, …, n tomar:
i 1
3-A:
uii  aii   lik uki ;
k 1
i 1
4-A:
uij  aij   lik ukj ;
k 1
i 1
1

a

l jk uki  ;

ji

uii 
k 1

Paso 4: Para el último término, un,n :
5-A:
l ji 
n 1
Tomar: unn  ann   lnk ukn
k 1
Paso 5: (Salida):
Desplegar L, U y L*U.
MÉTODO DE CHOLESKI
En el caso de que la matriz A, n x n, sea simétrica y definida positiva
(Xt . A . X > 0), entonces A admite una factorización en la forma A = L.Lt, siendo L
triangular inferior.
 l11 0

l21 l22
t
L.L  
 l31 l32

 l41 l42
0
0
l33
l43
0   l11 l21 l31 l41 
 a11



0   0 l22 l32 l42 
a21

 
 a31
0   0 0 l33 l43 
 


l44   0 0 0 l44 
 a41
a12
a13
a22
a23
a32
a33
a42
a43
a14 

a24 
a34 

a44 
Los 10 elementos de L se obtienen a partir de los elementos conocidos de A:
Algoritmo de Choleski
Entrada: Matriz A, L = zeros(n)
Salida: Matriz A, Matriz L, Matriz producto L*Lt.
Paso 1: Tomar l11  a11
Paso 2: Para j = 2,…, n tomar l j1 
a j1
l11
,
Paso 3: Para i = 2, … , n-1 seguir los pasos 4 y 5.
i 1
Paso 4: Tomar lii  aii   lik2
k 1
Paso 5:
Para j = i+1, … , n tomar:
i 1
1

l ji   a ji   l jk lik 
lii 
k 1

n 1
Paso 6: Tomar
lnn  ann   lnk2
k 1