# Hypergeometric function with a matrix argument

I am looking for the evaluation of a Hypergeometric function with a matrix argument as for example in Koev and Edelman or as showcased in this Wikipedia article.

From what I understand from Mathematica's documentation, it only accepts a scalar as the last argument.

• I didn't look at the links, but if your matrix is diagonalizable you can apply the scalar function to the eigenvalues and transform the result back afterwards. – Jens Aug 30 '14 at 19:17
• Is there a reason why a diagonalization approach doesn't work? Please provide some code with a minimal example. – Jens Aug 31 '14 at 5:33
• In the paper, the authors note that this function depends only on the eigenvalues of the matrix, so their algorithm is formulated only for diagonal matrices (and indeed, you may diagonalize your matrix and then assume that it started out this way -- but still, as others have said, Mathematica doesn't have this built-in). – Kellen Myers Aug 31 '14 at 16:41
• @Jens et al: please reopen the question. There's suprisingly much to be told. Thank you in advance. – Dr. Wolfgang Hintze Sep 4 '14 at 15:48
• @Dr.WolfgangHintze, it is no longer on hold. You can transfer your answer now. – RunnyKine Sep 4 '14 at 18:11

## 4 Answers

Matrix functions in MMA

First of all MMA does in general not support matrix arguments in its standard functions. Therefore there are special functions MatrixExp and MatrixPower available. But, as will be shown in this answer, it is possible to create user defined functions via infinite series, and it turns out that MMA is surprisingly well able in dealing with these matrix functions. The idea is simple: a function with a known series expansion can be transformed into a Matrix function letting z^k -> MatrixPower[z,k].

Hypergeometric matrix function in one variable

Let us start with the well known function $2F1(a,b,c;z)$ which in MMA is Hypergeometric2F1[a,b,c,z], and define

matrix2F1[a_, b_, c_, z_] :=
Sum[Pochhammer[a, k] Pochhammer[b, k]/
Pochhammer[c, k] MatrixPower[z, k]/k!, {k, 0, \[Infinity]}]

Here z is a matrix.

Example 1: Pauli matrix

z = PauliMatrix[1]

$\left( \begin{array}{cc} 0 & 1 \\ 1 & 0 \\ \end{array} \right)$

matrix2F1[a,b,c,z]

$\left( \begin{array}{cc} \frac{1}{2} (\text{Hypergeometric2F1}[a,b,c,-t]+\text{Hypergeometric2F1}[a,b,c,t]) & \frac{1}{2} (-\text{Hypergeometric2F1}[a,b,c,-t]+\text{Hypergeometric2F1}[a,b,c,t]) \\ \frac{1}{2} (-\text{Hypergeometric2F1}[a,b,c,-t]+\text{Hypergeometric2F1}[a,b,c,t]) & \frac{1}{2} (\text{Hypergeometric2F1}[a,b,c,-t]+\text{Hypergeometric2F1}[a,b,c,t]) \\ \end{array} \right)$

I found it very surprising that MMA computes the series and recognizes closed expressions without problems.

Example 2: general 2x2 matrix

z = {{p, q}, {r, s}}

$\left( \begin{array}{cc} p & q \\ r & s \\ \end{array} \right)$

matrix2F1[1, 2, 2, z t] // FullSimplify;
% // MatrixForm

$\left( \begin{array}{cc} \frac{1-s t}{1-t (s+q r t)+p t (-1+s t)} & -\frac{q t}{-1+t (p+s+q r t-p s t)} \\ -\frac{r t}{-1+t (p+s+q r t-p s t)} & \frac{1-p t}{1-t (s+q r t)+p t (-1+s t)} \\ \end{array} \right)$

Hypergeometric matrix function in two variables

For two variables there is a problem to be clarified first: two matrices in general do not commute, i.e. the result of a multiplication depends on the order. Hence $f(x,y)\neq f(y,x)$, in general. A relation between $f(x,y)$ and $f(y,x)$ become a Little simpler if x.y = A y.x + B with some scalars A and B. But let's us Abandon this question for a Moment and go to the hypergeometrics.

There are 4 different generalizations of the hypergeometric function, called Appell functions (http://mathworld.wolfram.com/AppellHypergeometricFunction.html).

Let's take the first one and define

matrixAppellF1[a_, b_, b1_, c_, x_, y_] :=
Sum[(Pochhammer[a, m + n] Pochhammer[b, m] Pochhammer[b1, n])/(
m! n! Pochhammer[c, m + n])
MatrixPower[x, m].MatrixPower[y, n], {m, 0, \[Infinity]}, {n,
0, \[Infinity]}]

Here we have taken the powers of x together to the left of the powers of y. This seems to be Kind of "natural", but there are other possibilities.

Example 3: AppelF1[a,b,b1,c,x,y]

First check result with scalars

AppellF1[1, 1, 1, 1, t , u ]

$\frac{1}{(1-t) (1-u)}$

x = PauliMatrix[1];
y = PauliMatrix[3];

matrixAppellF1[1, 1, 1, 1, t x, u y]

$\left( \begin{array}{cc} -\frac{1}{(-1+t) (1+t) (1-u)} & -\frac{t}{(-1+t) (1+t) (1+u)} \\ -\frac{t}{(-1+t) (1+t) (1-u)} & -\frac{1}{(-1+t) (1+t) (1+u)} \\ \end{array} \right)$

Summary

We have shown that MMA is very well suited to handle complicated functions with Matrix arguments.We simply have to replace in the power series Expansion of the usual function the power of the variable z^k with MatrixPower[z,k].

Surprisingly enough, MMA can do the infinite sums and provide closed expressions.

By the way: the there is no need for the user to think about eigenvalues or diagonalization.

Best regards, Wolfgang

• After correcting a missing ending in the first definition, I get stuck in line 3: matrix2F1[a,b,c,z] yields an error with 0^0 (using version 10). – Jens Sep 4 '14 at 19:06
• @Jens: sorry, I made a simple typing mistake in the basic formula for matrix2F1: I wrote z^k instead of MatrixPower[z,k]. I have corrected it. The results we ok. Might have been a Freudian error to make exactly the essential part wrong ;-). – Dr. Wolfgang Hintze Sep 5 '14 at 7:18
• Still waiting for your comments ... – Dr. Wolfgang Hintze Sep 5 '14 at 12:28
• Yes, now it all works as advertised (+1)! – Jens Sep 5 '14 at 16:57
• Wow, this is wonderful! Thanks so much guys! – Hirek Sep 8 '14 at 9:12

To expand on my comment to the question: there is in fact an "official" recommendation for how to approach this in the documentation of DiagonalizableMatrixQ, under "Applications". It uses the idea that for a function that can be expanded in a power series, we can always convert a matrix argument into set of a scalar arguments if the matrix is diagonalizable. The similarity transformation of the diagonalization can be pulled out of the power series, and the function only needs to be evaluated for the eigenvalues.

Since most built-in functions (including the hypergeometric ones) are listable, it will be efficient if we can use this listability to our advantage. With a matrix as an argument, listability is initially undesirable, but when we have the eigenvalues as a list it's a great feature.

So here is a more streamlined implementation of the suggestion in the documentation:

matrixEval[func_, arg_?DiagonalizableMatrixQ] := Module[{u, d, e},
{e, u} = Eigensystem[arg];
u = ConjugateTranspose[u];
d = func[e];
If[! ListQ[d],
d = func /@ e];
u.DiagonalMatrix[d].Inverse[u]
]

Using the Pauli x matrix as an example, I'll first show how the diagonalizing similarity transform wraps a general function:

Clear[f];
z = PauliMatrix[1];

matrixEval[f[#] &, z] // TraditionalForm

$$\left( \begin{array}{cc} \frac{f(-1)}{2}+\frac{f(1)}{2} & \frac{f(1)}{2}-\frac{f(-1)}{2} \\ \frac{f(1)}{2}-\frac{f(-1)}{2} & \frac{f(-1)}{2}+\frac{f(1)}{2} \\ \end{array} \right)$$

Now replace f by the desired function:

matrixEval[Hypergeometric2F1[a, b, c, #] &, z] // TraditionalForm

$$\left( \begin{array}{cc} \frac{\, _2F_1(a,b;c;-1)}{2}+\frac{\, _2F_1(a,b;c;1)}{2} & \frac{\, _2F_1(a,b;c;1)}{2}-\frac{\, _2F_1(a,b;c;-1)}{2} \\ \frac{\, _2F_1(a,b;c;1)}{2}-\frac{\, _2F_1(a,b;c;-1)}{2} & \frac{\, _2F_1(a,b;c;-1)}{2}+\frac{\, _2F_1(a,b;c;1)}{2} \\ \end{array} \right)$$

This shows the same structure as above because the diagonalization is the same.

Here is a test with a built-in function that actually has a matrix equivalent:

m = {{3, 2}, {1, 4}};

matrixEval[Exp, m] == MatrixExp[m]

(* ==> True *)

Now it remains to address the question of non-diagonalizable matrix arguments. This can be done in principle as follows:

Find the minimal polynomial of the matrix M. If it is of the form x^k then k is the degree of the matrix, meaning the power for which M^k == 0. This degree is then used in a series expansion of the function, and the matrix can be inserted into that by replacing Power with MatrixPower. Such a matrix is called nilpotent, and the good thing is that power series of such matrices always have at most k terms, so there is no need for infinite sums. I'll make this more concrete if the OP wants to work with such arguments. Probably the best way to actually do this in Mathematica is to use the JordanDecomposition.

So in summary, I would always choose the diagonalization approach rather than trying to rewrite the built in function using its infinite series definition because the latter requires you to first find the appropriate definition and then will also be slow. A series approach is, however, no problem when you have a nilpotent matrix argument, because you can simply apply Series to find the expression you need.

And of course, a power series approach is fine for any finite-dimensional matrix argument if the function in question has a finite number of terms in its series. Another case is where the series terminates but the matrix is infinite-dimensional, as in this answer. Then one can still get pretty good approximate results with a MatrixPower approach using truncated matrices.

• (at)Jens: please notice that I had already shown that the original question has a solution in MMA without caring about diagionalization. – Dr. Wolfgang Hintze Sep 5 '14 at 7:22
• @Dr.WolfgangHintze Sure, but my solution works without having to look up the series definition of the function. Going through the eigenvalues as I do here is also safer in cases where you have to worry about the radius of convergence of the function. For things like MatrixExp this isn't an issue because the function is absolutely convergent. But that's not the case generically. – Jens Sep 5 '14 at 17:02
• (at) Jens: part 1 of 2: (1) please remember that the question was how to treat hypergeometric functions of one an two variables with matrix arguments. I have made a natural proposal and solved exactly that question completely and surprisingly easy. (2) I point out again that I don't need to care about the complicated task of diagonalization, eigenvalues and so on. This is done "behind curtain" by MMA, as can be seen in my example of a general 2x2-matrix. Comment continued below. – Dr. Wolfgang Hintze Sep 6 '14 at 20:20
• (at) Jens: part 2 of 2: (3) The convergence question is also taken care of in my approach by allowing for a scalar parameter t. MMA would tell us during summing up the infinite series which range in t would lead to convergence or divergence. Furthermore, in the requested hypergeometric case this is not an issue. (4) the interesting questions of non-commutability I put up for the two variable case still is uncommented; also it seems no to be treated in the references quoted. – Dr. Wolfgang Hintze Sep 6 '14 at 20:21
• @Jens Thank you so much for your thoughtful response! Upon revisiting it, I believe there may be some miscommunication?! The coded up examples here build matrices of hypergeometric functions of each element in a matrix. What I am wondering about is that the hypergeometric function itself takes the matrix as argument and gives a scalar output. – Hirek Jan 30 '15 at 0:38

Mathematica does not support Hypergeometric functions of matrix arguments. Currently only the cited paper by Plamen and Koev at http://math.mit.edu/~plamen/files/hyper.pdf contains a link to code that implements it in MATLAB.

• The paper has pseudocode that is fairly easy to read and implement, so perhaps you can translate it to Mathematica? – rm -rf Aug 30 '14 at 15:58
• I might do that when I have enough time but I need to use my analytical results first for a paper I am writing. Where would I publish that though? – Hirek Aug 30 '14 at 18:50

You can use the new in M9 function MatrixFunction to do this. For instance:

MatrixFunction[Hypergeometric2F1[a, b, c, #]&, {{0, 1}, {1, 0}}] //TeXForm

$\left( \begin{array}{cc} \frac{\, _2F_1(a,b;c;-1)}{2}+\frac{\, _2F_1(a,b;c;1)}{2} & \frac{\, _2F_1(a,b;c;1)}{2}-\frac{\, _2F_1(a,b;c;-1)}{2} \\ \frac{\, _2F_1(a,b;c;1)}{2}-\frac{\, _2F_1(a,b;c;-1)}{2} & \frac{\, _2F_1(a,b;c;-1)}{2}+\frac{\, _2F_1(a,b;c;1)}{2} \\ \end{array} \right)$