Álgebra Linea Aplicada (Notas de Clase) PDF

Universidad Nacional de Colombia
Facultad de Ciencias
Escuela de Matemáticas
Invitación al álgebra lineal numérica
Carlos E. Mejía
Profesor Titular
Trabajo presentado en el marco del

Semestre Sabático Julio 2016 - Enero 2017
2
Índice general
I Algebra lineal 11
1. Matrices y vectores 13
1.1. Matrices . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 14
1.2. Multiplicación de matriz por vector . . . . . . . . . . . . . . . . . . . 17
1.3. Nota MATLAB . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 18
2. Matrices y determinantes 21
2.1. Producto de matrices . . . . . . . . . . . . . . . . . . . . . . . . . . . 21
2.2. Determinantes . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 24
2.3. Inversa de una matriz . . . . . . . . . . . . . . . . . . . . . . . . . . . 25
2.4. Nota MATLAB . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 26
3. Espacios vectoriales 27
3.1. Definición . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 27
3.2. Ejemplos de espacios vectoriales . . . . . . . . . . . . . . . . . . . . 28
3.3. Subespacios . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 28
3.4. Ejemplos de subespacios . . . . . . . . . . . . . . . . . . . . . . . . . 28
3.5. Independencia lineal . . . . . . . . . . . . . . . . . . . . . . . . . . . 29
3.6. Base de un espacio vectorial . . . . . . . . . . . . . . . . . . . . . . . 30
3.7. Sumas e intersecciones de subespacios . . . . . . . . . . . . . . . . . . 30
4. Sistemas lineales 33
4.1. Eliminación de Gauss . . . . . . . . . . . . . . . . . . . . . . . . . . . 33
4.2. Formulación matricial de la eliminación de Gauss . . . . . . . . . . . 36
4.3. Nota MATLAB . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 39
4.4. Invertibilidad . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 40
4.5. Nota MATLAB . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 41
3
4 ÍNDICE GENERAL
5. Forma escalonada reducida y rango 43

5.1. Operaciones elementales . . . . . . . . . . . . . . . . . . . . . . . . . 44
5.2. Matrices elementales . . . . . . . . . . . . . . . . . . . . . . . . . . . 45
5.3. Relación entre columnas y filas . . . . . . . . . . . . . . . . . . . . . 48
5.4. Forma normal del rango . . . . . . . . . . . . . . . . . . . . . . . . . 50
6. Matrices y transformaciones lineales 53

6.1. Definición . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 53
6.2. Representación matricial . . . . . . . . . . . . . . . . . . . . . . . . . 54
6.3. Composición y subespacios asociados . . . . . . . . . . . . . . . . . . 55
6.4. Ortogonalidad . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 56
6.5. Cuatro subespacios fundamentales . . . . . . . . . . . . . . . . . . . . 57
6.6. Valores y vectores propios . . . . . . . . . . . . . . . . . . . . . . . . 58
6.7. Nota MATLAB . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 59
7. Normas de vectores y matrices 61

7.1. Valor absoluto . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 61
7.2. Normas de vectores . . . . . . . . . . . . . . . . . . . . . . . . . . . . 62
7.3. Normas de matrices . . . . . . . . . . . . . . . . . . . . . . . . . . . . 64
7.4. Nota MATLAB . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 66
II Análisis matricial 67
8. Perturbaciones, errores y condicionamiento 69
8.1. Errores absoluto y relativo . . . . . . . . . . . . . . . . . . . . . . . . 70
8.2. Sistema lineal . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 71
8.3. Condicionamiento de normas matriciales . . . . . . . . . . . . . . . . 73
8.4. Suma y producto de matrices . . . . . . . . . . . . . . . . . . . . . . 73
8.5. Invertir una matriz . . . . . . . . . . . . . . . . . . . . . . . . . . . . 74
9. Condicionamiento y sistemas lineales 79

9.1. Lado derecho perturbado . . . . . . . . . . . . . . . . . . . . . . . . . 79
9.2. Caso general . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 80
9.3. Estabilidad de métodos directos . . . . . . . . . . . . . . . . . . . . . 81
10.Sistemas lineales especiales 83

10.1. Problemas con valores en la frontera . . . . . . . . . . . . . . . . . . 83
10.2. Ecuación del calor . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 87
ÍNDICE GENERAL 5
10.3. Ecuación de Laplace . . . . . . . . . . . . . . . . . . . . . . . . . . . 90

10.4. Nota MATLAB . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 94
11.Primeros algoritmos de factorización de matrices 95

11.1. Matriz diagonalizable . . . . . . . . . . . . . . . . . . . . . . . . . . . 95
11.2. Factorización LU . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 98
11.3. Factorización de Cholesky . . . . . . . . . . . . . . . . . . . . . . . . 100
12.Factorización QR 101
12.1. Matrices de Householder . . . . . . . . . . . . . . . . . . . . . . . . . 101
12.2. La factorización QR . . . . . . . . . . . . . . . . . . . . . . . . . . . 103
12.3. Valores propios . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 104
12.4. Nota MATLAB . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 107
13.Descomposición en valores singulares (SVD) 109

13.1. La factorización SVD . . . . . . . . . . . . . . . . . . . . . . . . . . . 110
13.2. Nota MATLAB . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 116
13.3. ¿Qué sigue? . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 116
6 ÍNDICE GENERAL
Prefacio
Estas notas de clase se escriben como apoyo para dos cursos: uno de la parte final
de un pregrado en matemáticas y otro de la parte inicial de una maestría en matemá-
ticas. En la Universidad Nacional de Colombia, sede Medellín, estos cursos se llaman
Algebra lineal aplicada y Aspectos numéricos del álgebra lineal respectivamente.
Del lector de este manuscrito se espera que haya tenido contacto con el álgebra
lineal y no necesariamente con el álgebra lineal numérica. Libros adecuados para ese
primer contacto son [S], [Lang], [M], [ND] y [P]. El libro de Noble y Daniel [ND] es
un libro introductorio de álgebra lineal numérica. El libro de Meyer [M] también es
un libro de álgebra lineal numérica pero es un poco más avanzado que el de Noble
y Daniel. Los otros tres libros son libros de texto para cursos de álgebra lineal de
segundo o tercer semestre de pregrado.
Desde ya advertimos que en estas notas no se citan todas las referencias bibliográ-
ficas que aparecen al final, dicha lista sirve de apoyo a las notas y es una invitación
a los lectores a conocer otros materiales sobre el tema.
En varias partes del manuscrito hay unas notas cortas llamadas Nota MATLAB.
Con ellas se busca motivar el uso de este software para el aprendizaje del álgebra li-
neal numérica. Recordemos que MATLAB es un nombre que surgió de la combinación
en una sola palabra de las tres primeras letras de las palabras Matrix Laboratory.
Para el autor, la mejor manera de aprender a usar a MATLAB es por medio de
ejemplos obtenidos de clase o libros o de algún otro documento físico o virtual. El
profesor Cleve Moler, inventor de MATLAB y fundador de la empresa Mathworks,
escribió dos excelentes libros que usan MATLAB y que se pueden comprar en SIAM
o descargarlos gratuitamente de https://www.mathworks.com/moler.html
Estas notas están divididas en dos partes. La primera, llamada Algebra lineal,
contiene material teórico fundamental. Allí hay definiciones y resultados relacionados
con vectores, matrices, transformaciones lineales, determinantes, sistemas lineales,
entre otros. La segunda, llamada Análisis matricial, incluye ejemplos y algorit-
mos seleccionados entre los muchos que hay en álgebra lineal numérica. La selección
se hizo con la idea de apoyar el trabajo interdisciplinario que tiene al álgebra li-
7
8 ÍNDICE GENERAL
neal numérica como uno de sus pilares. Por tanto el material tratado tiene que ver
con condicionamiento, estabilidad, factorizaciones de matrices, solución de sistemas
lineales por métodos directos y en menor medida, cálculo de valores propios.
En el prefacio de [S], Gilbert Strang escribe: Personalmente considero que mucha
más gente necesita álgebra lineal que cálculo. Hay muy buena justificación para esta
afirmación y de hecho el autor de estas notas tiene la misma creencia. Es que muchos
de los problemas de ingeniería y ciencias se modelan por medio de ecuaciones dife-
renciales ya sea ordinarias o parciales y lo usual para un investigador de ciencias o
ingeniería es resolver numéricamente estas ecuaciones. Eso significa llegar a un pro-
blema de álgebra lineal numérica. Más adelante en estas notas se explicará un poco
sobre este proceso.
Escribir este manuscrito es un anhelo que el autor tiene desde hace varios años.
El trabajo [AM] sobre métodos iterativos, escrito con Carlos D. Acosta, se puede
considerar predecesor de este libro pero trata temas más avanzados. Los trabajos
[Layton] y [vdG] son ejemplos inspiradores con mayor alcance que el presente trabajo
pero con el mismo estilo de distribución via web, siempre buscando llegar a una
amplia audiencia.
La escritura de este manuscrito se hizo durante el semestre sabático del que disfru-
tó el autor en el segundo semestre de 2016. Sea esta la oportunidad de agradecer a la
Universidad Nacional de Colombia por haberme concedido dicho semestre sabático.
Este libro está en construcción y asi se ofrece en la página web personal del autor
que se puede encontrar en la lista de docentes de la Escuela de Matemáticas de
la Universidad Nacional de Colombia, sede Medellín. Actualmente, en diciembre de
2016, las direcciones URL son:
Página web del autor http://www.unalmed.edu.co/~cemejia
Página en Google Sites: https://sites.google.com/unal.edu.co/cemejia/home
Página web de la Escuela de Matemáticas:
http://ciencias.medellin.unal.edu.co/escuelas/matematicas/
Nomenclatura
Escalares
Letras griegas minúsculas α, β, . . . o letras latinas o griegas minúsculas con sub-
índices xi , aij , λk , . . .
Vectores
Letras latinas minúsculas a, b, x, y, . . . o letras latinas mayúsculas con subíndice
Ai , A∗j , Ai∗ , . . . o la letra e minúscula con un subíndice por ejemplo ej para vectores
canónicos.
Matrices
Letras latinas, griegas o en estilo caligráfico mayúsculas A, B, Λ, V, . . . y ocasio-
nalmente, si no hay confusión, también se usan esas mismas letras con un subíndice.
Los elementos de la matriz A son aij o Aij , sus columnas son A∗j , sus filas son
Ai∗ . Los asteriscos se suprimen cuando no hay confusión, por ejemplo cuando en el
contexto solo se están utilizando las columnas de la matriz.
Los elementos del vector x son xi .
Principales superíndices y subíndices: i, j, k, l, m, n.
El número 0 se usa como escalar, vector y matriz.
La matriz identidad se denota I o si hay peligro de confusión, In .
LD: Linealmente dependiente.
LI: Linealmente independiente.
CL: Combinación lineal.
TL: Transformación lineal.
σ (A) : Espectro de la matriz A.
ρ (A) : Radio espectral de la matriz A.
9
10 ÍNDICE GENERAL
ej : Vector canónico. Es el vector columna cuya componente j − ésima es 1 y las

demás son 0.
det (A) : Determinante de la matriz A.
traza (A) : Traza de A.
AT : Traspuesta de la matriz A.
A∗ : Conjugada traspuesta de la matriz A.
condp (A) = kp (A) : Número de condición de la matriz A con respecto a la
p−norma, 1 ≤ p ≤ ∞. Se omite el subíndice si no hay confusión.
A−1 : Inversa de la matriz A.
rank (A) : Rango (rank) de A.
R (A) : Rango (range) de A. También se llama subespacio imagen de A.
ker (A) : Núcleo de A.
Si S es un subconjunto no vacío de un espacio vectorial, el subespacio generado
por S se denota Sg (S) o Sp (S) o Span (S) .
pλ (t) = det (A − λI) : Polinomio característico de A. Sus raices son los valores
propios de A.
Cp : Matriz compañera del polinomio p.
Matriz unitaria U ∈ Cn×n : Sus columnas constituyen una base ortonormal de
Cn .
Matriz ortogonal U ∈ Rn×n : Sus columnas constituyen una base ortonormal de
n
R .
Matriz de Hessenberg: H es de Hessenberg si hij = 0 para i > j + 1.
Matriz de Householder: U = I − 2ww T con w ∈ Rn y kwk2 = 1.
Parte I
Algebra lineal
11
Capítulo 1
Matrices y vectores
En el primer capítulo introducimos muchas de las notaciones y definiciones que

utilizaremos. Seguimos de cerca el libro [I] de Ilse Ipsen, profesora de North Carolina
State University. La profesora Ipsen regala una copia electrónica de su libro en su
página web http://www4.ncsu.edu/~ipsen.
Los vectores y las matrices son los principales objetos con los que trabajamos en
todo el libro. Entre los años 200 y 100 a. de C. hay evidencias de usos de matrices en
China y de forma menos clara en Babilonia. Los determinantes se anuncian en el libro
de Cardano Ars Magna (1645) y aparecen con más claridad hacia 1683 tanto en Japón
como en Europa. Leibniz es uno de los precursores para el concepto de determinante
y solo hacia 1812 con Cauchy la noción de determinante quedó establecida como la
tenemos ahora.
Aunque la eliminación gaussiana se vislumbra en una obra china titulada Nine

Chapters on the Mathematical Art del año 200 a. de C., es el propio Gauss quien la
usa por primera vez hacia 1809. El primero que usa el término matriz es Sylvester en
1850. La segunda mitad del siglo XIX es la época de avances en la teoría de matrices
y por tanto del álgebra lineal. Entre los matemáticos que más contribuyeron están
Sylvester, Cayley, Hamilton, Frobenius y Jordan.
Los datos históricos fueron tomados de

http://www-history.mcs.st-andrews.ac.uk/HistTopics/Matrices_and_determinants.html
13
14 CAPÍTULO 1. MATRICES Y VECTORES
1.1. Matrices
Llamamos matriz m × n a un rectángulo de números
 
a11 · · · a1n
 .. 
A =  ... . 
am1 · · · amn
con m filas y n columnas. Los elementos aij son escalares, es decir, elementos de
R el conjunto de los números reales o elementos de C, el conjunto de los números
complejos. Se escribe A ∈ Rm×n en el primer caso y A ∈ Cm×n en el segundo. Si
m = n, se dice que A es una matriz cuadrada de orden n.
Un vector fila es una matriz 1 × n y un vector columna es una matriz m × 1. El
conjunto de vectores columna con m componentes reales se denota Rm . Análogamente
se define Cm .
1.1.1. Submatrices
Una submatriz de una matriz A ∈ Rm×n es una matriz obtenida de A suprimiendo
filas y/o columnas. Más precisamente, dados dos conjunto de números enteros 1 ≤
i1 < i2 < . . . < ik ≤ m y 1 ≤ j1 < j2 < . . . < jl ≤ n entonces la matriz k × l
 
ai1 ,j1 ai1 ,j2 · · · ai1 ,jl
 ai2 ,j1 ai2 ,j2 · · · ai2 ,jl 
 
 .. .. .. 
 . . . 
aik ,j1 aik ,j2 · · · aik ,jl
es una submatriz de A. Se llama principal si es cuadrada y sus elementos diagonales

son elementos diagonales de A, es decir, k = l y ir = jr para r = 1, 2, . . . , k.
1.1.2. Matrices especiales

1. Matriz nula
La matriz nula m × n se denota 0 y todos sus elementos son cero. Si en una
matriz A alguno de sus elementos, digamos alj , es no nulo, entonces A 6= 0.
1.1. MATRICES 15
2. Matriz identidad
La matriz identidad de orden n, denotada In , es la matriz cuadrada
 
1
 .. 
 . 
1
con unos en la diagonal y ceros en todas las demás posiciones. En particular,
I1 = 1. El subíndice se suprime si no hay confusión y se habla simplemente de
la matriz identidad I.
3. Matriz diagonal
La matriz cuadrada A es diagonal si aij = 0 para todo i 6= j. En particular la
matriz nula es diagonal.
4. Matriz triangular superior (inferior)
La matriz cuadrada A es triangular superior (inferior) si aij = 0 para i > j
(i < j).
5. Matriz tridiagonal
La matriz cuadrada A es tridiagonal si aij = 0 para |i − j| > 1.
6. Matriz de permutación
Una matriz cuadrada es de permutación si contiene un solo uno en cada fila y
en cada columna. Aparte de esos unos, todos los demás elementos son cero. Se
puede pensar que es cualquier permutación de la matriz identidad.
7. Vectores canónicos
Las columnas de In se denominan vectores canónicos ei y por tanto
In = (e1 e2 · · · en ) .
1.1.3. Multiplicación por escalar

Sean A ∈ Rm×n y ρ ∈ R. El elemento (i, j) de la matriz ρA es
(ρA)ij = ρaij .
Si ρ = 0 entonces se obtiene la matriz nula.
La multiplicación por escalar es asociativa en el siguiente sentido
(αρ) A = α (ρA) .
1.1.4. Traza de una matriz cuadrada

Sea A ∈ Rn×n . Su traza es el número
n
X
traza (A) = aii .
i=1
Es decir, es la suma de los elementos en la diagonal de A.
1.1.5. Suma de matrices

Si A ∈ Rm×n y B ∈ Rm×n entonces
(A + B)ij = aij + bij .
La suma está bien definida en el sentido de la clausura: la suma de dos matrices de

Rm×n es una matriz de Rm×n . Cumple las siguientes condiciones:
1. Asociativa: Si A ∈ Rm×n , B ∈ Rm×n y C ∈ Rm×n , entonces
(A + B) + C = A + (B + C) .
2. Elemento idéntico: La matriz 0 ∈ Rm×n es tal que para toda matriz A ∈ Rm×n
se cumple
0 + A = A + 0 = A.
3. Conmutativa: Si A ∈ Rm×n y B ∈ Rm×n , entonces
A + B = B + A.
4. Distributivas: Si A ∈ Rm×n , B ∈ Rm×n , α ∈ R y β ∈ R entonces
(α + β) A = αA + βA y α (A + B) = αA + αB.
1.1.6. Combinación lineal

Dados n vectores V1 , V2 , . . . , Vn y n escalares α1 , α2 , . . . , αn , n ≥ 1, a la suma
n
X
αj V j
j=1
se le llama combinación lineal de los vectores Vi con coeficientes αi .

1.2. MULTIPLICACIÓN DE MATRIZ POR VECTOR 17
1.1.7. Producto interno o producto punto

 
w1
 w2 
 
Sean v = (v1 v2 · · · vn ) y w =  ..  . El producto interno de v y w es
 . 
wn
n
X
vw = vj wj .
j=1
Otras notaciones para el producto interno son v · w ó (v, w).
1.2. Multiplicación de matriz por vector

1.2.1. Matriz por vector columna
Sea A ∈ Rm×n con filas Ai∗ y columnas A∗i y sea v ∈ Rn con elementos vi .
Entonces    
A1∗ v1
 ..   .. 
A =  .  = (A∗1 A∗2 · · · A∗n ) , v= . .
Am∗ vn
Distinguimos dos formas de considerar el producto Av :
1. Av es un vector columna de productos internos: (Av)i = Ai∗ v.
n
X
2. Av es una combinación lineal de las columnas de A : Av = vi A∗i .
i=1
Ejemplo 1.1. La j−ésima columna de A se puede pensar como el producto A∗j =

Aej .
1.2.2. Vector fila multiplicado por matriz

Sean A como arriba y w = (w1 w2 · · · wm ) . En este caso también vemos el
producto de fila por matriz de dos maneras distintas:
1. wA es un vector fila de productos internos: wA = (wA∗1 wA∗2 · · · wA∗n ) .
m
X
2. wA es una combinación lineal de las filas de A : wA = vi Ai∗ .
i=1
1.2.3. Producto externo (outer product)

 
v1
 
Sean v =  ...  y w = (w1 w2 wm ). El producto vw es la matriz
vm
 
v1 w1 · · · v1 wn
 ..  .
vw =  ... . 
vm w1 · · · vm wn
Ejemplo 1.2. Una matriz de  Vandermonde

 de orden n con potencias de x se puede
1
 .. 
ver como el producto externo  .  (1 x x2 · · · xn−1 ) .
1
1.3. Nota MATLAB

Hay dos formas básicas de buscar ayuda en MATLAB. La primera es por medio
del cuadro de búsqueda que hay en la esquina derecha superior del escritorio y la
segunda es en la línea de comandos, por medio de la orden doc. Dos ensayos reco-
mendados son: En la ventana de búsqueda escribir gallery y en la línea de comandos
escribir doc gallery.
En general, la búsqueda en la línea de comandos por medio de la orden doc
es más efectiva. En el caso de gallery puede verse que doc lleva directamente a la
documentación sobre Gallery Test Matrices y en cambio si se usa la ventana de
búsqueda aparecen varios documentos entre los cuales el listado de primero es el que
nos interesa.
En estas notas sobre el uso de MATLAB, trabajamos matrices pequeñas, por
ejemplo matrices cuadradas de 8 filas y columnas, aun si se toman de Gallery. El
lector puede ensayar con tamaños mayores si desea.
a= gallery(’jordbloc’,8,4); % bloque de Jordan
c = gallery(’moler’,8); % matriz simetrica definida positiva
j = 3*a; % producto por escalar
f = a + c; % suma de matrices
v = trace(c); % traza
x = [4,1,-1,zeros(1,5)]; % vector fila
xa = x*a; % vector fila por matriz
1.3. NOTA MATLAB 19
y = xa’; % vector columna como traspuesta de vector fila

