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.