Wolfram Function Repository

Instant-use add-on functions for the Wolfram Language

Function Repository Resource:

 MatrixGeometricMean

Source Notebook

Compute the geometric mean of two matrices

Contributed by: Jan Mangaldan

ResourceFunction["MatrixGeometricMean"][a,b]

gives the geometric mean of the matrices a and b.

Details

The geometric mean of two matrices a and b is the unique positive definite matrix solution x of the equation , and is the generalization of the usual geometric mean for scalars.
The matrices a and b can be numerical or symbolic, but must be Hermitian and positive definite.
The geometric mean of two matrices can have complex-valued elements.

Examples

Basic Examples (2) 

The geometric mean of two exact 2×2 symmetric positive definite matrices:

In[1]:=
ResourceFunction["MatrixGeometricMean"][( { {2, 1}, {1, 2} } ), ( { {3, 2}, {2, 3} } )] // Simplify
Out[1]=

The geometric mean is also symmetric and positive definite:

In[2]:=
SymmetricMatrixQ[%] && PositiveDefiniteMatrixQ[%]
Out[2]=

Scope (2) 

Two symmetric positive definite matrices:

In[3]:=
ma = HilbertMatrix[4]; mb = Array[Min, {4, 4}];

Compute the geometric mean with machine arithmetic:

In[4]:=
ResourceFunction["MatrixGeometricMean"][N[ma], N[mb]]
Out[4]=

Compute the geometric mean with 24-digit precision arithmetic:

In[5]:=
ResourceFunction["MatrixGeometricMean"][N[ma, 24], N[mb, 24]]
Out[5]=

Compute the geometric mean of two random Hermitian positive definite matrices:

In[6]:=
ma = ConjugateTranspose[#] . # &[RandomComplex[1 + I, {4, 4}]]; mb = ConjugateTranspose[#] . # &[RandomComplex[1 + I, {4, 4}]]; ResourceFunction["MatrixGeometricMean"][ma, mb]
Out[7]=

Properties and Relations (6) 

MatrixGeometricMean of two 1×1 matrices is equivalent to GeometricMean:

In[8]:=
a = RandomReal[]; b = RandomReal[]; {ResourceFunction["MatrixGeometricMean"][{{a}}, {{b}}], GeometricMean[{a, b}]}
Out[9]=

The geometric mean is symmetric in its arguments:

In[10]:=
n = 3; a = Transpose[#] . # &[RandomReal[1, {n, n}]]; b = Transpose[#] . # &[RandomReal[1, {n, n}]]; ResourceFunction["MatrixGeometricMean"][a, b] == ResourceFunction["MatrixGeometricMean"][b, a]
Out[11]=

The geometric mean of a matrix and the identity is equivalent to the square root of the matrix:

In[12]:=
n = 3; a = Transpose[#] . # &[RandomReal[1, {n, n}]]; ResourceFunction["MatrixGeometricMean"][a, IdentityMatrix[n]] == MatrixPower[a, 1/2]
Out[13]=

If two matrices commute, their geometric mean is equivalent to the square root of their product:

In[14]:=
a = {{34, 12}, {12, 41}}; b = {{2100, 300}, {300, 2275}}; a . b == b . a
Out[15]=
In[16]:=
ResourceFunction["MatrixGeometricMean"][a, b] == MatrixPower[a . b, 1/2] // Simplify
Out[16]=

The geometric mean can be expressed in terms of MatrixPower:

In[17]:=
n = 3; a = Transpose[#] . # &[RandomReal[1, {n, n}]]; b = Transpose[#] . # &[RandomReal[1, {n, n}]]; eps = 10^(-8); mgm = ResourceFunction["MatrixGeometricMean"][a, b]; diff1 = Norm[mgm - MatrixPower[a . Inverse[b], 1/2] . b, Infinity]; diff2 = Norm[mgm - b . MatrixPower[Inverse[b] . a, 1/2], Infinity]; {diff1 < eps, diff2 < eps}
Out[18]=

The geometric mean can be expressed as an integral:

In[19]:=
n = 3; a = Transpose[#] . # &[RandomReal[1, {n, n}]]; b = Transpose[#] . # &[RandomReal[1, {n, n}]]; Norm[ResourceFunction["MatrixGeometricMean"][a, b] - 2/\[Pi] NIntegrate[ Inverse[(1 + t) Inverse[a] + (1 - t) Inverse[b]]/Sqrt[ 1 - t^2], {t, -1, 1}], \[Infinity]] // Chop
Out[20]=

Version History

  • 1.0.0 – 06 May 2022

Source Metadata

  • Citation:
    • Ando T., Li C.K., Mathias R., "Geometric means." Linear Algebra and its Applications, vol. 385, 305-334, 2004.
DOI: 10.1016/j.laa.2003.11.019
    • Bini D.A., Ianazzo B., "A note on computing matrix geometric means." Advances in Computational Mathematics, vol. 35, 175-192, 2011.
DOI: 10.1007/s10444-010-9165-0
    • Ianazzo B., "The geometric mean of two matrices from a computational viewpoint." arXiv:1201.0101 [math.NA], December 30, 2011.
    • Pusz W., Woronowicz S.L., "Functional calculus for sesquilinear forms and the purification map." Reports on Mathematical Physics, vol. 8, no. 2, 159-170, 1975.
DOI: 10.1016/0034-4877(75)90061-0

Related Resources

Author Notes

The matrix geometric mean for a pair of matrices is uniquely defined. The matrix geometric mean of three or more HPD matrices can be defined in many inequivalent ways, and all of them are nontrivial to implement robustly. Thus, that generalization is left for a future version.

License Information