GNU Scientific Library - Bugs: bug #52570, Inaccuracy of the Airy function...
You are not allowed to post comments on this tracker with your current authentication level.
bug #52570: Inaccuracy of the Airy function due to invocation of GSL's cosine function with large input parameters
Submitter: | Patrick Alken <psa> | ||
Submitted: | Fri 01 Dec 2017 05:32:49 PM UTC | ||
Category: | Runtime error | Severity: | 3 - Normal |
Operating System: | Status: | None | |
Assigned to: | None | Open/Closed: | Open |
Release: |
No files currently attached
Depends on the following items: None found
Items that depend on this one: None found
Carbon-Copy List
No changes have been made to this item
zhoulai =dot= fu =at= gmail =dot= com
Dear GSL developers,
This is a follow-up email of my previous report of regarding GSL’s Airy and
cosine functions (https://lists.gnu.org/archive/html/bug-gsl/2017-11/msg00011.html).
Here I would like to report another kind of unexpected results when using
the Airy function gsl_sf_Ai_e. Code to reproduce the issue:
#include <gsl/gsl_math.h>
#include <gsl/gsl_errno.h>
#include <gsl/gsl_sf_trig.h>
#include <gsl/gsl_sf_airy.h>
int main(){
gsl_sf_result res;
double x = -1.14e34;
gsl_mode_t mode = GSL_PREC_DOUBLE;
int status = gsl_sf_airy_Ai_e(x, mode, &res);
printf("status = %d, val = %g, err = %g\n", status,res.val,res.err);
}
Running the code gives: status = 0, val = -inf, err = inf
The calculation result here, val = -inf, is incorrect: Airy functions are
damped oscillatory for negative x; using Mathematica you get AiryAi
[-1.1e34]=-1.36e-9. Besides, the error estimate, inf, is obviously
imprecise, and the returned status, GSL_SUCCESS 0, may be somewhat
misleading.
Cause of the issue in my humble opinion: Airy_ai (-1.1e34) invokes
gsl_sf_cos_err_e (theta.val, theta.err, &cos_result) where theta.val =
-8.1145794715437919e+50, theta.err = 7.4985953321244595e+35 (according to
gdb). After this invocation, cos_result.val and cos_result.err become -inf
and inf respectively, which is clearly wrong, which stays away from its
[-1,1] bound. Since “it is known that the GSL gsl_sf_cos and gsl_sf_sin
functions fail for large inputs” (from https://lists.gnu.org/archive/
html/bug-gsl/2017-11/msg00012.html), maybe GSL's Airy functions, both the
first and the second kinds, should consider use libm’s cosine instead.
Attached is the same code as above.
Thanks.
Zhoulai Fu