Search code examples
pythonnumpyscipycomplex-numbers

Meijer G-function in Python and scipy


I'm in need of a Meijer G function in scipy. I read somewhere on the internet that due to its generality, the Meier G function is not supported as a special function in scipy, but everyone should write something up according to his personal use case.

My problem is that I have no experience whatsoever with complex integration. As LaTeX is forbidden here, here's what I'm trying to solve numerically:

enter image description here

(the first line being the general case, the second line my case that I'm trying to compute), with p(a), k, k2 given

As wikipedia states, there are three ways to get L:

  • L runs from −i∞ to +i∞ such that all poles of Γ(bj − s), j = 1, 2, ..., m, are on the right of the path, while all poles of Γ(1 − ak + s), k = 1, 2, ..., n, are on the left.
  • L is a loop beginning and ending at +∞, encircling all poles of Γ(bj − s), j = 1, 2, ..., m, exactly once in the negative direction, but not encircling any pole of Γ(1 − ak + s), k = 1, 2, ..., n.
  • L is a loop beginning and ending at −∞ and encircling all poles of Γ(1 − ak + s), k = 1, 2, ..., n, exactly once in the positive direction, but not encircling any pole of Γ(bj − s), j = 1, 2, ..., m.

How do I get L and solve the integral? The way I'm used to compute integrals over the reals is to

import numpy as np
myL = np.linspace(0, 1, 100)
densityL = myL[1] - myL[0]
myIntegral = (F(myL)*densityL).sum()

I'm not too much after efficiency, I'd prefer a simple and slow working example that I can use to understand the methodology.


Solution

  • For something as complicated, I really think you should avoid calculating the integral yourself, especially if you don't have experience with complex integration, and use a well tested existing implementation.

    Meijer G-function is implemented in mpmath and possibly in Sympy.