The eigenvectors and eigenvalues of matrix geometric mean This is a follow up question on from How to solve a non-homogeneous quadratic matrix equation?.
Given the matrix $G = A(A^{-1}M)^{1/2}=A^{1/2}(A^{-1/2}MA^{-1/2})^{1/2}A^{1/2}$, where $A=-H^{-1}$, with $H$ symmetric negative definite and $M$ symmetric positive definite. 
Is it possible to establish a clear relationship between the eigenvectors and eigenvalues of $G$, and those of $H$ and $M$?
Also, I have read that $G$ is on the midpoint of the geodesic linking $A$ and $M$, in the space of symmetric positive definite matrices, with distance between two matrices $A$ and $M$ defined as
$$ d(A,M) = \| \log(A^{-1/2}MA^{-1/2}) \|_{\mathrm{HS}},$$
where $ \| X \|_{\mathrm{HS}} = \mathrm{Trace} (X^TX)^{1/2}$
. Do you know if it would be possible to have a 3D graphic representation of this relationship for $2\times 2$ matrices, and if so, how one would go about it? 
Once again, many thanks. 
 A: As for visualization, I have seen SPD 2x2 matrices represented as ellipsoids (such as this one, which would be a 3D case), and the segment joining $A$ and $B$ in that Riemannian distance as a sequence of small ellipsoid slowly morphing one into the other (moving either in space or in time). It is a nice geometric visualization, but from a mathematical point of view it does not add too much insight.
I do not know of simple relations between the eigenvalues of $A$, $B$ and $A \#B$, apart from trivial cases where $A$ and $B$ have a common eigenvector. 
This might be because in some sense the "natural" operation on that matrix space is not orthogonal conjugation $A\mapsto QAQ^T$ but generic conjugation $A\mapsto MAM^T$, with $M$ not necessarily orthogonal. As you can read in Robert Bryant's answer, a suitable transformation of that kind can diagonalize two generic matrices $A$ and $B$ at the same time. So "what is preserved under orthogonal conjugation" might not be the best question to ask, in an abstract setting. In a concrete setting where you have the matrices and you want to do computations with them, well, I understand perfectly the logic of your question but I do not have a ready answer for it.
(Incidentally, if you are doing computations with matrix means, you might be interested in Bruno Iannazzo's Matlab toolbox
A: There are several known relations between eigenvalues of $G$ and $A, M$. 
However, as the following equality for $2\times 2$ matrices shows, these relations are as rich as the standard inequalities between eigenvalues of sums / products of positive matrices. 
More specifically, let $G$ denote the matrix geometric mean, then for $2\times 2$ matrices with unit determinant:
\begin{equation*}
  G = \frac{A+M}{\sqrt{\det(A+M)}},
\end{equation*}
so that $\lambda(G) \propto \lambda(A+M)$. For not unit determinants, one can suitably scale the matrices above to obtain a corresponding expression.

$\newcommand{\reals}{\mathbb{R}}\newcommand{\gm}{\sharp}$
For the general case, the following may be of interest. I'll use $A\gm_t M$ to denote the point on the geodesic between $A$ and $M$ given by
\begin{equation*}
   A\gm_t M := A^{1/2}(A^{-1/2}MA^{-1/2})^tA^{1/2},
\end{equation*}
so in particular the midpoint $G=A\gm_{1/2}M$ is the matrix geometric mean.

Theorem. Let $\prec_{\log}$ denote the log-majorisation order, i.e., for $x, y \in \reals_{++}^n$ ordered nonincreasingly, we say $x \prec_{\log} y$ if $\prod_{i=1}^{n-1}x_i \le \prod_{i=1}^{n-1}y_i$ and $\prod_{i=1}^nx_i=\prod_{i=1}^ny_i$. Then, for $A, M > 0$ and $t \in [0,1]$, we have the log-majorisation between the eigenvalues:
  \begin{equation*}
    \lambda(A\gm_t M) \prec_{\log} \lambda(A^{1-t}M^t) \prec_{\log} \lambda(A^{1-t})\lambda(M^t).
\end{equation*}

