05-ProblemasResueltos.Rmd

# Problemas Resueltos {#problems}

```{r echo = FALSE}
knitr::opts_chunk$set(echo = FALSE,
                      results = "asis",
                      message = FALSE,
                      warning = FALSE)

set.seed(4321)
``` 

En este capítulo planteamos varios problemas que recorren la mayoría de los contenidos de este libro y que intentan poner en práctica todas las técnicas de cada tema.

En cada problema resuelto, se hará referencia a la sección donde se ha dado la explicación que ayuda a la resolución del ejercicio.

## Espacios Vectoriales {#prob-ev}
```{r}
# Dimensión del espacio V
n <- 4#sample(4:5, size = 1)

# Definición de dos subespacios U y W
## Dimensiones
n_U <- sample(1:(n - 2), size = 1)
# n_W <- sample(1:(n - n_U), size = 1)
## Matrices aleatorias para sus cartesianas
A_U <- rMatrix(n = n_U, m = n, values = -1:3)
# A_W <- rMatrix(n = n_W, m = n, values = -1:3)

eqns_U <- write_system(A_U, b = matlab::zeros(c(nrow(A_U), 1))) %>%   str_flatten(", ")

# eqns_W <- write_system(A_W, b = matlab::zeros(c(nrow(A_W), 1))) %>%   str_flatten(", ")

vec0 <- matlab::zeros(c(n, 1))

# Hallar una base de U y W
genU <- solve_homogeneous(A_U)
bU <- linearly_independents(genU)
dim_U <- ncol(bU)

scalars <- rVector(n = dim_U)
genW <- cbind(bU %*% scalars, rMatrix(n, m = 1))
A_W <- parametric_to_cartesian(genW)$A %>% 
  remove_fraction()

# genW <- solve_homogeneous(A_W)
bW <- linearly_independents(genW)
dim_W <- ncol(bW)

# Calcular una base y la dimensión de:
# (a) U \cap W
A_Intersect <- rbind(A_U, A_W)
gen_Intersect <- solve_homogeneous(A_Intersect)
b_Intersect <- linearly_independents(gen_Intersect)
dim_Intersect <- ncol(b_Intersect)

# (b) U + W
gen_Sum <- cbind(bU, bW)
b_Sum <- linearly_independents(gen_Sum)
dim_Sum <- ncol(b_Sum)

# Comprobar el teorema de la dimensión
se_cumple <- dim_Sum == dim_U + dim_W - dim_Intersect

# Cambio de base
# Generar una base aleatoria y encontrar lo cambios a 
# y desde la canónica
repeat {
  
  B1 <- rMatrix(n = n, m = n * 4, values = -1:1)
  B1 <- linearly_independents(B1)

  if (ncol(B1) == n) break
  
}

# De B1 a C es trivial, es la propia B1.
# De C a B1 es su inversa.
P <- gauss_elimination(B1, eye(n), jordan = TRUE, diag1 = TRUE)

# Generar un vector aleatorio y expresarlo en las coordenadas
# de la base nueva
scalars <- matrix(sample(-2:4, size = n, replace = TRUE), ncol = 1)
v <- B1 %*% scalars
coord <- gauss_elimination(B1, v, jordan = FALSE, diag1 = FALSE)

unknowns <- c("x", "y", "z", "t", "u")[seq(n)]
generic_vector <- matrix(unknowns, ncol = 1)
```


__Ejercicio__

Consideremos como espacio vectorial $V = \mathbb{R}^{`r n`}$ y los subespacios siguientes:

```{r}
glue::glue("$$U = \\left\\{",
           "[to_latex(generic_vector)] \\in\\mathbb{R}^{[n]}:",
           "[str_flatten(write_system(A_U, b = zeros(n, 1), latex = TRUE, format = 'c'))]",
           "\\right\\}$$",
           .open = "[", .close = "]") %>% 
  cat()
glue::glue("$$W = \\left\\{",
           "[to_latex(generic_vector)] \\in\\mathbb{R}^{[n]}:",
           "[str_flatten(write_system(A_W, b = zeros(n, 1), latex = TRUE, format = 'c'))]",
           "\\right\\}$$",
           .open = "[", .close = "]") %>% 
  cat()

```

a. Hallar la dimensión y una base de $U$ y de $W$.
b. Hallar la dimensión y una base de $U\cap W$ y de $U+W$.
c. Comprobar el teorema de la dimensión.
d. Encontrar el subespacio suplementario de $U$.
e. Consideremos el vector $v$ dado por sus coordenadas en la base canónica \( `r glue::glue("{to_latex(v)}")` \). Consideremos también la base:

```{r}
glue::glue(
  "$$",
  "\\mathcal{B}_{1} = ",
  "\\left\\{[vectors_to_latex(B1)]\\right\\}",
  "$$",
  .open = "[", .close = "]"
) %>% 
  cat()
```

Encontrar las coordenadas de $v$ en la nueva base.

f. Encontrar las matrices de cambio de base entre $\mathcal{B}_1$ y $\mathcal{C}$ (en ambos sentidos), siendo $\mathcal{C}$ la base canónica.

__Solución__


