Quadratic variation of the Ornstein-Uhlenbeck process Let $(X_t)_{t\geq 0}$ be the zero-mean Ornstein-Uhlenbeck process such that $X_0 = 0$ almost surely, i.e.
$$X_t = \sigma e^{-\alpha t}\int_0^t e^{\alpha s}\,dB_s \quad \qquad (\triangle)$$ 
On the other hand, $(X_t)$ is the unique process that satisfies the SDE
$$dX_t = \alpha X_t\,dt + \sigma\,dB_t \quad X_0 = 0 \qquad (\square)$$
Since the SDE $(\square)$ satisfies the growth and the Lipschitz conditions, we know that the strong solution to this SDE, $(X_t)$, exists, is unique and continuous.
From the latter $[X,X]_t = [\int \alpha X_s\,ds + \int \sigma\,dB_s,\int \alpha X_s\,ds + \int \sigma\,dB_s]_t = [\int \sigma\,dB_s, \int \sigma\,dB_s]_t = \sigma^2 t$
Here I really needed continuity to make sure that the contribution of $\int \alpha X_s\,ds$ to the quadratic variation is $0$ since then this integral is of bounded variation and is continuous. 
Anyway, my question is how would one compute $[X,X]_t$ by just using $(\triangle)$?
I know that for semimartingales $M,N$
$$[G\cdot M,H\cdot N]_t = \int_{(0,t]}G_sH_s\,d[M,N]_t$$
I applied this result to my case as follows
\begin{align}[X,X]_t =& [\sigma e^{-\alpha t}\int_0^t e^{\alpha s}\,dB_s,\sigma e^{-\alpha t}\int_0^t e^{\alpha s}\,dB_s]\\ =& \sigma^2  [\int_0^t e^{-\alpha (t-s)}\,dB_s, \int_0^t e^{-\alpha (t-s)}\,dB_s]\\=& \sigma^2 \int_0^t e^{-2\alpha (t-s)}ds \neq \sigma^2 \end{align}
Clearly, I am doing something wrong in the last step. I feel like the $e^{-\alpha t}$ term in $X_t$ should be handled in a different way but I couldn't wrap my head around this.

In Formula
$$[G\cdot M,H\cdot N]_t = \int_{(0,t]}G_sH_s\,d[M,N]_t,\tag{1}$$
$G$ and $H$ are functions of $s$ only. However, in your case, $G=H=\sigma e^{-\alpha(t-s)}$, which are not functions of $s$ only, and then Formula$(1)$ is not applicable. 
