Matrices simétricas

Teorema espectral: si $M$ es simétrica, diagonaliza en una base ortonormal. Es decir, podemos escribir, para una matriz ortogonal $P$: $$ M = P^t\cdot D \cdot P. $$

Una matriz $P$ es ortogonal si $P^t\cdot P=P\cdot P^t=Id$.

Matrices definidas positivas

Una matriz simétrica es definida positiva si $$ \mathbf{x}^{t}M\mathbf{x}>0 $$ para cualquier vector $\mathbf{x}\neq \mathbf{0}$.

In [2]:
n = 5
M = diag(2*ones(1,n)) - diag(ones(1,n-1),1) - diag(ones(1,n-1),-1)

n =  5
M =

   2  -1  -0  -0  -0
  -1   2  -1  -0  -0
  -0  -1   2  -1  -0
  -0  -0  -1   2  -1
  -0  -0  -0  -1   2

El siguiente criterio ayuda a comprender las matrices definidas positivas:

Una matriz $M$ simétrica es definida positiva si y solo si todos sus autovalores son positivos.

In [3]:
eigs(M)

ans =

   3.73205
   3.00000
   2.00000
   1.00000
   0.26795
In [68]:
function esdf = es_definida_positiva(M)
    esdf = true;
    [nrows, ncols] = size(M);
    avals = eigs(M);
    for j = [1:nrows]
        l = avals(j);
        if l<=0
           esdf = false;
        end
    end
end
In [69]:
es_definida_positiva(M)

ans = 1
In [70]:
N = M - 0.29*eye(5);
es_definida_positiva(N)

ans = 0

La matriz inversa de $M$ es $$M^{-1} = P^t\cdot \left(\begin{array}{ccc} \frac{1}{d_1} & \dots & 0\\ \vdots & \ddots & \vdots\\ 0 & \dots & \frac{1}{d_n} \end{array}\right) \cdot P, $$ luego por lo tanto, $M^{-1}$ también es definida positiva.

In [39]:
for l=eigs(M)'
    1/l
end

ans =  0.26795
ans =  0.33333
ans =  0.50000
ans =  1
ans =  3.7321
In [41]:
Minv = inv(M);
eigs(Minv)

ans =

   3.73205
   1.00000
   0.50000
   0.33333
   0.26795

Ejercicio

Una matriz $M$ simétrica es definida positiva si y solo si todos sus menores principales tienen determinante positivo.

Dada una matriz $M$, el menor principal k-ésimo es la matriz formada por las primeras $k$ filas y las primeras $k$ columnas de $M$: $$ M_k=\left[ \begin{array}{cccc} a_{11} & a_{12} & \cdots & a_{1k}\\ a_{21} & a_{22} & \cdots & a_{2k}\\ \vdots & \vdots & & \vdots\\ a_{k1} & a_{k2} & \cdots & a_{kk} \end{array} \right] $$

  • Comprueba que tu función calcula correctamente si una matriz es definida positiva.
In [43]:
M = eye(5) - 0.1*ones(5);
es_definida_positiva(M)
M = eye(5) - 10*ones(5);
es_definida_positiva(M)

ans =  1
ans = 0
In [ ]:
In [ ]:

Raíces cuadradas de matrices

Si la matriz $M$ de dimensiones $n\times n$ es simétrica y definida positiva: $$M = P^t\cdot \left(\begin{array}{ccc} d_1 & \dots & 0\\ \vdots & \ddots & \vdots\\ 0 & \dots & d_n \end{array}\right) \cdot P, $$ y podemos definir una "matriz raíz cuadrada" de la forma: $$N=M^{1/2} = P^t\cdot \left(\begin{array}{ccc} \pm\sqrt{d_1} & \dots & 0\\ \vdots & \ddots & \vdots\\ 0 & \dots & \pm\sqrt{d_n} \end{array}\right) \cdot P. $$ Es fácil comprobar que $N$ es también simétrica y $N\cdot N = M$.

En general, una matriz cuadrada arbitraria tiene muchas raíces cuadradas, pero si $M$ es simétrica y definida positiva, exactamente una de sus raíces cuadradas es definida positiva.

In [53]:
M = [3 -1
    -1 3];
[P, D] = eig(M);
D_sqrt = diag(sqrt(diag(D)));
Pinv = inv(P);
M_sqrt = P*D_sqrt*Pinv

%Salvo error numerico, M_sqrt@M_sqrt es igual a M
M_sqrt*M_sqrt

M_sqrt =

   1.70711  -0.29289
  -0.29289   1.70711

ans =

   3.00000  -1.00000
  -1.00000   3.00000
In [62]:
% Curiosidad: todas las raíces cuadradas de M
M = [3 -1
    -1 3];
[P, D] = eig(M);
Pinv = inv(P);

for signo1 = [-1,1]
    for signo2 = [-1,1]
        fprintf('-----\n')
        signos = [signo1; signo2]
        D_sqrt = diag(signos.*sqrt(diag(D)));
        M_sqrt = P*D_sqrt*Pinv
        if es_definida_positiva(M_sqrt)
            fprintf('Es definida positiva\n')
        else
            fprintf('No es definida positiva\n')        
        end
        %Salvo error numerico, M_sqrt@M_sqrt es igual a M
        M_sqrt*M_sqrt
    end
end

-----
signos =

  -1
  -1

