I am trying to use the gsl library for calculating 1F1. I have some C code.
The following works and matches Mathematica's results for small real arguments z but not large arguments:
/* file 1F1.c */
#include <stdio.h>
#include <stdlib.h>
#include <gsl/gsl_sf.h>
#include <gsl/gsl_math.h>
double _1F1( double b,
double c,
double z)
{
double Fa = gsl_sf_hyperg_1F1(b, c, z);
return Fa;
}
int main()
{
double b = 100000.0;
double c = 14.0;
double z = 0.9;
double result;
result = _1F1(b, c, z);
printf("1F1(%.1f, %.1f, %.1f): %e\n", b, c, z, result);
/* result = 3.2042*10^236 in Mathematica*/
z = 10*z;
result = _1F1(b, c, z);
printf("1F1(%.1f, %.1f, %.1f): %e\n", b, c, z, result);
/* result = 9.52444189237289*10^794 in Mathematica*/
return 0;
}
I get:
$ gcc -std=c17 -Wall -pedantic 1F1.c -lm -lgsl -lcblas; ./a.out
1F1(100000.0, 14.0, 0.9): 3.204195e+236
1F1(100000.0, 14.0, 9.0): -nan
So,how do I get around and evaluate for large z, as Mathematica seems to be able to do?
I also tried using the Cephes library version, and got into the same problem.
Why does Mathematica not have this problem? They must be using some other formula for large arguments.