a. Para encontrar unas bases de $U$ y de $W$ podemos seguir las indicaciones sobre [cómo encontrar un sistema generador y una base a partir de la descripción de un subespacio](#base).

En el caso del subsepacio $U$, el sistema homogéneo dado por su ecuaciones cartesianas es:

```{r}
glue::glue("$$\\left\\{",
           "[str_flatten(write_system(A_U, b = zeros(n, 1), latex = TRUE))]",
           "\\right.$$",
           .open = "[", .close = "]") %>% 
  cat()
```

Este sistema, resuelto por Gauss, nos da las siguientes expresiones paramétricas:
```{r}
params <- c("\\alpha", "\\beta", "\\delta", "\\gamma")
glue::glue(
  "$$",
  "[to_latex(generic_vector)] = ",
  "[write_linear_combination(genU, vars = params)]",
  "$$",
  .open = "[", .close = "]"
) %>% 
  cat()
```

De estas paramétricas podemos [extraer un sistema generador](#base) y, [a partir de él, la base de $U$](#sistgenabase):
```{r}
glue::glue(
  "$$",
  "\\mathcal{B}_{U} = ",
  "\\left\\{[vectors_to_latex(bU)]\\right\\}",
  "$$",
  .open = "[", .close = "]"
) %>% 
  cat()
```

Y, así, $\mathrm{dim}(U) = `r ncol(bU)`$.

Repetimos todo el proceso para el subespacio $W$. Tomamos su sistema homogéneo asociado:

```{r}
glue::glue("$$\\left\\{",
           "[str_flatten(write_system(A_W, b = zeros(n, 1), latex = TRUE))]",
           "\\right.$$",
           .open = "[", .close = "]") %>% 
  cat()
```

Lo resolvemos por Gauss para obtener las ecuaciones paramétricas del conjunto solución:
```{r}
params <- c("\\alpha", "\\beta", "\\delta", "\\gamma")
PW <- gauss_elimination(A_W, matlab::zeros(nrow(A_W), 1))
glue_latex(
  "[glue_matrices(A_W, matlab::zeros(nrow(A_W), 1), latex = TRUE, fractions = TRUE)]",
  "\\sim",
  "[glue_matrices(PW$splits, latex = TRUE, fractions = TRUE)]",
  "\\Rightarrow\\ ",
  "[to_latex(generic_vector)] = ",
  "[write_linear_combination(genW, vars = params)]"
) %>% 
  cat()
```

Por tanto, una base de $W$ será:
```{r}
glue::glue(
  "$$",
  "\\mathcal{B}_{W} = ",
  "\\left\\{[vectors_to_latex(bW)]\\right\\}",
  "$$",
  .open = "[", .close = "]"
) %>% 
  cat()
```

Y, así, $\mathrm{dim}(W) = `r ncol(bW)`$.

b. El subespacio intersección $U\cap W$ es [el conjunto de vectores que verifican tanto las ecuaciones cartesianas de $U$ como las de $W$](#interseccion). Es decir, son el conjunto solución del sistema homogéneo:
```{r}
glue::glue("$$\\left\\{",
           "[str_flatten(write_system(A_Intersect, b = zeros(n, 1), latex = TRUE))]",
           "\\right.$$",
           .open = "[", .close = "]") %>% 
  cat()
```

```{r}

if (ncol(b_Intersect) > 0) {
  
params <- c("\\alpha", "\\beta", "\\delta", "\\gamma")
glue::glue(
  "Este sistema, resuelto por Gauss, nos da las siguientes expresiones paramétricas:\n",
  "$$",
  "[to_latex(generic_vector)] = ",
  "[write_linear_combination(gen_Intersect, vars = params)]",
  "$$",
  .open = "[", .close = "]"
) %>% 
  cat()
  
} else {
  
  glue::glue("Este sistema es compatible determinado, luego su única solución es el vector 0.") %>% 
    cat()
  
}
```

```{r}
if (ncol(b_Intersect) > 0) {
  
glue::glue(
  "Por tanto, una base de $U\\cap W$ será:\n",
  "$$",
  "\\mathcal{B}_{U\\cap W} = ",
  "\\left\\{[vectors_to_latex(b_Intersect)]\\right\\}",
  "$$",
  .open = "[", .close = "]"
) %>% 
  cat()
  
} else {
  
  glue::glue(
    "Luego:",
    "$$",
    "U\\cap W = \\{ 0 \\}",
    "$$",
    .open = "[", .close = "]"
  ) %>% 
    cat()
  
}

```

Y, así, $\mathrm{dim}(U\cap W) = `r ncol(b_Intersect)`$.

Para la suma $U+W$, [sabemos que $\mathcal{B}_U\cup\mathcal{B}_W$ es un sistema generador](#suma):

```{r}
glue::glue(
  "$$",
  "\\mathcal{B}_{U}\\cup\\mathcal{B}_{W} = ",
  "\\left\\{[vectors_to_latex(cbind(bU, bW))]\\right\\}",
  "$$",
  .open = "[", .close = "]"
) %>% 
  cat()
```

Basta con usar eliminación gaussiana en ese conjunto de vectores para [eliminar aquellos que sean linealmente dependientes](#sistgenabase).
```{r}

s <- gauss_elimination(t(cbind(bU, bW)))

cat("$$", to_latex(cbind(t(cbind(bU, bW))), fractions = TRUE), "\\sim", to_latex(cbind(s$U[s$ind_rows, ]), fractions = TRUE), "$$")
```

Aquellas filas que hayan acabado siendo enteras de ceros, se corresponden con vectores que originalmente eran dependientes linealmente de los demás.

Por tanto, eliminándolos, tenemos una base de $U+W$:
```{r}
glue::glue(
  "$$",
  "\\mathcal{B}_{U + W} = ",
  "\\left\\{[vectors_to_latex(b_Sum)]\\right\\}",
  "$$",
  .open = "[", .close = "]"
) %>% 
  cat()
```

Y, así, $\mathrm{dim}(U+W) = `r ncol(b_Sum)`$.

c. Podemos comprobar [el teorema de la dimensión](#th-dim):

$$\begin{array}{ccccccc}
\mathrm{dim}(U+W) & = & \mathrm{dim}(U) & + & \mathrm{dim}(W) & - & \mathrm{dim}(U\cap V) \\
`r dim_Sum` & = & `r dim_U` & + & `r dim_W` & - & `r dim_Intersect` \\
\end{array}$$

d. Calculemos [el espacio suplementario de $U$](#supl).

Tomamos una base de $U$, como la calculada antes, y le aplicamos Gauss para obtener un sistema de vectores en forma escalonada.

```{r}
PU <- gauss_elimination(t(bU), jordan = TRUE)
M <- PU$U
# Pivotes
  if (nrow(M) == 1) {

    pivots_idx <- which(M != 0)[1]
    params_idx <- setdiff(seq(ncol(M)),
                          pivots_idx)

  } else {

    pivots_idx <- apply(M, 1, function(r) which(r != 0)[1])
    pivots_idx_final <- pivots_idx[!is.na(pivots_idx)]

    params_idx <- setdiff(seq(ncol(M)), pivots_idx_final)

  }

# Canónica
C <- eye(n)

# Solo aquellos vectores de los que no tienen 1 
# en las posiciones de los pivotes
b_suplU <- matrix(C[, -pivots_idx], nrow = n)
```

```{r}
glue_latex(
  "[glue_matrices(t(bU), latex = TRUE, fractions = TRUE)]",
  "\\sim",
  "[glue_matrices(PU$splits, latex = TRUE, fractions = TRUE)]"
) %>% 
  cat()
```

Siguiendo las instrucciones de la [sección acerca del espacio suplementario](#supl), vamos a fijarnos en qué columnas están los pivotes: \(\{ `r str_flatten(pivots_idx, ", ")` \}\). 

De la base canónica $\mathcal{C}$ eliminamos aquellos vectores que tienen un 1 en las mismas posiciones que los pivotes. Lo que nos quede, ésa es la base del espacio suplementario de $U$:
```{r}
glue_latex(
  "\\mathcal{B}_{\\mathrm{supl}} =",
  "\\left\\{",
  "[vectors_to_latex(b_suplU)]",
  "\\right\\}"
) %>% cat()
```

```{r}
params <- letters[seq(ncol(B1))]
str <- c("$", str_flatten(params, ", "), "$") %>% str_flatten()
```

e. Unos escalares `r str` serán [las coordenadas de $v$ en la base $\mathcal{B}_1$ si se verifica](#coord):

\small
```{r}
glue::glue(
  "$$",
  "[to_latex(v)] = [write_linear_combination(B1, vars = params)]",
  "$$",
  .open = "[", .close = "]"
  
) %>% 
  cat()
```

\normalsize
lo que equivale al siguiente sistema de ecuaciones lineales, que se puede resolver usando métodos Gaussianos, y nos da:
```{r}
glue::glue(
  "$$\\left\\{",
  "[str_flatten(write_system(B1, v, latex = TRUE, vars = params))]",
  "\\right.",
  "\\Rightarrow",
  " v = [to_latex(matrix(params, ncol = 1))]_{\\mathcal{B}_1} = [to_latex(scalars)]_{\\mathcal{B}_1}",
  "$$",
  .open = "[", .close = "]"
) %>% 
  cat()
```

f. El procedimiento para construir la [matriz de cambio de base](#cambiobase) de una base $\mathcal{B}$ a otra $\mathcal{B}'$ es el siguiente:

- Calcular las coordenadas de cada vector en $\mathcal{B}$ en la base $\mathcal{B}'$. Eso supone la resolución de $n$ sistemas de ecuaciones lineales, y el procedimiento es igual que el visto en la sección anterior.
- Construir la matriz poniendo, en forma de columna, las coordenadas calculadas: la primera columna se corresponde con las calculadas para el primer vector de $\mathcal{B}$, la segunda columna, para las del segundo, etc.

Otro aspecto teórico importante es que la matriz de $\mathcal{B}$ a $\mathcal{B}'$ es la inversa de la del cambio de base de $\mathcal{B}'$ a $\mathcal{B}$. Si una de ellas es fácil de calcular, la otra se puede construir usando, por ejemplo, el método de Gauss-Jordan para el cálculo de inversas.

Hay un caso fácil: el cambio de base desde $\mathcal{B}_1$ a la base canónica $\mathcal{C}$ es la propia matriz cuyas columnas son los vectores de $\mathcal{B}_1$:
```{r}
glue::glue(
  "$$P_{\\mathcal{B}_1\\to\\mathcal{C}} = ",
  "[to_latex(B1)]",
  "$$",
  .open = "[", .close = "]"
) %>% 
  cat()
```

La matriz del cambio de base opuesto la vamos a calcular [buscando la matriz inversa](#cambiobase) de la de $\mathcal{B}_1$ a $\mathcal{C}$.

Para ello, se puede usar la eliminación de Gauss-Jordan:
```{r}
glue::glue(
  "$$",
  "[glue_matrices(B1, eye(n), latex = TRUE)]",
  "\\sim",
  "[glue_matrices(P$splits, fractions = TRUE, latex = TRUE)]",
  "$$",
  .open = "[", .close = "]"
) %>% 
  cat()
```

Luego
```{r}
glue::glue(
  "$$P_{\\mathcal{C}\\to\\mathcal{B}_1} = ",
  "[to_latex(P$splits[[2]],fractions = TRUE)]",
  "$$",
  .open = "[", .close = "]"
) %>% 
  cat()
```


## Aplicaciones lineales {#prob-apli}

```{r, results = "asis"}
unknowns <- c("x", "y", "z", "t")
n <- dimV <- sample(3:4, size = 1)
m <- dimW <- sample(2:3, size = 1)

A <- rMatrix(n = dimW, m = dimV)

str0 <- to_latex(matrix(unknowns[seq(dimV)], ncol = 1))

str <- write_system(A, latex = TRUE, format = "c")

repeat {
  
  B1 <- rMatrix(n = n, m = n * 4, values = -2:1)
  B1 <- linearly_independents(B1)

  if (ncol(B1) == n) break
  
}

repeat {
  
  B2 <- rMatrix(n = m, m = m * 4, values = 0:2)
  B2 <- linearly_independents(B2)

  if (ncol(B2) == m) break
  
}

gen_vector_V <- matrix(get_unknowns(m = n, 
                                    latex = TRUE), 
                       ncol = 1)
```


__Ejercicio__

Consideremos la aplicación $f:V=\mathbb{R}^{`r dimV`}\to W=\mathbb{R}^{`r dimW`}$ dada por:

```{r results="asis"}
cat("$$f", str0, " = \\left(", str, "\\right)$$")
```

a. Encontrar la matriz asociada a $f$ en las bases canónicas y en las bases:
```{r}
glue_latex(
  "\\mathcal{B}_{1} = \\{v_i\\} = ",
  "\\left\\{[vectors_to_latex(B1)]\\right\\}"
) %>% 
  cat()

glue_latex(
  "\\mathcal{B}_{2} = \\{w_i\\} = ",
  "\\left\\{[vectors_to_latex(B2)]\\right\\}"
) %>% 
  cat()
```

b. Dar una base y la dimensión del núcleo y de la imagen de $f$.

c. Identificar si $f$ es inyectiva o sobreyectiva.

d. Comprobar el teorema de la dimensión de rango y nulidad.

e. Calcular una base de $f(U)$, donde $U$ es el subespacio de $V$ dado por:
```{r}
U <- rMatrix(n = n, m = 2)
AU <- parametric_to_cartesian(U)$A
glue_latex(
  "U = \\left\\{",
  "[to_latex(generic_vector_V)]\\in\\mathbb{R}^{[n]}:",
  "[write_system(AU, matlab::zeros(nrow(AU), 1), latex = TRUE, fractions = TRUE, format = 'c')]",
  "\\right\\}"
) %>% 
  cat()
```


__Solución__


a. Para hallar la [matriz asociada](#matriz-apli) a la aplicación $f$ en las bases canónicas, basta con mirar los coeficientes de las variables $x,y,\ldots$ en la expresión de $f$ y ponerlos por columnas. De esta forma:
```{r}
glue_latex(
  "A =",
  "[glue_matrices(A, latex = TRUE, fractions = TRUE)]"
) %>% 
  cat()
```

Para calcular la matriz $A'$ asociada a $f$ [con respecto a las nuevas bases](#cambiobase-apli), tenemos dos opciones:

1. Calcular la imagen de los vectores de $\mathcal{B}_1$ y encontrar sus coordenadas con respecto a la base $\mathcal{B}_2$, y ponerlas por columnas.
2. Calcular las matrices $P$ y $Q$, de [cambio de base](#cambiobase) desde la canónica en $V$ y en $W$ a $\mathcal{B}_1$ y $\mathcal{B}_2$, respectivamente, y usarlas para calcular $A' = Q^{-1}\ A\ P$.

Ambos caminos son completamente equivalentes y, de hecho, uno se deduce del otro. Por tanto, en cada situación debemos ver cuál de los dos nos puede resultar más sencillo.

En este caso, optamos por el [camino 1](#matriz-apli).

Llamemos $Y = (f(v_1)|f(v_2)|\ldots|f(v_n))$ a la matriz que tiene por columnas las imágenes de los vectores de $\mathcal{B}_1$ mediante $f$:
```{r}
Y <- A %*% B1
glue_latex(
  "Y = ",
  "[glue_matrices(Y, latex = TRUE, fractions = TRUE)]"
) %>% 
  cat()
```
y $B_2 = (w_1|w_2|\ldots|w_m)$ a la matriz que tiene por columnas los vectores de $\mathcal{B}_2$.

Entonces $A' = B_2^{-1}\ Y$, que podemos calcular mediante Gauss-Jordan, partiendo de $(B_2|Y)$ y llegando a $(I_m|A')$:
```{r}
Ap <- gauss_elimination(B2, Y, jordan = TRUE, diag1 = TRUE)
glue_latex(
  "[glue_matrices(B2, Y, latex = TRUE, fractions = TRUE)]",
  "\\sim",
  "[glue_matrices(Ap$splits, latex = TRUE, fractions = TRUE)]"
) %>% 
  cat()
```

Luego 
```{r}
glue_latex(
  "A' =",
  "[glue_matrices(Ap$splits[[2]], latex = TRUE, fractions = TRUE)]"
) %>% 
  cat()
```

b. Comenzamos por estudiar el [núcleo](#nucleo) de $f$.

Sabemos que el núcleo $\mathrm{Ker}\ f$ es el conjunto de vectores $v\in V$ tal que $f(v) = 0$, y que eso coincide con el conjunto de soluciones del sistema de ecuaciones homogéneo cuya matriz de coeficientes es $A$.

Todo esto se resume en que las ecuaciones cartesianas de $\mathrm{Ker}\ f$ son aquellas cuya matriz de coeficientes es $A$:
```{r}
glue_latex(
  "\\mathrm{Ker}\\ f =",
  "\\left\\{",
  "[to_latex(gen_vector_V)]\\in V:",
  "[write_system(A, matlab::zeros(m, 1), latex = TRUE, format = 'c')]",
  "\\right\\}"
) %>% 
  cat()
```

Resolvemos este sistema por Gauss para [encontrar las ecuaciones paramétricas](#base) y, a partir de ahí, una [base](#sisgenabase) del núcleo:
```{r}
zeros_vec <- matlab::zeros(m, 1)
gauss_Ker <- gauss_elimination(A, zeros_vec, jordan = TRUE)
glue_latex(
  "[glue_matrices(A, zeros_vec, latex = TRUE, fractions = TRUE)]",
  "\\sim",
  "[glue_matrices(gauss_Ker$splits, latex = TRUE, fractions = TRUE)]"
) %>% 
  cat()
```
Luego nos quedan las siguientes ecuaciones paramétricas para $\mathrm{Ker}\ f$:
```{r}
genKer <- solve_homogeneous(A)
params <- c("\\alpha", "\\beta", "\\gamma", "\\delta")
glue_latex(
  "[to_latex(gen_vector_V)] =",
  "[write_linear_combination(genKer, vars = params)]"
) %>% 
  cat()
```

Entonces la base del núcleo de $f$ es:
```{r}
glue_latex(
  "\\mathcal{B}_{\\mathrm{Ker}\\ f} = ",
  "\\left\\{",
  "[vectors_to_latex(genKer)]",
  "\\right\\}"
) %>% 
  cat()
```
y \(\mathrm{dim}(\mathrm{Ker}\ f) = `r ncol(genKer)`\).

Con respecto a la imagen de $f$, hemos de recordar que [$\mathrm{Im}\ f$ es el subespacio de $W$ generado por los vectores que forman las columnas de la matriz asociada $A$](#imagen).

Por tanto, un [sistema generador](#gen) de $\mathrm{Im}\ f$ será:
```{r}
glue_latex(
  "\\mathcal{G} =",
  "\\left\\{",
  "[vectors_to_latex(A)]",
  "\\right\\}"
) %>% 
  cat()
```

Desde este sistema generador, [eliminamos los vectores linealmente dependientes](#sisgenabase) para obtener una base de $\mathrm{Im}\ f$. Usamos Gauss con los vectores de $\mathcal{G}$ por _filas_, y eliminamos aquellos que, al final, se hayan convertido en filas de ceros:
```{r}
gauss_Im <- gauss_elimination(t(A))
glue_latex(
  "[glue_matrices(t(A), latex = TRUE, fractions = TRUE)]",
  "\\sim",
  "[glue_matrices(gauss_Im$splits[[1]], latex = TRUE, fractions = TRUE)]"
) %>% 
  cat()
```

Luego
```{r}
bImf <- linearly_independents(A)
glue_latex(
  "\\mathcal{B}_{\\mathrm{Im}\\ f} =",
  "\\left\\{[vectors_to_latex(bImf)]\\right\\}"
) %>% 
  cat()
```
y \(\mathrm{dim}(\mathrm{Im}\ f) = `r ncol(bImf)`\).

c. Para determinar si $f$ es [inyectiva](#inyectividad) o [sobreyectiva](#sobrey), nos fijaremos en las dimensiones de su núcleo y de su imagen:

- $f$ será inyectiva si y solo si se cumple $\mathrm{dim}(\mathrm{Ker}\ f) = 0$. Sabemos que \(\mathrm{dim}(\mathrm{Ker}\ f) = `r ncol(genKer)`\), luego $f$ `r ifelse(ncol(genKer) == 0, "", "no")` es inyectiva.
- $f$ será sobreyectiva si y solo si se cumple que \(\mathrm{dim}(\mathrm{Im}\ f) = \mathrm{dim}(W) = \mathrm{dim}(\mathbb{R}^{`r m`})\), y esto `r ifelse(ncol(bImf) == m, "", "no")` se cumple, como hemos visto en el apartado anterior, luego $f$ `r ifelse(ncol(bImf) == m, "", "no")` es sobreyectiva.

d. El [teorema de la dimensión para núcleo e imagen](#th-dim-nucleo) dice que $\mathrm{dim}(V) = \mathrm{dim}(\mathrm{Ker}\ f) + \mathrm{dim}(\mathrm{Im}\ f)$. Lo podemos comprobar fácilmente:
$$\begin{array}{ccccc}
\mathrm{dim}(\mathbb{R}^{`r n`}) & = & \mathrm{dim}(\mathrm{Ker}\ f) & + & \mathrm{dim}(\mathrm{Im}\ f) \\
`r n` & = & `r ncol(genKer)` & + & `r ncol(bImf)` \\
\end{array}$$

e. Dado el subespacio $U$, debemos encontrar un sistema de vectores $\mathcal{G}$ que [lo genere](#gen), calcular su imagen $f(\mathcal{G})$, y eso será un sistema generador de $f(U)$. A partir de ahí, podemos [extraer una base](#sisgenabase) de $f(U)$.

Para calcular un sistema generador de $U$, pasamos sus [ecuaciones cartesianas](#base) a paramétricas, usando métodos Gaussianos:
```{r}
U <- solve_homogeneous(AU)
glue_latex(
  "\\left\\{[write_system(AU, zero_vector(AU), fractions = TRUE, latex = TRUE)]\\right.",
  "\\Rightarrow",
  "[to_latex(gen_vector_V)] = [write_linear_combination(U, vars = params)]"
) %>% 
  cat()
```

De aquí podemos determinar una base de $U$, si eliminamos sus vectores linealmente dependientes:
```{r}
bU <- linearly_independents(U)
glue_latex(
  "\\mathcal{B}_U = \\left\\{",
  "[vectors_to_latex(bU)]",
  "\\right\\}"
) %>% 
  cat()
```

Ahora procedemos a calcular $f(\mathcal{B}_U)$:
```{r}
fbU <- A %*% bU
glue_latex(
  "f(\\mathcal{B}_U) = \\left\\{",
  "[vectors_to_latex(fbU)]",
  "\\right\\}"
) %>% 
  cat()
```
que es un sistema generador de $f(U)$. Eliminamos los vectores linealmente dependientes (en este caso es sencillo) para obtener una base de $f(U)$:
```{r}
glue_latex(
  "\\mathcal{B}_{f(U)} =",
  "\\left\\{[vectors_to_latex(fbU)]\\right\\}"
) %>% 
  cat()
```
y así \(\mathrm{dim}(f(U)) = `r ncol(fbU)`\).


## Diagonalización {#prob-diag}

```{r}
# Dimensión del problema
n <- dim <- 3

# Autovalores para crear una matriz
lambdas <- sample(setdiff(-1:2, 0), 
                  size = dim - 1, 
                  replace = TRUE) %>% 
  sort()
lambdas <- c(lambdas, lambdas[dim - 1])
D <- diag(lambdas)

# Matriz de paso aleatoria
# Vamos a exigir que si los 3 autovalores son diferentes,
# entonces P sea regular.
# Si el determinante es no nulo, será diagonalizable
repeat {
  
  P <- rMatrix(n = dim, values = -2:3)
  
  if (det(P) != 0) break
  
}

# Calculamos A = inv(P)AP
s <- gauss_elimination(P, D %*% P, 
                       diag1 = TRUE, 
                       jordan = TRUE)
A <- s$splits[[2]] %>% remove_fraction(byrow = FALSE)

# Polinomio característico
p <- charpoly(A)
p[abs(p) < 1.e-7] <- 0
# Raíces del polinomio característico
# Autovalores con multiplicidad
lambdas2 <- polyroots(p)
p_latex <- poly2latex(p, var = "\\lambda")

multiplicidad_algebraica <- lambdas2$mult
autovalores <- round(lambdas2$root)

o <- order(abs(autovalores))
autovalores <- autovalores[o]
multiplicidad_algebraica <- multiplicidad_algebraica[o]

D <- diag(rep(autovalores, times = multiplicidad_algebraica))

# Por cada autovalor, sacamos una base del subsepacio asociado
gen <- list()
for (L in autovalores) {

  B <- A - L * eye(dim)
  B[abs(B) < 1.e-7] <- 0
  gen <- append(gen, list(solve_homogeneous(B)))
  
}

b_eigenvectors <- do.call(cbind, gen)
# foo <- lapply(gen, vectors_to_latex) %>% lapply(cat)
# glue_matrices(gen) %>% cat()

# Miramos las multiplicidades geométricas
multiplicidad_geometrica <- sapply(gen, ncol)

# Si coinciden las multiplicidades algebraica y geométrica
# es diagonalizable
diagonalizable <- all(multiplicidad_algebraica == multiplicidad_geometrica)

if (diagonalizable) {
  
  P2 <- do.call(cbind, gen)
  

}
```

__Ejercicio__

Consideremos el siguiente endomorfismo \(f:\mathbb{R}^{`r n`}\to\mathbb{R}^{`r n`}\), dado por 
```{r}
glue_latex(
  "f[to_latex(generic_vector(n))] = ",
  "\\left([write_system(A, latex = TRUE, format = 'c')]\\right)"
) %>% 
  cat()
```

a. Hallar todos los autovalores de $f$.
b. Hallar los subespacios de autovectores de $f$.
c. Comprobar si $f$ es diagonalizable y dar una base $\mathcal{B}$ de \(\mathbb{R}^{`r n`}\) tal que la matriz asociada a $f$ en dicha base sea diagonal.
d. Expresar y calcular, como suma de potencias de grado menor que \(`r n`\), tanto $A^{-1}$ como $A^4$, siendo $A$ la matriz asociada a $f$ en la base canónica.

__Solución__

Como paso previo a todos los cálculos, debemos determinar [la matriz asociada](#matriz-apli) a $f$ en las bases canónicas:
```{r}
glue_latex(
  "A = [glue_matrices(A, latex = TRUE, fractions = TRUE)]"
) %>% 
  cat()
```

a. Para [calcular los autovalores](#eigenvalue) de $f$, o equivalentemente, los de $A$, debemos construir el _polinomio característico_:

```{r}
A_lambdaI <- to_fraction(A, latex = TRUE)
for (i in seq(nrow(A))) {
  
  if (A[i, i] != 0) {
    
    A_lambdaI[i, i] <- paste0(A_lambdaI[i, i], "-\\lambda")
    
  } else {
    
    A_lambdaI[i, i] <- "-\\lambda"
    
  }
  
}
glue_latex(
  "\\begin{array}{rcl}",
  "p(\\lambda) & = & \\mathrm{det}(A - \\lambda\\ I) = \\\\",
  " & = & \\mathrm{det}\\left(",
  "[glue_matrices(A, latex = TRUE, fractions = TRUE)] - \\lambda [glue_matrices(eye(n), latex = TRUE)]",
  "\\right) = \\\\",
  " & = & \\mathrm{det}[glue_matrices(A_lambdaI, latex = TRUE)] = \\\\",
  " & = & [p_latex]",
  "\\end{array}"
) %>% 
  cat()
```

Factorizamos el polinomio para hallar sus raíces (los autovalores) más fácilmente:
```{r}
terms <- -autovalores %>% to_fraction(latex = TRUE)
terms[autovalores < 0] <- paste0("+", terms[autovalores < 0])
terms <- paste0("(\\lambda ", terms, ")")
powers <- paste0("^{", multiplicidad_algebraica, "}")
powers[multiplicidad_algebraica == 1] <- "" 
factorization <- glue::glue("{terms}{powers}") %>% 
  stringr::str_flatten("\\cdot ")
glue_latex(
  "p(\\lambda) = [factorization]"
) %>% 
  cat()
```

Por tanto, los autovalores de $f$ son \(\{`r str_flatten(autovalores, ", ")`\}\).

b. Para calcular el [subespacio asociado a un autovector](#eigenspace) $\lambda$, debemos hallar la solución del sistema homogéneo $(A-\lambda\ I)x = 0$:

```{r}
for (av_idx in seq_along(autovalores)) {
  
  av <- autovalores[av_idx]
  
  cat("- Para el autovalor \\(\\lambda = ", av, "\\):\n")
  
  B <- A - av * eye(n)
  B[abs(B) < 1.e-7] <- 0
  
  av_str <- -av %>% to_fraction(latex = TRUE)
  if (av < 0) {
    
    av_str <- paste0("+ ", av_str)
    
  }
  
  glue_latex(
    "\\begin{array}{rcl}",
    "(A [av_str] I)[to_latex(generic_vector(n))] = 0 &",
    "\\Leftrightarrow & ",
    "[glue_matrices(B, latex = TRUE, fractions = TRUE)][to_latex(generic_vector(n))] = 0 \\\\",
    "& \\Leftrightarrow & ",
    "\\left\\{[write_system(B, zero_vector(n), latex = TRUE, fractions = TRUE)]\\right.",
    "\\end{array}"
  ) %>% 
    cat()
  
  cat("\n\n")
  
  cat("Resolvemos este sistema por Gauss-Jordan, encontrando la forma paramétrica de su conjunto solución, que es \\(U_{", av, "}\\):\n")
  
  params <- c("\\alpha", "\\beta", "\\gamma", "\\delta")
  s <- gauss_elimination(B, zero_vector(n), jordan = TRUE)
  
  glue_latex(
    "[glue_matrices(B, zero_vector(n), latex = TRUE, fractions = TRUE)]",
    "\\sim",
    "[glue_matrices(s$splits, latex = TRUE, fractions = TRUE)]",
    "\\Rightarrow",
    "[to_latex(generic_vector(n))] = [write_linear_combination(gen[[av_idx]], vars = params)]"
  ) %>% 
    cat()
  
  cat("\n\n")
  
  cat("De aquí que una base de \\(U_{", av, "}\\) sea:\n")
  
  glue_latex(
    "\\mathcal{B}_{U_{[av]}} = ",
    "\\left\\{[vectors_to_latex(gen[[av_idx]])]\\right\\}"
  ) %>% 
    cat()
  
  cat("\n")
  
}
```

c. Estudiamos si $f$ es [diagonalizable](#diagonalizable). Para ello, tenemos en cuenta la multiplicidad algebraica y la geométrica de cada autovalor, como presentamos en la siguiente tabla:

```{r}
# L <- as.matrix(lambdas2[o, ])
L <- matrix(c(autovalores, multiplicidad_algebraica, multiplicidad_geometrica), ncol = 3)

multi_table2 <- rbind(c("\\text{Autovalor}", "\\text{Mult. Algebraica}", "\\text{Mult. Geométrica}"), L)
```

$$ `r glue_matrices(multi_table2, 
                  latex = TRUE, 
                  ldeco = "", rdeco = "")` $$
                  
Con esto, queda demostrado que el endomorfismo $f$ __`r ifelse(diagonalizable, "sí", "no")`__ es diagonalizable.

```{r}
if (diagonalizable) {
  
  cat("Además, la base $\\mathcal{B}$ de autovectores es:\n")
  
  subspaces <- paste0("\\mathcal{B}_{", autovalores, "}") %>% 
    str_flatten("\\cup\\ ")
  
  glue_latex(
    "\\mathcal{B} = [subspaces] = \\left\\{[vectors_to_latex(b_eigenvectors)]\\right\\}"
  ) %>% cat()
  
  cat("\n\n")
  
  cat("La matriz de cambio de base de $\\mathcal{B}$ a $\\mathcal{C}$ es la que tiene los elementos de $\\mathcal{B}$ por columnas:\n")
  
  glue_latex(
    "P = [glue_matrices(b_eigenvectors, latex = TRUE, fractions = TRUE)]"
  ) %>% 
    cat()
  
  cat("\n\n")
  
  cat("Y la matriz diagonal, asociada a $f$ en la nueva base, es la que tiene los autovalores en la diagonal, que se calcula como:\n")
  
  glue_latex(
    "\\begin{array}{rcl}",
    "D & = & P^{-1}\\ A\\ P = \\\\",
    " & = & [glue_matrices(b_eigenvectors, latex = TRUE, fractions = TRUE)]^{-1}[glue_matrices(A, latex = TRUE, fractions = TRUE)][glue_matrices(b_eigenvectors, latex = TRUE, fractions = TRUE)]\\\\",
    " & = & [glue_matrices(D, latex = TRUE, fractions = TRUE)]",
    "\\end{array}"
  ) %>% 
    cat()
  
}
```

d. Ahora aplicaremos el [Teorema de Cayley-Hamilton](#th-cayley) para poder expresar la inversas de $A$ y $A^4$ como suma de potencias de grado menor que \(`r n`\).

El teorema me dice que la matriz $A$ anula su polinomio característico, es decir:
```{r}
pA <- poly2latex(p, var = "A", is_matrix = TRUE)
glue_latex(
  "p(A) = [pA] = 0"
) %>% 
  cat()
```

Despejando la identidad, nos queda:
```{r}
p_inv <- -p[-length(p)] / (p[length(p)])
p_dec <- p
p_dec[length(p)] <- 0
p_indep <- p[length(p)]

glue_latex(
  "\\begin{array}{rcl}",
  "\\displaystyle I & = & [to_fraction(1 / p_indep, latex = TRUE)] \\left([poly2latex(p_dec, var = 'A', is_matrix = TRUE)]\\right) = \\\\",
  " & = & \\left([poly2latex(p_inv, var = 'A', is_matrix = TRUE)]\\right)\\ A",
  "\\end{array}"
) %>% 
  cat()
```

Luego tenemos
```{r}

matrices <- list(eye(n))
    
for (i in seq(n - 1)) {
  tmp <- A %*% matrices[[1]]
  tmp[abs(tmp) < 1.e-7] <- 0
  
  matrices <- append(list(tmp), 
                     matrices)
  
}

inversa <- reduce2(.x = matrices, 
                   .y = p_inv, 
                   .f = function(ac, x, y) 
                     ac + x * y, 
                   .init = 0)

p_inv_str <- to_fraction(p_inv, latex = TRUE)
p_inv_str[p_inv > 0] <- paste0("+", p_inv_str[p_inv > 0])
str <- c()
for (i in seq(n)) {
  
  str <- c(str,
           glue::glue(
             "{p_inv_str[i]}{glue_matrices(matrices[[i]], latex = TRUE, fractions = TRUE)}"
           )
  )
  
}

str <- str_flatten(str)

glue_latex(
  "\\begin{array}{rcl}",
  "A^{-1} & = & [poly2latex(p_inv, var = 'A', is_matrix = TRUE)] =\\\\",
  " & = & [str] =\\\\",
  " & = & [glue_matrices(inversa, latex = TRUE, fractions = TRUE)]",
  "\\end{array}"
) %>% cat()
```

Calculemos ahora $A^4$. Tenemos que partir de despejar el término \(A^{`r n`}\) de la igualdad que nos proporciona el Teorema de Cayley-Hamilton:
```{r}
p_power <- -p[-1] / p[1]
glue_latex(
  "A^{[n]} = [poly2latex(p_power, var = 'A', is_matrix = TRUE)]"
) %>% 
  cat()
```

Para $A^4$, tenemos:
```{r}
p_4 <- c(p_power, 0)
glue_latex(
  "A^4 = A\\cdot A^3 = A \\left([poly2latex(p_power, var = 'A', is_matrix = TRUE)]\\right) = [poly2latex(p_4, var = 'A', is_matrix = TRUE)]"
) %>% 
  cat()
```

Si queremos dejarlo como suma de potencias de $A$ de grado menor que \(`r n`\), debemos sustituir aquí \(A^{`r n`}\) por la expresión que hemos despejado hace un momento.

```{r}
p_4_part <- p_4[-1]

p_4_red <- p_4[1] * p_power + p_4_part

p_4_red_str <- to_fraction(p_4_red, latex = TRUE)
pos <- setdiff(which(p_4_red > 0), 1)
p_4_red_str[pos] <- paste0("+", p_4_red_str[pos])
str2 <- c()
for (i in seq(n)) {
  
  if (p_4_red[i] != 0) {
    
    str2 <- c(str2,
              glue::glue(
                "{p_4_red_str[i]}{glue_matrices(matrices[[i]], latex = TRUE, fractions = TRUE)}"
              )
    )
    
  }
  
}

str2 <- str_flatten(str2)

potencia <- reduce2(.x = matrices, 
                   .y = p_4_red, 
                   .f = function(ac, x, y) 
                     ac + x * y, 
                   .init = 0)


str <- c(
  glue::glue(
    "{to_fraction(p_4[1], latex = TRUE)} \\left( {poly2latex(p_power, var = 'A', is_matrix = TRUE)} \\right)"
  ),
  ifelse(p_4[2] > 0, "+", ""),
  glue::glue(
    "{poly2latex(p_4_part, var = 'A', is_matrix = TRUE)}"
  )
) %>% str_flatten()

glue_latex(
  "\\begin{array}{rcl}",
  "A^4 & = & [str] = \\\\",
  " & = & [poly2latex(p_4_red, var = 'A', is_matrix = TRUE)] = \\\\",
  " & = & [str2] = \\\\",
  " & = & [glue_matrices(potencia, latex = TRUE, fractions = TRUE)]",  
  "\\end{array}"
  ) %>% 
  cat()

cat("\n")
```

## Espacios Euclídeos {#prob-espeuclideo}

__Ejercicio__

```{r}
set.seed(4321)
n <- 4
u <- rVector(n)
v <- rVector(n)

repeat {
  
  bU <- rMatrix(n = n, m = 2) %>% 
    linearly_independents()
  
  if (ncol(bU) == 2) break
  
}


repeat {
  
  n <- 4
  lambdas <- sample(-1:0, size = n, replace = TRUE) %>% 
    sort()
  D <- diag(lambdas)
  repeat {
    
    B <- rMatrix(n = n, m = n * 4, values = 0:1) %>% 
      linearly_independents()
    
    if (ncol(B) == n) break
    
  }
  c(B, s) %<-% gram_schmidt(B)
  
  B <- B %>% t() %>% 
    remove_fraction() %>% t()
  
  A <- t(B) %*% D %*% B
  A[abs(A) < 1.e-7] <- 0
  
  p <- charpoly(A)
  p[abs(p) < 1.e-7] <- 0
  
  L <- polyroots(p)
  
  roots <- L$root
  
  if (max(roots - floor(roots)) == 0) break
  
}
p_latex <- poly2latex(p, var = "\\lambda")
```


a. Consideremos los vectores 
```{r}
glue_latex(
  "u = [to_latex(u, fractions = TRUE)]\\quad\\quad v= [to_latex(v, fractions = TRUE)]"
) %>% 
  cat()
```
Hallar el ángulo que forman $u$ y $v$, la norma de cada uno de ellos, $d(u,v)$, y normalizarlos.

b. Encontrar el complemento ortogonal del subespacio generado por el siguiente sistema de vectores:
```{r}
glue_latex(
  "\\mathcal{B} = \\left\\{",
  "[vectors_to_latex(bU)]",
  "\\right\\}"
) %>% 
  cat()
```
c. Diagonalizar ortogonalmente la matriz
```{r}
glue_latex(
  "A = [glue_matrices(A, latex = TRUE, fractions = TRUE)]"
) %>% 
  cat()
```
hallando una base ortonormal de autovectores y expresando de forma explícita la matriz diagonal en función de la matriz $A$ y de la matriz de cambio de base.

__Solución__

a. Vamos a calcular tanto el [producto escalar](#espeuc) de $u$ y $v$ como sus [normas](#norm):
```{r}
prod_str <- "\\cdot"
u_str <- to_fraction(u, latex = TRUE)
u_str[u < 0] <- paste0("(", u_str[u < 0], ")")
v_str <- to_fraction(v, latex = TRUE)
v_str[v < 0] <- paste0("(", v_str[v < 0], ")")
dot_prod <- glue::glue("{u_str}{prod_str}{v_str}") %>% 
  str_flatten("+")
glue_latex(
  "\\langle u, v \\rangle = \\langle [to_latex(u)], [to_latex(v)]\\rangle = [dot_prod] = [sum(u * v)]" 
) %>% 
  cat()

u2_str <- glue::glue("{u_str}^2") %>% str_flatten("+")
u_vector <- Vector$new(u)
norm_u <- u_vector$norm()
glue_latex(
  "\\|u\\| = \\sqrt{\\langle u, u \\rangle} = \\sqrt{[u2_str]} = [norm_u$to_latex()]"
) %>% 
  cat()

v2_str <- glue::glue("{v_str}^2") %>% str_flatten("+")
v_vector <- Vector$new(v)
norm_v <- v_vector$norm()
glue_latex(
  "\\|v\\| = \\sqrt{\\langle v, v \\rangle} = \\sqrt{[v2_str]} = [norm_v$to_latex()]"
) %>% 
  cat()

```

Entonces
```{r}
uv <- Number$new(sum(u * v))
norm_u_inv <- Number$new(norm_u)
norm_u_inv$inverse()
norm_u_inv$simplify()
norm_u_inv$rationalize()
norm_v_inv <- Number$new(norm_v)
norm_v_inv$inverse()
norm_v_inv$simplify()
norm_v_inv$rationalize()
norm_prod <- Number$new(norm_u_inv)$prod(norm_v_inv)
cos_theta <- uv$prod(norm_prod)
cos_theta$simplify()
cos_theta$rationalize()
glue_latex(
  "\\mathrm{ang}(u, v) = \\arccos\\frac{\\langle u, v \\rangle}{\\|u\\| \\|v\\|} = \\arccos\\frac{[sum(u * v)]}{[norm_u$to_latex()][norm_v$to_latex()]} = \\arccos [cos_theta$to_latex()]"
) %>% 
  cat()
```

Pasamos a normalizar cada vector:
```{r}
u_latex <- u_vector$to_latex()
u_norm <- u_vector$prod(norm_u_inv)
v_latex <- v_vector$to_latex()
v_norm <- v_vector$prod(norm_v_inv)
glue_latex(
  "\\frac{u}{\\|u\\|} = [norm_u_inv$to_latex()][u_latex] = [u_norm$to_latex()],\\quad\\quad \\frac{v}{\\|v\\|} = [norm_v_inv$to_latex()][v_latex] = [v_norm$to_latex()]"
) %>% 
  cat()
```

Por último, calculamos la [distancia](#norm):
```{r}
u_v <- u - v
u_v_vector <- Vector$new(u_v)
glue_latex(
  "d(u,v) = \\|u-v\\| = \\left\\|[u_latex] - [v_latex]\\right\\| = \\left\\|[u_v_vector$to_latex()]\\right\\| = [u_v_vector$norm()$to_latex()]"
) %>% 
  cat()
```

b. Para calcular el [complemento ortogonal](#compl-ortho) del subespacio $U$, partimos de la base $\mathcal{B}$, y buscamos aquellos vectores \(v = `r to_latex(generic_vector(n))`\) de $V$ cuyo producto escalar por los vectores de la base sea 0.

Concretamente, usando el producto escalar usual en \(\mathbb{R}^{`r n`}\), esto equivale a buscar las soluciones del sistema de ecuaciones lineales homogéneo: 
```{r}
glue_latex(
  "[write_system(t(bU), zero_vector(t(bU)), fractions = TRUE, latex = TRUE)]"
) %>% 
  cat()
```
donde los coeficientes de cada fila son las coordenadas del vector correspondiente de $\mathcal{B}$. De esta forma, este sistema representa las ecuaciones cartesianas del espacio $U^{\perp}$.

Podemos resolver el sistema usando Gauss-Jordan, y despejando para encontrar las ecuaciones paramétricas del conjunto solución:
```{r}
s <- gauss_elimination(t(bU), diag1 = FALSE, jordan = TRUE)
params <- c("\\alpha", "\\beta", "\\delta")
bUt <- solve_homogeneous(t(bU))
glue_latex(
  "[glue_matrices(t(bU), latex = TRUE, fractions = TRUE)]",
  "\\sim",
  "[glue_matrices(s$splits, latex = TRUE, fractions = TRUE)]",
  "\\Rightarrow",
  "[to_latex(generic_vector(n))] = [write_linear_combination(bUt, vars = params)]"
) %>% 
  cat()
```
de donde podemos deducir una [base](#base) para $U^{\perp}$:
```{r}
bUt <- bUt %>% linearly_independents()
glue_latex(
  "\\mathcal{B}' = \\left\\{",
  "[vectors_to_latex(bUt)]",
  "\\right\\}"
) %>% 
  cat()
```

c. Comenzamos por hallar el [polinomio característico](#eigenvalue) de $A$ y calculando sus [autovalores](#eigenvalue).

```{r}
p_latex <- poly2latex(p, var = "\\lambda")

multiplicidad_algebraica <- L$mult
autovalores <- round(L$root)

o <- order(abs(autovalores))
autovalores <- autovalores[o]
multiplicidad_algebraica <- multiplicidad_algebraica[o]

D <- diag(rep(autovalores, times = multiplicidad_algebraica))

# Por cada autovalor, sacamos una base del subsepacio asociado
gen <- list()
basis <- list()
for (av in autovalores) {
  
  B <- A - av * eye(n)
  B[abs(B) < 1.e-7] <- 0
  H <- solve_homogeneous(B) %>% 
    linearly_independents()
  
  # cat("H\n")
  # glue_latex(
  #   "[glue_matrices(H, latex = TRUE, fractions = TRUE)]"
  # ) %>% 
  #   cat()
  # cat("\n")

  gen <- append(gen, list(H))
  
  # Ortogonalizamos con GS
  Bi <- lapply(seq(ncol(H)),
               function(i) Vector$new(H[, i]))
  
  basis <- append(basis, GS(Bi))
  
}

b_eigenvectors <- do.call(cbind, gen)

normas <- lapply(basis, function(v) v$norm()) %>% 
  lapply(function(s) s$simplify()) %>% 
  lapply(function(s) s$rationalize())

normas_str <- normas %>% 
  sapply(function(s) s$to_latex())

normas_str <- glue::glue("\\frac{1}{[normas_str]}", 
                         .open = "[", .close = "]")  

normas <- normas %>% 
  lapply(function(s) s$inverse())

base_str <- basis %>% 
  sapply(function(v) v$to_latex())

Bpp <- glue::glue("{normas_str}{base_str}") %>% 
  str_flatten(", ")

base2 <- lapply(seq_along(basis), function(i) basis[[i]]$prod(normas[[i]])) %>% 
  lapply(function(s) s$simplify()) %>% 
  lapply(function(s) s$rationalize())

base2_str <- base2 %>% 
  sapply(function(v) v$to_latex()) %>% 
  str_flatten(", ")


A_lambdaI <- to_fraction(A, latex = TRUE)
for (i in seq(nrow(A))) {
  
  if (A[i, i] != 0) {
    
    A_lambdaI[i, i] <- paste0(A_lambdaI[i, i], "-\\lambda")
    
  } else {
    
    A_lambdaI[i, i] <- "-\\lambda"
    
  }
  
}
glue_latex(
  "\\begin{array}{rcl}",
  "p(\\lambda) & = & \\mathrm{det}(A - \\lambda\\ I) = \\\\",
  " & = & \\mathrm{det}\\left(",
  "[glue_matrices(A, latex = TRUE, fractions = TRUE)] - \\lambda [glue_matrices(eye(n), latex = TRUE)]",
  "\\right) = \\\\",
  " & = & \\mathrm{det}[glue_matrices(A_lambdaI, latex = TRUE)] \\\\",
  " & = & [p_latex]",
  "\\end{array}"
) %>% 
  cat()
```

Igualándolo a 0 y resolviendo la ecuación que queda, llegamos a que los autovalores son: \(\lambda = `r str_flatten(autovalores, ", ")`\), ya que:
```{r}
terms <- -autovalores %>% to_fraction(latex = TRUE)
terms[autovalores < 0] <- paste0("+", terms[autovalores < 0])
terms[autovalores == 0] <- ""
terms <- paste0("\\lambda ", terms)
terms[autovalores != 0] <- paste0("(", terms[autovalores != 0], ")")
powers <- paste0("^{", multiplicidad_algebraica, "}")
powers[multiplicidad_algebraica == 1] <- "" 
factorization <- glue::glue("{terms}{powers}") %>% 
  stringr::str_flatten("\\cdot ")
glue_latex(
  "p(\\lambda) = [factorization]"
) %>% 
  cat()
```

Recorremos ahora cada autovalor, encontrando una base del subespacio propio asociado, y ortonormalizándola usando el [método de Gram-Schmidt](#gram-schmidt).

```{r}
all_bases <- list()
for (av_idx in seq_along(autovalores)) {
  
  av <- autovalores[av_idx]
  
  cat("\n- Para el autovalor \\(\\lambda = ", av, "\\):\n")
  
  B <- A - av * eye(n)
  B[abs(B) < 1.e-7] <- 0
  
  av_str <- -av %>% to_fraction(latex = TRUE)
  if (av < 0) {
    
    av_str <- paste0("+ ", av_str)
    
  }
  if (av == 0) {
    
    av_str <- "- 0"
    
  }
  
  glue_latex(
    "\\begin{array}{rcl}",
    "(A [av_str] I)[to_latex(generic_vector(n))] = 0",
    "& \\Leftrightarrow & ",
    "[glue_matrices(B, latex = TRUE, fractions = TRUE)][to_latex(generic_vector(n))] = 0 \\Leftrightarrow \\\\",
    "& \\Leftrightarrow & ",
    "\\left\\{[write_system(B, zero_vector(n), latex = TRUE, fractions = TRUE)]\\right. \\\\",
    "\\end{array}"
  ) %>% 
    cat()
  
  cat("\n\n")
  
  cat("Resolvemos este sistema por Gauss-Jordan, encontrando la forma paramétrica de \\(U_{", av, "}\\):\n")
  
  params <- c("\\alpha", "\\beta", "\\gamma", "\\delta")
  s <- gauss_elimination(B, zero_vector(n), jordan = TRUE)
  
  glue_latex(
    "[glue_matrices(B, zero_vector(n), latex = TRUE, fractions = TRUE)]",
    "\\sim",
    "[glue_matrices(s$splits, latex = TRUE, fractions = TRUE)]"    ) %>% 
    cat()
  
  glue_latex(
    "\\Rightarrow",
    "[to_latex(generic_vector(n))] = [write_linear_combination(gen[[av_idx]], vars = params)]"
  ) %>% 
    cat()
  
  cat("\n\n")
  
  cat("Luego una base de \\(U_{", av, "}\\) es:\n")
  
  glue_latex(
    "\\mathcal{B}_{U_{[av]}} = ",
    "\\left\\{[vectors_to_latex(gen[[av_idx]])]\\right\\}"
  ) %>% 
    cat()
  
  cat("\n")
  
  if (ncol(gen[[av_idx]]) > 1) {
    
      cat("El siguiente paso es _ortonormalizar_ esta base, mediante el método de Gram-Schmidt. Llamemos $v_i$ a los vectores de esta base hallada.\n")

    B <- gen[[av_idx]]
    
    cat("Llamamos\n")
    
    v1 <- matrix(B[, 1], ncol = 1)
    Bp <- v1
    glue_latex(
      "u_1 = v_1 = [to_latex(v1, fractions = TRUE)]"
    ) %>% 
      cat()

    cat("\n")
    
    cat("Procedemos con el resto de vectores:\n")
    cat("\n")
    for (i in seq(2, ncol(B))) {
      
      Ui <- Bp
      vi <- matrix(B[, i], ncol = 1)
      
      uivi <- (matrix(rep(vi, i - 1), ncol = i - 1) * Ui) %>% 
        colSums()
      uiui <- colSums(Ui * Ui)
      
      ii <- seq(i - 1)
      
      proy_form <- glue::glue("\\frac{\\langle v_[i], u_[ii] \\rangle}{\\langle u_[ii], u_[ii] \\rangle} u_[ii]", .open = "[", .close = "]") %>% 
        str_flatten("+")
      
      ui_str <- sapply(seq(i - 1), function(j) glue_matrices(matrix(Ui[, j], ncol = 1), latex = TRUE, fractions = TRUE))
      
      proy_val <- glue::glue("\\frac{[to_fraction(uivi, latex = TRUE)]}{[to_fraction(uiui, latex = TRUE)]} [ui_str]", .open = "[", .close = "]") %>% 
        str_flatten("+")
      
      coefs <- matrix(rep(uivi/uiui, nrow(Ui)), nrow = nrow(Ui), byrow = TRUE)
      proyec <- matrix(rowSums(coefs * Ui), ncol = 1)
      ui_new <- vi - proyec
      Bp <- cbind(Bp, ui_new)
      
      previous_ui <- glue::glue("u_{ii}") %>% str_flatten(", ")
      
      cat("Construimos:\n")
      glue_latex(
        "\\mathcal{U}_[i-1] = \\mathcal{L}(\\{[previous_ui]\\}) = ",
        "\\mathcal{L}(\\left\\{[vectors_to_latex(Ui)]\\right\\})"
      ) %>% str_replace_all("\n", " ") %>% 
        cat()
      cat("\n")
      
      cat("Entonces\n")
      
      glue_latex(
        "\\mathrm{proy}_{\\mathcal{U}_[i-1]}(v_[i])  = [proy_form]  =  [proy_val]  = [glue_matrices(proyec, fractions = TRUE, latex = TRUE)]"
      ) %>% 
        str_replace_all(pattern = "\n", replacement = " ") %>% 
        cat()
      
      cat("\n")
      
      cat("Y llamamos\n")
      
      glue_latex(
        "u_[i]  = v_[i] - \\mathrm{proy}_{\\mathcal{U}_[i - 1]}(v_[i]) = ",
        "[to_latex(vi, fractions = TRUE)] - [glue_matrices(proyec, fractions = TRUE, latex = TRUE)] = ",
        " [glue_matrices(ui_new, fractions = TRUE, latex = TRUE)]"
      ) %>% str_replace_all(pattern = "\n", replacement = " ") %>% cat()
      
      cat("\n")
      
    }    
    
    cat("\nY así hemos llegado a un sistema ortogonal formado por los vectores $u_i$.\n")
    
  } else {
    
    cat("\nComo solo tenemos un vector en la base, ya forma un sistema ortogonal.\n")
    
    Bp <- gen[[av_idx]]
    
  }
  
  cat("Debemos normalizar ahora los vectores de este sistema ortogonal que tenemos, multiplicando cada uno por el inverso de su norma, llegando a la siguiente base ortonormal de $U_{", av, "}$.\n")
  
  this_basis <- lapply(seq(ncol(Bp)), 
                       function(i) Vector$new(Bp[, i]))
  
normas <- lapply(this_basis, function(v) v$norm()) %>% 
  lapply(function(s) s$simplify()) %>% 
  lapply(function(s) s$rationalize())

normas_str <- normas %>% 
  sapply(function(s) s$to_latex())

normas_str <- glue::glue("\\frac{1}{[normas_str]}", 
                         .open = "[", .close = "]")  

normas <- normas %>% 
  lapply(function(s) s$inverse())

base_str <- this_basis %>% 
  sapply(function(v) v$to_latex())

Bpp <- glue::glue("{normas_str}{base_str}") %>% 
  str_flatten(", ")

this_base2 <- lapply(seq_along(this_basis), function(i) this_basis[[i]]$prod(normas[[i]])) %>% 
  lapply(function(s) s$simplify()) %>% 
  lapply(function(s) s$rationalize())

all_bases <- append(all_bases, this_base2)

base2_str <- this_base2 %>% 
  sapply(function(v) v$to_latex()) %>% 
  str_flatten(", ")

if (length(this_base2) > 2) {
  
  glue_latex(
    "\\begin{array}{rcl}",
    "\\mathcal{B}'_{U_{[av]}} & = & \\left\\{",
    "[Bpp]\\right\\} = \\\\",
    " & = & \\left\\{[base2_str]\\right\\} \\\\",
    "\\end{array}"
  ) %>% 
    cat()
  
} else {
  
  glue_latex(
    "\\mathcal{B}'_{U_{[av]}}  =  \\left\\{",
    "[Bpp]\\right\\} = ",
    " \\left\\{[base2_str]\\right\\}"
  ) %>% 
    cat()
  
}
  
  cat("\n\n")
  
}
```
Definimos
```{r}
str <- glue::glue("\\mathcal{B}'_{U_{[autovalores]}}", 
           .open = "[", .close = "]") %>% 
  str_flatten(" \\cup ")
base2_str <- all_bases %>% 
  sapply(function(v) v$to_latex()) %>% 
  str_flatten(", ")

glue_latex(
  "\\mathcal{B} = [str] = \\left\\{[base2_str]\\right\\}"
) %>% 
  cat()
```
Es una base ortonormal de \(\mathbb{R}^{`r n`}\), formada únicamente por autovectores de la matriz $A$.

Además, la matriz $P$ del [cambio de base](#cambiobase) de esta base $\mathcal{B}$ a la canónica, formada por los vectores de $\mathcal{B}$ puestos por columnas, es _ortogonal_ y nos proporciona la relación entre $A$ y la matriz diagonal $D$ de los autovalores:
\small
```{r}
P <- vectors2matrix(all_bases)
glue_latex(
  "\\begin{array}{rcl}",
  "D & \\!\\! = \\!\\! & P^{\\text{t}}\\ A\\ P = \\\\",
  " & \\!\\! = \\! & \\!\\!\\! [glue_matrices(t(P), latex = TRUE)]\\!\\!\\! [glue_matrices(A, latex = TRUE, fractions = TRUE)]\\!\\!\\! [glue_matrices(P, latex = TRUE)] \\\\",
  " &\\!\\! = \\! & \\!\\! [glue_matrices(D, latex = TRUE, fractions = TRUE)]",
  "\\end{array}"
) %>% 
  cat()
```
\normalsize