M_sqrt =

  -1.70711   0.29289
   0.29289  -1.70711

No es definida positiva
ans =

   3.00000  -1.00000
  -1.00000   3.00000

-----
signos =

  -1
   1

M_sqrt =

   0.29289  -1.70711
  -1.70711   0.29289

No es definida positiva
ans =

   3.00000  -1.00000
  -1.00000   3.00000

-----
signos =

   1
  -1

M_sqrt =

  -0.29289   1.70711
   1.70711  -0.29289

No es definida positiva
ans =

   3.00000  -1.00000
  -1.00000   3.00000

-----
signos =

   1
   1

M_sqrt =

   1.70711  -0.29289
  -0.29289   1.70711

Es definida positiva
ans =

   3.00000  -1.00000
  -1.00000   3.00000

Hay otra factorización de matrices simétricas y definidas positivas que también es muy interesante:

Descomposición LDLT: Una matriz simétrica y definida positiva $S$ admite una descomposición: $$ S = L\cdot D \cdot L^t $$ donde $L$ es una matriz triangular inferior con unos en la digonal.

Descomposición de Choleski: Una matriz simétrica y definida positiva $S$ admite una descomposición: $$ S = L\cdot L^t $$ donde $L$ es una matriz triangular inferior.

En este segundo caso, no estamos hablando de la raíz cuadrada, ya que $S$ no es $L\cdot L$, sino $L\cdot L^t$.

In [63]:
n = 5
M = diag(2*ones(1,n)) - diag(ones(1,n-1),1) - diag(ones(1,n-1),-1)

n =  5
M =

   2  -1  -0  -0  -0
  -1   2  -1  -0  -0
  -0  -1   2  -1  -0
  -0  -0  -1   2  -1
  -0  -0  -0  -1   2

Observación: esta matriz es tridiagonal, pero su inversa tiene todas sus entradas no nulas...

In [64]:
inv(M)

ans =

   0.83333   0.66667   0.50000   0.33333   0.16667
   0.66667   1.33333   1.00000   0.66667   0.33333
   0.50000   1.00000   1.50000   1.00000   0.50000
   0.33333   0.66667   1.00000   1.33333   0.66667
   0.16667   0.33333   0.50000   0.66667   0.83333

Ejercicio

  • Para la matriz simétrica que definimos antes (para un valor de $n$ de por ejemplo 100): $$A_5= \left(\begin{array}{rrrrr} 2.0 & -1.0 & 0.0 & 0.0 & 0.0 \\ -1.0 & 2.0 & -1.0 & 0.0 & 0.0 \\ 0.0 & -1.0 & 2.0 & -1.0 & 0.0 \\ 0.0 & 0.0 & -1.0 & 2.0 & -1.0 \\ 0.0 & 0.0 & 0.0 & -1.0 & 2.0 \end{array}\right) $$ encuentra el comando necesario y calcula la factorización de Cholesky.
  • Resuelve el sistema de ecuaciones $A\cdot \mathbf{x}= \mathbf{b}$, para $n=100$, y para el vector $\mathbf{b}$ que tiene entrada $i$ igual a $i(n-i)$, usando la factorización de Cholesky. Comprueba que obtienes el mismo resultado usando A\b, usando la factorización PLU, y usando la factorización de Cholesky.
  • (opcional) Mide los tiempos necesarios para resolver el sistema usando cada método.
In [ ]:
In [ ]:
In [ ]:

Matrices banda

Una matriz $n\times n$ es una matriz banda si existen dos enteros $p,q$ con $1\leq p,q< n$ con la propiedad que $a_{ij}=0$ cuando $p\leq j-i$ o $q\leq i-j$.

El ancho de banda de una matriz banda se define como $w=p+q-1$.

Por ejemplo, si $p=q=2$, la matriz es tridiagonal, y el ancho de banda es 3: tridiagonal

Ejercicio

Comprueba empíricamente, con varios ejemplos, los siguientes resultados:

  • Si una matriz banda $A$ con enteros $p$ y $q$ es de diagonal estrictamente dominante, entonces admite una descomposición LU (sin permutación). Nos preguntamo si:
    • $L$, además de ser triangular inferior, también es matriz banda, con el mismo $p$ de la matriz $A$ (pero $q=1$ al ser triangular inferior)
    • $U$, además de ser triangular superior, también es matriz banda, con el mismo $q$ de la matriz $A$ (pero $p=1$ al ser triangular superior)
  • Si $A$ es una matriz simétrica, definida positiva, y banda (necesariamente con $p=q$, por ser simétrica), entonces la matriz $L$ de su descomposición de Cholesky también es banda (¿pero cuánto valen $p$ y $q$?).
In [ ]:
In [ ]:

Ejercicio

Busque ejemplos,o contraejemplos, a las siguientes preguntas:

  • Si una matriz es triangular inferior, su traspuesta es triangular superior.
  • Si una matriz es triangular inferior, su inversa es triangular inferior.
  • Si una matriz es banda, ¿su inversa es banda?
  • Si una matriz es simétrica, ¿su inversa es simétrica?
  • Si una matriz es definida positiva, ¿su inversa es definida positiva?
  • Si una matriz es simétrica y definida positiva, ¿coincide necesariamente su factorización LU con su factorización de Cholesky?
In [ ]:
In [ ]:
In [ ]:
In [ ]:
In [ ]:
In [ ]: