04-EspaciosEuclideos.Rmd

# Espacios Euclídeos {#espeuc}

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

set.seed(1234)
``` 

En este capítulo, damos un paso para acercarnos a los conceptos geométricos clásicos de distancias y ángulos, abstrayéndolos en el contexto de los espacios vectoriales.

Para ello, precisaremos de la noción de _producto escalar_: Una aplicación $g:V\times V\to \mathbb{R}$, donde $V$ es un espacio vectorial sobre $\mathbb{R}$, es un __producto escalar__ si verifica:

- Es simétrica: $g(u, v) = g(v, u)$ para todo $u,v\in V$.
- Es bilineal: $g(\alpha u + \beta v, w) = \alpha g(u, w) + \beta g(v, w)$ para todos los $\alpha,\beta\in\mathbb{R}$, $u,v,w\in V$. Por simetría, se tiene lo mismo para la linealidad en la segunda componente.
- Es definida positiva: $g(v, v) > 0$ para todo $v\in V\setminus\{0\}$, $g(v, v) = 0$ si, y sólo si $v = 0$.

En resumen, se suele definir un producto escalar como una _forma bilineal, simétrica y definida positiva_.

__Notación__: El producto escalar de dos vectores $u$ y $v$ suele recibir diversas notaciones, aunque las más usadas son $u\cdot v$ y $\langle u, v \rangle$.

A un espacio vectorial $V$ donde podemos definir un producto vectorial lo llamamos __espacio euclídeo__, y lo solemos denotar mediante el par $\left(V, \langle\cdot,\cdot\rangle\right)$.

---

__Ejemplo__

En $\mathbb{R}^n$, disponemos de un _producto escalar usual_, de forma que si $u, v\in\mathbb{R}^n$, entonces
$$\langle u,v \rangle = u_1v_1+u_2v_2+\ldots+u_nv_n$$
siendo $u_i, v_i$, $i=1,\ldots,n$ las coordenadas de $u$ y de $v$, respectivamente, en una base dada (la canónica como caso particular). A este producto vectorial es al que nos referiremos cuando usemos la notación $u\cdot v$.

Podemos generalizar este producto escalar, definiendo unos pesos $\alpha_i\in\mathbb{R}^{+}$ para cada componente:
$$\langle u,v \rangle = \alpha_1u_1v_1+\alpha_2u_2v_2+\ldots+\alpha_nu_nv_n$$

Podemos construir productos escalares de la siguiente forma: Tomemos una matriz $A\in\mathcal{M}_n(\mathbb{R})$ simétrica y definida positiva. Entonces la aplicación $\langle\cdot,\cdot\rangle:\mathbb{R}^n\times\mathbb{R}^n\to\mathbb{R}$ definida por
$$\langle u, v \rangle = u^{\text{t}}\ A\ v$$
es un producto escalar. A esta matriz $A$ se le denomina _matriz de Gram_, y todo producto escalar tiene asociada una. Los dos ejemplos anteriores se correspondían con tomar como $A$ la matriz identidad o una matriz diagonal con elementos positivos $\alpha_i\in\mathbb{R}^{+}$ en la diagonal.

---

__¿Qué preguntas vamos a responder en este capítulo?__

- [¿Qué son las normas y distancias?](#norm)
- [¿Qué significa el concepto de ortogonalidad?](#ortho)
- [¿Qué es el complemento ortogonal de un subespacio?](#compl-ortho)
- [¿Cómo calculamos la proyección de un vector sobre un subespacio?](#proyec)
- [¿Qué es y qué utilidad tiene una base ortonormal?](#ortho-basis)
- [¿Cómo construimos una base ortonormal?](#gram-schmidt)
- [¿Qué es una matriz ortogonal y una aplicación ortogonal?](#ortho-matrix)
- [¿Qué es la diagonalización ortogonal?](#diag-ortho)
- [Qué matrices son diagonalizables ortogonalmente](#diag-ortho2)

## Qué son las normas y distancias {#norm}
\sectionmark{Normas y distancias}

```{r}
definition("Una __norma__ vectorial sobre un espacio $V$ (sobre un cuerpo $\\mathbb{K}$) es un operador que sirve para _medir_ la longitud de un vector (es decir, su distancia al origen), y se define formalmente como una aplicación $\\|\\cdot\\|:V\\to\\mathbb{R}$ que verifica:\n\n- Para todo $v\\in V$, $\\|v\\| \\ge 0$, y $\\|v\\| = 0$ si, y solo si, $v = 0$.\n- Para todo $c\\in\\mathbb{K}$, $v\\in V$, se tiene $\\|cv\\| = |c|\\cdot \\|v\\|$.\n- Para todo $u,v\\in V$, $\\|u+v\\|\\le\\|u\\|+\\|v\\|$. Esta propiedad se llama _desigualdad triangular_ (en un triángulo, la suma de la longitud de dos de sus lados --$\\|u\\|$ y $\\|v\\|$-- es siempre menor que la longitud de su otro lado --$\\|u+v\\|$--).", "Norma")
```

Cualquier operador $\|\cdot\|$ que verifique lo anterior es una norma vectorial. Sin embargo, podemos construir normas a partir de un producto escalar:

```{r}
theorem("Si $\\langle\\cdot,\\cdot\\rangle$ es un producto escalar en el espacio $V$, entonces la aplicación $\\|\\cdot\\|:V\\to\\mathbb{R}$ dada por $\\|v\\| = \\sqrt{\\langle v, v\\rangle}$ es una __norma derivada del producto escalar__ sobre $V$.")
```

Al definir la norma vectorial a partir de un producto escalar, además, tenemos la siguiente propiedad que será de gran interés en la [siguiente sección](#ortho):

```{r}
theorem("Si $\\|\\cdot\\|$ es una norma derivada del producto escalar $\\langle\\cdot,\\cdot\\rangle$, entonces
$$|\\langle u, v \\rangle| \\le \\|u\\|\\cdot\\|v\\|$$
para todo $u,v\\in V$.", "Desigualdad de Cauchy-Schwarz")
```

Un concepto implícitamente relacionado con el de norma vectorial es el de _distancia_. 

```{r}
definition("Una __distancia__ (también llamada _métrica_) es un operador $d:V\\times V\\to R$ que verifica las siguientes propiedades:\n\n- $d(u, v) = 0$ si, y solo si, $u = v$.\n- Es simétrica: $d(u, v) = d(v, u)$.\n- Cumple la desigualdad triangular: $d(u, v) \\le d(u, w) + d(w, v)$.", "Distancia")
```

Además, de esta definición se deduce que $d(u, v) \ge 0$ para todo $u,v\in V$.

Podemos deducir una distancia a partir de una norma, para así completar las relaciones existentes entre los conceptos aquí explicados:

```{r}
theorem("Si $\\|\\cdot\\|$ es una norma sobre un espacio vectorial $V$, entonces la aplicación $d:V\\times V\\to \\mathbb{R}$ definida por
$$d(u, v) = \\|u-v\\|$$
para todo $u,v\\in V$ es una _distancia_.")
```

__Nota__: Existen muchas más normas de las que se deducen a partir de un producto escalar, al igual que más _distancias_.

---

__Ejemplo__

Podemos deducir la expresión de la norma (a veces llamada _módulo_) de un vector $v\in\mathbb{R}^n$.

Si $v=\left(\begin{array}{c}v_1\\ v_2\\\vdots\\ v_n\end{array}\right)$, entonces
$$\|v\| = \sqrt{v\cdot v} = \sqrt{v_1^2+v_2^2+\ldots+v_n^2}$$

Y así, además, deducimos la expresión para la distancia euclídea usual en $\mathbb{R}^n$:
$$d(u,v) = \|u - v\| = \sqrt{(u_1-v_1)^2+\ldots+(u_n-v_n)^2}$$

Como ejemplo, podemos calcular la norma del vector
```{r}
v <- rVector(n = 4)
glue_latex(
  "v = [to_latex(v)]"
) %>% 
  cat()
v_str <- as.character(v)
v_str[v < 0] <- paste0("(", v_str[v < 0], ")")
glue_latex(
  "\\|v\\| = \\sqrt{[str_flatten(glue::glue('{v_str}^2'), '+')]} = \\sqrt{[sum(v^2)]}"
) %>% 
  cat()
```


---

## Qué significa el concepto de ortogonalidad {#ortho}
\sectionmark{Ortogonalidad}

Si partimos de la _desigualdad de Cauchy-Schwarz_ que [hemos visto antes](#norm), podemos llegar a que 
$$-1 \le \frac{\langle u, v\rangle}{\|u\|\cdot\|v\|} \le 1$$
para todo $u,v\in V$.

Por otro lado, la función coseno restringida a $[0, \pi]$, es decir, $\cos:[0,\pi]\to[-1,1]$, es biyectiva. Eso quiere decir que existe un único $\theta\in[0, \pi]$ tal que $\cos\theta = \frac{\langle u, v\rangle}{\|u\|\cdot\|v\|}$. 

Esto nos permite la siguiente definición:
```{r}
definition("Llamamos __ángulo entre los vectores $u,v\\in V$__ al único $\\theta\\in[0,\\pi]$ tal que $\\cos\\theta = \\frac{\\langle u, v\\rangle}{\\|u\\|\\cdot\\|v\\|}$.", "Ángulo entre vectores")
```

```{r}
definition("Diremos que dos vectores $u$ y $v$ de $V$ son __ortogonales__ (o __perpendiculares__) si el ángulo que forman es $\\frac{\\pi}{2}$, es decir, si $\\langle u, v \\rangle = 0$.", "Vectores ortogonales")
```

La noción de _ortogonalidad_ en el espacio euclídeo se corresponde con la perpendicularidad de vectores que todos tenemos de forma intuitiva. Lo podemos extender también a ortogonalidad con respecto a un espacio vectorial: Un vector $v\in V$ es __ortogonal al subespacio $U$__ de $V$ si, y sólo si, es ortogonal a todos y cada uno de los vectores de $U$, es decir, $\langle v, u \rangle = 0$ para todo $u\in U$.

El siguiente resultado nos proporciona una ayuda a la hora de determinar la ortogonalidad entre subespacios vectoriales.

```{r}
theorem("Sea $U$ un subespacio de $V$, y sea $\\mathcal{B} = \\{u_1,\\ldots,u_k\\}$ una base de $U$. Un vector $v\\in V$ es ortogonal a $U$ si, y solo si, es ortogonal a cada $u_i$: $\\langle v, u_i \\rangle = 0$ para todo $i=1,\\ldots,k$.")
```

Un __sistema ortogonal__ es un conjunto de vectores $\{v_1,\ldots,v_m\}\subset V$ donde los $v_i$ son ortogonales dos a dos, es decir, $\langle v_i, v_j \rangle = 0$ para todo $i\ne j$.

__Algunas propiedades y caracterizaciones de los vectores y sistemas ortogonales__

```{r}
theorem("Dos vectores $u$ y $w$ son ortogonales si, y solo si, $\\|u+w\\|^2 = \\|u\\|^2 + \\|v\\|^2$", "Generalización del teorema de Pitágoras")
```

```{r}
theorem("Si un sistema de vectores es ortogonal, entonces es linealmente independiente. Por tanto, en un espacio $V$ de dimensión finita $n$, cualquier conjunto de $n$ vectores ortogonales entre sí es una base.")
```

Este último resultado es interesante, pues nos constata que cualquier sistema ortogonal con tantos vectores como la dimensión del espacio vectorial es una base.

---

__Ejemplo__

```{r}
n <- 4
# Construimos un espacio U dando su base
# Nos aseguramos que tenga dimensión 2
repeat {
  
  bU <- rMatrix(n, m = 2)
  bU <- linearly_independents(bU)
  if (ncol(bU) == 2) break
  
}
# Vamos a construir sus ecuaciones cartesianas:
AU <- parametric_to_cartesian(bU)$A %>% remove_fraction()
bU <- solve_homogeneous(AU)

# Vamos a buscar un vector ortogonal a U
# Encontramos la base del complemento ortogonal
bUt <- solve_homogeneous(t(bU)) %>%  
  linearly_independents() %>% t() %>% 
  remove_fraction() %>% t()

# Vector v aleatorio en el ortogonal
scalars <- rVector(n = ncol(bUt))
v <- bUt %*% scalars
w <- rVector(n = n)

# Construyamos ahora una base aleatoria y la 
# ortonormalizamos usando GS
repeat {
  
  B <- rMatrix(n = n, m = n * 4, values = 0:1) %>% 
    linearly_independents()
  
  if (ncol(B) == n) break
  
}
c(B_ortho, s) %<-% gram_schmidt(B)

B_ortho[] <- B_ortho %>% t() %>% 
  remove_fraction() %>% t() %>% 
  round()

```

Consideremos el espacio $U$ dado por:
```{r}
glue_latex(
  "U = \\left\\{",
  "[to_latex(generic_vector(n))]\\in\\mathbb{R}^{[n]}:",
  "[write_system(AU, zero_vector(AU), fractions = TRUE, latex = TRUE, format = 'c')]",
  "\\right\\}"
) %>% 
  cat()
```
y los vectores
```{r}
glue_latex(
  "v = [to_latex(v)]\\quad\\quad w = [to_latex(w)]"
) %>% 
  cat()
```

Calculemos el ángulo entre $v$ y $w$:
```{r}
v_str <- as.character(v)
v_str[v < 0] <- paste0("(", v_str[v < 0], ")")
w_str <- as.character(w)
w_str[w < 0] <- paste0("(", w_str[w < 0], ")")
glue_latex(
  "\\theta = \\arccos \\frac{v\\cdot w}{\\|v\\|\\cdot\\|w\\|}"
) %>% 
  cat()
prod_str <- "\\cdot"
glue_latex(
  "v\\cdot w = [glue::glue('{v_str}{prod_str}{w_str}') %>% str_flatten('+')] = [sum(v * w)]"
) %>% 
  cat()
glue_latex(
  "\\|v\\| = \\sqrt{[glue::glue('{v_str}^2') %>% str_flatten('+')]} = \\sqrt{[sum(v^2)]}"
) %>% 
  cat()
glue_latex(
  "\\|w\\| = \\sqrt{[glue::glue('{w_str}^2') %>% str_flatten('+')]} = \\sqrt{[sum(w^2)]}"
) %>% 
  cat()

vv <- Vector$new(v)
vw <- Vector$new(w)
dot_prod <- vv$prod(vw)
nv <- vv$norm()
nv_orig <- vv$norm()
nw <- vw$norm()
nvnw <- nv$prod(nw)
nvnw$inverse()
nvnw$prod(dot_prod)
nvnw$simplify()
nvnw$rationalize()
```

Luego
```{r}
glue_latex(
  "\\theta = \\arccos \\frac{[dot_prod$to_latex()]}{[nv_orig$to_latex()][nw$to_latex()]} = \\arccos [nvnw$to_latex()]"
) %>% 
  cat()
```

Vamos a comprobar ahora que $v$ es ortogonal al subespacio $U$. Para ello, a partir de las ecuaciones cartesianas de $U$, vamos a [hallar una base](#base) suya, resolviendo el sistema por Gauss-Jordan y pasando por las ecuaciones paramétricas:
```{r}
s <- gauss_elimination(AU, jordan = TRUE, diag1 = TRUE)
params <- c("\\alpha", "\\beta", "\\delta", "\\gamma")
glue_latex(
  "[glue_matrices(AU, latex = TRUE, fractions = TRUE)]",
  "\\sim",
  "[glue_matrices(s$splits, latex = TRUE, fractions = TRUE)]",
  "\\Rightarrow",
  "[to_latex(generic_vector(n))] = [write_linear_combination(bU, vars = params)]"
) %>% 
  cat()
```

De aquí que la base de $U$ sea:
```{r}
glue_latex(
  "\\mathcal{B}_U = \\left\\{u_i\\right\\} = \\left\\{",
  "[vectors_to_latex(bU)]",
  "\\right\\}"
) %>% 
  cat()
```

Como hemos visto antes, para comprobar que $v$ es ortogonal a $U$, nos basta con comprobar que lo es a cada vector de $\mathcal{B}_U$, hallando su producto escalar:
```{r}
for (i in seq(ncol(bU))) {
  
  ui <- matrix(bU[, i], ncol = 1)
  ui_str <- as.character(ui)
  ui_str[ui < 0] <- paste0("(", ui_str[ui < 0], ")")
  glue_latex(
    "\\langle v, u_[i] \\rangle = ",
    "\\langle [to_latex(v)], [to_latex(ui)] \\rangle = ",
    "[glue::glue('{v_str}{prod_str}{ui_str}') %>% str_flatten('+')] = ",
    "[sum(v * ui)]"
  ) %>% 
    cat()
  
}
```

Consideremos ahora el siguiente sistema de vectores:
```{r}
glue_latex(
  "S = \\{s_i\\} = \\left\\{",
  "[vectors_to_latex(B_ortho)]",
  "\\right\\}"
) %>% 
  cat()
```

Podemos comprobar fácilmente que es un sistema ortogonal, hallando los productos escalares de los vectores de $S$ dos a dos:
```{r}
for (i in seq(ncol(B_ortho) - 1)) {
  
  vi <- matrix(B_ortho[, i], ncol = 1)
  vi_str <- as.character(vi)
  vi_str[vi < 0] <- paste0("(", vi_str[vi < 0], ")")
  
  for (j in seq(i + 1, ncol(B_ortho))) {
    
    vj <- matrix(B_ortho[, j], ncol = 1)
    vj_str <- as.character(vj)
    vj_str[vj < 0] <- paste0("(", vj_str[vj < 0], ")")

    glue_latex(
      "\\langle s_[i], s_[j] \\rangle = ",
      "\\langle [to_latex(vi)], [to_latex(vj)] \\rangle = ",
      "[glue::glue('{vi_str}{prod_str}{vj_str}') %>% str_flatten('+')] = ",
      "[sum(vi * vj)]"
    ) %>% 
      cat()
      
  }
  
}
```

Tenemos, por tanto, un sistema de `r ncol(B_ortho)` vectores ortogonales en \(\mathbb{R}^{`r n`}\), y, como hemos comentado antes, son linealmente independientes, luego forman una base del espacio vectorial.



## Qué es el complemento ortogonal de un subespacio {#compl-ortho}
\sectionmark{Complemento Ortogonal}

Consideremos un espacio vectorial $V$ con producto escalar $\langle \cdot,\cdot \rangle$, y sea $U$ un subespacio de $V$.

```{r}
definition("Llamamos __complemento ortogonal de $U$__ al conjunto de vectores que son ortogonales a todos los de $U$:
$$U^{\\perp} = \\{v\\in V: \\langle v, u \\rangle = 0\\,\\,\\text{para todo }u\\in U\\}$$", "Complemento ortogonal")
```


Tenemos las siguientes propiedades del complemento ortogonal de un subespacio:

```{r}
theorem("\n- $U^{\\perp}$ es un subespacio vectorial de $V$.\n- $U\\cap U^{\\perp} = \\{0\\}$.\n- Si $V$ es de dimensión finita, $\\mathrm{dim}(U) + \\mathrm{dim}(U^{\\perp}) = \\mathrm{dim}(V)$.", "Propiedades del complemento ortogonal")
```

Como consecuencia, este resultado nos dice que __$V$ se puede descomponer como la suma directa de $U$ y de su complemento ortogonal__: $U\oplus U^{\perp} = V$.

__¿Cómo calculamos el complemento ortogonal de un subespacio?__

Hemos de hallar los vectores que son ortogonales a un subespacio $U$. Como hemos [visto antes](#ortho), un vector $v$ es ortogonal a $U$ si, y solo si, es ortogonal a todos los vectores de una base de $U$.

Así pues, partimos de que tenemos una base $\mathcal{B} = \{u_1,\ldots,u_m\}$ de $U$ y tomemos un vector genérico $v\in V$. Para que $v\in U^{\perp}$, debe cumplirse que $\langle v, u_i\rangle = 0$ para todo $i=1,\ldots,m$.

Cada una de las restricciones $\langle v, u_i\rangle = 0$ nos proporciona una ecuación lineal que debe verificar un vector para poder pertenecer a $U^{\perp}$.

Tomamos, por tanto, todas las ecuaciones lineales que se deducen de las restricciones anteriores: ese sistema de ecuaciones nos representa las ecuaciones cartesianas del espacio $U^{\perp}$.

En el caso concreto de un espacio $V=\mathbb{R}^n$ con el _producto escalar euclídeo_, podemos ver qué aspecto tienen las ecuaciones del tipo $\langle v, u_i \rangle = 0$. Supongamos que los vectores $v$ y $u_i$ tienen por coordenadas, en la base canónica, las siguientes:
$$v = \left(
\begin{array}{c}
x_1\\ x_2\\ \vdots \\ x_n\\
\end{array}\right), \quad 
u_i = \left(
\begin{array}{c}
\alpha_{i,1}\\ \alpha_{i,2}\\ \vdots \\ \alpha_{i,n}\\
\end{array}\right)$$

Entonces
$$
\langle v, u_i \rangle = 0 \Leftrightarrow 
\langle \left(
\begin{array}{c}
x_1\\ x_2\\ \vdots \\ x_n\\
\end{array}\right) , 
\left(
\begin{array}{c}
\alpha_{i,1}\\ \alpha_{i,2}\\ \vdots \\ \alpha_{i,n}\\
\end{array}\right) \rangle = 0 \Leftrightarrow \alpha_{i,1}x_1+\alpha_{i,2}x_2+\ldots+\alpha_{i,n}x_n = 0
$$

Si nos fijamos, nos queda una ecuación lineal cuyas incógnitas son las coordenadas canónicas de $v$ y cuyos coeficientes son las coordenadas de $u_i$. Si repetimos el proceso para todo $i=1,\ldots,m$, tenemos el siguiente sistema de ecuaciones:
$$Av = 0\quad\text{donde}\quad A = \left(
\begin{array}{cccc}
\alpha_{1,1} & \alpha_{1,2} & \ldots & \alpha_{1,n} \\
\alpha_{2,1} & \alpha_{2,2} & \ldots & \alpha_{2,n} \\
\vdots & \vdots & \ddots & \vdots \\
\alpha_{m,1} & \alpha_{m,2} & \ldots & \alpha_{m,n} \\
\end{array}
\right)$$

Es decir, las _ecuaciones cartesianas_ de $U^{\perp}$ son las que tienen por matriz de coeficientes a las coordenadas de los vectores de la base de $U$ puestos _por filas_. De hecho, si llamamos $B$ a la matriz que resulta de poner por columnas los vectores de $\mathcal{B}$, entonces $A = B^{\mathrm{t}}$.

A partir de esas ecuaciones cartesianas, ya podríamos [encontrar una base](#sisgena base) del subespacio $U^{\perp}$.

---

__Ejemplo__

```{r}
bW <- rMatrix(n = n, m = 2)
bWt <- solve_homogeneous(t(bW))
```

Consideremos el espacio $U$ generado por la base siguiente:
```{r}
glue_latex(
  "\\mathcal{B} = \\{u_i\\} = \\left\\{",
  "[vectors_to_latex(bW)]",
  "\\right\\}"
) %>% 
  cat()
```

Para hallar su complemento ortogonal, suponemos un vector genérico 
```{r}
glue_latex(
  "v = [to_latex(generic_vector(n))]\\in \\mathbb{R}^{[n]}"
) %>% 
  cat()
```
y hacemos su producto escalar por todos los vectores de $\mathcal{B}$, imponiendo que valga 0:
```{r}
for (i in seq(ncol(bW))) {
  
  ui <- matrix(bW[, i])
  
  glue_latex(
    "\\langle v, u_[i] \\rangle = 0",
    "\\Leftrightarrow",
    "\\langle [to_latex(generic_vector(n))], [to_latex(ui)] \\rangle = 0",
    "\\Leftrightarrow",
    "\\ [write_system(t(ui), 0, latex = TRUE, format = 'c')]"
  ) %>% 
    cat()
}
```

Hemos llegado, por tanto, al siguiente sistema de ecuaciones que representan las ecuaciones cartesianas de $U^{\perp}$:
```{r}
glue_latex(
  "\\left\\{[write_system(t(bW), zero_vector(t(bW)), latex = TRUE)]\\right."
) %>% 
  cat()
```
Nos podemos dar cuenta de que la matriz del sistema es la que resulta de poner los vectores de la base $\mathcal{B}$ _por filas_.

Como siempre, podemos resolver el sistema por Gauss-Jordan para [deducir la forma paramétrica](#param2cartesian) de la solución, y, por tanto, un [sistema generador](#sisgen) de $U^{\perp}$:
```{r}
params <- c("\\alpha", "\\beta", "\\gamma", "\\delta")
s <- gauss_elimination(t(bW), jordan = TRUE, diag1 = TRUE)
glue_latex(
  "[glue_matrices(t(bW), latex = TRUE, fractions = TRUE)]",
  "\\sim",
  "[glue_matrices(s$splits, latex = TRUE, fractions = TRUE)]",
  "\\Rightarrow",
  "[to_latex(generic_vector(n))] = [write_linear_combination(bWt, vars = params)]"
) %>% 
  cat()
```

De aquí que una base del complemento ortogonal $U^{\perp}$ sea:
```{r}
glue_latex(
  "\\mathcal{B}' = \\left\\{",
  "[vectors_to_latex(bWt)]",
  "\\right\\}"
) %>% 
  cat()
```

Podemos comprobar fácilmente que los vectores de $\mathcal{B}'$ son ortogonales a los de $\mathcal{B}$.

---

## Cómo calculamos la proyección de un vector sobre un subespacio {#proyec}
\sectionmark{Proyecciones}

```{r}
definition("Consideremos un espacio euclídeo $\\left(V, \\langle\\cdot,\\cdot\\rangle\\right)$. Consideremos un subespacio $U$, del que conocemos una [base _ortogonal_](#ortho) $\\mathcal{B}_U = \\{u_1,\\ldots,u_m\\}$.  Definimos la __proyección _ortogonal_ de un vector $v$ sobre el subespacio $U$__, como:\n$$\\mathrm{proy}_U(v) = \\frac{\\langle v, u_1\\rangle}{\\langle u_1, u_1 \\rangle}u_1 + \\frac{\\langle v, u_2\\rangle}{\\langle u_2, u_2 \\rangle}u_2 + \\ldots + \\frac{\\langle v, u_m\\rangle}{\\langle u_m, u_m \\rangle}u_m$$", "Proyección ortogonal")
```

Es decir, la proyección $\mathrm{proy}_U(v)$ es un vector del subespacio $U$, ya que es combinación lineal de los elementos de su base, y además esta combinación lineal tiene como coeficientes el producto escalar de $v$ por cada elemento de la base, dividido por la norma de $u_i$ al cuadrado ([recordemos que $\|x\| = \sqrt{\langle x,x \rangle}$](#norm)).

Evidentemente, si los elementos de la base de $U$ tienen norma 1, el denominador desaparece.

La proyección de $v$ sobre el subespacio $U$ tiene la propiedad siguiente:

```{r}
theorem("Si $u$ es un vector del subespacio $U$, y $v\\in V$, entonces $d(u, v) \\ge d(\\mathrm{proy}_U(v), v)$, y la única forma para que se dé la igualdad es que $u = \\mathrm{proy}_U(v)$.")
```

Esto significa que la proyección ortogonal de $v$ sobre $U$ es el vector de $U$ más cercano (según la distancia $d$) a $v$.

Como caso particular, podríamos definir la __proyección ortogonal de un vector $v$ sobre otro vector $u$__ como
$$\mathrm{proy}_u(v) = \frac{\langle v, u \rangle}{\langle u, u\rangle}u$$
(realmente es la proyección del vector $v$ sobre el [subespacio generado](#sisgen) por el vector $u$).

Podemos entonces concluir que, si la base ortogonal de $U$ es $\mathcal{B}_U = \{u_1,\ldots,u_m\}$, entonces
$$\mathrm{proy}_U(v) = \mathrm{proy}_{u_1}(v)+\ldots+\mathrm{proy}_{u_m}(v)$$

A partir del hecho de que la base de $U$ es _ortogonal_, podemos deducir que:

```{r}
theorem("Si $U$ es un subespacio de $V$, $v\\in V$, y una base ortogonal de $U$ es $\\mathcal{B}_U = \\{u_1,\\ldots,u_m\\}$, entonces el vector $v - \\mathrm{proy}_U(v)$ es ortogonal a $U$, es decir, $v - \\mathrm{proy}_U(v) \\in U^{\\perp}$ está en el [complemento ortogonal](#compl-ortho) de $U$.")
```

Esto, junto con que $\mathrm{proy}_U(v)\in U$, nos justifican el siguiente resultado, llamado __teorema de la descomposición ortogonal__:

```{r}
theorem("En un espacio euclídeo $(V, \\langle\\cdot,\\cdot\\rangle)$, dado un subespacio $U$, y $v\\in V$, existe una única descomposición de $v$ como suma de dos vectores ortogonales, $v = v_1 + v_2$, tales que $v_1\\in U$ y $v_2\\in U^{\\perp}$. Concretamente, $v_1 = \\mathrm{proy}_U(v)$ y $v_2 = v - v_1 = \\mathrm{proy}_{U^{\\perp}}(v)$.", "Descomposición ortogonal")
```

Podemos siempre descomponer un vector de $V$ como suma de dos vectores, uno de $U$ y otro de $U^{\perp}$, que coinciden con las proyecciones sobre ambos subespacios.

---

__Ejemplo__

```{r}
repeat {

  bU <- rMatrix(n = n, m = 2) %>% 
    linearly_independents() %>% 
    t() %>% remove_fraction() %>% t()

  if (ncol(bU) == 2) break  
  
}
c(bU, s) %<-% gram_schmidt(bU)
bU <- (bU / s) %>% t() %>% 
  remove_fraction() %>% t()
```

Consideremos el subespacio $U$ cuya base es
```{r}
glue_latex(
  "\\mathcal{B} = \\{u_i\\} = \\left\\{",
  "[vectors_to_latex(bU)]",
  "\\right\\}"
) %>% 
  cat()
```

Notemos que $\mathcal{B}$ es una base _ortogonal_. 

Consideremos el vector
```{r}
v <- rVector(n = n)
glue_latex(
  "v = [glue_matrices(v, latex = TRUE, fractions = TRUE)]"
) %>% 
  cat()
```

Buscamos la _descomposición ortogonal_ de $v$ con respecto al subespacio $U$. Sabemos que podemos escribir $v = v_1  + v_2$, donde $v_1 = \mathrm{proy}_U(v)$ y $v_2 = v - v_1 \in U^{\perp}$.
```{r}
ii <- seq(ncol(bU))
str <- glue::glue("\\mathrm{proy}_{u_[ii]}(v)",
                  .open = "[", .close = "]") %>% 
  str_flatten("+")
str2 <- glue::glue("\\frac{\\langle v, u_[ii]\\rangle}{\\langle u_[ii], u_[ii] \\rangle} u_[ii]",
                   .open = "[", .close = "]") %>% 
  str_flatten("+")
glue_latex(
  "v_1 = \\mathrm{proy}_U(v) = [str] = [str2]"
) %>% 
  cat()
```

Si calculamos esas cantidades
```{r}
v_str <- as.character(v)
v_str[v < 0] <- paste0("(", v_str[v < 0], ")")
prod_str <- "\\cdot"
uv <- rep(0, ncol(bU))
uu <- rep(0, ncol(bU))
for (i in seq(ncol(bU))) {
  
  ui <- matrix(bU[, i], ncol = 1)
  ui_str <- as.character(ui)
  ui_str[ui < 0] <- paste0("(", ui_str[ui < 0], ")")
  
  uv[i] <- sum(v * ui)

  glue_latex(
    "\\langle v, u_[i] \\rangle = \\langle [to_latex(v)], ",
    "[to_latex(ui)] \\rangle = ",
    "[glue::glue('{v_str}{prod_str}{ui_str}') %>% str_flatten('+')] = [sum(v * ui)]"
  ) %>% 
    cat()
  
}

for (i in seq(ncol(bU))) {
  
  ui <- matrix(bU[, i], ncol = 1)
  ui_str <- as.character(ui)
  ui_str[ui < 0] <- paste0("(", ui_str[ui < 0], ")")
  
  uu[i] <- sum(ui * ui)

  glue_latex(
    "\\langle u_[i], u_[i] \\rangle = \\langle [to_latex(ui)], ",
    "[to_latex(ui)] \\rangle = ",
    "[glue::glue('{ui_str}^2') %>% str_flatten('+')] = [sum(ui * ui)]"
  ) %>% 
    cat()
}
```
podemos sustituirlas en la expresión de más arriba y llegamos a
```{r}
bU_cols <- lapply(ii, function(i) to_latex(matrix(bU[, i], ncol = 1)))
coef <- to_fraction(uv / uu, latex = TRUE)
coef[uv < 0] <- paste0("(", coef[uv < 0], ")")
str <- glue::glue("{coef}{prod_str}{bU_cols}") %>% 
  str_flatten("+")
C <- matrix(rep(uv / uu, nrow(bU)), nrow = nrow(bU), byrow = TRUE)
proy <- matrix(rowSums(C * bU), ncol = 1)
glue_latex(
  "v_1 = [str] = [glue_matrices(proy, latex = TRUE, fractions = TRUE)]"
) %>% 
  cat()
```

Por tanto,
```{r}
glue_latex(
  "v_2 = [glue_matrices(v, latex = TRUE, fractions = TRUE)] - [glue_matrices(proy, latex = TRUE, fractions = TRUE)] = ",
  "[glue_matrices(v - proy, latex = TRUE, fractions = TRUE)]"
) %>% 
  cat()
```

Por lo comentado anteriormente, sabemos que $v_2 = \mathrm{proy}_{U^{\perp}}(v)\in U^{\perp}$, lo cual se podría comprobar haciendo el producto escalar $\langle v_2, u_i \rangle$ y viendo que su valor es 0 para todo $i$.

---

## Qué es y qué utilidad tiene una base ortonormal {#ortho-basis}
\sectionmark{Bases ortonormales}

```{r}
definition("Un vector $v\\in V$ se dice __unitario__ si $\\|v\\| = 1$.", "Vector unitario")
```

Es fácil crear vectores unitarios: dado $u\in V$, si definimos $v = \frac{1}{\|u\|}u$, este vector tiene norma 1. A este proceso se le llama _normalizar_.

Un sistema de vectores $\{v_1,\ldots,v_k\}$ se llama __sistema ortonormal__ si es un [sistema ortogonal](#ortho) y cada $v_i$ es un vector unitario ($i=1,\ldots, k$).

En relación a lo comentado de las propiedades de los sistemas ortogonales, en nuestro caso podemos decir que un sistema ortonormal $\{v_1,\ldots,v_n\}$, donde $\mathrm{dim}(V) = n$, es una __base ortonormal__.

El interés de una base ortonormal es que es sencillo calcular las [coordenadas](#coord) de cualquier vector con respecto a dicha base:

```{r}
theorem("Sea $\\mathcal{B} = \\{v_1,\\ldots,v_n\\}$ una base _ortonormal_ del espacio vectorial $V$. Entonces, para cada $v\\in V$, tenemos:\n$$v = \\langle v, v_1 \\rangle v_1 + \\ldots + \\langle v, v_n \\rangle v_n$$\nluego las coordenadas de $v$ en la base $\\mathcal{B}$ las podemos calcular como:\n$$[v]_{\\mathcal{B}} = \\left(\\begin{array}{c}\\langle v, v_1\\rangle \\\\\\langle v, v_2 \\rangle \\\\\\vdots \\\\\\langle v, v_n \\rangle \\\\\\end{array}\\right)_{\\mathcal{B}}$$", "Coordenadas en base ortonormal")
```


---

__Ejemplo__

Consideremos ahora el sistema de vectores ortogonales del [ejemplo anterior](#ortho):
```{r}
glue_latex(
  "S = \\{s_i\\} = \\left\\{",
  "[vectors_to_latex(B_ortho)]",
  "\\right\\}"
) %>% 
  cat()
```

A partir de $S$, formemos un sistema ortonormal $S'$, donde cada vector de $S'$ es un vector de $S$, que se ha normalizado.
```{r}
base <- lapply(seq(ncol(B_ortho)), function(i) Vector$new(matrix(B_ortho[, i], ncol = 1)))
base_norm <- lapply(base, 
                    function(v) {
                      norma <- v$norm()
                      v$prod(norma$inverse())
                    })

base <- lapply(base, 
                    function(v) v$simplify())

base_norm <- lapply(base_norm, 
                    function(v) v$rationalize())

str_vectors <- lapply(base_norm, 
       function(v) v$to_latex()) %>% 
  str_flatten(", ")
glue_latex(
  "S' = \\left\\{",
  "[str_vectors]",
  "\\right\\}"
) %>% 
  cat()
```

$S'$ forma una base _ortonormal_ de \(\mathbb{R}^{`r n`}\). 
Si consideramos el vector 
```{r}
glue_latex(
  "v = [glue_matrices(v, latex = TRUE, fractions = TRUE)]"
) %>% 
  cat()
```
podemos hallar sus coordenadas en la base $S'$, ya que cada coordenada de $v$ será el producto escalar de $v$ por el correspondiente vector de $S'$:
```{r}

v <- Vector$new(v)
producto <- list()
for (i in seq_along(base)) {
  
  si <- base[[i]]
  p <- v$prod(si)
  p$simplify()
  p$rationalize()
  producto[[i]] <- p
  
  glue_latex(
    "\\langle v, s'_[i] \\rangle = ",
    "\\langle [v$to_latex()], [si$to_latex()] \\rangle = ",
    "[p$to_latex()]"
  ) %>% 
    cat()
}
```

Entonces, podremos escribir:
```{r}
coord <- Vector$new(producto)
glue_latex(
  "v = [coord$to_latex()]_{S'}"
) %>% 
  cat()
```

---


## Cómo construimos una base ortonormal {#gram-schmidt}
\sectionmark{Método de Gram-Schmidt}

En las secciones anteriores hemos usado bases ortogonales y ortonormales, pero cabe preguntarnos si la existencia de tales bases está asegurada. La respuesta nos la da el siguiente resultado:

```{r}
theorem("Sea $V$ un espacio vectorial euclídeo de _dimensión finita_. Entonces $V$ tiene una base ortonormal.")
```

La demostración de este resultado teórico realmente nos presenta un método constructivo que permite hallar una base ortonormal del espacio $V$ a partir de una base cualquiera.

__Método de Ortonormalización de Gram-Schmidt__

Consideremos una base $\mathcal{B} = \{v_1,\ldots,v_n\}$ en $V$. Vamos a utilizar el método de Gram-Schmidt para construir una base ortonormal a partir de ella.

El proceso se divide en dos fases:

- En la primera fase, construiremos una base ortogonal $\mathcal{B}'$.
- En un segundo paso, cada vector de $\mathcal{B}'$ se normaliza, de forma que obtendremos una base ortonormal.

_Primera fase_

El proceso del primer paso es incremental.

Llamamos $u_1 = v_1$, el primer vector de la base $\mathcal{B}$.

Ahora consideramos el segundo vector $v_2$. Por un [resultado teórico](#proyec) sobre las proyecciones que vimos antes, si llamamos $\mathcal{U}_1 = \mathcal{L}(\{u_1\})$, entonces $v_2 - \mathrm{proy}_{\mathcal{U}_1}(v_2)$ es un vector ortogonal a $U_1$, en particular es ortogonal a $u_1$. 

Llamemos entonces $u_2 = v_2 - \mathrm{proy}_{\mathcal{U}_1}(v_2)$.
Entonces, el sistema $\{u_1, u_2\}$ es ortogonal.

Llamamos ahora $\mathcal{U}_2 = \mathcal{L}(\{u_1,u_2\})$. Entonces, por el mismo motivo que antes, si llamamos $u_3 = v_3 - \mathrm{proy}_{\mathcal{U}_2}(v_3)$, tenemos que $u_3$ es ortogonal a $\mathcal{U}_2$, luego $\{u_1, u_2, u_3\}$ es un sistema ortogonal.

De esta forma, sucesivamente, en el paso $m$, vamos construyendo $\mathcal{U}_m = \mathcal{L}(\{u_1,\ldots,u_m\})$ y consideramos $u_{m+1} = v_{m+1} - \mathrm{proy}_{\mathcal{U}_m}(v_{m+1})$, que es un vector ortogonal a $U_m$ y, por tanto, $\{u_1,\ldots,u_{m+1}\}$ es un sistema ortogonal.

Una vez se haya completado el proceso, se habrá construido el conjunto $\mathcal{B}' = \{u_1,\ldots,u_n\}$. Como estamos en un espacio $V$ de dimensión $n$, y tenemos $n$ vectores ortogonales, luego [linealmente independientes](#ortho), éstos forman una base. Luego $\mathcal{B}'$ es una base ortogonal.

_Segunda fase_

La segunda fase es más sencilla, ya que implica únicamente normalizar cada vector de $\mathcal{B}'$.

De esta forma, la base ortonormal es la formada por 
$$\left\{\frac{u_1}{\|u_1\|}, \ldots, \frac{u_n}{\|u_n\|}\right\}$$

---

__Ejemplo__

Consideremos la base $\mathcal{B}$ en \(\mathbb{R}^{`r n`}\) dada por
```{r}
repeat {
  
  B <- rMatrix(n = n, m = 4 * n, values = -1:1) %>% 
    linearly_independents()
  
  if (ncol(B) == n) break
  
}

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

Vamos a usar el método de Gram-Schmidt para encontrar una base ortonormal, a partir de $\mathcal{B}$.

_Fase 1_

Llamamos 
```{r}
v1 <- matrix(B[, 1], ncol = 1)
Bp <- v1
glue_latex(
  "u_1 = v_1 = [to_latex(v1)]"
) %>% 
  cat()
```

A partir de este momento, procedemos de forma recursiva:
```{r}
cat("\n")
for (i in seq(2, n)) {
  
  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(
    "\\begin{array}{rcl}",
    "\\mathrm{proy}_{\\mathcal{U}_[i-1]}(v_[i]) & = & [proy_form] \\\\",
    " & = & [proy_val] \\\\",
    " & = & [glue_matrices(proyec, fractions = TRUE, latex = TRUE)]\\\\",
    "\\end{array}"
  ) %>% 
    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)] - [glue_matrices(proyec, fractions = TRUE, latex = TRUE)] = ",
    " [glue_matrices(ui_new, fractions = TRUE, latex = TRUE)]"
  ) %>% str_replace_all(pattern = "\n", replacement = " ") %>% cat()
  
  cat("\n")
  
}
```

Por tanto, la base ortogonal a la que llegamos es
```{r}
glue_latex(
  "\\mathcal{B}' = \\{u_i\\} = \\left\\{",
  "[vectors_to_latex(Bp)]",
  "\\right\\}"
) %>% 
  cat()
```

_Fase 2_

Ahora debemos multiplicar cada vector $u_i$ en $\mathcal{B}'$ por el inverso de su norma, y llegamos a los vectores $\frac{1}{\|u_i\|}u_i$, que nos forman la base ortonormal siguiente:
```{r}
base <- lapply(seq(ncol(Bp)), function(v) Vector$new(Bp[, v]))
normas <- lapply(base, 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 <- base %>% 
  sapply(function(v) v$to_latex())

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

base2 <- lapply(seq_along(base), function(i) base[[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(", ")

glue_latex(
  "\\begin{array}{rcl}",
  "\\mathcal{B}'' & = & \\left\\{",
  "[Bpp]",
  "\\right\\} = \\\\",
  " & = & \\left\\{",
  "[base2_str]",
  "\\right\\} \\\\",
  "\\end{array}"
) %>% 
  cat()
```


---


## Qué es una matriz ortogonal y una aplicación ortogonal {#ortho-matrix}
\sectionmark{Matrices y aplicaciones ortogonales}

```{r}
definition("Una matriz cuadrada $Q$ se llama __ortogonal__ si $Q^{\\text{t}}Q = I$ o, equivalentemente, $Q^{-1} = Q^{\\text{t}}$.", "Matriz ortogonal")
```

Es decir, una matriz ortogonal es aquella cuya inversa es su propia traspuesta.

```{r}
theorem("Una matriz $A\\in\\mathcal{M}_n(\\mathbb{R})$ es ortogonal si, y sólo si, sus columnas forman una base _ortonormal_ de $\\mathbb{R}^n$.")
```

```{r}
theorem("Sea $A$ una matriz ortogonal. Entonces $\\mathrm{det}(A) = \\pm 1$.")
```

```{r}
definition("Un endomorfismo $f:V\\to V$ se llama __ortogonal__ si respeta el producto escalar, es decir, si $\\langle\\cdot,\\cdot\\rangle$ es el producto escalar en $V$, entonces $\\langle f(u), f(v) \\rangle = \\langle u, v \\rangle$ para todo $u,v\\in V$.", "Endomorfismo ortogonal")
```

Como la aplicación ortogonal conserva el producto escalar, también conserva todas las cantidades que hemos definido a partir de él: [la norma, la distancia](#norm) y [el ángulo](#ortho) entre vectores de $V$: $\|f(v)\| = \|v\|$, $d(f(u), f(v)) = d(u,v)$ y $\mathrm{ang}(f(u),f(v)) = \mathrm{ang}(u, v)$ para todo $u,v\in V$.

```{r}
theorem("Una aplicación $f:V\\to V$ es ortogonal si, y sólo si, la imagen de toda base ortonormal es de nuevo una base ortonormal de $V$.")
```

¿Qué relación existe entre matrices ortogonales y aplicaciones ortogonales?

```{r}
theorem("Sea $f:V\\to V$ una aplicación lineal y sea $A$ su matriz asociada respecto a una base ortonormal de $V$. Entonces $f$ es ortogonal si, y solo si, la matriz $A$ es ortogonal.")
```

---

__Ejemplo__

Como hemos comentado, si tomamos una base ortonormal y ponemos sus vectores por columnas, nos forman una matriz _ortogonal_. Por tanto, podemos considerar la base $\mathcal{B}''$ del [ejemplo anterior](#gram-schmidt), y construir la matriz ortogonal asociada:
```{r}
base2_str <- base2 %>% 
  sapply(function(v) v$to_latex(ldeco = "", rdeco = "")) %>% 
  str_flatten(" ")
matriz <- glue::glue("\\left({base2_str}\\right)")
glue_latex(
  "Q = [matriz]"
) %>% 
  cat()
```

Es fácil comprobar que, en este caso, $Q^{\text{t}}\ Q = I$.

Hay otros ejemplos de matrices ortogonales. En particular, una tipología concreta de matrices ortogonales en $\mathbb{R}^2$ es muy utilizada en la práctica. Llamemos:
$$R(\alpha) = \left(
\begin{array}{cc}
\cos(\alpha) & \sin(\alpha) \\
-\sin(\alpha) & \cos(\alpha) \\
\end{array}\right)$$

Son las denominadas __matrices de rotación__. Están asociadas a la aplicación lineal que representa un giro o rotación de $\alpha$ radianes en sentido positivo (_antihorario_) alrededor del origen de coordenadas.

Así 
$$
\left(
\begin{array}{c}
x' \\ y' \\
\end{array}\right) = R(\alpha)\left(
\begin{array}{c}
x \\ y \\
\end{array}\right)$$
son las coordenadas del vector girado, en función de las coordenadas del vector original.

Estas matrices de rotación se pueden extender a dimensiones superiores. Por ejemplo, en $\mathbb{R}^3$, la rotación de $\alpha$ radianes alrededor del _eje $X$_ se representa mediante la matriz ortogonal
$$R(\alpha) = \left(
\begin{array}{ccc}
1 & 0 & 0 \\
0 & \cos(\alpha) & \sin(\alpha) \\
0 & -\sin(\alpha) & \cos(\alpha) \\
\end{array}\right)
$$

---

## Qué es la diagonalización ortogonal {#diag-ortho}
\sectionmark{Diagonalización ortogonal}

Recordemos [la definición de endomorfismo o de matriz diagonalizable](#diagonalizable). Una aplicación lineal $f:V\to V$ es diagonalizable si y solo si existe una base de $V$ formada por autovectores de $f$, lo cual equivale a que la representación matricial de $f$ en dicha base sea una matriz diagonal $D$. En representación matricial, esto se traduce en que una matriz $A\in\mathcal{M}_n(\mathbb{R})$ es diagonalizable si y solo si, existe una matriz regular $P$ tal que $D = P^{-1}\ A\ P$. Las columnas de $P$ están formadas por los autovectores de la base de la definición.

Veamos cómo influye el hecho de tener bases ortonormales en un espacio euclídeo en la diagonalización, lo que da lugar a la denominada _diagonalización ortogonal_.

```{r}
definition("Diremos que una aplicación lineal $f:V\\to V$ es __diagonalizable ortogonalmente__ si y solo si existe una base ortonormal en $V$ formada por autovectores de $f$. ", "Aplicación diagonalizable ortogonalmente")
```

```{r}
definition("En matrices, una matriz $A$ se dice __diagonalizable ortogonalmente__ si existe una matriz $P$ _ortogonal_ tal que $D = P^{\\text{t}}\\ A\\ P$, donde $D$, la matriz asociada a $f$ en la base de autovectores, es diagonal con autovalores en la diagonal principal.", "Matriz diagonalizable ortogonalmente")
```

---

__Ejemplo__

```{r}
set.seed(3412)
repeat {
  
  n <- 4
  lambdas <- sample(-1:1, 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")

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()
  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(", ")

```

Consideremos la matriz $A$ dada por 
```{r}
glue_latex(
  "A = [glue_matrices(A, fractions = TRUE, latex = TRUE)]"
) %>% 
  cat()
```
asociada a un endomorfismo \(f:\mathbb{R}^{`r n`} \to \mathbb{R}^{`r n`}\) en la base canónica.

Si tomamos la base
```{r}
glue_latex(
  "\\mathcal{B} = \\left\\{",
  "[base2_str]",
  "\\right\\}"
) %>% 
  cat()
```
resulta ser una [base ortonormal](#ortho-basis) formada por autovectores de $f$ (en [la siguiente sección](#diag-ortho2) detallaremos el proceso completo).

Por tanto, esta matriz (y el endomorfismo $f$) es diagonalizable ortogonalmente. De hecho, si consideramos como matriz $P$ la del [cambio de base](#cambiobase) de $\mathcal{B}$ a la canónica, podemos calcular la expresión diagonal, teniendo en cuenta que $P$ es una matriz ortogonal:

\small
```{r}
P_str <- vectors2matrix(base2)
glue_latex(
  "\\begin{array}{rcl}",
  "D &\\!\\! = &\\!\\! P^{\\text{t}}\\ A\\ P = \\\\",
  "  &\\!\\! = &\\!\\! [glue_matrices(P_str, latex = TRUE)]^{\\text{t}} [glue_matrices(A, latex = TRUE, fractions = TRUE)] [glue_matrices(P_str, latex = TRUE)] \\\\",
  "  &\\!\\! = &\\!\\! [glue_matrices(D, latex = TRUE, fractions = TRUE)]",
  "\\end{array}"
) %>% 
  cat()
```
\normalsize
donde podemos comprobar que los elementos de la diagonal de $D$ son autovalores de $A$.

---

## Qué matrices son diagonalizables ortogonalmente {#diag-ortho2}
\sectionmark{Diagonalización de matrices simétricas}

Supongamos una matriz $A\in\mathcal{M}_n(\mathbb{R})$ que sea [diagonalizable ortogonalmente](#diag-ortho), es decir, tal que existe $P$ ortogonal con $P^{\text{t}}\ A\ P = D$ una matriz diagonal.

Como $D$ es diagonal, en particular tenemos que $D^{\text{t}} = D$, luego: 
$$P^{\text{t}}\ A\ P = D = D^{\text{t}} = (P^{\text{t}}\ A\ P)^{\text{t}} = P^{\text{t}}\ A^{\text{t}}\ \left(P^{\text{t}}\right)^{\text{t}} = P^{\text{t}}\ A^{\text{t}}\ P$$
Luego, multiplicando por $P$ a la izquierda en los extremos de la igualdad, y por $P^{\text{t}}$ en la derecha, recordando que $P$ es una [matriz ortogonal](#ortho-matrix), llegamos a que $A = A^{\text{t}}$.

Esto significa que:

```{r}
theorem("Si una matriz $A\\in\\mathcal{M}_n(\\mathbb{R})$ es _diagonalizable ortogonalmente_, entonces es simétrica.")
```

Por tanto, una matriz _no simétrica_ podrá ser [diagonalizable](#diagonalizable) (dependiendo de si cumple los criterios ya estudiados), pero _nunca_ puede ser _diagonalizable ortogonalmente_.

Podemos estudiar entonces algunas de las propiedades de las matrices simétricas en lo que se refiere a diagonalización:

```{r}
theorem("Sea $A\\in\\mathcal{M}_n(\\mathbb{R})$ una matriz simétrica. Entonces:\n\n- Todos los autovalores de $A$ son reales.\n- Si $\\lambda_1$ y $\\lambda_2$ son dos autovalores distintos, y sus subespacios asociados son $U_{\\lambda_1}$ y $U_{\\lambda_2}$, respectivamente, entonces $U_{\\lambda_1}$ es ortogonal a $U_{\\lambda_2}$.\n- La multiplicidad geométrica de cada autovalor de $A$ coincide con su multiplicidad algebraica.")
```

En consecuencia

```{r}
theorem("Una matriz $A\\in\\mathcal{M}_n(\\mathbb{R})$ es diagonalizable ortogonalmente si y solo si es simétrica.")
```


__¿Cómo diagonalizamos ortogonalmente una matriz simétrica?__

El procedimiento es idéntico al que se conoce para [diagonalizar](#diagonalizable) una matriz cualquiera, con las siguientes puntualizaciones:

- Al ser una matriz simétrica, tenemos asegurada que es diagonalizable (incluso ortogonalmente).
- Al calcular la base del subespacio asociado a un autovalor $\lambda$, podemos hacer un paso extra, usando el [método de Gram-Schmidt](#gram-schmidt), para tener dicha base en forma ortonormal.
- Al finalizar, al igual que hacíamos en la digonalización _estándar_, podemos unir todas las bases de los distintos subespacios asociados a los autovalores y tenemos una base ortonormal del espacio vectorial $V$.

De forma esquemática, incluyendo todo lo que conocemos de diagonalización, queda:

1. Calcular el [polinomio característico](#eigenvalue) $p(\lambda)$ de $A$, y sus raíces, los autovalores de $A$. Sean $\lambda_1,\ldots,\lambda_k$ los autovalores.
2. Para cada autovalor $\lambda$:
    1. Determinar una [base](#base) $\mathcal{B}_{\lambda}$ del [subespacio propio](#eigenspace) $U_{\lambda}$. 
    2. Aplicar [el método de Gram-Schmidt](#gram-schmidt) a $\mathcal{B}_{\lambda}$ para obtener una [base _ortonormal](#ortho-basis) de $U_{\lambda}$ a la que llamaremos $\mathcal{B}'_{\lambda}$.
3. Llamamos $\mathcal{B} = \mathcal{B}'_{\lambda_1}\cup\ldots\cup\mathcal{B}'_{\lambda_k}$. Entonces $\mathcal{B}$ es una base ortonormal de $V$. Además, la matriz $P$ de [cambio de base](#cambiobase) de $\mathcal{B}$ a la canónica es una matriz [ortogonal](#ortho-matrix), tiene por columnas las coordenadas de los vectores de $\mathcal{B}$, y verifica que $D = P^{\text{t}}\ A\ P$ es una matriz diagonal con los autovalores en la diagonal principal (con multiplicidades incluidas).

---

__Ejemplo__

Consideremos la matriz
```{r}
glue_latex(
  "A = [glue_matrices(A, latex = TRUE, fractions = TRUE)]"
) %>% 
  cat()
```

Como es simétrica, es diagonalizable ortogonalmente. Vamos a seguir los pasos vistos para encontrar su forma diagonal y la base ortonormal de autovectores.

_Paso 1: Encontrar los autovalores de $A$_

Para calcular los autovalores de $A$, comenzamos por construir su [polinomio característico](#eigenvalue):
```{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()
```

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()
```

_Paso 2: Encontrar una base ortonormal de cada subespacio propio_

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

```{r}
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, "}\\) será:\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")
    
  }
  
  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(gen[[av_idx]])), 
                       function(i) Vector$new(gen[[av_idx]][, 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())

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

  glue_latex(
    "\\mathcal{B}'_{U_{[av]}} = \\left\\{",
    "[Bpp]\\right\\} = \\left\\{[base2_str]\\right\\}"
  ) %>% 
    cat()
  
  cat("\n\n")
  
}
```

_Paso 3: Hallar una base ortonormal formada por autovectores_

Definimos
```{r}
str <- glue::glue("\\mathcal{B}'_{U_{[autovalores]}}", 
           .open = "[", .close = "]") %>% 
  str_flatten(" \\cup ")
base2_str <- base2 %>% 
  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(base2)
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

Recordemos que en la diagonal de $D$ se encuentran los autovalores (en el mismo orden en que se han puesto los autovectores en las columnas de $P$), con sus respectivas multiplicidades.

---