3
$\begingroup$

I need to evaluate the hypergeometric function $_2F_1$ with $|z| > 1$ as in Wolfram Language with GSL but the GSL documentation says the $_2F_1$ needs $|z| < 1$.

Is there any way I can use GSL and $_2F_1$ with $|z| > 1$? Any reference about the subject?

$\endgroup$
1

3 Answers 3

2
$\begingroup$

The reason GSL doesn't do this is that you need an analytic continuation. That's what Mathematica does. For example, see here and here. A very similar question was asked here, but they propose a roll-your-own solution.

Edit: Useful paper reference here per my link above (thanks for pulling it out @JoseM)

$\endgroup$
2
1
$\begingroup$

Please, consider this as a toy or educational code... If you want to do it yourself consider the use of a variable precision library like GMP.

In my case z > 2, so reading the paper on the forum suggested by @spencer-bryngelson, Computation of Hypergeometric Functions you get Eq.4.20 (use table 13) or you can use DLMF 15.8.2 as suggested by @a-rural-reader (for me this's more clear because you can see a discontinuity if $$a - b \in \mathbb{Z}$$ see Hypergeometric2F1 too)

#include <stdio.h>
#include <stdlib.h>

#include <gsl_sf.h>
#include <gsl_math.h>

double _2F1(double a,
            double b,
            double c,
            double z)
{
  /* a - b !€ Z */
  double bd = b - 1e-15;

  double A = (gsl_sf_gamma(bd - a)*gsl_sf_gamma(c)*pow(-z, -a))/(gsl_sf_gamma(b)*gsl_sf_gamma(c - a));
  double Fa = gsl_sf_hyperg_2F1(a, a - c + 1.0, a - bd + 1.0, 1.0/z);
  double B = (gsl_sf_gamma(a - bd)*gsl_sf_gamma(c)*pow(-z, -b))/(gsl_sf_gamma(a)*gsl_sf_gamma(c - bd));
  double Fb = gsl_sf_hyperg_2F1(bd, bd - c + 1.0, bd - a + 1.0, 1.0/z);

  return A*Fa + B*Fb;
}

int main()
{
  /* Mathematica
  Hypergeometric2F1[2., 3., 4., 5.0]
  = 0.156542 + 0.150796i
  */

  double a = 2.0;
  double b = 3.0;
  double c = 4.0;
  double z = 5.0;

  double result = _2F1(a, b, c, z);
  printf("2F1(%.1f, %.1f, %.1f, %.1f): %f\n", a, b, c, z, result);
  /* result = 0.15625*/

  /*
  Hypergeometric2F1[2., 3., 4., 10.0]
  Out[30]= 0.03985 + 0.0188496 I
  */
  z = 10.0;

  result = _2F1(a, b, c, z);
  printf("2F1(%.1f, %.1f, %.1f, %.1f): %f\n", a, b, c, z, result);
  /* result = 0.0341796875*/

  return 0;
}

As you can see the results are far from accurate.

If you don't want to do it yourself, I think Arb seems a good choice.

$\endgroup$
1
  • $\begingroup$ Nice follow up, José. $\endgroup$
    – user20857
    Dec 1, 2021 at 19:17
0
$\begingroup$

You can use a transformation formula. Here is how I do in R:

library(gsl)
Gauss2F1 <- function(a,b,c,x){
    if(x>=0 & x<1){
        hyperg_2F1(a,b,c,x)
    }else{
        hyperg_2F1(c-a,b,c,1-1/(1-x))/(1-x)^b
    }
}
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.