linear algebra - The arithmetic-geometric mean for symmetric positive definite matrices

A while back, I wanted to see if the notion of the arithmetic-geometric mean could be extended to a pair of symmetric positive definite matrices. (I considered positive definite matrices only since the notion of the matrix square root is a bit intricate for other kinds of matrices.)

I expected that some complications would arise since, unlike scalar multiplication, matrix multiplication is noncommutative. Another complication would be that the product of two symmetric matrices need not be symmetric (though the positive definiteness is retained, so one can still speak of a principal matrix square root).

By analogy with the scalar AGM, I considered the iteration

$$\mathbf A_0=\mathbf A \; ; \; \mathbf B_0=\mathbf B$$ $$\mathbf A_{i+1}=\frac12(\mathbf A_i+\mathbf B_i) \; ; \; \mathbf B_{i+1}=\sqrt{\mathbf A_i \mathbf B_i}$$

I cranked up a short Mathematica routine:

matAGM[u_, v_] := First[FixedPoint[
      {Apply[Plus, #]/2, MatrixPower[Apply[Dot, #], 1/2]} &, {u, v}]] /;
      MatrixQ[u, InexactNumberQ] && MatrixQ[v, InexactNumberQ]

and decided to try it out on randomly generated SPD matrices.

(A numerical note: Mathematica uses the numerically stable Schur decomposition in computing matrix functions like the matrix square root.)

I found that for all of the randomly generated pairs of SPD matrices I tried, the process was convergent (though the rate of convergence is apparently not as fast as the scalar AGM). As expected, the order matters: matAGM[A, B] and matAGM[B, A] are usually not equal (and both results are unsymmetric) unless A and B commute (for the special case of diagonal A and B, the result is the diagonal matrix whose entries are the arithmetic-geometric means of the corresponding entries of the pair.)

I now have three questions:

1. How do I prove or disprove that this process converges for any pair of SPD matrices? If it is convergent, what is the rate of convergence?




2. Is there any relationship between matAGM[A, B] and matAGM[B, A] if the two matrices A and B do not commute?




3. Is there any relationship between this matrix arithmetic-geometric mean and the usual scalar arithmetic-geometric mean? Would, say, arithmetic-geometric means of the eigenvalues of the two matrices have anything to do with this?

---

(added 8/12/2011)

More digging around has me convinced that I should indeed be considering the formulation of the geometric mean by Pusz and Woronowicz:

$$\mathbf A\#\mathbf B=\mathbf A^{1/2}(\mathbf A^{-1/2}\mathbf B\mathbf A^{-1/2})^{1/2}\mathbf A^{1/2}$$

as more natural; the proof of convergence is then simplified, as shown in the article Willie linked to. However, I'm still wondering why the original "unnatural" formulation seems to be convergent (or else, I'd like to see a pair of SPD matrices that cause trouble for the unnatural iteration). I am also interested in how elliptic integrals might crop up in here, just as they did for the scalar version of the AGM.