z = c*y; % producto de matriz por vector columna
w = z*x; % producto externo
s = x*z; % producto interno
Capítulo 2
Matrices y determinantes
En el segundo capítulo sobre matrices y vectores definimos varias operaciones con

matrices y algunas clases especiales de matrices. También hacemos referencia a los
determinantes que en el resto del libro serán utilizados como herramienta. El material
es elemental y por tanto aparece prácticamente en todos los libros del tema. Para
escribir estas notas se consultaron primeros capítulos de [I], [L], [Layton] y [Ortega].
2.1. Producto de matrices

Sean A ∈ Rm×n y B ∈ Rn×r . El número de columnas de A debe ser igual al
número de filas de B. Su producto es la matriz D = AB ∈ Rm×r cuyo elemento (i, j)
es
Xn
dij = aik bkj .
k=1
Pensando en las filas y columnas de las matrices A y B, este producto

 se puede
A1∗
 A2∗ 
 
presentar de varias formas: Sean A = (A∗1 A∗2 · · · A∗n ) =  ..  y B =
 . 
Am∗
 
B1∗
 B2∗ 
 
(B∗1 B∗2 · · · B∗r ) =  .. , es decir, la columna j − ésima de A es A∗j , la
 . 
Bn∗
fila i − ésima de A es Ai∗ y análogamente para B.
21
22 CAPÍTULO 2. MATRICES Y DETERMINANTES
1. AB = (AB∗1 AB∗2 ··· AB∗r ) .

 
A1∗ B
 A2∗ B 
 
2. AB =  .. .
 . 
Am∗ B
3. (AB)ij = Ai∗ B∗j . Matriz de productos internos.

n
X
4. AB = A∗i Bi∗ . Suma de productos externos.
i=1
2.1.1. Propiedades del producto de matrices

1. Im A = AIn = A. Identidad.
2. A (BC) = (AB) C. Asociatividad.
3. A (B + C) = AB + AC y (A + B) C = AC + BC. Distributividades.
Ejemplo 2.1. El producto no es conmutativo:

0 1 0 0 0 2 0 0 0 1 0 0
= pero = .
0 0 0 2 0 0 0 2 0 0 0 0
Nota 2.2. El producto de A consigo misma n veces se denota An .
Definición 2.3. Una matriz cuadrada es
1. Una involución si A2 = I.
2. Idempotente o proyector si A2 = A.
3. Nilpotente si existe un entero positivo s tal que As = 0.
Definición 2.4. Sea A ∈ Rm×n . Su traspuesta, denotada AT , es la matriz cuya

componente (i, j) es AT ij = aji . Es decir, las filas (columnas) de A son las columnas
(filas) de AT .
Definición 2.5. Sea A ∈ Cm×n . Su traspuesta conjugada, denotada A∗ , es la matriz

cuya componente (i, j) es (A∗ )ij = aji el conjugado de aji .
2.1. PRODUCTO DE MATRICES 23
2.1.2. Propiedades de las traspuestas

1. Para matrices reales, la conjugada traspuesta es lo mismo que la traspuesta.
T
2. AT = A, (A∗ )∗ = A.
3. (βA)T = βAT , (βA)∗ = βA∗ .

4. (A + B)T = AT + B T , (A + B)∗ = A∗ + B ∗ .
5. (AB)T = B T AT , (AB)∗ = B ∗ A∗ .
2.1.3. Propiedades de productos de vectores

   
x1 y1
   
Sean x =  ...  , y =  ...  vectores de Cn . Por definición de producto
xn yn
Xn n
X
interno, x∗ y = xi yi y xT y = xi yi .
i=1 i=1
   
x1 y1
   
Definición 2.6. Decimos que los vectores x =  ...  y y =  ...  de Cn son
xn yn
∗
ortogonales si x y = 0. Para vectores reales la condición es la misma, pero se prefiere
la siguiente escritura: xT y = 0.

x1
Si z = x1 + ix2 ∈ C entonces x = ∈ R2 y además, i2 = −1. El valor
x2
p √
absoluto o módulo del número complejo z es |z| = x21 + x22 = xT x. Para los
vectores x y y de Cn también se cumplen los siguientes enunciados:
1. y ∗x = (x∗ y).
2. y T x = xT y.
3. x∗ x = 0 si y solo si x = 0.
4. Si x es real entonces xT x = 0 si y solo si x = 0.
n
X
5. In = ei eTi .
i=1
Definición 2.7. Una matriz A ∈ Cn×n es

1. Simétrica si AT = A.
2. Hermitiana si A∗ = A.
3. Antisimétrica (skew-symmetric) si AT = −A.
4. Antihermitiana (skew-hermitian) si A∗ = −A.
Proposición 2.8. Si A ∈ Cm×n entonces AAT y AT A son simétricas pero AA∗ y
A∗ A son hermitianas. Además, A + AT es simétrica y A + A∗ es hermitiana.
Definición 2.9. Una matriz A ∈ Rn×n es simétrica definida positiva si es si-
métrica y para todo vector no nulo x ∈ Rn , xT Ax > 0. Una matriz A ∈ Cn×n es
hermitiana definida positiva si es hermitiana y para todo vector no nulo x ∈ Cn ,
x∗ Ax > 0. Se dice que son semidefinidas positivas si en lugar del signo > se usa
el signo ≥ .
2.2. Determinantes
El determinante de una matriz cuadrada A de Rn×n o Cn×n se define asi:
X
det (A) = εσ aσ(1)1 aσ(2)2 · · · aσ(n)n ,
σ∈Sn
donde εσ = 1 si σ es par y εσ = −1 si σ es impar y además denotamos Sn al grupo

de permutaciones de los números {1, 2, . . . , n} .
Esta definición se puede encontrar en [C], [Ortega] y [H]. Esta última fuente
contiene un tratamiento completo sobre el tema.
Otras definiciones para el determinante de una matriz se basan en expansión en
cofactores. Más precisamente, para la matriz A denotamos Aij a la matriz obtenida
de A suprimiendo la fila i y la columna j. Al número (−1)i+j det (Aij ) se le denomina
cofactor del elemento aij . Las definiciones basadas en cofactores son:
n
X
det (A) = aij (−1)i+j det (Aij ) para cada i (2.1)
j=1
Xn
det (A) = aij (−1)i+j det (Aij ) para cada j.
i=1
Las propiedades que presentamos a continuación, son las que queremos destacar
sobre determinantes.
2.3. INVERSA DE UNA MATRIZ 25
Teorema 2.10. Propiedades de los determinantes
1. Si dos de las filas (columnas) de A son iguales, entonces det (A) = 0.
2. Si A tiene una fila (columna) nula, entonces det (A) = 0.
3. Si se intercambian dos filas (columnas) de A, entonces el determinante de la

nueva matriz es − det (A) .
4. Si se multiplica una fila (columna) de A por un escalar β el determinante de

la nueva matriz es β det (A) .

5. det AT = det (A) , det (A∗ ) = det (A).
6. Si un múltiplo escalar de una fila (columna) de A se suma a otra fila (columna),

el determinante se conserva.
7. Si A y B son matrices cuadradas del mismo tamaño, det (AB) = det (A) det (B) .
n
Y
8. Si A es triangular superior o triangular inferior entonces det (A) = aii .
i=1
2.3. Inversa de una matriz

Una matriz cuadrada A, real o compleja, tiene inversa (también se dice es no
singular o es invertible) si existe una matriz del mismo tamaño denotada A−1 tal que
AA−1 = A−1 A = I.
Si A no tiene inversa se dice que es singular o no invertible.

Proposición 2.11. La inversa es única.
Demostración. Sea A una matriz cuadrada y sean B y C matrices del mismo tamaño
tales que
AB = BA = I y AC = CA = I.
Entonces
B = BI = B (AC) = (BA) C = IC = C.
Proposición 2.12. Sean A ∈ Rn×n , x, b ∈ Rn .

1. Suponga que x 6= 0. Se cumple:
Si Ax = 0 entonces A es singular.
Si A es no singular entonces Ax 6= 0.
2. Supongamos Ax = b. Si b 6= 0 entonces x 6= 0.
Esta es la primera aparición de la ecuación Ax = b en estas notas. Una ecuación de

esta forma se denomina indistintamente sistema de ecuaciones lineales, sistema
lineal o ecuación lineal. Generalmente se exige que A sea una matriz invertible
pues tal condición garantiza una única solución para el siguiente problema: Dados
A no singular y un vector compatible b, hallar x tal que Ax = b. La existencia y
unicidad de solución para el sistema Ax = b se consideran posteriormente.
Proposición 2.13. Si A es invertible entonces AT y A∗ son invertibles. Sus inversas

son −1 T ∗
AT = A−1 y (A∗ )−1 = A−1 .
Proposición 2.14. Una matriz idempotente es la identidad o es singular.
Demostración. Supongamos que A es idempotente, es decir, A2 = A. Entonces 0 =

A (A − I) . Si A = I, se cumple la igualdad. Si A 6= I, hay una columna no nula en
A−I. Si la llamamos x, entonces Ax = 0 y por proposición anterior, A es singular.
2.4. Nota MATLAB

a = gallery(’jordbloc’,8,4); % bloque de Jordan
c = gallery(’moler’,8); % matriz simetrica definida positiva
a1 = inv(a); % inversa
g = a’; % traspuesta
p = a*c; % producto de matrices
q = a.*c; % producto componente a componente
h = det(a); % determinante
Capítulo 3
Espacios vectoriales
En este capítulo nos concentramos en espacios vectoriales de dimensión finita y

lo que se incluye es una breve revisión de conceptos básicos. Para estudiar espacios
vectoriales en un contexto más amplio se recomienda [H]. Las principales referencias
que consultamos para escribir este capítulo son [HK], [I], [Lang], [L] y [M].
Probablemente el primer libro que se refiere a las propiedades que cumple un
espacio vectorial lo publicó Peano en 1888. En ese libro también aparecen por primera
vez las operaciones unión e intersección para conjuntos y el signo ∈ .
Para saber más de esta historia recomendamos
http://www-history.mcs.st-andrews.ac.uk/HistTopics/Abstract_linear_spaces.html
Advertencia: A los elementos de espacios vectoriales se les llama vectores.
3.1. Definición
Un espacio vectorial sobre el campo de los números reales es un conjunto V junto
con dos operaciones
+ : V × V −→ V y · : R×V −→ V
tales que
1. (V, +) es un grupo abeliano, es decir:
La suma es asociativa: v + (w + z) = (v + w) + z para todo v, w, z ∈ V.

La suma tiene un elemento idéntico: Existe un único vector 0 ∈ V tal que
v + 0 = 0 + v = v para todo v ∈ V.
27
28 CAPÍTULO 3. ESPACIOS VECTORIALES
Para cada vector v ∈ V existe un único vector, su inverso aditivo −v, tal
que v + (−v) = (−v) + v = 0.
La suma es conmutativa: v + w = w + v para v, w ∈ V.
2. (αβ) v = α (βv) para α, β ∈ R y v ∈ V.
3. (α + β) v = αv + βv para α, β ∈ R y v ∈ V.
4. α (v + w) = αv + αw, para α ∈ R y v, w ∈ V.
5. 1v = v para v ∈ V.
3.2. Ejemplos de espacios vectoriales

1. Rn es un espacio vectorial sobre R con la suma y el producto por escalar usuales.
2. Rm×n es un espacio vectorial sobre R con la suma y el producto por escalar

usuales.
3. V es el conjunto de las funciones f de [0, 1] en R con la suma y el producto

por escalar usuales.
4. El conjunto de los polinomios con coeficientes reales con la suma y el producto

por escalar usuales.
3.3. Subespacios
Sea V un espacio vectorial sobre R y sea W ⊆ V , W 6= φ. El subconjunto W
es un subespacio de V si con las operaciones heredadas de V también es un espacio
vectorial. De forma equivalente puede decirse que W es un subespacio de V si y solo
si αv + βw ∈ W para todo v, w ∈ W, α, β ∈ R.
Esta definición significa que lo único que debe revisarse es la clausura de las
operaciones en el subconjunto.
3.4. Ejemplos de subespacios

1. Para cada espacio vectorial V los subconjuntos {0} y V son los subespacios
triviales de V .
3.5. INDEPENDENCIA LINEAL 29
2. W = {A ∈ Rn×n : A es simétrica} es un subespacio de Rn×n .
3. Supongamos que Y es el conjunto de las funciones continuas f de [0, 1] en R

con la suma y el producto por escalar usuales. Entonces Y es un subespacio
del espacio V descrito en 3.2.
4. Sea S un subconjunto no vacío de un espacio vectorial V. El conjunto S puede

ser finito o infinito. El conjunto de todas las combinaciones lineales (ver 1.1.6)
con elementos de S es un subespacio de V llamado subespacio generado por S
y se denota Sg (S) ó Span (S) ó Sp (S) .
3.5. Independencia lineal

3.5.1. Conjunto linealmente dependiente (LD)
Un subconjunto S no vacío de un espacio vectorial V se llama linealmente
dependiente si existen vectores V1 , V2 , . . . , Vk ∈ S y escalares α1 , α2 , . . . , αk ∈ R, no
k
X
todos nulos, tales que αj Vj = 0. Esto es lo que se llama una combinación lineal
j=1
no trivial para el vector 0.
3.5.2. Conjunto linealmente independiente (LI)

Un subconjunto S no vacío de un espacio vectorial V se llama linealmente
independiente si no es linealmente dependiente. En otras palabras, en cualquier
k
X
combinación lineal del cero αj Vj = 0 con vectores V1 , V2 , . . . , Vk ∈ S y escalares
j=1
α1 , α2 , . . . , αk ∈ R, obligatoriamente se cumple que todos los escalares son nulos:
α1 = α2 = . . . = αk = 0.
3.5.3. Ejemplos
1. Si S es un subconjunto de un espacio vectorial V y 0 ∈ S, entonces S es un
subconjunto linealmente dependiente.
2. En Rn el conjunto B = {e1 , e2 , . . . , en } de los vectores canónicos (ver 1.1.2,

numeral 7) de Rn es linealmente independiente.

1 1
3. S = , es un subconjunto linealmente independiente de R2 .
0 2
3.6. Base de un espacio vectorial

Un subconjunto B de un espacio vectorial V es una base de V si
1. B es un subconjunto linealmente independiente de V.
2. Sg (B) = V.
Nota 3.1. En el ejemplo 2 de la lista 3.5.3 describimos la base canónica de Rn .

Nota 3.2. Existen espacios vectoriales en los que hay bases infinitas. Por ejemplo, el
conjunto B = {1, x, x2 , . . .} es una base del espacio de los polinomios con coeficientes
reales con la suma y el producto por escalar usuales.
Nota 3.3. En lo que sigue, todos los espacios vectoriales con los que se trabaja son
finito dimensionales, a no ser que se diga lo contrario.
Teorema 3.4. Sea V un espacio vectorial con bases B1 y B2 de n y m elementos
respectivamente. Entonces n = m.
Demostración. Ver Teorema 3, cap. II de [Lang].
Definición 3.5. Sea V un espacio vectorial que tiene una base con n elementos. Al
número n se le llama dimensión de V y se le denota dim (V ) .
Nota 3.6. Se dice que el espacio vectorial trivial {0} tiene dimensión 0.
Teorema 3.7. Sea V un espacio vectorial de dimensión n. Sea W un subespacio
distinto de {0} de V. Entonces dim (W ) ≤ n.
Demostración. Ver Teorema 6, cap. II de [Lang].
3.7. Sumas e intersecciones de subespacios

Sean W y Y subespacios de un espacio vectorial V.
3.7.1. Suma
La suma de W y Y es W + Y = {w + y : w ∈ W, y ∈ Y }
3.7. SUMAS E INTERSECCIONES DE SUBESPACIOS 31
3.7.2. Intersección
La intersección de W y Y es W ∩ Y = {x : x ∈ W y x ∈ Y }
Teorema 3.8. La suma y la intersección de subespacios de un espacio vectorial V

son subespacios de V.
Nota 3.9. La unión de subespacios de V no necesariamente es un subespacio de V.
Definición 3.10. Sean W y Y subespacios de un espacio vectorial V y sea Z =

W + Y. Si se cumple que W ∩ Y = {0} se dice que Z es la suma directa de W y Y
y se escribe Z = W ⊕ Y.
Teorema 3.11. Sean W y Y subespacios de un espacio vectorial de dimensión finita

V . Supongamos que Z = W ⊕ Y. Entonces
1. Cada z ∈ Z se escribe de forma única como z = w + y con w ∈ W y y ∈ Y.
2. dim (Z) = dim (W ) + dim (Y ) .
Demostración. Ver Theorem 2.26 en [L] y Teorema 6 del capítulo 2 de [HK].

Capítulo 4
Sistemas lineales
En este capítulo consideramos la solución de sistemas lineales e introducimos el

método directo por excelencia que es la eliminación de Gauss. Este es el algoritmo
básico del álgebra lineal y sirve para resolver exactamente el sistema lineal Ax = b
para A matriz no singular. Es más utilizado cuando la matriz A es densa, es decir,
la mayoría de sus elementos son no nulos.
Para matrices en las que la mayoría de sus elementos son ceros, que son llamadas
ralas o poco densas (sparse), es más común utilizar métodos iterativos. Con
estos métodos se pretende aproximar la solución del sistema Ax = b, no se pretende
encontrarla exactamente. Son especialmente útiles cuando las matrices son grandes
y ralas. Recomendamos [Layton] para apreciar el efecto llamado Maldición de la
dimensión (Curse of dimensionality) que es el que obliga a recurrir a estrategias
iterativas y [AM] para una introducción a estos métodos.
En [SV] hay un recuento histórico sobre métodos iterativos escrito por dos de los
principales investigadores del tema. Se reportan avances en el tema realizados durante
el siglo XX. Para una historia de los métodos directos y en particular de la eliminación
de Gauss, recomendamos la monografía que Carl D. Meyer ofrece en su página web
http://meyer.math.ncsu.edu/Meyer/PS_Files/GaussianEliminationHistory.pdf
Para escribir este capítulo se consultaron [I], [L], [M] y [Ortega].
4.1. Eliminación de Gauss

Sean A ∈ Rn×n , x, b ∈ Rn . Interesa conocer condiciones necesarias y suficientes
para existencia y unicidad de soluciones para el sistema
Ax = b. (4.1)
33
34 CAPÍTULO 4. SISTEMAS LINEALES
Escrito por filas, el sistema (4.1) es

n
X
aij xj = bi , i = 1, . . . n (4.2)
j=1
Introducimos la eliminación de Gauss por medio de la demostración de un teorema

central a la manera clásica de [Ortega].
Teorema 4.1. Si det (A) 6= 0 entonces el sistema (4.1) tiene una única solución.
Demostración. Empezamos demostrando que existe una solución. Supongamos que
a11 6= 0. Multiplicamos la primera ecuación de (4.2) por a−1
11 a21 y la sustraemos de
la segunda ecuación. Enseguida multiplicamos la primera ecuación por a−1 11 a31 y la
sustraemos de la tercera ecuación. Continuamos en esta forma hasta llegar al sistema
a11 x1 + a12 x2 + · · · + a1n xn = b1

(1) (1) (1)
a22 x2 + · · · + a2n xn = b2
.. .. .. (4.3)
. . .
(1) (1) (1)
an2 x2 + · · · + ann xn = bn
donde los (1)′ s en los superíndices indican que pudieron haber cambiado. El resultado
es que se eliminó la incógnita x1 de las últimas n−1 ecuaciones. El sistema modificado
(4.3) y el sistema (4.2) tienen el mismo vector x como solución. Enseguida suponemos
(1)
que a22 6= 0 y repetimos el mismo proceso para el sistema de n − 1 ecuaciones con
n − 1 incógnitas para eliminar la incógnita x2 de las últimas n − 2 ecuaciones. La
repetición de este procedimiento lleva a un sistema con matriz triangular que tiene
la siguiente forma:
a11 x1 +a12 x2 · · · +a1n xn = b1

(1) (1) (1)
a22 x2 · · · +a2n xn = b2
.. .. .. (4.4)
. . .
(n−1) (n−1)
ann xn = bn
El sistema modificado (4.4) y el sistema (4.2) tienen el mismo vector x como solución.
(n−1)
Finalmente, suponemos que ann 6= 0. De la última ecuación en (4.4) despejamos
para obtener el valor de xn . Conocido xn , de la penúltima ecuación, encontramos
el valor de xn−1 . Continuamos este procedimiento de sustitución regresiva (back-
substitution) y obtenemos el vector x que es solución del sistema (4.4) y por tanto
del sistema original (4.2). Como en cada uno de los pasos de sustitución regresiva la
4.1. ELIMINACIÓN DE GAUSS 35
solución obtenida es única, entonces encontramos una única solución para el sistema
(1) (n−1)
(4.2) siempre que los elementos a11 , a22 , . . . , ann ,llamados pivotes, sean no nulos.
Concluimos la demostración utilizando la hipótesis det (A) 6= 0 para garantizar que
siempre se pueden econtrar pivotes no nulos. Supongamos que a11 = 0. Buscamos en
la primera columna un elemento no nulo, digamos ak1 e intercambiamos las filas 1 y
k en los dos lados de la igualdad (4.2). El elemento no nulo de la primera columna
existe porque si la columna fuera nula entonces det (A) = 0. El sistema con las filas
intercambiadas tiene la misma solución x que el sistema original excepto el orden
de las componentes en el vector x. Más precisamente, si pensamos en (4.2) como un
sistema de n ecuaciones, el sistema con dos ecuaciones intercambiadas tiene las mis-
mas soluciones que el sistema original. Para considerar el análogo al sistema reducido
(1)
(4.3), hacemos el mismo razonamiento. Si a22 = 0, escogemos un elemento no nulo,
(1)
digamos ak2 con k > 2, e intercambiamos las filas 2 y k. La existencia del elemento no
nulo en las filas siguientes a la segunda de la segunda columna de (4.3) se justifica de
nuevo por la hipótesis det (A) 6= 0. Si la columna 1 de la submatriz (n − 1) × (n − 1)
de A, consistente en las filas 2, . . . , n y las columnas 2, . . . , n, fuera nula, entonces
el determinante de la submatriz sería 0 y por la expansión en cofactores (2.1) para
la columna 1 de la matriz en (4.3), se tendría det (A) = 0. Debe tenerse en cuenta
que si hay intercambios de filas el determinante cambia a lo sumo en el signo. Este
(i−1)
argumento se puede continuar y concluir que los pivotes aii , i = 3, . . . , n son no
nulos o se pueden hacer no nulos por intercambio de filas.
Enseguida relacionamos vectores canónicos con inversas de matrices y sistemas
lineales.
Consideremos de nuevo el sistema (4.1) y supongamos que det (A) 6= 0. Ahora
definamos n sistemas lineales
Ax(i) = ei , i = 1, . . . , n. (4.5)
Basados en el teorema (4.1), estos

sistemas tienen únicas soluciones x(i) para cada
i. La matriz C = x(1) x(2) . . . x(n) que tiene estas soluciones como columnas, cumple
AC = I. (4.6)
Por numeral 7 de Teorema 2.10, 1 = det (I) = det (A) det (C) y como det (A) 6= 0
entonces det (C) 6= 0.
Lo que hacemos ahora es plantear sistemas de la forma (4.5) para la matriz C
que tiene determinante no nulo. Llegamos a una matriz X tal que
CX = I. (4.7)
Al premultiplicar (4.7) por A y usar asociatividad, encontramos que X = A. Recor-

dando la definición y la unicidad de la inversa de una matriz, concluimos que C =
A−1 . De manera que hemos demostrado
Teorema 4.2. Sea A ∈ Rn×n . Si det (A) 6= 0 entonces A es no singular.
Nota 4.3. Si Ax = b y det (A) 6= 0, entonces el sistema tiene una única solución
dada por
x = A−1 b.
Aunque de importancia teórica, esta expresión no proporciona un método práctico
para encontrar la solución de una ecuación lineal.
Teorema 4.4. Sea A ∈ Rn×n . Si A es no singular, entonces det (A) 6= 0.
Demostración. Suponemos que existe la inversa (única) A−1 de A la cual cumple
AA−1 = A−1 A = I.
Tomando determinante, 1 = det (I) = det (AA−1 ) = det (A) det (A−1 ) . Ninguno de
estos factores puede ser nulo, en particular, det (A) 6= 0.
4.2. Formulación matricial de la eliminación de Gauss

En esta sección introducimos las factorizaciones A = LU y P A = LU con L
triangular inferior con unos en la diagonal, U triangular superior y P matriz de per-
mutación. Estas factorizaciones sirven para formular matricialmente la eliminación
de Gauss.
Empezamos por el caso más sencillo de la eliminación de Gauss, es decir, todos
los divisores en el proceso son no nulos y no se necesitan intercambios de filas. Los
cómputos que llevan al sistema (4.3) se pueden escribir en términos de matrices como
L(1) Ax = L(1) b,
donde  
1
 −l21 1 
 
 −l31 0 1 
L(1) =  (4.8)
 .. .. .. 
 . . . 
−ln1 0 · · · 0 1
4.2. FORMULACIÓN MATRICIAL DE LA ELIMINACIÓN DE GAUSS 37
ai1
con li1 = , i = 2, . . . , n.
a11
El paso k se realiza por premultiplicación por la matriz
 
1
 . 
 0 .. 
 
 0 1 
(k) 
L =  .. .. .. , (4.9)
. . −l . 
 k+1,k 
 . 
 .
. 
0 0 −ln,k 0 1
(k−1)
aik
donde lik = (k−1)
, i = k + 1, . . . , n.
akk
El sistema triangular (4.4) se puede escribir
b
LAx b con L
= Lb, b = L(n−1) · · · L(2) L(1) . (4.10)
La matriz L b es triangular inferior con unos en su diagonal. Sea U la matriz

triangular superior del sistema (4.4). El trabajo anterior indica que LAb = U. La
matriz Lb tiene determinante 1, por tanto es no singular. Su inversa L = L
b−1 también
es triangular inferior con unos en su diagonal y cumple que
A = LU. (4.11)
En 1.1.2, numeral 6, definimos matriz de permutación y ese es precisamente el

concepto que requerimos ahora. Si P es una matriz de permutación, se puede pensar
de P como la matriz que se obtiene de la identidad I por permutación de dos de sus
filas, digamos las filas j y k. El producto P A es la matriz que se obtiene de A por
intercambio de sus filas j y k.
Volvemos a la eliminación de Gauss y su formulación matricial. Si un divisor es
cero, se necesita un intercambio de filas de A. Supongamos que a11 = 0 y que la
primera fila se debe intercambiar con la fila k. La nueva matriz es P (1) A, donde
P (1) es una matriz de permutación. La matriz L(1) de (4.8) con los elementos de A
sustituidos por los de P (1) A, premultiplica a P (1) A o sea L(1) P (1) A.
(1)
Para el segundo paso, si a22 = 0, se requiere una permutación P (2) y el resultado
es P (2) L(1) P (1) A. Finalmente se obtiene la expresión similar a (4.10) dada por
b
LAx b
= Lb, b = L(n−1) P (n−1) · · · L(2) P (2) L(1) P (1)
L (4.12)
En esta expresión P (i) = I si no se requirió intercambio en dicho paso. Las

matrices L(i) son triangulares inferiores con unos en la diagonal pero Lb no tiene que
b es invertible por ser producto de
ser triangular inferior. Lo que si se sabe es que L
matrices invertibles y entonces se puede escribir
b−1 U
A=L (4.13)
b−1 no tiene que ser triangular inferior.

que es análogo a (4.11) aunque L
Teorema 4.5. Si A es una matriz invertible de Rn×n entonces existe una matriz de
permutación P tal que
P A = LU (4.14)
donde L es triangular inferior con unos en la diagonal y U es triangular inferior.
Demostración. Supongamos que en el paso i se requiere intercambio con la fila k. Si

a la matriz original A se le hubieran intercambiado las filas i y k antes de empezar
la eliminación, en el paso i habríamos encontrado el elemento diagonal no nulo que
requerimos. La estrategia es entonces la siguiente: Realicemos la eliminación de Gauss
completa en la matriz A y registremos los intercambios de fila que se necesitaron.
Agrupemos dichos intercambios en la matriz de permutación P. El resultado de la
eliminación con intercambios es el mismo que el de la eliminación sin intercambios
pero realizado a la matriz P A. Por tanto
P A = LU.
 
4 4 4 4
 2 2 2 2 
Ejemplo 4.6. Sea A =  
 1 0 1 0  . En la eliminación de Gauss se necesita
1 1 1 1
intercambiar filas 2 y 3 para poder finalizarla. Sea
 
1 0 0 0
 0 0 1 0 
P =  0 1 0 0 .

0 0 0 1
4.3. NOTA MATLAB 39
El producto P A es  
4 4 4 4
 1 0 1 0 
PA = 
 2

2 2 2 
1 1 1 1
y la eliminación de Gauss de P A no requiere intercambios y usa los multiplicadores
l21 = − 41 , l31 = − 12 y l41 = − 41 . Dicha eliminación lleva a la matriz
 
4 4 4 4
 0 −1 0 −1 
U =



Por tanto  
1
 1/4 1 
L(1) =
 1/2 0 1


1/4 0 0 1
es tal que P A = L(1) U.
Nota 4.7. Se acostumbra omitir los elementos nulos de las matrices siempre que no
haya confusión.
4.3. Nota MATLAB

Tomado de la ayuda de MATLAB al invocar doc lu
a = [1, 2, 3; 4, 5, 6; 7, 8, 0];
[l, u, p] = lu (a) ;
Las matrices que resultan son:

     
1 7 8 0 0 0 1
l =  1/7 1 , u= 6/7 3 , p= 1 0 0 
4/7 1/2 1 9/2 0 1 0
Nótese que se hubiera podido realizar eliminación de Gauss sin intercambio de filas.
Aquí aparecen intercambios debido a la estrategia de pivoteo parcial.
4.4. Invertibilidad
Terminamos esta sección con algunas equivalencias para invertibilidad de matri-
ces.
Teorema 4.8. Sean A ∈ Rn×n . Las siguientes condiciones son equivalentes:
1. det (A) 6= 0.
2. Hay una matriz n × n denotada A−1 tal que AA−1 = A−1 A = I.

3. Ax = b tiene una única solución para todo b ∈ Rn .
4. Ax = 0 tiene solamente la solución trivial x = 0.
5. Las columnas (filas) de A son linealmente independientes.
Demostración. Ya se demostró la equivalencia de 1. y 2. y que 1. implica 3. Además,

4. es un caso particular de 3. Veamos la equivalencia de 4. y 5. para columnas.
Escribimos A por columnas en la forma
A = (A∗1 A∗2 A∗n ) .
La ecuación Ax = 0 se escribe
n
X
0 = Ax = xi A∗i .
i=1
Si se cumple 4. entonces las columnas de A son linealmente independientes. Recípro-

camente, si las columnas de A son LI entonces Ax = 0 tiene una única solución que
es la trivial x = 0.
La condición para las filas es consecuencia de la igualdad det (A) = det AT . Fina-
lizamos la demostración probando que 4. implica 1. Supongamos que se cumple 4.
pero det (A) = 0. Por ser nulo el determinante, el proceso de eliminación de Gauss
lleva a una matriz reducida de la forma
 
∗ ··· ··· ∗ ··· ··· ∗
 .. .. .. 
 . . . 
 
b  ∗ ∗ · · · · · · ∗  A(1) a B
A= = (4.15)
 0 ∗ ··· ∗  0 0 C
 .. . . .. 
 . . . 
0 ∗ ··· ∗
4.5. NOTA MATLAB 41
Los asteriscos indican elementos que pueden ser no nulos y en la columna k − ésima
el elemento diagonal y todos los que hay debajo son cero. Nótese que si k = 1, la
primera columna de A es nula y asi Ae1 = 0, lo cual contradice 4. Por consiguiente
k ≥ 2.
La matriz A(1) es triangular superior, invertible, de tamaño (k − 1) × (k − 1) , a es
un vector de tamaño (k − 1) , B es (k − 1) × (n − k) y C es (n − k + 1) × (n − k) .
El sistema A(1) x(1) = a tiene solución x(1) , por tanto
 
(1) x(1) (1) (1)
A a B   A x −a
−1 = = 0.
0 0 C 0
0
Esto significa que hay un vector no nulo x b tal que Ab bx = 0. Pero, de acuerdo con
(4.12),
Ab = L(k−1) P (k−1) · · · L(1) P (1) A
donde las matrices L(i) son las matrices triangulares inferiores y las matrices P (i) son
las matrices de permutación que surgen en la eliminación de Gauss. Tenemos:
bx = L(k−1) P (k−1) · · · L(1) P (1) Ab
0 = Ab x.
Como las matrices L(i) y P (i) son invertibles, entonces Ab

x = 0. Como x
b 6= 0, esto es
una contradicción.
Nota 4.9. Hay otras equivalencias para la no singularidad de una matriz y también
hay otros métodos de demostración, por ejemplo aquellos basados en valores propios.
4.5. Nota MATLAB

a = [4,4,4,4;2,2,2,2;1,0,1,0;1,1,1,1];
[l,u,p] = lu(a);
Este código MATLAB permite encontrar las matrices que aparecen arriba. Debe
tenerse en cuenta que la rutina lu de MATLAB aplica por defecto la estrategia de
pivoteo parcial que consiste en que en la j − ésima columna, de la posición (j, j)
hacia abajo, se escoge como pivote el primer elemento no nulo que tiene el mayor
valor absoluto de esa parte de la columna.
Capítulo 5
Forma escalonada reducida y rango
Este capítulo trata de eliminación de Gauss para matrices rectangulares y para

matrices cuadradas que no tienen que ser invertibles. Empezamos con la forma esca-
lonada por filas y la definición de las tres operaciones elementales fila de las cuales dos
se han usado ya antes. Después llegamos a la forma escalonada reducida por filas a
la definición de rango que es una de las que más se usarán en el resto del manuscrito.
Para escribir este capítulo consultamos principalmente [M], [L] y [Ortega].
Consideremos A ∈ Rm×n . En el caso rectangular también podemos aplicar eli-
minación de Gauss con intercambio de filas pero se puede presentar la situación
siguiente: el elemento diagonal es cero y todos los elementos debajo de él en la mis-
ma columna también son nulos. En tal caso, se debe pasar a la siguiente columna en
la misma fila y continuar el proceso. Se llega a la siguiente eliminación modificada
de Gauss:
Suponga que U es la matriz que se obtiene después de completar i − 1 pasos de
la eliminación. Para ejecutar el paso i − ésimo, se siguen las siguientes instrucciones:
1. Mirando de izquierda a derecha en U, localizar la primera columna que contiene

un elemento no nulo en la fila i o debajo de ella. La llamamos U∗j .
2. El pivote para el paso i es el elemento de la posición (i, j) .
3. De ser necesario, se intercambia la fila i con una fila inferior para traer un
número distinto de cero a la posición (i, j) y después, se anulan los elementos
de las posiciones inferiores a la del pivote (i, j) .
4. Si la i − ésima fila de U y todas las que están debajo son nulas, el proceso está
completo.
43
44 CAPÍTULO 5. FORMA ESCALONADA REDUCIDA Y RANGO
Definición 5.1. Se dice que una matriz A ∈ Rm×n , con filas Ai∗ y columnas A∗j , está
en forma escalonada por filas (row echelon form), si se cumplen las dos condiciones
siguientes:
1. Si Ai∗ es una fila nula, entonces todas las filas debajo de ella también son nulas.
2. Si el primer elemento no nulo de Ai∗ está en la posición j, entonces todos los

elementos debajo de la posición i en las columnas A∗1 , A∗2 , . . . , A∗j son nulos.
Las matrices en forma escalonada por filas tienen la siguiente forma:

 
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
 ∗ ∗ ∗ ∗ ∗ ∗ 
 
 ∗ ∗ ∗ ∗ ∗ 
 
 ∗ ∗ 
 
 
5.1. Operaciones elementales

Sea A ∈ Rm×n . Las siguientes operaciones sobre las filas Ai∗ de A se denomi-
nan operaciones elementales. Se pueden definir las operaciones análogas sobre las
columnas de la matriz.
1. E1: Intercambio de las filas i y j.
2. E2: Reemplazar la fila i por un múltiplo escalar no nulo de ella misma.
3. E3: Reemplazar la fila j por Aj∗ + βAi∗ con β 6= 0.
Nota 5.2. Hasta ahora hemos utilizado operaciones E1 y E3 en la eliminación de

Gauss y en la construcción de formas escalonadas por filas. Para llegar a formas
escalonadas reducidas por filas, necesitamos también operaciones de tipo E2.
Definición 5.3. Una matriz A ∈ Rm×n está en forma escalonada reducida por filas
(reduced row echelon form) si se cumplen las siguientes condiciones:
1. A está en forma escalonada por filas.
2. El primer elemento de cada fila no nula es un uno.

5.2. MATRICES ELEMENTALES 45
3. Todos los elementos encima de cada pivote en su misma columna son ceros.
Las matrices en forma escalonada reducida por filas tienen la siguiente forma:
 
1 ∗ 0 0 ∗ ∗ 0 ∗
 1 0 ∗ ∗ 0 ∗ 
 
 1 ∗ ∗ 0 ∗ 
 
 1 ∗ 
 
 
 
1 2 3 3
Ejemplo 5.4. Sea A =  2 4 6 9  . Usando solamente operaciones elemen-
2 6 7 6
tales de tipo E1 y E3, realizamos el siguiente proceso de eliminación que no tiene
estrategias de pivoteo:
     
1 2 3 3 1 2 3 3 1 2 3 3
 2 4 6 9  −→  0 0 0 3  −→  0 2 1 0  .
2 6 7 6 0 2 1 0 0 0 0 3
Para obtener la forma escalonada reducida, hacemos lo siguiente:

     
1 2 3 3 1 2 3 3 1 0 2 0
 0 2 1 0  −→  0 1 1/2 0  −→  0 1 1/2 0  .
0 0 0 3 0 0 0 3 0 0 0 1
Aquí usamos operaciones elementales de tipos E2 y E3.
Nota 5.5. Para una matriz invertible A su forma escalonada reducida es la matriz
identidad. El proceso de eliminación en este caso se denomina de Gauss Jordan. Se
caracteriza por la utilización de operaciones E2 para la obtención de los unos de la
diagonal.
5.2. Matrices elementales

Las matrices de la forma
I − uv T , (5.1)
con u y v vectores de Rn tales que v T u 6= 1, se denominan matrices elementales.

Siempre son invertibles, sus inversas tienen la forma
−1 uv T
I − uv T =I− (5.2)
vT u − 1
y por tanto también son matrices elementales.
La relación de estas matrices con las operaciones elementales definidas en (5.1)
es la siguiente: El resultado de realizar una operación elemental de alguno de los tres
tipos en la matriz identidad, es una matriz elemental. Para ilustrarlo, proponemos
lo que ocurre con matrices 3 × 3.
1. Al realizar una operación E1 en I, consistente en intercambiar filas 1 y 2, se

obtiene  
0 1 0
S (1) =  1 0 0 
0 0 1
y S (1) = I − vv T , donde v = e1 − e2 .
2. Al realizar una operación E2 en I, consistente en multiplicar su segunda fila
por un escalar β, se obtiene
 
1 0 0
S (2) =  0 β 0 
0 0 1
y S (2) = I − (1 − β) e2 eT2 .
3. Al realizar una operación E3 en I, consistente en tomar β veces la fila 1 y
sumarla a la fila 3, se obtiene
 
1 0 0
S (3) =  0 1 0 
β 0 1
y S (3) = I + βe3 eT1 .
Proposición 5.6. Propiedades de las matrices elementales
1. Si una matriz A se pre-multiplica por una matriz elemental Q la matriz re-

sultante es la matriz A con sus filas modificadas de acuerdo con la operación
elemental asociada con Q.
5.2. MATRICES ELEMENTALES 47
2. Si una matriz A se post-multiplica por una matriz elemental Q, la matriz resul-

tante es A con sus columnas modificadas de acuerdo con la operación elemental
asociada a Q.
Teorema 5.7. La matriz cuadrada A es invertible si y solo si es un producto de

matrices elementales.
Demostración. Si A es invertible, entonces por medio de eliminación de Gauss Jordan
(Nota 5.5) en la que se usan solamente operaciones elementales fila se puede llegar
a la matriz identidad. Si las matrices elementales asociadas con esas operaciones
elementales se denotan Qi , entonces
I = Qk Qk−1 · · · Q2 Q1 A.
Por tanto A = (Qk Qk−1 · · · Q2 Q1 )−1 = Q−1 −1 −1
1 Q2 · · · Qk que es un producto de ma-
trices elementales.
Definición 5.8. Equivalencia
Si una matriz B se puede obtener de una matriz A por medio de operaciones ele-
mentales fila y columna se escribe A ∼ B y se dice que A es equivalente a B. Esto
lo podemos reescribir asi: A ∼ B si y solo si hay matrices invertibles Q1 y Q2 tales
que Q1 AQ2 = B.
Nota 5.9. 1. Si todas las operaciones elementales se hacen en las filas entonces
Q2 = I y si todas las operaciones elementales se hacen en las columnas entonces
Q1 = I.
2. Como las matrices Q1 y Q2 son invertibles entonces A ∼ B si y solo si B ∼ A.
3. También se da A ∼ A.
4. Finalmente, también se da la transitividad A ∼ B y B ∼ C =⇒ A ∼ C.
5. Los tres numerales anteriores indican que ∼ es una relación de equivalencia.
Proposición 5.10. Para cada matriz A ∈ Rm×n hay una única forma escalonada
reducida por filas. Se denota EA .
Demostración. La demostración puede leerse en el Example 3.9.2 de [M].
Nota 5.11. La unicidad de la forma escalonada reducida de A implica la unicidad de
las posiciones para los pivotes y por tanto también quedan unívocamente determina-
dos el número de pivotes y el número de columnas de A que contienen las posiciones
de los pivotes.
Definición 5.12. A las columnas de A (EA ) que contienen las posiciones de los
pivotes se les llama columnas básicas de A (EA ). A las demás columnas de A
(EA ) se les llama columnas no básicas.
Definición 5.13. Para cada matriz A ∈ Rm×n el número de filas no nulas en su

forma escalonada reducida se denomina rango (rank) de A y se denota rank (A) .
Nota 5.14. El rango es un protagonista en estas notas. Utilizaremos las siguientes

igualdades:
rank (A) = número de pivotes
= número de f ilas no nulas en EA
= número de columnas básicas de A.
5.3. Relación entre columnas y filas

Proposición 5.15. Si A ∼ B por filas entonces las relaciones lineales entre columnas
de A también se dan entre las columnas de B. Es decir,
n
X n
X
B∗k = αj B∗j si y solo si A∗k = αj A∗j .
j=1 j=1
En particular, las relaciones entre las columnas de A y EA son las mismas de manera
que las columnas no básicas de A son combinaciones lineales de las columnas básicas
de A.
Si A ∼ B por columnas, entonces las relaciones lineales entre filas de A también se
dan entre las filas de B.
Demostración. Si A ∼ B por filas entonces QA = B para una matriz invertible Q.

La j − ésima columna de B es
B∗j = (QA)∗j = QA∗j .
Luego, si
n
X
A∗k = αj A∗j
j=1
entonces n
X
QA∗k = αj QA∗j
j=1
5.3. RELACIÓN ENTRE COLUMNAS Y FILAS 49
que corresponde a
n
X
B∗k = αj B∗j .
j=1
n
X
Recíprocamente, si B∗k = αj B∗j entonces la multiplicación por Q−1 lleva a A∗k =
j=1
n
X
αj A∗j .
j=1
Nota 5.16. Para A y EA la situación es: Cada columna no básica E∗k de EA es

una combinación lineal de las j columnas básicas que hay a la izquierda de E∗k .
Los coeficientes que se necesitan son las primeras j componentes de E∗k . Por la
proposición anterior también es verdadero lo siguiente:
Cada columna no básica A∗k de A es una combinación lineal de las j columnas básicas
que hay a la izquierda de A∗k . Los coeficientes que se necesitan son las primeras j
componentes de E∗k .
Corolario 5.17. Sea A ∈ Rm×n . Las columnas básicas de A son un conjunto LI.
Demostración. Escribimos a A por columnas
A = (A1 A2 . . . An )
y suponemos que el conjunto de columnas básicas de A es {Ar1 , Ar2 , . . . , Ark } . El
conjunto correspondiente en EA es {Qr1 , Qr2 , . . . , Qrk } .
Demostración. Sea
k
X
0= βi Ari
i=1
una CL del cero con las columnas básicas de A. Por proposición anterior, 0 =
Xk
βi Qri . Por el desfase de los unos en EA , esto obliga βi = 0 para todo i.
i=1
Nótese que si la matriz A ∈ Rm×n se escribe por columnas, es decir, A =

(A1 A2 An ) y T es el conjunto de las columnas básicas de A, entonces Sg (T ) =
Sg ({A1 , A2 , . . . , An }) .
Definición 5.18. Sea A = (A1 A2 An ) ∈ Rm×n . El subespacio de Rm dado por
Sg ({A1 , A2 , . . . , An }) se denomina rango (range) de A o espacio columna de A o
imagen de A y se denota R (A) . Por corolario 5.17, dim R (A) = rank (A) .
 
X n r 1
 
Nota 5.19. Si w ∈ R (A) entonces w = rj Aj = Ar donde r =  ...  . Esto
j=1 rn
explica el nombre imagen de la definición anterior.
5.4. Forma normal del rango

Utilizando operaciones elementales fila y columna, se puede demostrar lo siguien-
te:
Proposición 5.20. Si A ∈ Rm×n con rank (A) = r, entonces

Ir 0
A ∼ Nr = .
0 0
Nr se denomina forma normal del rango (rank normal form). Es la matriz más
simplificada a la que se puede llegar por medio de operaciones elementales.
Demostración. Se sabe que A ∼ EA por filas y esto significa que existe una matriz no
singular P tal que P A = EA . Las r columnas básicas de EA son vectores canónicos.
Si se intercambian columnas hasta que estos vectores canónicos queden en la posición
superior izquierda, es porque hay una matriz no singular Q1 tal que

Ir M
P AQ1 = EA Q1 = .
0 0
Al multiplicar ambos lados por

Ir −M
Q2 =
0 I
llegamos a
Ir M Ir −M Ir 0
P AQ1 Q2 = = .
0 0 0 I 0 0
Es decir, A ∼ Nr .
Terminamos este capítulo con dos proposiciones que facilitan el trabajo con el
rango.
5.4. FORMA NORMAL DEL RANGO 51
Proposición 5.21. Si A y B son matrices m × n entonces:
1. A ∼ B si y solo si rank (A) = rank (B) . En particular, multiplicación por

matriz invertible genera una matriz con el mismo rango (rank).
2. A ∼ B por filas si y solo si EA = EB .
3. A ∼ B por columnas si y solo si EAT = EBT .
Demostración. Ver [M] pag. 138.
Proposición 5.22. La trasposición no cambia el rango. Es decir, si A es una matriz

m × n entonces

rank (A) = rank AT y rank (A) = rank (A∗ ) .
Demostración. Ver [M] pag. 139.

Capítulo 6
Matrices y transformaciones lineales
Las transformaciones lineales, también llamadas funciones lineales, operadores li-

neales y homomorfismos, son las funciones entre espacios vectoriales que respetan la
suma y el producto por escalar. Veremos que las matrices se pueden identificar como
transformaciones lineales aunque el estudio de matrices se puede hacer sin recurrir a
dicha identificación (ver por ejemplo, [I].) Nuestra aproximación está basada en [L] y
también consultamos [M] y [Ortega]. Los principales temas que tratamos son: defini-
ción, representación matricial de una transformación lineal y los cuatro subespacios
fundamentales de una transformación lineal y por tanto de una matriz.
Antes de empezar el desarrollo del tema, advertimos que siempre que mencione-
mos un espacio vectorial será de dimensión finita y el campo asociado será el de los
números reales.
6.1. Definición
Sean V y W espacios vectoriales. Entonces T : V −→ W es una transformación
lineal (TL) de V a W si se cumple
T (βx + ρy) = βT (x) + ρT (y) para todo β, ρ ∈ R y todo x, y ∈ V.
Ejemplo 6.1. 1. La función 0 : V −→ W tal que para todo v ∈ V, 0 (v) = 0 es

una TL.
2. La función idéntica i : V −→ V tal que para todo v ∈ V, i (v) = v es una TL.
3. Sea A ∈ Rm×n y sea T : Rn −→ Rm definida asi: T (v) = Av es una transfor-

mación lineal pues T (βx + ρy) = A (βx + ρy) = βAx+ρAy = βT (x)+ρT (y) .
53
54 CAPÍTULO 6. MATRICES Y TRANSFORMACIONES LINEALES
Cuando no hay confusión, se omite el paréntesis en T (v) y se escribe simplemente

T v.
6.2. Representación matricial

Sea T : Rn −→ Rm una TL y supongamos que B1 = {V1 , V2 , . . . , Vn } y B2 =
{W1 , W2 , . . . , Wn } son bases de Rn y Rm respectivamente. Definimos la representación
matricial de T en las bases B1 y B2 como la matriz A ∈ Rm×n
 
a11 · · · a1n
 .. 
A =  ... .  (6.1)
am1 · · · amn
m
X
cuyas componentes cumplen que T Vi = aki Wk para cada i = 1, . . . , n. Es decir, la
k=1
columna i − ésima de A está conformada por los (únicos) coeficientes que permiten
escribir T Vi como combinación lineal de los elementos de la base B2 de Rm . Definamos
las matrices de las bases asi:
C = (V1 V2 . . . Vn ) , F = (W1 W2 . . . Wm )
y denotemos A∗i a la columna i − ésima de A, es decir,
 
a1i
 
A∗i =  ...  y A = (A∗1 A∗2 . . . A∗n ) .
ami
Puede verse que T Vi = F A∗i .
Por supuesto la matriz A depende de las bases B1 y B2 . Además, dado un vector
v ∈ Rn , con escritura (única) en términos de la base B1 dada por
Xn
v= ρj Vj
j=1
entonces T v queda determinado por el valor

 de T en los elementos de la base ya que
X n ρ1
 .. 
Tv = ρj T Vj . Denotemos x =  .  . Entonces v = Cx y
j=1 ρn
n
X n
X n
X
T Cx = T v = ρj T Vj = ρj F A∗j = F ρj A∗j = F Ax. (6.2)
j=1 j=1 j=1
6.3. COMPOSICIÓN Y SUBESPACIOS ASOCIADOS 55
Como esto es válido para todo x pues es válido para todo v, desearíamos concluir que
T C es igual a F A, pero debemos ser cuidadosos pues no disponemos de una definición
para T C. Si las bases B1 y B2 son las bases canónicas respectivas, entonces C y F
son matrices identidad y la igualdad que nos gustaría obtener es T = A pero de
nuevo, falta significado para una igualdad como ésta.
La forma correcta de considerar a las matrices como transformaciones lineales es
a través de la transformación lineal T definida en el numeral 3 del ejemplo 6.1 como
T v = Av. A esta TL también la llamamos A y se habla indistintamente de la matriz
A o de la transformación lineal A. Afortunadamente pensar en A como TL y como
matriz casi nunca genera confusión. A manera de resumen proponemos la siguiente
equivalencia:
A ∈ Rm×n ⇐⇒ A : Rn −→ Rm es la TL que a cada x ∈ Rn asigna Ax. (6.3)
6.3. Composición y subespacios asociados

6.3.1. Composición de transformaciones lineales
Sean A ∈ Rm×n , B ∈ Rn×p con sus transformaciones lineales asociadas
A : Rn −→ Rm B : Rp −→ Rn
y .
x 7−→ Ax x 7−→ Bx
La composición de A y B, que esquemáticamente se acostumbra escribir
B A
p
R −→ R −→ Rm
n
x 7−→ Bx 7−→ ABx,
es la TL
C
p
R −→ Rm
x 7−→ Cx = ABx.
La matriz de la TL C es AB.
6.3.2. Subespacios asociados

Sea A : Rn −→ Rm la transformación lineal asociada con la matriz A ∈ Rm×n .
Definición 6.2. El rango de A, denotado R (A) , es
R (A) = {w ∈ Rm : w = Av para algún v ∈ Rn } . (6.4)
Esta definición coincide con la definición 5.18 pues los vectores de la forma Av son
las combinaciones lineales de las columnas de A.
Definición 6.3. El espacio nulo de A, denotado ker (A), es
ker (A) = {v ∈ Rn : Av = 0} .
Al espacio nulo de A también se le llama núcleo (kernel) de A.
Teorema 6.4. El rango de A es un subespacio de Rm y el espacio nulo de A es un

subespacio de Rn .
Demostración. Sean w, z ∈ R (A) y sean x, y ∈ Rn tales que Ax = w y Ay = z.
Entonces αw + βz = αAx + βAy = A (αx + βy) . Luego αw + βz ∈ R (A) .
Para el espacio nulo, sean z, w ∈ ker (A) . Entonces Az = 0 y Aw = 0. El elemento
αz + βw es tal que A (αz + βw) = αAz + βAw = 0 + 0 = 0. Luego αz + βw ∈
ker (A) .
6.4. Ortogonalidad
Sea R = {V1 , V2 , . . . , Vk } un subconjunto de vectores no nulos de Rn . Se dice que
R es ortogonal si ViT Vj = 0 para i 6= j y ortonormal si ViT Vj = δij , donde δij es
el delta de Kronecker, es decir,

1 si i = j,
δij =
0 si i 6= j.
Definición 6.5. Sea S un subespacio de Rn . El complemento ortogonal de S,

denotado S ⊥ , se define asi:

S ⊥ = v ∈ Rn : v T s = 0 para todo s ∈ S .
Teorema 6.6. Sean R, S subespacios de Rn . Entonces
1. El complemento ortogonal de S es un subespacio de Rn .
2. S ⊕ S ⊥ = Rn .
6.5. CUATRO SUBESPACIOS FUNDAMENTALES 57
⊥
3. S ⊥ = S.
4. R es un subespacio de S si y solo si S ⊥ es un subespacio de R⊥ .
5. (R + S)⊥ = R⊥ ∩ S ⊥ .
6. (R ∩ S)⊥ = R⊥ + S ⊥ .
Teorema 6.7. Sea A : Rn −→ Rm . Entonces

1. ker (A)⊥ = R AT .

2. R (A)⊥ = ker AT .
Demostración. Demostremos la parte 2. x ∈ R (A)⊥ ⇐⇒ (Ay)

T
x = 0 para todo y
T T T T
⇐⇒ y A x = 0 para todo y ⇐⇒ A x = 0 ⇐⇒ x ∈ ker A .
6.5. Cuatro subespacios fundamentales

Los cuatro subespacios fundamentales de la transformación lineal A : Rn −→ Rm
de rango rank (A) = r, son: R (A) , R (A)⊥ , ker (A) y ker (A)⊥ . Por teorema 6.7,
teorema 6.6 y proposición 5.22,

r = dim R (A) = dim R AT = dim ker (A)⊥
y
n − r = dim R (A)⊥ = dim ker AT .
Teorema 6.8. Sea A : Rn −→ Rm . Entonces
n = dim ker (A) + rank (A) = dim ker (A) + dim R (A) .
Demostración. Por
teorema 6.6, n = dim ker (A) + dim
ker (A)⊥ . Por teorema 6.7,
⊥
ker (A) = R AT y por proposición 5.22, dim R AT = dim R (A) = rank (A) .
Definición 6.9. Sea A : Rn −→ Rm . Decimos que A es sobreyectiva o sobre si

R (A) = Rm . Decimos que A es inyectiva o 1-1 si Ax = Ay =⇒ x = y para todo
par de vectores x, y en Rn .
Una formulación equivalente de inyectividad es: A es inyectiva si y solo si ker (A) =

{0} y por tanto dim ker (A) = 0.
Teorema 6.10. Sean A ∈ Rm×n y B ∈ Rn×p . Entonces
1. R (AB) ⊆ R (A) y si B tiene filas LI entonces R (AB) = R (A) .
2. ker (B) ⊆ ker (AB) y si A tiene columnas LI entonces ker (B) = ker (AB) .
Demostración. 1. Si y ∈ R (AB) , y = ABx y asi y ∈ R (A) . Si B tiene filas LI,

B T tiene columnas LI.

rank (B) = rank B T = n
y por tanto R (B) = Rn . Si y ∈ R (A) entonces y = Ax para algún x ∈ Rn =

R (B) o sea y ∈ R (AB) .
2. Si x ∈ ker (B) , Bx = 0, luego ABx = 0 y asi x ∈ ker (AB) . Si A tiene

columnas LI, R (A) = Rn , rank (A) = n y ker (A) = {0}. De manera que si
ABx = 0 es porque Bx = 0.
Finalizamos este capítulo con las definiciones de valores y vectores propios que
serán importantes posteriormente.
6.6. Valores y vectores propios

Para una matriz cuadrada A, los escalares λ y los vectores x 6= 0 que satisfacen
Ax = λx
se denominan valores y vectores propios de A respectivamente. Cada pareja (λ, x)

se denomina par propio o pareja propia de A. El conjunto de todos los valores propios
de A se denomina Espectro de A y se denota σ (A) . Se cumple:
1.
λ ∈ σ (A) ⇐⇒ A − λI es singular ⇐⇒ det (A − λI) = 0. (6.5)
2.
{x 6= 0 : x ∈ ker (A − λI)}
es el conjunto de los vectores propios asociados con λ. Al subespacio ker (A − λI)
se le llama subespacio propio de A asociado con el valor propio λ.
6.7. NOTA MATLAB 59
En general el espectro de A contiene números complejos aun en el caso en que A

sea real. Es que la búsqueda de valores propios se hace resolviendo la ecuación
pλ (t) = det (A − λI) = 0 (6.6)
que es un polinomio en λ con coeficientes reales si A es real y complejos si A es

compleja.
El polinomio (6.6) se denomina polinomio característico de A.
Finalmente, el radio espectral de A, denotado ρ (A) , se define como
ρ (A) = máx {|λ| : λ ∈ σ (A)} , (6.7)
es decir, es el máximo valor absoluto o módulo en el espectro de A.
6.7. Nota MATLAB

a = gallery(’fiedler’,8,1:5);
[v,d] = eig(a); % valores y vectores propios
Capítulo 7
Normas de vectores y matrices
En álgebra lineal numérica es de especial importancia tener una idea del tamaño
de una perturbación o de un error. De hecho uno de los asuntos que se estudian es
el condicionamiento, es decir, si en un problema los datos se perturban, ¿qué tan
grande es la diferencia entre la solución exacta del problema y la solución calculada?
Si esta diferencia se parece a la perturbación en los datos, se habla de buen condicio-
namiento y si es muy grande en comparación con la perturbación, se habla de mal
condicionamiento.
Para medir vectores y matrices se usan normas. En este capítulo definimos y
describimos las principales normas. Para escribirlo consultamos principalmente [I],
[M] y [vdG].
7.1. Valor absoluto

Si β ∈ R su valor absoluto es

β, si β ≥ 0
|β| =
−β, si β < 0.
Si β ∈ C su valor absoluto o módulo se define asi: Sea β = b1 + ib2 . Entonces,

q q
2 2
|β| = b1 + b2 = ββ
donde β es el conjugado de β, es decir, β = b1 − ib2 .

El valor absoluto de números complejos o reales tiene las siguientes propiedades:
1. β 6= 0 =⇒ |β| > 0. El valor absoluto es definido positivo.
61
62 CAPÍTULO 7. NORMAS DE VECTORES Y MATRICES
2. |αβ| = |α| |β| . El valor absoluto es homogéneo.
3. |α + β| ≤ |α| + |β| . El valor absoluto satisface la desigualdad triangular.
7.2. Normas de vectores

La noción de norma es la extensión del valor absoluto para vectores.
7.2.1. Definición
Sea N : Cn −→ R. Entonces N es una norma vectorial si para todo x, y ∈ Cn ,
1. x 6= 0 =⇒ N (x) > 0. La norma es definida positiva.
2. N (βx) = |β| N (x) . La norma es homogénea.
3. N (x + y) ≤ N (x) + N (y) . La norma satisface la desigualdad triangular.
7.2.2. Observaciones
1. La norma vectorial en Rn se define de la misma manera.
2. Si N es una norma en Cn o en Rn , entonces N (0) = 0.
La segunda observación se deriva de la propiedad 2 de la definición de norma

vectorial 7.2.1 tomando como escalar al número 0. Debido a esto, es común ver que
la propiedad 1 de la definición 7.2.1 se reescriba asi:
N (x) ≥ 0 y N (x) = 0 ⇐⇒ x = 0.
La notación más común para normas es k·k .
7.2.3. Norma euclideana o norma 2 o 2-norma
n
! 21
X 1
kxk2 = |xi |2 = (x∗ x) 2 (7.1)
i=1
Para demostrar que (7.1) es una norma vectorial se utiliza la desigualdad de Cauchy-
Schwarz.
7.2. NORMAS DE VECTORES 63
Teorema 7.1. (Desigualdad de Cauchy-Schwarz) Si x, y ∈ Cn entonces
|x∗ y| ≤ kxk2 kyk2 .
Demostración. Se puede consultar por ejemplo en [vdG] Theorem 2.4.
Teorema 7.2. La 2-norma es una norma vectorial
Demostración. Sean x, y ∈ Cn y β ∈ C.
1. Supongamos x 6= 0. Entonces alguna de sus componentes es no nula, digamos

n
! 21
X 1
xj 6= 0. kxk2 = |xi |2 ≥ |xj |2 2 = |xj | > 0.
i=1
n
! 21 n
! 21 n
! 21
X X X
2. kβxk2 = |βxi |2 = |β|2 |xi |2 = |β|2
|xi |2 = |β| kxk2 .
i=1 i=1 i=1
3. kx + yk22 = (x + y)∗ (x + y) = x∗ x + y ∗x + x∗ y + y ∗ y
≤ kxk22 + 2 kxk2 kyk2 + kyk22 = (kxk2 + kyk2 )2 .
7.2.4. La norma 1 o 1-norma

n
X
kxk1 = |xi | .
i=1
7.2.5. La norma ∞ o norma uniforme
kxk∞ = máx |xi | .

i
7.2.6. La norma p o p − norma
n
! p1
X
kxkp = |xi |p .
i=1
Las normas 1 y 2 son casos particulares de la norma p. Demostrar que en efecto

son normas es un ejercicio que vale la pena intentar. Para demostrar que la norma p
satisface la desigualdad triangular se utiliza una desigualdad denominada de Hölder.
1 1
Teorema 7.3. (Desigualdad de Hölder) Sean x, y ∈ Cn y p
+ q
= 1 para reales p, q
mayores que 1. Entonces
|x∗ y| ≤ kxkp kykq .
Demostración. Sugerimos ver el ejercicio 5.1.12 de [M].
7.3. Normas de matrices

