Let us consider the linear equation system $$\bf Ax = b$$
We can formulate it's normal equations:
$${\bf A}^T{\bf Ax=A}^T{\bf b}$$
but these are often harder to solve, because ${\bf A}^T{\bf A}$ has worse condition number than $\bf A$. Some time ago I stumbled upon this approach to formulate a linear equation system without normal equations:
$$\begin{bmatrix}\bf I&\bf A\\{\bf A}^T&\bf 0\end{bmatrix}\begin{bmatrix}\bf r\\\bf x\end{bmatrix}=\begin{bmatrix}\bf b\\\bf 0\end{bmatrix}$$
Here we can see the constructed matrix is symmetric so we can avoid the normal equations, but at the expense of blowing up the space to dim(r)+dim(x).
How can we prove / realize that this choice will be better/faster?
This does not seem to work in general.
The $\bf A$ matrix can for example be such that ${\bf B} = \begin{bmatrix}\bf I&\bf A\\{\bf A}^T&\bf 0\end{bmatrix}$ has non-positive eigenvalues. In this case Conjugate Gradient does not converge as it demands positive definiteness. We can circumvent by ${\bf B}^2{\bf x = Bb}$ but then the equation will meet same troubles as we tried avoiding with dodging normal equations in the first place.
However, if we instead run Biconjugate Gradient in Matlab, it seems to converge nicely.
I haven't read about the differences about this biorthogonalization thing or why it could work better.