I am trying to understand why the following theorem holds only for $\sigma_1$ and not for $\sigma_n$ ($n \in \mathbb{N} \ \backslash \ \{ 1 \}$). I would like to know if my attempt to answer my own question is right.
Theorem 8.38. Let $L$ be a self-adjoint operator satisfying $(8.5)$ and $(8.6)$. Then the minimum eigenvalue is simple and has a positive eigenfunction.
Proof. If $u$ is an eigenfunction of $\sigma_1$, then it follows from the formula $(8.92)$ that $|u|$ is one also. But then Harnack inequality, Theorem $8.21$, we must have $|u|$ positive (a.e.) in $\Omega$ and hence $\sigma_1$ has a positive eigenfunction. This argument also shows that the eigenfunctions of $\sigma_1$ are either positive or negative and hence it is impossible that two of them are orthogonal, whence $V_1$ must be one-dimensional and $\sigma_1$ simple. $\square$
Relevant informations are the same done in this topic. Another relevant information is
we may characterize the higher eigenvalues of $L$ through the formula
$$(8.94) \ \sigma_m = \inf \ \{ J(u) \ | \ u \not\equiv 0, (u,v) = 0, \forall v \in \{ V_1, \cdots, V_{m-1} \} \}$$
$V_n$ is the eigenspace associated to the eigenvalue $\sigma_n$ ($n \in \mathbb{N}$) and $(\cdot,\cdot)$ is the inner product on $H = W_0^{1,2}(\Omega)$. This was proved on pages $16$ and $17$ of this lecture notes (it is the proof of the lemma $29$).
Before explain what I tried to do to answer my own question, I would like to explain the proof given by Gilbarg and Trudinger.
First, the first line of the proof is justified in the topic linked with relevant informations and the lemma $29$ of the lecture notes linked. The positivity of $|u|$ is immediate from the Harnack's inequality. The impossibility of two eigenfunctions of $\sigma_1$ be orthogonal follows from the observation that if $u$ and $v$ are two of these eigenfunctions, then $u = u^+$ and $v = v^+$ without loss of generality (if $v$ is negative, then we work with $-v^-$) and
$$(u,v) = (u^+,v^+) = \int_{\Omega} u^+v^+ dx + \sum_{i=1}^n \int_{\Omega} \frac{u^+}{\partial x_i} \frac{v^+}{\partial x_i} dx \geq \int_{\Omega} u^+v^+ dx > 0,$$
where $\sum_{i=1}^n \int_{\Omega} \frac{u^+}{\partial x_i} \frac{v^+}{\partial x_i} dx \geq 0$ because this is an inner product (perhaps, the trivial part to show that this defines an inner product is show that it is positive definitive, but this follows from the fact that the seminorm $|u|^2 := \sum_{i=1}^n \int_{\Omega} \left| \frac{u}{\partial x_i} \right|^2 dx$ is a norm on $H = W_0^{1,2}(\Omega)$ from the Poincaré's inequality).
This impossibility ensures that $V_1$ is one dimensional because if $u$ and $v$ are two eigenfunctions of $\sigma_1$, then
$$v = \frac{(u,v)}{(u,u)} u + v^{\perp},$$
where $v^{\perp}$ denotes the normal component of $v$ with respect to $u$, which implies that $v^{\perp}$ is an eigenfunction of $\sigma_1$, hence $v^{\perp}$ is not orthogonal to $u$, a contradiction, therefore $V_1$ is one-dimensional.
Thus, I think the particularity of the eigenfunctions of $\sigma_1$ are used to show that $|u|$ is an eigenfunction of $\sigma_1$ in the topic linked above. Although it is true that $J(u) = J(|u|)$, it is not true that if $(u,v) = 0, \forall v \in \{ V_1, \cdots, V_{m-1} \}$, then $(|u|,v) = 0$ also because
$$(u^+ - u^-, v) = (u,v) = 0 \Longrightarrow (u^+,v) = (u^-,v),$$
but
$$(u^+ + u^-, v) = (u^+,v) + (u^-,v) = 2 (u^+,v) > 0,$$
where the inequality above holds because $v$ is an eigenfunction, hence $v \not\equiv 0$, and $u^+ \not\equiv 0$, otherwise, $u^-$ would be also, hence $u$ would not be an eigenfunction of $\sigma_1$. Thus, we can not apply the lemma $29$ of the lecture notes linked for $\sigma_m$ when $m > 1$).
Thanks in advance!