La definición de una norma matricial se enuncia de forma análoga a la de una
norma vectorial pero requiere de una condición adicional.
7.3.1. Definición
Sea N : Cm×n −→ R. Entonces N es una norma matricial si para todo A, B ∈
Cm×n ,
1. A 6= 0 =⇒ N (A) > 0. La norma es definida positiva.
2. N (βA) = |β| N (A) . La norma es homogénea.
3. N (A + B) ≤ N (A) + N (B) . La norma satisface la desigualdad triangular.
4. N (AB) ≤ N (A) N (B) . Propiedad submultiplicativa.
Asi como antes, la notación usual para norma matricial es k·k . Además, tal como
lo discutimos después de la definición 7.2.1, la propiedad 1. se puede reescribir asi:
kAk ≥ 0 y kAk = 0 si y solo si A = 0.
Para matrices reales la definición es análoga. La primera norma que definimos es

una generalización de la 2-norma vectorial.
7.3. NORMAS DE MATRICES 65
7.3.2. Norma de Frobenius

Si A ∈ Rm×n , su norma de Frobenius se define asi:
! 21
X 1
kAkF = |aij |2 = traza AT A 2 ,
i,j
donde la traza de A fue definida en 1.1.4. Puede observarse que esta norma se pue-
de pensar como la norma 2 del vector obtenido al hacer un vector fila (columna)
consistente en todas las filas (columnas) de A, una enseguida de la otra.
7.3.3. Normas matriciales inducidas

Para normas vectoriales definidas en Rm y Rn se puede definir una norma matri-
cial en Rm×n de la siguiente manera:
kAk = máx kAxk , (7.2)

kxk=1
donde A ∈ Rm×n y x ∈ Rn .
La existencia de este máximo se debe a que toda norma vectorial es una función
continua y el disco kxk = 1 de Rn es un conjunto cerrado y acotado. Más detalles se
pueden encontrar en la sección 5.1 de [M].
De la definición (7.2) se puede deducir que para cualquier norma inducida se
cumple la siguiente condición de compatibilidad
kAxk ≤ kAk kxk . (7.3)
Entre las normas inducidas se destacan las que corresponden a las normas vec-
toriales 1, ∞ y 2. Todas ellas tienen expresiones alternativas a (7.2) para su cálculo
que se resumen asi:
Teorema 7.4. Para A ∈ Rm×n
1. m
X
kAk1 = máx kAxk1 = máx |aij | .
kxk1 =1 j
i=1
2. m
X
kAk∞ = máx kAxk∞ = máx |aij | .
kxk∞ =1 i
j=1
3. p
kAk2 = máx kAxk2 = ρ (AT A).
kxk2 =1
Demostración. Para demostraciones de estas igualdades recomendamos la sección

5.2 de [M].
7.4. Nota MATLAB

a = rand(8); % matriz cuadrada escogida de forma aleatoria
u = abs(3); % valor absoluto de numero real
k = complex(3,-4); % numero complejo 3 - 4i
t = abs(k); % valor absoluto o modulo
b = a(:,1); % vector con la primera columna de a
r1 = norm(a,inf); % norma inf de a
r2 = norm(b,1); % norma 1 de b
r3 = norm(a); % norma 2 de a (es la norma por defecto en MATLAB)
r4 = norm(a,’fro’); % Norma de Frobenius
Parte II
Análisis matricial
67
Capítulo 8
Perturbaciones, errores y
condicionamiento
Cuando se consulta en Google la frase Garbage in garbage out (abreviada GIGO)

aparecen más de 10 millones de resultados. El primero de ellos es el siguiente:
garbage in, garbage out
phrase of garbage
1.
used to express the idea that in computing and other spheres,
incorrect or poor quality input will always produce faulty output.
La frase garbage in, garbage out se puso de moda alrededor de 1960 y se refería
a todo tipo de errores, incluyendo errores en el programa de computador. En este
capítulo deseamos tratar ideas preliminares acerca de cómputos con datos inexactos
(perturbados.) Suponemos que los cómputos se hacen de forma exacta o con progra-
mas de computador correctamente escritos y que los errores (perturbaciones) en los
datos no se pueden evitar, ya sea porque vienen de algún proceso de aproximación
previo o porque se agregaron artificialmente para imitar datos obtenidos en la vida
diaria por medio de instrumentos de medición.
Distinguiremos dos categorías de matrices:
1. Matrices bien condicionadas: Algoritmos que se aplican a estas matrices
con datos perturbados siempre producen soluciones que son cercanas a las
soluciones exactas y las aproximan correctamente.
2. Matrices mal condicionadas: Cuando los diversos algoritmos se aplican a
a estas matrices con datos perturbados, pueden producir soluciones que son
lejanas a las soluciones exactas y que son inservibles como aproximaciones.
69
70 CAPÍTULO 8. PERTURBACIONES, ERRORES Y CONDICIONAMIENTO
Insistimos en que hay ocasiones en que las matrices mal condicionadas no causan
problemas en los cómputos y que la frase Garbage in garbage out no se puede aplicar
de forma general en este caso pero si sirve como motivación y para que estemos
prevenidos.
Un concepto del análisis numérico que se aplica mucho en álgebra lineal numé-
rica es el de estabilidad. Se usa para calificar algoritmos asi como el concepto de
condicionamiento se aplica a matrices.
Trefethen y Bau en [TB] pag. 104 dicen:
A stable algorithm gives nearly the right answer to nearly the right question (Un
algoritmo estable da una respuesta casi igual a la correcta cuando la pregunta es casi
la misma que la original.)
Por su parte, Ipsen en [I] pag. 54 dice:
An algorithm is (very informally) numerically stable in exact arithmetic if each
step in the algorithm is not worse conditioned than the original problem. (Un al-
goritmo es, hablando informalmente, estable numéricamente en aritmética exacta si
cada paso del algoritmo no es más mal condicionado que el problema original.)
No intentaremos dar definiciones para estabilidad pero advertimos que hay varias.
Nos quedaremos con las ideas expuestas en los dos libros citados arriba y confiamos
que sean suficientes.
En la redacción de este capítulo utilizamos principalmente [I], [vdG], [TB] y [Da].
En este capítulo consideramos la propagación de errores en suma, multiplicación e
inversión de matrices y en sistemas de ecuaciones. En los siguientes capítulos regre-
saremos a este tema varias veces.
8.1. Errores absoluto y relativo

Definición 8.1. Supongamos que el número real x∗ es una aproximación del número
real x. Entonces:
1. El error absoluto en x∗ es |x − x∗ | .
|x − x∗ |
2. Si x 6= 0, el error relativo en x∗ es .
|x|
Ejemplo 8.2. Supongamos x = 10. El error absoluto que resulta de aproximar x con
el número x∗ = 10,2 es 0,2 y el error relativo es 0,02.
Supongamos ahora que x = 0,1. El error absoluto que resulta de aproximar x con el
número x∗ = 0,2 es 0,1 y el error relativo es 1.
8.2. SISTEMA LINEAL 71
8.2. Sistema lineal

Las definiciones de error se extienden a vectores y matrices pero en lugar de valor
absoluto se usan normas. Para observar propagación de errores en datos, el caso más
sencillo de estudiar es el de un sistema lineal de dos ecuaciones con dos incógnitas
que corresponde a dos rectas en el plano. Si las rectas son casi paralelas se presenta
mal condicionamiento como lo muestra el ejemplo de la siguiente nota.
Nota 8.3. MATLAB

Pequeñas perturbaciones pueden causar grandes efectos en las soluciones de los pro-
blemas modificados con respecto a la solución exacta. El siguiente código trabaja con
un sistema lineal con una matriz 2 × 2 moderadamente mal condicionada.
a = [1 1/7;0.95 0.14];
ca = cond(a); % numero condicional de a
b = [1;0]; % lado derecho
x = a\b; % solucion exacta
a1 = [1 1/7;1 0.14]; % matriz perturbada
x1 = a1\b; % Solucion con a perturbada
ea1 = norm(x1-x); % error absoluto en x1
er1 = norm(x1-x)/norm(x); % error relativo en x1
b2 = [1;-0.2]; % lado derecho perturbado
x2 = a\b2; % solucion con lado derecho perturbado
ea2 = norm(x2-x); % error absoluto en x2
er2 = norm(x2-x)/norm(x); % error relativo en x2
los_b = [b,b2]; % los dos vectores b
los_x = [x x1 x2]; % las tres soluciones: x exacta, x1, x2
los_ea = [ea1 ea2]; % los dos errores absolutos
los_er = [er1 er2]; % los dos errores relativos
y1 = @(x) 7*(1-x); % recta 1
y2 = @(x) -(0.95/0.14)*x; % recta 2
xx = linspace(-100,100); % dominio
plot(xx,y1(xx),’r’,xx,y2(xx),’k–’) % grafica
grid on
hold on
plot(x(1),x(2),’kd’,’MarkerSize’,12,’MarkerFaceColor’,[0.5,0.5,0.5],...
’MarkerEdgeColor’,’r’)
plot(x1(1),x1(2),’ks’,’MarkerSize’,12,’MarkerFaceColor’,[0.5,0.5,0.5],...
’MarkerEdgeColor’,’b’)
Dos rectas casi paralelas

800
recta 1
600 recta 2
x
400 x1
x2
200
−200
−400
−600
−800
−100 −50 0 50 100
Figura 8.1: Mal condicionamiento
plot(x2(1),x2(2),’ko’,’MarkerSize’,12,’MarkerFaceColor’,[0.5,0.5,0.5],...
’MarkerEdgeColor’,’g’)
legend(’recta 1’,’recta 2’,’x’,’x1’,’x2’)
title(’Dos rectas casi paralelas’)
%%%

1 1
los_b = ,
0 −0,2
32,6667 −49 39,3333
los_x = , ,
−221,6667 350 −268,3333
los_ea = [577,4705 47,1405]
los_er = [2,5773 0,2104]
De estos resultados se puede concluir que perturbar esta matriz puede llevar a un
sistema con solución que en nada aproxima a la solución del sistema original. Una
perturbación del vector del lado derecho generó un sistema con solución menos lejana
a la solución del sistema original.
En el código aparece el número de condición de una matriz. Lo definimos a con-
tinuación.
8.3. CONDICIONAMIENTO DE NORMAS MATRICIALES 73
Definición 8.4. Para una norma matricial k·k dada, el número de condición de una
matriz invertible A se define como cond (A) = κ (A) = kAk kA−1 k .
Veremos que este número surge naturalmente en estudios de condicionamiento y
que 1 ≤ κ (A) < ∞. Mientras mayor es κ (A) más mal condicionada es la matriz.
8.3. Condicionamiento de normas matriciales

Las normas matriciales están bien condicionadas en el siguiente sentido.
Proposición 8.5. Si A, E ∈ Rm×n entonces
|kA + Ek − kAk| ≤ kEk .
Demostración. Por desigualdad triangular,
kA + Ek ≤ kAk + kEk y kA + Ek − kAk ≤ kEk .
Ahora, A = (A + E) − E y
kAk ≤ kA + Ek + kEk
luego
− kEk ≤ kA + Ek − kAk .
De la proposición anterior también podemos decir que prueba la continuidad de

la norma como función de las componentes de la matriz.
8.4. Suma y producto de matrices

Teorema 8.6. (Propagación de perturbaciones en la suma de matrices) Sean U, V, Uε , Vε ∈
Rm×n tales que U, V, U + V son no nulas y Uε , Vε se toman como perturbaciones de
U y V respectivamente. Entonces
k(Uε + Vε ) − (U + V )kp kUkp + kV kp
≤ máx {εU , εV } (8.1)
kU + V kp kU + V kp
donde
kUε − Ukp kVε − V kp
εU = y εV = . (8.2)
kUkp kV kp
Demostración. k(Uε + Vε ) − (U + V )kp ≤ kUε − Ukp + kVε − V kp

= kUkp εU + kV kp εV ≤ kUkp + kV kp máx {εU , εV } .
Teorema 8.7. (Propagación de perturbaciones en la multiplicación de matrices)

Sean U, Uε ∈ Rm×n y V, Vε ∈ Rn×q tales que U, V, UV son no nulas y como antes,
Uε , Vε se toman como perturbaciones de U y V respectivamente. Entonces
kUε Vε − UV kp kUkp kV kp
≤ (εU + εV + εU εV ) (8.3)
kUV kp kUV kp
donde εU y εV están dados por (8.2).
Demostración. Sean Uε = U+E y Vε = V +F. Entonces Uε Vε −UV = (U + E) (V + F )−
UV
= UF + EV + EF. Tomando normas, contando con desigualdad triangular y pro-
piedad submultiplicativa, llegamos a kUε Vε − UV kp ≤ kUkp kF kp + kEkp kV kp +
kEkp kF kp =
kF k kEk kEk kF k

kUkp kV kp kV kp + kU kp + kU kp kV kp .
p p p p
Nota 8.8. Observamos:

1. En lo sucesivo, si no hay confusión, se suprime el subíndice p de la norma.
kU k +kV k kU k kV k
2. Los números kUp+V k p y kUpV k p se denominan números de condición para la
p p
suma y la multiplicación respectivamente.
3. Los resultados (8.1) y (8.3) indican que el error que se genera al utilizar apro-
ximaciones de las matrices, es del mismo orden de magnitud que las perturba-
ciones en las matrices. Por extensión a lo establecido anteriormente, se dice
que las operaciones suma y producto son bien condicionadas o son estables
con respecto a perturbaciones en los datos.
8.5. Invertir una matriz

Ahora interesa saber cómo se propagan los errores cuando se invierte una matriz.
Teorema 8.9. (Inversa de identidad perturbada) Si E ∈ Rn×n y kEk < 1 entonces
I + E es no singular y
1 1
≤ (I + E)−1 ≤ .
1 + kEk 1 − kEk
8.5. INVERTIR UNA MATRIZ 75

Si además kEk ≤ 1
2
entonces (I + E)−1 ≤ 2.
Demostración.
Supongamos kEk < 1 y trabajemos por absurdo. Si I+E es singular, hay un vector
x 6= 0 tal que (I + E) x = 0. De aquí se deduce kxk ≤ kExk ≤ kEk kxk < kxk . Esto
es una contradicción. Para el primer estimado hacemos lo siguiente:
I = (I + E) (I + E)−1 (8.4)
y
1 = kIk = (I + E) (I + E)−1 ≤ (1 + kEk) (I + E)−1 .
Para la segunda desigualdad, de (8.4) obtenemos
I = (I + E)−1 + E (I + E)−1
y por tanto
1 = kIk ≥ (I + E)−1 − E (I + E)−1 . (8.5)

Por otro lado, E (I + E)−1 ≤ kEk (I + E)−1 de donde

(I + E)−1 − E (I + E)−1 ≥ (I + E)−1 − kEk (I + E)−1

= (I + E)−1 (1 − kEk) > 0.
Volvemos a (8.5),

1 = kIk ≥ (I + E)−1 − E (I + E)−1

= (I + E)−1 − E (I + E)−1

= (I + E)−1 (1 − kEk) .
1 1 1
Finalmente, si kEk ≤ 2
entonces 1 − kEk ≥ 2
y por tanto 1−kEk
≤ 2.
Para matrices generales se tiene:
Teorema 8.10. Sea A ∈ Rn×n no singular y sea E ∈ Rn×n tal que kA−1 Ek < 1.
Entonces, A + E es invertible y
kA−1 k
(A + E)−1 ≤ .
1 − kA−1 Ek

Si además, kA−1 k kEk ≤ 1
2
entonces (A + E)−1 ≤ 2 kA−1 k .
Demostración. Escribimos A + E = A (I + A−1 E) . Como kA−1 Ek < 1, el teorema

8.9 dice que I + A−1 E es invertible y por tanto A + E también es invertible. Su
inversa es −1 −1 −1
(A + E)−1 = A I + A−1 E = I + A−1 E A .
La submultiplicatividad y el teorema 8.9 permiten concluir el resultado deseado. Para
la afirmación final, basta tener en cuenta que kA−1 Ek ≤ kA−1 k kEk ≤ 21 .
Ya estamos listos para probar un teorema sobre el error que se comete cuando en
lugar de invertir una matriz A, se invierte una aproximación de A dada por A + E.
Teorema 8.11. (Inversa de una matriz perturbada) Sea A ∈ Rn×n no singular y sea
E ∈ Rn×n tal que kA−1 Ek < 1. Entonces,
kA−1 Ek
(A + E)−1 − A−1 ≤ A−1 .
1 − kA−1 Ek
Además, si kA−1 k kEk ≤ 12 , entonces

(A + E)−1 − A−1 kEk
≤ 2κ (A)
kA−1 k kAk
donde κ (A) = kAk kA−1 k .
Demostración. Llamemos F = A−1 E.
(A + E)−1 − A−1 = (I + F )−1 A−1 − A−1

= (I + F )−1 − I A−1
= − (I + F )−1 F A−1
usando la igualdad I = (I + F )−1 (I + F ) . Por teorema 8.9,

(A + E)−1 − A−1 ≤ (I + F )−1 kF k A−1
kF k
≤ A−1 .
1 − kF k
Para la segunda parte, si kA−1 k kEk ≤ 12 ,

(A + E)−1 − A−1 ≤ 2 kF k A−1
8.5. INVERTIR UNA MATRIZ 77
kEk
pues (I + F )−1 ≤ 2. Pero kF k ≤ kA−1 k kEk = kAk kA−1 k y asi
kAk

(A + E)−1 − A−1 −1 kEk
≤ 2 kAk A . (8.6)
kA−1 k kAk
kEk
Nota 8.12. Notemos que el número es el error relativo en A + E. Por tanto el
kAk
estimado (8.6) dice que una cota para el error relativo en (A + E)−1 es un múltiplo
del error relativo en A+E. ¿Qué múltiplo? 2κ (A) . Por consiguiente, mientras mayor
sea κ (A) mayor puede ser el error relativo en (A + E)−1 . Con razón al número κ (A)
se le llama número de condición para la matriz A.
Capítulo 9
Condicionamiento y sistemas lineales
En este capítulo nos interesa estimar el efecto que tienen perturbaciones en la

matriz y/o el lado derecho de un sistema lineal sobre la solución del sistema. La
estimación depende del número de condición de la matriz y sirve solamente para
advertir de una posible discrepancia significativa en la solución del nuevo problema.
Sean A ∈ Rn×n no singular y b ∈ Rn . Por teorema 4.8 el sistema
Ax = b (9.1)
tiene una única solución. Ahora, nos interesa considerar el sistema perturbado
(A + E) z = b + r (9.2)
y estimar el error en la nueva solución z. En primer lugar estudiamos el caso en el que

A no ha sido perturbada y después consideramos el caso general. Para la escritura
tuvimos en cuenta principalmente los libros [I] y [A].
9.1. Lado derecho perturbado

Consideremos el sistema
Ay = b + r (9.3)
que interpretamos asi: b + r es una aproximación de b. Si hacemos e = y − x y
restamos (9.3)-(9.1), obtenemos la ecuación para el error que es
Ae = r. (9.4)
79
80 CAPÍTULO 9. CONDICIONAMIENTO Y SISTEMAS LINEALES
Se dice que el error e satisface la misma ecuación lineal pero con el residual r
kek krk
como lado derecho. Nótese que y son los errores relativos en y y b + r
kxk kbk
respectivamente.
De (9.4), krk = kAek ≤ kAk kek y como e = A−1 r, entonces también kek ≤
kA−1 k krk . La primera desigualdad la dividimos por kAk kxk, la segunda desigualdad
la dividimos por kxk y obtenemos
krk kek kA−1 k krk
≤ ≤ .
kAk kxk kxk kxk
Enseguida usamos las desigualdades

kbk ≤ kAk kxk y kxk ≤ A−1 kbk
y llegamos a
1 krk kek −1 krk
≤ ≤ kAk A . (9.5)
kAk kA−1 k kbk kxk kbk
Con la notación establecida κ (A) = kAk kA−1 k , vemos claramente que a mayor
tamaño del número de condición κ (A) mayor es el rango en el que puede estar el
error relativo en y. Además, la desigualdad (9.5) también dice que el error relativo
en y está acotado por dos múltiplos constantes del error relativo en b + r.
9.2. Caso general

El caso general (9.2) se desarrolla de forma similar. La suposición que hacemos
es kA−1 k kEk < 1 que implica kA−1 Ek < 1, lo que garantiza que A + E es invertible.
El resultado que se da en este caso es

kek κ (A) krk kEk
≤ + . (9.6)
kxk kEk kbk kAk
1 − κ (A)
kAk
De teorema 8.10 sabemos

kA−1 k kA−1 k
(A + E)−1 ≤ ≤ .
1 − kA−1 Ek 1 − kA−1 k kEk
Ahora, de (9.2) si hacemos e = z − x, tenemos Ax + Ex + (A + E) e = b + r. Como
Ax = b, entonces
e = (A + E)−1 (r − Ex) .
9.3. ESTABILIDAD DE MÉTODOS DIRECTOS 81
Luego

kek = (A + E)−1 (r − Ex)

≤ (A + E)−1 kr − Exk
kA−1 k
≤ (krk + kEk kxk)
1 − kA−1 k kEk

kAk kA−1 k krk kEk
= + kxk .
−1
kEk kAk kAk
1 − kAk kA k
kAk
Luego

kek κ (A) krk kEk
≤ +
kxk kEk kAk kxk kAk
1 − κ (A)
kAk

