Off[General::spell1]
If[$VersionNumber == 5. , 
            HessenbergDecomposition[A_] = 
              Developer`HessenbergDecomposition[A]];

dim = 3;

Print["Transpose[A] es la transpuesta de A"]

\(\(\(
\)\(Print["Una matriz rectangular real"]

    Mrect = Table[Random[], {dim}, {dim + 3}]; 

    \(\(Mrect // MatrixForm\) // 
        Print;\)



    Print["Una matriz cuadrada con rango <dim "]

    \(MSing = Table[Random[], {dim - 1}, {dim}];\)

    \(MSing = Transpose[MSing] . MSing;\)

    MSing // MatrixForm

    Det[MSing]




    Print["Una matriz cuadrada real"]

    \(Mreal = Table[Random[], {dim}, {dim}];\)

    Mreal // MatrixForm


    Print["Una matriz cuadrada compleja"]

    \(Mcompl = Table[Random[Complex], {dim}, {dim}];\)

    Mcompl // 
      MatrixForm



    Print["Una matriz triangular superior"]

    \(MtrS = 
        Table[If[i > j, 0. , Random[]], {i, dim}, {j, 
            dim}];\)

    MtrS // MatrixForm


    Print["Una matriz simetrica"]

    \(MSim = Mreal + Transpose[Mreal];\)

    MSim // 
      MatrixForm



    Print["Una matriz antisimetrica"]

    \(MASim = Mreal - Transpose[Mreal];\)

    MASim // 
      MatrixForm



    Print["Una matriz hermitica"]

    \(MHer = Mcompl + Conjugate[Transpose[Mcompl]];\)

    MHer // 
      MatrixForm



    Print["Una matriz antihermitica"]

    \(MAHer = Mcompl - Conjugate[Transpose[Mcompl]];\)

    MAHer // MatrixForm

    Print["U otra matriz hermitica..."]

    \(MHer2 = I\ MAHer;\)

    MHer // 
      MatrixForm



    Print["Una matriz positiva"]

    \(MPosit = MHer . Conjugate[Transpose[MHer]];\)

    MPosit // 
      MatrixForm



    \)\)\)

\(\(\( (*Funciones\ invariantes\ sobre\ matrices\ cuadradas*) \)\(
\)\(
\)\(\(Print["Determinante"];\)

    Det[Mreal]

    \(Print["Traza"];\)

    Tr[Mreal]\)\)\)

(*Para trabajar con espacios filas y columnas*)

Print["Reduccion Gauss Jordan por filas de Mrect"]; 
RowReduce[Mrect] // MatrixForm
Print["Espacio nulo de Mrect"]; 
NS = NullSpace[Mrect]; NS // MatrixForm
Print["Mrect actuando sobre los generadores del espacio nulo (las \
filas de NS )"]; 
Mrect . Transpose[NS] // Chop // 
      MatrixForm

Inversas y pseudo inversas

Print["La pseudoinversa de Mrect"]
PInvMrect = PseudoInverse[Mrect]; 
PInvMrect // MatrixForm
Print["De un lado actua como la inversa"]
Mrect . PInvMrect // Chop // MatrixForm
Print["Pero del otro..."]
PInvMrect . Mrect // Chop // 
      MatrixForm

Print["La inversa de una matriz cuadrada"]
Minv = Inverse[Mreal]; 
Minv . Mreal // MatrixForm
Mreal . Minv // 
      MatrixForm


Print["La inversa a derecha de Mrect"]
MrectinvDer = 
      Transpose[Mrect] . 
        Inverse[Mrect . Transpose[Mrect]]; 
Mrect . MrectinvDer // 
      MatrixForm

Print["La inversa a izquierda de Mrect... Problemas: \
Transpose[Mrect].Mrect es singular...."]
MrectinvIz = 
      Inverse[Transpose[Mrect] . Mrect] . 
        Transpose[Mrect]; 
MrectinvIz . Mrect // 
      MatrixForm

Formas canonicas

{S, diag} = JordanDecomposition[Mreal]; 
Print["Descomposicion de Jordan: matriz diagonal"]; 
diag // MatrixForm
Print["Descomposicion de Jordan: matriz de similaridad"]; 
S // MatrixForm
Print["Descomposicion de Jordan: verificacion"]; 
Mreal -  S . diag . Inverse[S] // MatrixForm

Print["Caso singular:"]
{S, Jor}=JordanDecomposition[(({{0., 1., 0.}, {0., 0., 0.}, {0., 0., 0.}}))];
Print["Bloque de Jordan..."]
Jor // MatrixForm 
Print["Pero dado que minimas \
fluctuaciones numéricas nos sacan del caso singular, el resultado
de JordanDecomposition suele ser una diagonalizacion:"]
{S, Jor}=JordanDecomposition[Mreal.(({{0., 1., 0.}, {0., 0., 0.}, {0., 0., 0.}})).Inverse[Mreal]];
Jor // MatrixForm

Print["Caso hermitico:"]; 
S =. ; diag =. ; 
{S, diag} = JordanDecomposition[MHer]; 
Print["Los autovectores son ortogonales"]
Table[
          S[[i]] . Conjugate[S[[j]]], {i, dim}, {j, 
            dim}] // Chop // MatrixForm

Print["Descomposicion en valores singulares"]

{U, S, V} = 
      SingularValueDecomposition[Mrect]; 
U // MatrixForm
S // MatrixForm
V // MatrixForm
Print["Verificacion..."]
U . S . Transpose[V] - Mrect // MatrixForm

Print["Descomposicion de Schur"]; 
{Q, T} = SchurDecomposition[Mreal]; 
Q // MatrixForm
T // MatrixForm
Print["Verificacion"]; 
Q . T . Conjugate[Transpose[Q]] - Mreal // 
      MatrixForm

Print["Descomposicion QR"]; 
{Q, R} = QRDecomposition[Mreal]; 
Q // MatrixForm
R // MatrixForm
Print["Verificacion"]; 
Conjugate[Transpose[Q]] . R - Mreal // 
      MatrixForm

\(\(\(
\)\(Print["Hessemberg Decomposition"]

    \({P, H} = HessenbergDecomposition[Mreal];\)

    P // MatrixForm

    Print["H tiene la `forma de Hessenberg` (triangular superior + primera fila subdiagonal)"]\ 

    H // MatrixForm

    Print["P es unitaria"]

    P . Conjugate[Transpose[P]] // MatrixForm

    Print["Verificacion"]

    P . H . Conjugate[Transpose[P]] - Mreal // 
      MatrixForm






    \(Print["Hessemberg Decomposition sobre matrices simetricas"];\)

    P =. ; H =. ;

    \({P, H} = HessenbergDecomposition[MHer];\)

    P // MatrixForm

    Print["Si MHer es Hermitica o simétrica, H resulta tridiagonal"]\ 

    H // MatrixForm

    Print["P es unitaria"]

    P . Conjugate[Transpose[P]] // MatrixForm

    Print["Verificacion"]

    P . H . Conjugate[Transpose[P]] - MHer // 
      MatrixForm



    \)\)\)

$VersionNumber == 5.

Print["Norma de matrices"]


Print["Norm[Mrect] da la 'norma' definida por el mayor valor singular \
de Mrect"]
Norm[Mrect]
SingularValueDecomposition[Mrect][[2]][[1, 
      1]]

\(Print["Norm[Mrect,Frobenius] da la 'norma Eucliudea'"]
\)
Norm[Mrect, "Frobenius"]
Sqrt[
      Tr[Mrect . 
          Transpose[
            Conjugate[Mrect]]]]

Print["Otra forma de verlo: v es un vector cuyas componentes son las \
de la matriz"]
v = Flatten[Mrect]; 
(Sqrt[Conjugate[v] . v])

Print["Exponencial de una matriz"]

Mexp = MatrixExp[Mreal]; 
Mexp // MatrixForm

Print["Autovalores de la matriz original"]; 
Eigenvalues[
      Mreal]


Print["Autovalores de la matriz original exponenciados"]; 
Eigenvalues[Mreal] // Exp

Print["Autovalores de la matriz exponenciada"]; 
Eigenvalues[
      Mexp]


Print["Evaluacion como serie de potencias..."]; 
Mexp1 = IdentityMatrix[dim] + Mreal; 
Mexp2 = Mexp1 + Mreal . Mreal/2!; 
Mexp3 = 
      Mexp2 + Mreal . Mreal . Mreal/3!; 
Mexp4 = 
      Mexp3 + Mreal . Mreal . Mreal . 
            Mreal/4!; 
Print["Mejora al aumentar el orden del desarrollo..."]
G1 = {Mexp1 - Mexp // Norm, 
Mexp2 - Mexp // 
            Norm, 
Mexp3 - Mexp // 
            Norm, 
Mexp4 - Mexp // Norm} // 
        ListPlot[#1, PlotJoined -> True, 
            PlotStyle -> 
              Hue[ .1]] &

Print["Lo mismo, pero deshaciendonos de la traza..."]
lac = Tr[Mreal]/dim;
Mtmp = Mreal - lac IdentityMatrix[dim];
Print["Tr[Mtmp]=", 
        Tr[Mtmp]];

Mexp1 = IdentityMatrix[dim] + Mtmp;
Mexp2 = Mexp1 + Mtmp . Mtmp/2!;
Mexp3 = 
        Mexp2 + Mtmp . Mtmp . Mtmp/3!;
Mexp4 = 
        Mexp3 + Mtmp . Mtmp . Mtmp . 
              Mtmp/4!;
Print["Mejora al aumentar el orden del desarrollo..."]
G2 = {Exp[lac] Mexp1 - Mexp // 
            Norm, 
Exp[lac] Mexp2 - Mexp // 
            Norm, 
Exp[lac] Mexp3 - Mexp // 
            Norm, 
Exp[lac] Mexp4 - Mexp // Norm} // 
        ListPlot[#1, PlotJoined -> True, 
            PlotStyle -> 
              Hue[ .4]] &


Show[G1, 
      G2]