κ (A) krk kEk
≤ + .
kEk kbk kAk
1 − κ (A)
kAk
De esta forma queda estimado el error relativo en z por medio de los errores relativos
en b + r y en A + E. El coeficiente depende de κ (A) y de la relación kA−1 k kEk < 1
kEk kEk
que corresponde a 1 − κ (A) > 0. Si = 0,001 y κ (A) = 900 entonces
kAk kAk
κ (A)
= 9000. De manera que existe la posibilidad de tener un rango amplio
kEk
1 − κ (A)
kAk
para el error relativo en z.
9.3. Estabilidad de métodos directos

La mayoría de los métodos directos para resolver sistemas lineales se basan en
factorizaciones como veremos en los capítulos finales de estas notas: LU, LLT , QR
son algunas de estas factorizaciones. El siguiente teorema analiza la propagación del
error cuando la matriz del sistema está factorizada. Se parte de la siguiente base:
Para el sistema Ax = b si A = S1 S2 , con S1 y S2 invertibles, entonces el sistema
se resuelve por medio de dos sistemas lineales
S1 y = b
S2 x = y.
82 CAPÍTULO 9. CONDICIONAMIENTO Y SISTEMAS LINEALES
El ejemplo que ya conocemos en estas notas es: A = LU. El sistema se resuelve

asi:
Ly = b
Ux = y.
Teorema 9.1. Sean A ∈ Rn×n no singular, b ∈ Rn no nulo y suponga que Ax = b.

Supongamos
kEk
A + E = S1 S2 , ε A =
kAk
kr1 k
S1 y = b + r 1 , ε 1 =
kbk
kr2 k
S2 z = y + r 2 , ε 2 = .
kyk
1
Si kA−1 k kEk ≤ 2
entonces
kz − xk
≤ 2κ (A) (εA + ε1 + ε)
kxk
donde −1 −1
S S1
ε= 2 ε (1 + ε1 ) .
(A + E)−1 2
Demostración. Ver Fact 3.14 de [I].

Capítulo 10
Sistemas lineales especiales
La estructura de la matriz A en un sistema Ax = b sugiere frecuentemente la

manera de resolver el sistema de una forma eficiente. El ejemplo más simple es cuan-
do la matriz A es triangular, ya sea superior o inferior. En esos casos la eliminación
de Gauss se reduce a una sustitución regresiva (back substitution) o progresiva (for-
ward substitution) respectivamente. Cuando la matriz del sistema A se obtiene de
la discretización de ecuaciones diferenciales por diferencias finitas o por elementos
finitos, es muy común que A sea de una de estas clases: diagonalmente dominante,
tridiagonal, tridiagonal por bloques o simétrica definida positiva. Para cada una de
ellas hay estrategias recomendadas para la solución del sistema Ax = b. Además, si
las ecuaciones diferenciales son no lineales, es de esperarse que se deba recurrir a un
método iterativo, por ejemplo, el método de Newton y que haya que resolver varios
sistemas lineales con matrices de algunas de las formas mencionadas arriba.
En este capítulo presentamos algunas matrices que vienen de discretización de
ecuaciones diferenciales. Para escribirlo utilizamos como referencias los libros [Da],
[M], [KC] y [O] además de las notas de clase [Me]. No es casualidad que dos de estos
libros y las notas sean de análisis numérico, es que estamos en un punto de gran
interacción entre varias disciplinas.
10.1. Problemas con valores en la frontera

La primera fuente de matrices especiales que estudiamos corresponde a ecuaciones
diferenciales ordinarias. Nos interesa la solución numérica del siguiente problema:
y ′′ = f (t, y, y ′) (10.1)
y (a) = α, y (b) = β
83
84 CAPÍTULO 10. SISTEMAS LINEALES ESPECIALES
Con unas condiciones de borde un poco más generales y bajo ciertas hipótesis
sobre f y sus derivadas parciales ft , fy y fy′ , Herbert B. Keller probó un teorema
de existencia y unicidad para el problema (10.1) que puede consultarse en la sección
8.9 de [KC].
Hay un caso particular que es de interés, es el caso lineal dado por
y ′′ = u (t) + v (t) y + w (t) y ′ (10.2)
y (a) = α, y (b) = β
Para discretizar las ecuaciones en los problemas (10.1) y (10.2), utilizamos dife-
rencias finitas, es decir, sustitución de derivadas por diferencias a partir del Teorema
de Taylor, que es probablemente el teorema del cálculo de mayor importancia en
análisis numérico.
El intervalo [a, b] lo dividimos en n + 1 subintervalos [ti , ti+1 ] , i = 0, . . . n, de
longitud uniforme h > 0. Aquí ti = a + ih, 0 ≤ i ≤ n + 1. La relación entre h y n es
b−a
h= .
n+1
La variable continua y la sustituimos por la variable discreta Y, por tanto, Yi es la
aproximación de y (ti ) . En los puntos interiores del dominio discreto, o sea de t1 a
tn , las derivadas las sustituimos por las siguientes diferencias:
Yi−1 − 2Yi + Yi+1
y ′′ se cambia por
h2
Yi+1 − Yi−1
y ′ se cambia por
2h
Además, las condiciones de borde las imponemos asi: Y0 = α, Yn+1 = β. Las
versiones discretas de (10.1) y (10.2) son:
Y0 = α (10.3)

Yi−1 − 2Yi + Yi+1 Yi+1 − Yi−1
2
= f ti , Yi , para i = 1, . . . n
h 2h
Yn+1 = β
y
Y0 = α (10.4)
Yi−1 − 2Yi + Yi+1 Yi+1 − Yi−1
2
= ui + vi Yi + wi para i = 1, . . . n
h 2h
Yn+1 = β,
10.1. PROBLEMAS CON VALORES EN LA FRONTERA 85
donde ui = u (ti ) , vi = v (ti ) y wi = w (ti ) .

Antes de continuar introducimos un concepto de álgebra lineal.
Definición 10.1. Una matriz A ∈ Rn×n se dice que es estrictamente diagonalmente
dominante (EDD) si X
|aii | > |aij | , i = 1, . . . n. (10.5)
j6=i
Proposición 10.2. Si A ∈ Rn×n es una matriz EDD entonces A es invertible.

Demostración. Supongamos que A es EDD y que existe un vector x 6= 0 tal que
Ax = 0. Sea
|xk | = máx |xj | .
1≤j≤n
n
X
Entonces |xk | > 0 y como Ax = 0 entonces aij xj = 0 para i = 1, . . . n. Luego
j=1
para k se tiene
X X

|akk | |xk | = akj xj ≤ |xk | |akj |

j6=k j6=k
y esto contradice (10.5).
10.1.1. Caso lineal

Continuamos con los problemas discretos (10.3) y (10.4) pero empezamos con el
más sencillo que es (10.4). Lo reordenamos
Y0 = α (10.6)

h h
−1 − wi Yi−1 + 2 + h2 vi Yi + −1 + wi Yi+1 = −h2 ui , i = 1, . . . , n
2 2
Yn+1 = β,
y usamos las siguientes notaciones:
h
ai = −1 − wi para i = 0, 1, . . . , n − 1
2
di = 2 + h2 vi para i = 1, . . . , n
h
ci = −1 + wi para i = 1, . . . , n
2
2
bi = −h ui para i = 1, . . . , n.
El sistema (10.6) se escribe en forma de sistema tridiagonal

    
d 1 c1 Y1 b1 − a0 α
 a1 d2 c2   Y2   b2 
    . 
 a d c   ..   . 
 2 3 3  .   . 
 .. .. ..    =  ..  (10.7)
 . . .    . 
    
 an−2 dn−1 cn−1     bn−1 
an−1 dn Yn bn − cn β

h
Si h es suficientemente pequeño de manera que wi < 1 para todo i y además
2
vi > 0 para todo i, entonces la matriz del sistema es EDD ya que

h h
2 + h2 vi > 2 = −1 − wi + −1 + wi
2 2

h h
= 1 + wi + 1 − wi .
2 2
Proposición 10.3. El sistema (10.7) tiene una única solución.
Demostración. Como la matriz del sistema es EDD, entonces es no singular por
proposición 10.2 y la conclusión se sigue de teorema 4.8.
10.1.2. Caso no lineal

Estamos interesados en la solución numérica del problema (10.3) que repetimos
aquí por comodidad
Y0 = α (10.8)

Yi−1 − 2Yi + Yi+1 Yi+1 − Yi−1
2
= f ti , Yi , para i = 1, . . . n
h 2h
Yn+1 = β.
Como la función f puede ser no lineal, planteamos una iteración de Newton ([O]
chapter 8) para resolver numéricamente el sistema (10.8) que reordenamos asi:

2 Yi+1 − Yi−1
Fi (Y ) = −Yi−1 + 2Yi − Yi+1 + h f ti , Yi , , i = 1, . . . , n,
2h
con Y0 = α y Yn+1 = β.
10.2. ECUACIÓN DEL CALOR 87
∂Fi
La matriz jacobiana de F es F ′ y el elemento (i, j) de F ′ es ∂Y j
. Por la caracte-
rística especial de Fi que depende solamente de Yi−1 , Yi y Yi+1 , entonces la matriz
jacobiana es tridiagonal y está dada por
 ∂F ∂F 
1 1
∂Y1 ∂Y2
 ∂F2 ∂F2 ∂F2 
 ∂Y1 ∂Y2 ∂Y3 
 ∂F3 ∂F3 ∂F3 
 ∂Y2 ∂Y3 ∂Y4 
F (Y ) = 
′
 .. .. .. 

 . . . 
 ∂Fn−1 ∂Fn−1 ∂Fn−1 
 ∂Yn−2 ∂Yn−1 ∂Yn 
∂Fn ∂Fn
∂Yn−1 ∂Yn
Los elementos genéricos son:

∂Fi h ∂f Yi+1 − Yi−1
= −1 − ti , Yi ,
∂Yi−1 2 ∂y ′ 2h

∂Fi 2 ∂f Yi+1 − Yi−1
=2+h ti , Yi ,
∂Yi ∂y 2h

∂Fi h ∂f Yi+1 − Yi−1
= −1 + ti , Yi ,
∂Yi+1 2 ∂y ′ 2h
en donde se sigue teniendo en cuenta que i = 1, . . . , n, Y0 = α y Yn+1 = β.
Recordemos que el método de Newton, en la iteración k − ésima, resuelve el
sistema lineal
F ′ Y (k) Z = −F Y (k)
y el nuevo iterado es Y (k+1) = Y (k) + Z. Es decir, en cada iteración para la solución
numérica del problema con valores en la frontera (10.1), se debe resolver un sistema
tridiagonal para el cual la eliminación de Gauss es mucho más sencilla.
10.2. Ecuación del calor

Esta sección es una parte del capítulo 7 de [Me]. La ecuación diferencial parcial
de tipo parabólico más sencilla es la ecuación unidimensional
∂u ∂2u
= 2 , a < x < b. (10.9)
∂t ∂x
La variable espacial x se toma en un intervalo [a, b] y la variable temporal t se toma
en el intervalo semi-infinito [0, ∞). Las condiciones de borde se definen para x = a y
x = b y para todo t. La condición inicial se define para x ∈ (a, b) y t = 0. La ecuación

(10.9) tiene validez en (a, b) × (0, ∞) . Esta ecuación da cuenta de procesos difusivos
lineales, por ejemplo, la conducción de calor. Debido a ésto, muy frecuentemente nos
referimos a la variable dependiente como temperatura, aunque este modelo se puede
utilizar en muchos otros procesos.
Consideremos discretizaciones uniformes en espacio y tiempo dadas por
xi = a + ih, i = 0, 1, ... y tj = jk, j = 0, 1, ...
donde h y k son los tamaños de paso en las direcciones de espacio y tiempo respec-
tivamente. Ambos son reales positivos. Es conveniente definir también
k
r= .
h2
10.2.1. Diferencias finitas

En todas las aproximaciones de diferencias finitas utilizamos la variable discreta
Vij para representar la aproximación de u en el punto de la malla (xi , tj ) que se
b−a
calcula por medio del método numérico. Sea h = y supongamos que tenemos
n+1
condiciones de borde de tipo Dirichlet para x = a = x0 y x = b = xn+1 , es decir,
los elementos V0j y Vn+1
j
son conocidos para todo j. Además, definimos la variable
j j T
discreta vectorial V = V1 , V2 , ..., Vnj para representar el vector de aproximaciones
j
en el nivel temporal j. Nótese que una condición inicial significa que el vector V 0 es
conocido.
10.2.2. Métodos más comunes

Distinguimos varias expresiones de diferencias finitas que aproximan (10.9).Veamos
algunas:
Método explícito La temperatura Vij+1 se escribe en términos de temperaturas

en niveles anteriores en la escala temporal.
j
Vij+1 − Vij Vi+1 − 2Vij + Vi−1
j
=
k h2
que conduce a
Vij+1 = Vij + r Vi−1
j
− 2Vij + Vi+1
j
. (10.10)
10.2. ECUACIÓN DEL CALOR 89
Sea A = trid (r, 1 − 2r, r) la matriz tridiagonal de orden n con elementos diagonales
1 − 2r y elementos super y sub diagonales iguales a r.
 
1 − 2r r
 r 1 − 2r r 
 
 .. .. .. 
A= . . . 
 
 r 1 − 2r r 
r 1 − 2r
La igualdad (10.10) se puede escribir
V j+1 = AV j (10.11)
lo que indica que
V j = Aj V 0 . (10.12)
Los iterados consecutivos se consiguen por aplicación de potencias de A a la condición
inicial.
Los métodos explícitos, para que sean útiles, requieren de una condición sobre los
parámetros de discretización. En este caso la condición es
1
r≤ . (10.13)
2
La utilidad no es un término matemático, en realidad lo que se busca es que el método
numérico sea estable y a métodos con restricciones se les llama condicionalmente
estables. Detalles sobre este tema pueden verse en [KC] y [Me].
Método implícito La temperatura Vij+1 se escribe en términos de temperaturas

en el mismo nivel y en niveles anteriores de la escala temporal,
Vij+1 − Vij V j+1 − 2Vij+1 + Vi−1
j+1
= i+1 ,
k h2
que se puede escribir
(1 + 2r) Vij+1 − rVi−1
j+1 j+1
− rVi+1 = Vij . (10.14)
Sea A = trid (−r, 1 + 2r, −r) la matriz tridiagonal de orden n con elementos diago-
nales 1 + 2r y elementos super y sub diagonales iguales a −r.
 
1 + 2r −r
 −r 1 + 2r −r 
 
 .. .. .. 
A= . . . 
 
 −r 1 + 2r −r 
−r 1 + 2r
La igualdad (10.14) se puede escribir
AV j+1 = V j . (10.15)
Es decir, en cada paso temporal, se debe resolver un sistema de ecuaciones lineales

con matriz simétrica, tridiagonal y diagonalmente dominante (por tanto no singular).
Agregamos finalmente que para este método implícito no hay una restricción
sobre los parámetros de discretización del estilo de (10.13). Este método es incondi-
cionalmente estable.
Métodos mixtos dependientes de un parámetro α ∈ [0, 1] La forma general

es
Vij+1 − Vij V j+1 − 2Vij+1 + Vi−1
j+1 j
Vi+1 − 2Vij + Vi−1
j
= α i+1 + (1 − α) .
k h2 h2
Los dos métodos anteriores son casos particulares de éste, basta hacer α = 0 y 1
respectivamente. El más importante de los métodos mixtos es el que corresponde a
1
α = que se llama de Crank-Nicolson. Está dado por la igualdad
2
Vij+1 − Vij 1 j+1

= 2 Vi+1 − 2Vij+1 + Vi−1
j+1 j
+ Vi+1 − 2Vij + Vi−1
j
k 2h
que también podemos escribir
j+1
− rVi−1 + 2 (1 + r) Vij+1 − rVi+1
j+1 j
= rVi−1 + 2 (1 − r) Vij + rVi+1
j
. (10.16)
Sean A = trid (−r, 2 (1 + r) − r) y B = trid (r, 2 (1 − r) , r) matrices de orden n. La

ecuación (10.16) se escribe como el sistema
AV j+1 = BV j ,
con V j+1 como incógnita. Este es el sistema que debe resolverse para cada iteración
temporal.
Anotamos finalmente que este método también es incondicionalmente estable.
10.3. Ecuación de Laplace

En dos variables, la ecuación de Laplace es
∆u = uxx + uyy = 0 (10.17)

10.3. ECUACIÓN DE LAPLACE 91
para u una función de (x, y) . Se define en regiones de R2 y no depende de la variable

temporal. Con base en la ecuación de Laplace, nos interesa presentar la solución
numérica de un sencillo problema conocido como problema de Dirichlet, que aparece
frecuentemente en ciencias e ingeniería.
Sea Ω un subconjunto abierto de R2 con frontera ∂Ω. Queremos resolver numéri-
camente un caso particular del siguiente problema de Dirichlet:
∆u = uxx + uyy = 0 en Ω (10.18)

u (x, y) = g (x, y) en ∂Ω
u es continua en Ω.
El caso particular de interés es utilizar el cuadrado unitario como dominio.
Ω = {(x, y) : 0 < x < 1, 0 < y < 1} .
Para trabajar por diferencias finitas debemos hacer una malla de puntos en todo
el cuadrado. Presentamos únicamente el caso en el que la malla es uniforme y el
parámetro de discretización es el mismo en las dos direcciones, x y y.
1
(xi , yj ) = (ih, jh) , i, j = 0, . . . , n + 1, h = .
n+1
Para el problema discreto utilizamos la variable Vij para aproximar a u (xi , yj ) . Nó-
tese que
V0j = g (0, yj ) , Vn+1

j
= g (1, yj ) , Vi0 = g (xi , 0) y Vin+1 = g (xi , 1) (10.19)
son conocidos, vienen de valores en la frontera.

Siguiendo con la técnica de diferencias finitas, el problema discreto que sirve para
resolver (10.18) es
j
Vi−1 − 2Vij + Vi+1
j
Vij−1 − 2Vij + Vij+1
+ =0
h2 h2
para los puntos interiores de la malla que corresponden a i, j = 1, . . . , n. Reordenan-
do, se llega a
4Vij − Vi−1
j j
− Vi+1 − Vij−1 − Vij+1 = 0. (10.20)
Consideremos el caso particular n = 4 ilustrado en la figura adjunta.
Deseamos escribir (10.20) en la forma Av = b con A una matriz 16 × 16, v el
vector que corresponde a los Vij de los nodos interiores de la malla y b el vector de
lado derecho que aparece debido a los valores conocidos en la frontera.
1 = y5
y4
y3
y2
y1
0 = y0 x1 x2 x3 x4 x5 = 1
0 = x0
Figura 10.1: Malla de 16 nodos interiores
Hagamos
v1 = V11 , v2 = V21 , v3 = V31 , v4 = V41 ,

v5 = V12 , v6 = V22 , v7 = V32 , v8 = V42 ,
v9 = V13 , v10 = V23 , v11 = V33 , v12 = V43 ,
v13 = V14 , v14 = V24 , v15 = V34 , v16 = V44 .
Con esta notación y basados en (10.20) y (10.19), el sistema lineal que se debe
10.3. ECUACIÓN DE LAPLACE 93
resolver es:
4v1 − v2 − v5 = V01 + V10 (10.21)

4v2 − v1 − v3 − v6 = V20
4v3 − v2 − v4 − v7 = V30
4v4 − v3 − v8 = V40 + V51
4v5 − v1 − v6 − v9 = V02
4v6 − v2 − v7 − v10 − v5 =0
4v7 − v3 − v8 − v11 − v6 =0
4v8 − v4 − v7 − v12 = V52
4v9 − v5 − v10 − v13 = V03
4v10 − v6 − v11 − v14 − v9 =0
4v11 − v7 − v12 − v15 − v10 =0
4v12 − v8 − v16 − v11 = V53
4v13 − v9 − v14 = V15 + V04
4v14 − v10 − v15 − v13 = V25
4v15 − v11 − v16 − v14 = V35
4v16 − v12 − v15 = V54 + V45
La matriz de este sistema es tridiagonal por bloques dada por
 
T −I
 −I T −I 
A=

 (10.22)
−I T −I 
−I T
donde
   
4 −1 −1
 −1 4 −1   −1 
T =

 y −I =
 
.

−1 4 −1 −1
−1 4 −1
El vector del lado derecho es

 
V01 + V10
 V20 
 
 V30 
 
 V40 + V51 
 
 V02 
 
 0 
 
 0 
 
 V52 
b=



 V03 
 0 
 
 0 
 
 V53 
 
 V1 + V04
5 
 
 V25 
 
 V35 
V54 + V45
10.4. Nota MATLAB

1. Para resolver problemas con valores en la frontera MATLAB ofrece las rutinas
bvp4c y bvp5c.
2. Al escribir doc gallery en la línea de comandos, aparece el sistema de ayuda

para una galería de matrices que se sugieren como ejemplo o como modelo. La
matriz (10.22) se conoce como poisson en esta galería. Las matrices tridiago-
nales aparecen en la galería como tridiag.
Capítulo 11
Primeros algoritmos de factorización

de matrices
La semejanza de matrices es una herramienta teórica de gran importancia, basta

pensar en el estudio de matrices diagonalizables, en la forma de Jordan y en el teore-
ma de triangularización de Schur, temas que se pueden consultar en muchos libros,
por ejemplo en [M] y [O]. La semejanza es un caso particular de factorización pero
hay otros casos de factorización de matrices y varios de ellos son de interés en álgebra
lineal numérica. Ya conocemos el algoritmo de la factorización LU. En este capítulo
la volvemos a considerar y también introducimos el algoritmo de Cholesky. Además,
explicamos la forma de resolver sistemas lineales por medio de una factorización
de la matriz del sistema. Las principales fuentes que consultamos son [A], [I], [Da],
[JR], [O] y [M]. Para que la exposición sea lo más completa posible mencionamos
brevemente los conceptos de diagonalización y triangularización de matrices.
11.1. Matriz diagonalizable

Definición 11.1. Dos matrices A y B de tamaño n × n se dice que son semejantes
si existe una matriz no singular P tal que P −1 AP = B.
Teorema 11.2. Si A y B de tamaño n×n son semejantes, entonces sus espectros son
iguales o sea σ (A) = σ (B) o en otras palabras, tienen los mismos valores propios.
Demostración. Basta demostrar que tienen el mismo polinomio característico. Su-
95
96CAPÍTULO 11. PRIMEROS ALGORITMOS DE FACTORIZACIÓN DE MATRICES
pongamos B = P −1AP. Entonces

det (A − λI) = det P −1 P det (A − λI)
= det P −1 det (A − λI) det P

= det P −1 (A − λI) P

= det P −1 AP − λI = det (B − λI) .
Definición 11.3. Una matriz cuadrada A se dice que es diagonalizable si es se-

mejante a una matriz diagonal. Por teorema 11.2, los valores propios de A aparecen
en la diagonal de la matriz diagonal semejante.
Definición 11.4. Un conjunto completo de vectores propios para una matriz cua-
drada A es cualquier conjunto de n vectores propios que es LI.
No todas
las matrices
tienen un
conjunto
completo
de vectores
propios.
Por ejem-
1 1 1 1 1
plo, A = es tal que A = , de manera que 1, es un
0 1 0 0 0
1 x1
par propio de A. Llamemos z = . No hay un vector x = tal que {z, x}
0 x2
es LI y ambos son vectores
propios
de
A. El
único valor propio de A es 1, de manera
x1 + x2 x1
que Ax = x significa = , por tanto x2 = 0 y resulta que x es
x2 x2
un múltiplo escalar de z por lo cual no puede ser que {z, x} sea LI. A las matrices
que no tienen un conjunto completo de vectores propios se les llama defectivas o
deficientes.
Definición 11.5. Se dice que una matriz U ∈ Cn×n es unitaria si sus columnas
constituyen una base ortonormal de Cn .
Definición 11.6. Se dice que una matriz U ∈ Rn×n es ortogonal si sus columnas
constituyen una base ortonormal de Rn .
Proposición 11.7. 1. U ∈ Cn×n es unitaria si y solo si U ∗ U = UU ∗ = I si y
solo si U −1 = U ∗ .
2. U ∈ Rn×n es ortogonal si y solo si U T U = UU T = I si y solo si U −1 = U.
Teorema 11.8. Una matriz cuadrada A tiene un conjunto completo de vectores
propios si y solo si es diagonalizable. Es decir, A tiene un conjunto completo de
vectores propios si y solo si existe una matriz invertible P tal que P −1 AP = D donde
D es una matriz diagonal.
11.1. MATRIZ DIAGONALIZABLE 97
Demostración. La idea principal es que si (λi , Vi ) , i = 1, . . . , n es un conjunto de

pares propios de A y el conjunto {V1 , V2 , . . . , Vn } es LI, entonces la matriz P =
(V1 V2 · · · Vn ) es no singular y P −1 AP = D con
 
λ1
 λ2 
 
D= .. .
 . 
λn
Teorema 11.9. (Schur) Toda matriz cuadrada real o compleja es unitariamente

semejante a una matriz triangular superior. Es decir, si A es n × n existen una
matriz unitaria U y una matriz triangular superior T tales que U ∗ AU = T. Además,
los elementos diagonales de T son los valores propios de A.
Demostración. Ver [M] chapter 7 o theorem 3.6 de [JR].
Teorema 11.10. (Ejes principales) Toda matriz cuadrada hermitiana (simétrica) A

es unitariamente (ortogonalmente) semejante a una matriz diagonal. Los elementos
en la diagonal de la matriz diagonal son los valores propios de A.
Demostración. Consultar theorem 7.4 de [A].
Nota 11.11. 1. Si U T AU = D es la semejanza ortogonal para la matriz simétrica

A, entonces las columnas de U conforman una base ortonormal de vectores
propios de A.
2. Toda matriz cuadrada es triangularizable, es decir, es semejante a una matriz

triangular superior. Sin embargo, no toda matriz cuadrada es diagonalizable.
3. Hay un teorema de Schur para matrices reales que usa matrices ortogonales para
la semejanza. Sin embargo, la matriz T no tiene que ser triangular superior en
este caso, pues en la diagonal pueden aparecer bloques 2 × 2 correspondientes
a los valores propios complejos de la matriz A. Para los detalles ver theorem
9.45 de [Da].
4. Sea A una matriz cuadrada. No puede haber un algoritmo basado en

funciones básicas de los elementos de la matriz, para lograr que en un
número finito de pasos se encuentre una matriz P tal que P −1 AP = T con
T triangular superior. Funciones básicas son: suma y resta, multiplicación y
división, además de raíces k − ésimas. Si tal procedimiento existiera estaría-

mos encontrando de forma sencilla los valores propios de la matriz A, lo que
significa que tendríamos un procedimiento basado en funciones elementales so-
bre los coeficientes de A para factorizar polinomios con coeficientes reales de
cualquier grado (polinomios característicos.) Esto no se puede hacer, pues el
famoso teorema de Abel dice: Si n ≥ 5 y p (t) es un polinomio de grado n con
coeficientes reales, entonces en general no hay manera de expresar las raices
de p (t) en términos de funciones básicas de los coeficientes de p (t) .
5. Sea p (t) un polinomio de grado n con coeficientes reales. Existe una matriz
llamada Matriz Compañera Cp de p (t) y es tal que el polinomio característico de
Cp es precisamente p (t) . Es decir, cada matriz tiene un polinomio característico
asociado y cada polinomio es el polinomio característico de una matriz.
6. Repetimos la frase definitiva de la página 192 en [TB]: Any eigenvalue solver

must be iterative. El objetivo de un algoritmo para aproximar valores propios
es producir una sucesión de números que convergen rápidamente a los valores
propios.
Un resultado útil que relaciona norma y radio espectral de una matriz es el

siguiente:
Proposición 11.12. Sea A una matriz cuadrada. Si λ es un valor propio de A,

entonces |λ| ≤ kAkp para 1 ≤ p ≤ ∞. En particular, ρ (A) ≤ kAkp .
Demostración. Sea (λ, x) un par propio de A. Omitimos el subíndice de la norma

por simplificar la escritura. Entonces x 6= 0 y podemos tomar kxk = 1.
|λ| = |λ| kxk = kλxk = kAxk ≤ kAk kxk = kAk .
11.2. Factorización LU
El algoritmo para la factorización LU con pivoteo parcial es el siguiente:
11.2.1. Algoritmo: Factorización LU con pivoteo parcial

Entrada: Matriz A ∈ Rn×n
11.2. FACTORIZACIÓN LU 99
Salida: Matriz de permutación P, matriz triangular inferior con 1′ s en la diagonal

L y matriz triangular superior U tales que P A = LU.
Pasos:
1. Si n = 1, entonces P = 1, L = 1 y U = A.
2. Si n > 1, escoja una matriz de permutación Pn tal que

α a
Pn A =
d An−1
donde α es el primer elemento de la primera columna, de arriba hacia abajo,
que tiene la máxima magnitud en dicha columna; esto indica que |α| ≥ kdk∞ .
Factorizamos
1 0 α a
Pn A =
l In−1 0 S
donde l = dα−1 y S = An−1 − la.
3. Calcular Pn−1 S = Ln−1 Un−1 donde Pn−1 es una matriz de permutación, Ln−1
es una matriz triangular inferior con 1′ s en la diagonal y Un−1 es triangular
superior.

1 0 1 0 α a
4. P = Pn , L = ,U= .
0 Pn−1 Pn−1 l Ln−1 0 Un−1
La eliminación de Gauss con pivoteo parcial sigue el siguiente algoritmo
11.2.2. Algoritmo: Eliminación de Gauss con pivoteo parcial

Entrada: Matriz invertible A ∈ Rn×n , vector b ∈ Rn .
Salida: Solución de Ax = b.
1. Factorizar P A = LU por algoritmo 11.2.1.

2. Resolver el sistema Ly = P b.
3. Resolver el sistema Ux = y.
Nota 11.13. El primer sistema a resolver tiene matriz triangular inferior con unos
en la diagonal, el cual se resuelve con una sustitución progresiva. El segundo sistema
tiene matriz triangular superior y por tanto se resuelve con una sustitución regresiva.
En los dos casos se trata de sistemas simplificados.
11.3. Factorización de Cholesky

En el caso en que A sea hermitiana definida positiva o simétrica definida positiva,
la eliminación de Gauss es mucho más sencilla pues existe la llamada factorización
de Cholesky que además es única.
11.3.1. Algoritmo: Factorización de Cholesky

Entrada: A ∈ Rn×n simétrica definida positiva
Salida: Matriz triangular inferior L con elementos diagonales positivos tal que
A = LLT .
√
1. Si n = 1, L = A.
2. Si n > 1, se realizan las siguientes particiones y factorizaciones:

1 1 1
α aT α2 0 1 0 α 2 α− 2 aT
A= = 1 ,
a An−1 aα− 2 In−1 0 S 0 In−1
donde S = An−1 − aα−1 aT .
3. Calcule S = Ln−1 LTn−1 , donde Ln−1 es triangular inferior con elementos diago-
nales positivos.
1
a2 0
4. L = 1 .
aα− 2 Ln−1
11.3.2. Nota MATLAB

L = chol(A,’lower’) calcula la matriz L tal que A = L ∗ LT . En MATLAB
se chequea la calidad del resultado invocando la siguiente norma: err = norm(A -
L*L’,’fro’). El número err debe ser del orden de 10−14 . La matriz A debe ser
simétrica definida positiva.
Capítulo 12
Factorización QR
Toda matriz real invertible A tiene una única factorización A = QR donde Q es

ortogonal y R es triangular superior con elementos diagonales positivos. En general,
para una matriz A ∈ Cm×n con m ≥ n, hay una matriz unitaria Q ∈ Cm×m y una
matriz
triangular
superior R ∈ Cn×n con elementos diagonales no negativos tales que
R
A=Q .
0
En el caso de un sistema lineal Ax = b con A invertible, la factorización A = QR
permite una fácil estrategia de solución para el sistema, a saber:
Ax = b
QRx = b
Rx = QT b
y el último sistema es triangular superior y por tanto se resuelve por sustitución

regresiva.
Para desarrollar este tema es conveniente tener presentes otros temas relaciona-
dos como los reflectores de Householder, rotaciones de Givens y proceso de Gram-
Schmidt. Fuentes para estos materiales son [A], [I], [TB] y [M].
12.1. Matrices de Householder

Definición 12.1. Una matriz H se llama matriz de Hessenberg si hij = 0 para
i > j + 1.
Consideremos de nuevo las notas que escribimos después del teorema de Schur
11.9. En la práctica, ¿cómo hacer para aproximar los valores propios de una matriz A?
101
102 CAPÍTULO 12. FACTORIZACIÓN QR
Se hace en dos pasos: en el primer paso, por medio de transformaciones de semejanza

sobre una matriz A podemos llegar a una matriz de Hessenberg o sea, es posible
encontrar una matriz invertible S tal que S −1 AS = H con H matriz de Hessenberg.
El primer paso es un método directo, que enseguida vemos cómo se realiza. El segundo
paso es un proceso iterativo: se construye una sucesión de matrices de Hessenberg
que converge a una matriz triangular superior.
Empezamos considerando una clase especial de matrices ortogonales llamadas
matrices o reflectores de Householder.
Definición 12.2. Si w ∈ Rn con kwk2 = 1, definimos la matriz de Householder

correspondiente a w asi
U = I − 2ww T . (12.1)
Proposición 12.3. Si U es una matriz de Householder entonces
U T = U, U 2 = I.
Es decir, U es simétrica y ortogonal.
Las matrices de Householder se usan para transformar un vector no nulo en uno

que tenga ceros en la mayoría de sus componentes.
Sea b 6= 0 un vector de Rn . Queremos encontrar una matriz de Householder U
tal que las componentes de Ub de la (r
+ 1) a la n sean cero para algún 1 ≤ r < n.
c
Escribimos m = n − r + 1 y b = con c ∈ Rr−1 , d ∈ Rm .
d
Escogemos
Or−1 Ir−1 0
w= yU= (12.2)
v 0 Im − 2vv T
con v ∈ Rm y kvk2 = 1. Pedimos

T
α
Im − 2vv d= (12.3)
0m−1
para algún α. Por ser Im − 2vv T ortogonal, |α| = kdk2 = S. Definimos
p = v T d.
De (12.3),
α
d − 2pv = . (12.4)
0m−1
12.2. LA FACTORIZACIÓN QR 103
A Q R
= 0
0
m×n m×m m×n
Figura 12.1: QR
Multiplicamos miembro a miembro por v T y usamos kvk2 = 1

T T T α
v d − 2pv v = −p y v = αv1 .
0m−1
Por tanto p = −αv1 . Lo sustituímos en la primera componente del vector en (12.4)
d1 + 2αv12 = α. (12.5)
Tomamos el signo de α de manera que sign (α) = −sign (d1 ) . De (12.5) se despeja
v12 y de aquí se obtiene v1 , cualquiera de los dos valores se puede usar. Enseguida se
consigue p = −αv1 y con p se sigue de (12.4) que
dj
vj = , j = 2, . . . , m.
2p
Ya queda precisado el vector w y por tanto la matriz ortogonal U.
12.2. La factorización QR
Sea A una matriz cuadrada real. Veamos que por medio de matrices de Househol-
der se pueden obtener Q ortogonal y R triangular superior tales que A = QR.
Sea
Pr = I − 2w (r) w (r)T , r = 1, . . . , n − 1
donde w (r) se escoge como en (12.2). Escribimos a A en términos de sus columnas
A = (A∗1 A∗2 · · · A∗n ) .
P1 A = (P1 A∗1 P1 A∗2 · · · P1 A∗n ) .

Escogemos b = A∗1 , P1 y w (1) con el procedimiento visto antes para que P1 A contenga
ceros debajo de la diagonal en su primera columna. Escogemos P2 de forma similar,
en tal forma que P2 P1 A tiene ceros en su segunda columna debajo de la diagonal.
En este caso b es la segunda columna de P1 A.
Cuando concluye este procedimiento, queda construída la matriz triangular su-
perior
R = Pn−1 · · · P1 A.
Si en el paso r de la construcción todos los elementos debajo de la diagonal de la
columna r son ceros, entonces se toma Pr = I. Finalmente se define
QT = Pn−1 · · · P1
y se obtiene A = QR.
12.3. Valores propios

Regresamos a los dos pasos descritos al principio de la sección 12.1 para calcular
los valores propios de una matriz. El paso 1 consiste en utilizar matrices de Hou-
seholder para reducir a A a una matriz semejante que es de Hessenberg. El paso 2
es un procedimiento iterativo que pasamos a explicar enseguida.
Sean Q ortogonal y R triangular superior tales que A = QR. Sea A1 = A.
Definimos la sucesión de matrices Am , Qm y Rm por
Am = Qm Rm , Am+1 = Rm Qm , m = 1, 2, . . . (12.6)
donde Am = Qm Rm es la factorización QR de Am .
Nótese que Am+1 = QTm Am Qm , lo que implica que todas las matrices de la sucesión
son ortogonalmente semejantes a A. En realidad
Am+1 = QTm · · · QT1 A1 Q1 · · · Qm . (12.7)

12.3. VALORES PROPIOS 105
Sean
Pm = Q1 · · · Qm y Um = Rm · · · R1 . (12.8)
La matriz Pm es ortogonal y la matriz Um es triangular superior.
De (12.7) obtenemos
Am+1 = PmT A1 Pm , m ≥ 1.
Enseguida usamos (12.7) con m en lugar de m − 1 para llegar a
Q1 · · · Qm−1 Am = A1 Q1 · · · Qm−1 .
Ahora,
Pm Um = Q1 · · · Qm Rm · · · R1
= Q1 · · · Qm−1 Am Rm−1 · · · R1
= A1 Q1 · · · Qm−1 Rm−1 · · · R1
= A1 Pm−1 Um−1 .
Como A1 = P1 U1 = Q1 R1 , iteración de la igualdad Pm Um = A1 Pm−1 Um−1 lleva a
Pm Um = Am
1 , m ≥ 1.
Este proceso iterativo es mucho más sencillo si la matriz a la que se le aplica ya

ha pasado por el paso 1 y por tanto es una matriz de Hessenberg.
Teorema 12.4. La sucesión Am definida en (12.6) converge a una matriz triangular

superior que contiene los valores propios de A en su diagonal. Si A es simétrica, la
sucesión converge a una matriz diagonal.
Demostración. Ver theorem 9.6 de [A].
Enseguida aplicamos teorema 9.1 para demostrar que el algoritmo QR es numé-

ricamente estable. En la demostración siguiente se requiere que la norma matricial
sea k·k2 . Los principales resultados sobre k·k2 que se usan aparecen en el siguiente
lema.
Lema 12.5. Sea A ∈ Rn×n . Entonces

1. AT 2 = kAk2 .
2. Si U es una matriz ortogonal n × n, entonces kUk2 = 1 y kUAk2 = kAk2 =

kUAkF . (Para norma Frobenius basta pensar en su definición como traza AT A .)
3. La norma de Frobenius es compatible con la norma 2, es decir,

kAxk2 ≤ kAkF kxk2


A1∗ x
 
(se piensa en Ax como vector de productos internos  ...  y se aplica
An∗ x
desigualdad de Cauchy-Schwarz.)
Demostración.
El espectro de AAT es el mismo de AT A pues si AAT y = µy entonces

2
AT A AT y = µAT y. Análogo en la otra dirección, por tanto AT 2 = ρ AAT =

ρ AT A = kAk22 .

Ahora, si U es unitaria, kUk22 = ρ U T U = ρ (I) = 1. Además, kUAxk22 = (UAx)T (UAx) =
xT AT U T UAx = xT AT Ax = kAxk22 , por tanto kUAk2 = kAk2 .
Teorema 12.6. Sean A ∈ Rn×n no singular, b ∈ Rn no nulo y suponga que Ax = b.
Supongamos
kEk2
A + E = QR, εA =
kAk2
kr1 k2
Qy = b + r1 , ε1 =
kbk2
kr2 k2
Rz = y + r2 , ε2 = .
kyk2
1
Si kA−1 k kEk ≤ 2
entonces
kz − xk2
≤ 2κ (A) (εA + ε1 + ε2 (1 + ε1 ))
kxk2
Demostración. Se aplica teorema 9.1 con S1 = Q, S2 = R. Haciendo las sustituciones
correspondientes, el factor que aparece en ese teorema acompañando a ε2 (1 + ε1 ) es
kR−1 k2 kQ−1 k2 2 T

. Por ser Q ortogonal, kQk 2 = ρ Q Q = ρ (I) = 1. Además,
k(A+E) k2
−1

(A + E)−1 = (QR)−1 = R−1 Q−1 = R−1 .
2 2 2 2
La última igualdad es de nuevo debido a la ortogonalidad de Q pues

−1 −1
R Q = Q−T R−T = QR−T = R−T = R−1 .
2 2 2 2 2
12.4. NOTA MATLAB 107
12.4. Nota MATLAB

[Q,R] = qr(A) donde A es m×n, produce R, una matriz m×n que es trapezoidal
superior y una matriz Q que es m × m unitaria tales que A = Q∗ R.La matriz
R
T
trapezoidal superior es triangular superior si m = n y es de la forma con T
0
triangular superior si m > n.
Capítulo 13
Descomposición en valores singulares

(SVD)
Esta es la factorización principal para matrices rectangulares, no solo porque se

usa en la definición de pseudoinversa y en la aproximación de mínimos cuadrados,
sino también porque es una de las factorizaciones que más información sobre la matriz
aporta, por ejemplo rango (rank), subespacio imagen (range) y núcleo. Pero estas
razones son insuficientes de acuerdo con Carla D. Martin y Mason A. Porter quienes
en [MP] exhiben una lista impresionante de aplicaciones de la SVD que va desde
ingeniería y ciencias hasta política partidista en Estados Unidos.
La SVD fue descubierta en la primera mitad del siglo XIX y de su desarrollo
temprano hay un recuento [S4] hecho por un especialista, el profesor G.W. Stewart.
Para escribir estas notas utilizamos principalmente los libros [Da], [I], [M], [JR] y
[TB].
Los principales temas que tratamos en el capítulo son la demostración del teorema
de la factorización SVD, las aproximaciones de rango 1, relación con el rango (rank)
y relación con los valores propios.
Iniciamos con un lema sobre matrices simétricas.
Lema 13.1. 1. Si A es real simétrica entonces sus valores propios son reales.
2. Si A es real simétrica definida (semidefinida) positiva entonces sus valores pro-
pios son positivos (no negativos).
Demostración. Supongamos Ax = λx con x 6= 0. Sea x el conjugado de x. Entonces
Ax = Ax y λx = λx. Como Ax = λx, entonces λ, x es un par propio de A. Además,
X X
xT x = xT x = xi xi = |xi |2 > 0 pues x 6= 0.
i i
109
110 CAPÍTULO 13. DESCOMPOSICIÓN EN VALORES SINGULARES (SVD)
A U > VT
= 0
0
m×n m×m m×n n×n
Figura 13.1: SVD
Concluimos que xT x = xT x es un número real positivo. Enseguida consideramos la

ecuación Ax = λx. Al premultiplicarla con xT llegamos a
xT Ax = xT λx = λxT x.
Por otro lado, de Ax = λx obtenemos
xT Ax = xT λx = λxT x = λxT x.
Finalmente usamos la simetría de A.
xT Ax = (Ax)T x = xT AT x = xT Ax
lo cual indica λxT x = λxT x. Cancelando, λ = λ, es decir λ ∈ R.
La positividad se desprende de Ax = λx con x 6= 0 pues
0 < xT Ax = λxT x = λ kxk22
lo cual implica λ > 0.
13.1. La factorización SVD

Teorema 13.2. (Descomposición en valores singulares SVD) Sea A ∈ Rm×n . En-
tonces existen matrices ortogonales U ∈ Rm×m , V ∈ Rn×n tales que
A = UΣV T ,
13.1. LA FACTORIZACIÓN SVD 111

D
donde Σ = si m ≥ n o Σ = F 0 si m < n, con
0
   
σ1 σ1
 σ2   σ2 
   
D= ..  oF = .. 
 .   . 
σn σm
y σ1 ≥ σ2 ≥ . . . ≥ σp ≥ 0 donde p = mı́n {m, n} . La última parte se acostumbra

resumir asi: Σ = diag (σ1 , σ2 , . . . , σp ) ∈ Rm×n . A pesar de haber hecho el enunciado
con esta generalidad, solo tratamos el caso m ≥ n.
Demostración. La matriz AT A es simétrica semidefinida positiva. Digamos que sus

valores propios, que son reales no negativos, son λi = σi2 , i = 1, . . . , n. Llamemos Vi
a sus vectores propios asociados. Supongamos además que σ1 ≥ σ2 ≥ . . . ≥ σr > 0 y
σr+1 = . . . = σn = 0. Hagamos Ψ1 = (V1 V2 . . . Vr ) , Ψ2 = (Vr+1 Vr+2 . . . Vn ) y V =
(Ψ1 Ψ2 ) . Con base en la primera observación después del teorema de ejes principales
11.10, las columnas de V son una base ortonormal de vectores propios de AT A y esto
significa
ViT AT AVi = σi2 y ViT AT AVj = 0 para i 6= j. (13.1)
Definimos los vectores

1
Ui = AVi para i = 1, . . . , r. (13.2)
σi
Forman un conjunto ortogonal pues
1 1 1
UiT Uj = (AVi )T AVj = V T AT AVj = δij
σi σj σi σj i

1 si i = j
donde δij es la delta de Kronecker definida asi: δij = .
0 si i 6= j
Sea Λ1 = (U1 U2 . . . Ur ) y escojamos Λ2 = (Ur+1 Ur+2 . . . Um ) con la condiciones
UjT A = 0 y kUj k2 = 1
que automáticamente permiten que U = (Λ1 Λ2 ) sea una matriz ortogonal cuyas
columnas son una base ortonormal de Rm . Finalmente,

 1 T T 
σ
V1 A
   11 V T AT 
 σ2 2 
U1T  .. 
 UT   . 
 2   
T  1
U AV =  ..  A (V1 V2 . . . Vn ) =  σr VrT AT  A (V1 V2 . . . Vn )
 .  
 UT 
UmT  r+1 
 .. 
 . 
T
Um
 1 2

σ
σ1 1
 1 2
σ 
 σ2 2 
 .. 
 . 
 
= 1 2
σ =Σ
 σr r 
 
 
 
Nota 13.3. 1. Los elementos de la diagonal σi se denominan valores singulares

de A.
2. Los valores propios de AT A y AAT son los mismos, es decir, {σi }ni=1 . Esto se
debe a que
AAT = UΣV T V ΣT U T = UΣΣT U T .
Pero ΣΣT = diag (σ12 , σ22 , . . . , σn2 ), es decir, AAT es unitariamente semejante a
la misma matriz diagonal que AT A.
3. Las columnas de U se denominan vectores singulares a izquierda. Se acostum-

bran escribir U = (U1 , U2 , . . . , Um ) .
4. Las columnas de V se llaman vectores singulares a derecha. Se acostumbra

escribir V = (V1 , V2 , . . . , Vn ) .
5. Como trabajamos solamente el caso m ≥ n, entonces hay n valores singulares

para A (pero puede haber repeticiones en este conteo, también puede haber
ceros.)
6. σ1 = máxi σi y σn = mı́ni σi .
7. Los valores singulares son únicos pero los vectores singulares no lo son. Sin
embargo, por abuso del lenguaje, se acostumbra hablar de la descomposición en
valores singulares de A.
Proposición 13.4. (SVD reducida) Si A = UΣV T es la SV D de A ∈ Rm×n (m ≥
n), entonces se puede obtener también la siguiente factorización
A = UΣ1 V T
donde U = (U1 , U2 , . . . , Un ) está compuesto por las primeras n columnas de U, por

tanto es una matriz de columnas ortonormales y Σ1 = diag (σ1 , σ2 , . . . , σn ) ∈ Rn×n .
Proposición 13.5. Sean σ1 ≥ σ2 ≥ . . . ≥ σn los valores singulares de A ∈ Rm×n
(m ≥ n). Entonces
1. kAk2 = σ1 .
2. σn = mı́n kAxk2 .
kxk2 =1
n
! 21
X
3. kAkF = σi2 .
i
1
4. Cuando A es n × n no singular, kA−1 k2 = σn
.
σ1
5. Cuando A es n × n no singular, κ2 (A) = kAk2 kA−1 k2 = σn
.
6. rank(A) = r, el número de valores singulares no nulos.

Demostración. 1. Se usa lema 12.5: kAk2 = UΣV T 2 = kΣk2 = máx σi = σ1 .

D
2. A = U V T . Sea z tal que kzk2 = 1 y kAzk2 = mı́n kAxk2 . Sea
0 kxk2 =1
T
y = V z.
X 1
T
mı́n kAxk2 = kAzk2 = UΣ1 V z 2 = kΣ1 yk2 = 2
σi |yi | 2 2
kxk2 =1
≥ σn kyk2 = σn .
Para la otra desigualdad,

σn = kΣ1 en k2 = U T AV en 2 = kAV en k2 ≥ mı́n kAxk2 .
kxk2 =1
X 12
T
3. kAkF = UΣV F = kΣkF = σi2 .
4. Como A es no singular entonces Σ es no singular, es decir, todos los valores

singulares son positivos. A−1 = V Σ−1 U T y el resultado sigue de parte 1.
5. Consecuencia de partes 1 y 2.
6. Como A = UΣV T con U, V no singulares, por definición 5.8 A ∼ Σ y por

proposición 5.21 tienen el mismo rango (rank). Además rank (Σ) = número de
valores singulares no nulos.
Los valores singulares están bien condicionados en el siguiente sentido.

Proposición 13.6. Sean A, A + E ∈ Rm×n , m ≥ n, con valores singulares σ1 ≥
· · · ≥ σn y η1 ≥ · · · ≥ ηn respectivamente. Entonces
|σ1 − η1 | ≤ kEk2 y |σn − ηn | ≤ kEk2 .
Demostración. Para σ1 y η1 se sigue de parte 1 de la proposición anterior y de

la proposición 8.5. Para σn y ηn pensamos en A como transformación lineal y en
especial como función continua definida en Rn . Tomamos y ∈ Rn tal que kyk2 = 1 y
kAyk2 = σn .
ηn = mı́n k(A + E) xk2 ≤ k(A + E) yk2 ≤ kAyk2 + kEyk2

kxk2 =1
= σn + kEyk2 ≤ σn + kEk2 .
Para la otra desigualdad, tomamos y tal que k(A + E) yk2 = ηn y kyk2 = 1.
σn = mı́n kAxk2 ≤ kAyk2 = k(A + E) y − Eyk2

kxk2 =1
≤ k(A + E) yk2 + kEyk2 = ηn + kEyk2 ≤ ηn + kEk2 .
Hay unas matrices de rango (rank) 1 estrechamente relacionadas con la SVD. Se

trata de los productos exteriores.

Proposición 13.7. Sean u ∈ Rm y v ∈ Rn no nulos. Entonces rank uv T = 1.
Demostración. Recordemos que uv T = (uv1 uv2 . . . uvn ) o sea todas las columnas
son múltiplos del vector u.
Proposición 13.8. Si A ∈ Rm×n con m ≥ n y rank (A) = r tiene la SVD dada en

teorema 13.2 entonces r
X
A= σj Uj VjT .
j=1
Demostración. Sean U = (Ur Um−r ) y V = (Vr Vn−r ) donde Ur , Um−r contienen las
primeras r columnas y las últimas m − r columnas de U respectivamente. Análogo
para V . Resulta que
 
V1T r
  X
A = Ur Σr VrT = (σ1 U1 · · · σr Ur )  ...  = σj Uj VjT , (13.3)
VrT j=1
donde Σr = diag (σ1 , σ2 , . . . , σr ) .

Nota 13.9. La fórmula (13.3) es una de las más aplicadas en lo que respecta a SVD.
Para una aplicación básica al procesamiento de imágenes sugerimos consultar sección
8.4 de [ND]. Los autores de este libro aprecian tanto esta aplicación que la usaron
para la carátula.
Teorema 13.10. Si A ∈ Rm×n con m ≥ n y rank (A) = r tiene la SVD A = UΣV T
dada en teorema 13.2 y U = (Ur Um−r ) , V = (Vr Vn−r ) donde Ur , Vr contienen las
primeras r columnas de U y V respectivamente, entonces:
1. Las columnas de Ur forman una base ortonormal para R (A) .
2. Las columnas de Vn−r forman una base ortonormal para ker (A) .

3. Las columnas de Vr forman una base ortonormal para R AT .

4. Las columnas de Um−r forman una base ortonormal para ker AT .
Demostración. 1. Sin pérdida de generalidad suponemos A 6= 0. De (13.3) A =

Ur Σr VrT y Σr es diagonal sin ceros en la diagonal.
 
σ1 V1T
 
Σr VrT =  ... 
σr VrT
y como V es unitaria entonces {V1 , . . . , Vr } es un conjunto X

LI. Afirmamos que
el conjunto {σ1 V1 , . . . , σr Vr } también es LI pues si 0 = αi σi Vi entonces
αi σi = 0 para todo i y como todos los σi son no nulos, entonces son los αi los
que deben anularse. Esto implica que Σr VrT tiene filas LI y por teorema 6.10,
R (A) = R (Ur ) .
2. Sin pérdida de generalidad suponemos ker (A) 6= {0} . De (13.3) A = Ur Σr VrT
y Σr es diagonal sin ceros en la diagonal.

Ur Σr = σ1 U1 · · · σr Ur

y estas columnas son LI, luego ker (A) = ker VrT . Por teorema 6.7,

ker (A) = ker VrT = R (Vr )⊥ .
El subespacio ortogonal de R (Vr ) es R (Vn−r ) y asi queda demostrado el teo-
rema.
3. Aplicar parte 1 a AT .
4. Aplicar parte 2 a AT .
13.2. Nota MATLAB

[U,S,V] = svd(A) sirve para efectuar la SVD de A. Las matrices U y V, que
cumplen A = USV T , son ortogonales (unitarias)
con m y n filas respectivamente. Si
D
A es m × n y m ≥ n, S es de la forma donde D es diagonal n × n.
0
s = svds(A,k)sirve para hallar los k valores singulares más grandes.
13.3. ¿Qué sigue?

Esta corta sección sirve para cerrar el capítulo sobre SVD y, al menos por ahora,
el manuscrito completo. La última sección del artículo [MP] tiene el siguiente título:
Everywhere you go, always take the SVD with you. En mi caso, esto es lo
que yo hago con el álgebra lineal numérica, pues me hace falta en mi trabajo
investigativo.
El objetivo es que muchos lectores de estas notas eventualmente se vuelvan usua-
rios o investigadores del álgebra lineal numérica y deban mantenerla a la mano.
Este manuscrito es una invitación al álgebra lineal numérica, no es un tratado ni
un libro de texto. Sugerencias para mejorarlo, son bienvenidas.
Bibliografía
[AM] Acosta, C.D. y C.E. Mejía, Cuaderno de Investiga-

ción: Métodos Iterativos para Sistemas Lineales, Universi-
dad Nacional de Colombia, Manizales, 2005. Disponible en
http://medellin.unal.edu.co/~cemejia/doc/Cuaderno_Cdacosta_Cemejia.pdf
[A] Atkinson, K.E. An introduction to numerical analysis, 2a. ed. John Wiley
and sons, 1989.
[C] Ciarlet, P.G. Introduction to numerical linear algebra and optimisation,

Cambridge University Press, 1989.
[Da] Datta, B.N. Numerical linear algebra and applications, 2a. ed. SIAM, 2010.
[D] Demmel, J.W., Applied Numerical Linear Algebra, SIAM, 1997.
[GO] Golub, G. y J.M. Ortega, Scientific computing, an introduction with parallel

computing, Academic Press, 1993.
[GV] Golub, G. y C.F. Van Loan, Matrix computations, 4a. ed. The Johns Hop-
kins University Press, 2013.
[HH] Higham, D.J. y N.J. Higham, Matlab Guide, 2a. ed. SIAM, 2005.
[HS] Hirsch, M.W. y S. Smale, Differential equations, dynamical systems and

linear algebra, Academic Press, 1974.
[HK] Hoffman, K. y R. Kunze, Algebra lineal, Prentice Hall Hispanoamericana,

1973.
[HJ] Horn, R.A. y C.R. Johnson, Matrix Analysis, Cambridge University Press,
1990.
[H] Hungerford, T.W. Algebra, Springer, 1987.
117
118 BIBLIOGRAFÍA
[I] Ipsen, I.C.F. Numerical Matrix Analysis, linear systems and least squares,
SIAM, 2009.
[JR] Johnson, L.W. y R.D. Riess, Numerical Analysis, 2a. ed. Addison-Wesley,
1982.
[KC] Kincaid, D. y W. Cheney, Numerical analysis, mathematics of scientific

computing, 3a. ed. Brooks/Cole Thomson Learning, 2002.
[Lang] Lang, S. Algebra lineal, Fondo Educativo Interamericano, 1974.
[L] Laub, A.J. Matrix analysis for scientists and engineers, SIAM, 2005.
[Layton] Layton, W. y M. Sussman, Numerical Linear Algebra, 2014, disponible en

http://www.lulu.com
[MP] Martin, C.D. y M.A. Porter, The extraordinary SVD. arXiv:1103.2338v5,

2012.
[Me] Mejía, C.E. Invitación al análisis numérico, 2002, Notas de clase disponibles
en http://medellin.unal.edu.co/~cemejia/doc/vinculos.html
[M] Meyer, C.D. Matrix analysis and applied linear algebra, SIAM, 2000.
[Mu] Murty, K.G. Computational and algorithmic linear algebra and n-

dimensional geometry, World Scientific Publishing, 2015.
[ND] Noble, B. y J.W. Daniel, Applied Linear Algebra, 3a. ed. Prentice-Hall,
1988.
[Ortega] Ortega, J.M. Matrix theory, a second course, Springer, 1987.
[O] Ortega, J.M. Numerical Analysis, a second course, SIAM, 1990.
[P] Poole, D. Algebra lineal, una introducción moderna, 2a. ed. Cengage Lear-
ning, 2007.
[SV] Saad, Y. y H.A. van der Vorst, Iterative solution of linear systems in
the 20th century, Journal of Computational and Applied Mathematics 123
(2000) 1-33.
[Se] Sewell, G. Computational methods of linear algebra, 2a. ed. Wiley-

Interscience, 2005.
BIBLIOGRAFÍA 119
[S1] Stewart, G.W., Matrix Algorithms Volume I Basic Decompositions, SIAM,

1998.
[S2] Stewart, G.W. Afternotes on Numerical Analysis, SIAM, 1996.
[S3] Stewart, G.W. Afternotes goes to Graduate School, SIAM, 1998.
[S4] Stewart, G.W. On the early history of the singular value decomposition,
SIAM Review 35 (1993) 551-566.
[S] Strang, G. Algebra lineal y sus aplicaciones, 4a. ed. Thomson, 2007.
[TB] Trefethen, L.N. y D. Bau III, Numerical linear algebra, SIAM, 1997.
[vdG] van de Geijn, R.A. Linear Algebra: Foundations to Frontiers. A Co-

llection of Notes on Numerical Linear Algebra, 2014. Disponible en
http://www.ulaff.net
[W] Watkins, D.A. Fundamentals of matrix computations, 2a. ed. Wiley-

Interscience, 2002.

Álgebra Linea Aplicada (Notas de Clase) PDF

Cargado por

Información del documentohacer clic para expandir la información del documento

Copyright:

Formatos disponibles

Álgebra Linea Aplicada (Notas de Clase) PDF

Cargado por

Información del documento

Título original

Derechos de autor

Formatos disponibles

Compartir este documento

Compartir o incrustar documentos

Opciones para compartir

¿Le pareció útil este documento?

¿Este contenido es inapropiado?

Copyright:

Formatos disponibles

Álgebra Linea Aplicada (Notas de Clase) PDF

Cargado por

Copyright:

Formatos disponibles

Universidad Nacional de Colombia

Invitación al álgebra lineal numérica

Trabajo presentado en el marco del

5. Forma escalonada reducida y rango 43

6. Matrices y transformaciones lineales 53

7. Normas de vectores y matrices 61

9. Condicionamiento y sistemas lineales 79

10.Sistemas lineales especiales 83

10.3. Ecuación de Laplace . . . . . . . . . . . . . . . . . . . . . . . . . . . 90

11.Primeros algoritmos de factorización de matrices 95

13.Descomposición en valores singulares (SVD) 109

ej : Vector canónico. Es el vector columna cuya componente j − ésima es 1 y las

En el primer capítulo introducimos muchas de las notaciones y definiciones que

Aunque la eliminación gaussiana se vislumbra en una obra china titulada Nine

Los datos históricos fueron tomados de

es una submatriz de A. Se llama principal si es cuadrada y sus elementos diagonales

1.1.2. Matrices especiales

1.1.3. Multiplicación por escalar

1.1.4. Traza de una matriz cuadrada

Es decir, es la suma de los elementos en la diagonal de A.

1.1.5. Suma de matrices

(A + B)ij = aij + bij .

La suma está bien definida en el sentido de la clausura: la suma de dos matrices de

1. Asociativa: Si A ∈ Rm×n , B ∈ Rm×n y C ∈ Rm×n , entonces

3. Conmutativa: Si A ∈ Rm×n y B ∈ Rm×n , entonces

4. Distributivas: Si A ∈ Rm×n , B ∈ Rm×n , α ∈ R y β ∈ R entonces

1.1.6. Combinación lineal

se le llama combinación lineal de los vectores Vi con coeficientes αi .

1.1.7. Producto interno o producto punto

Otras notaciones para el producto interno son v · w ó (v, w).

1.2. Multiplicación de matriz por vector

Ejemplo 1.1. La j−ésima columna de A se puede pensar como el producto A∗j =

1.2.2. Vector fila multiplicado por matriz

1.2.3. Producto externo (outer product)

Ejemplo 1.2. Una matriz de  Vandermonde

1.3. Nota MATLAB

y = xa’; % vector columna como traspuesta de vector fila

En el segundo capítulo sobre matrices y vectores definimos varias operaciones con

2.1. Producto de matrices

Pensando en las filas y columnas de las matrices A y B, este producto

1. AB = (AB∗1 AB∗2 ··· AB∗r ) .

3. (AB)ij = Ai∗ B∗j . Matriz de productos internos.

2.1.1. Propiedades del producto de matrices

2. A (BC) = (AB) C. Asociatividad.

Ejemplo 2.1. El producto no es conmutativo:

Nota 2.2. El producto de A consigo misma n veces se denota An .

Definición 2.3. Una matriz cuadrada es

3. Nilpotente si existe un entero positivo s tal que As = 0.

Definición 2.4. Sea A ∈ Rm×n . Su traspuesta, denotada AT , es la matriz cuya

Definición 2.5. Sea A ∈ Cm×n . Su traspuesta conjugada, denotada A∗ , es la matriz

2.1.2. Propiedades de las traspuestas

3. (βA)T = βAT , (βA)∗ = βA∗ .

2.1.3. Propiedades de productos de vectores

Definición 2.7. Una matriz A ∈ Cn×n es

donde εσ = 1 si σ es par y εσ = −1 si σ es impar y además denotamos Sn al grupo

Teorema 2.10. Propiedades de los determinantes

1. Si dos de las filas (columnas) de A son iguales, entonces det (A) = 0.

2. Si A tiene una fila (columna) nula, entonces det (A) = 0.