bugGNU Scientific Library - Bugs: bug #53876, gsl_sf_hyperg_2F1_renorm missing...

 
 

You are not allowed to post comments on this tracker with your current authentication level.

bug #53876: gsl_sf_hyperg_2F1_renorm missing factor

Submitter:  Charles Hussong <chussong>
Submitted:  Fri 11 May 2018 02:08:22 AM UTC
   
 
Category:  Accuracy problem Severity:  3 - Normal
Priority:  5 - Normal Operating System: 
Status:  None Assigned to:  None
Open/Closed:  Open Release: 

Discussion

Wed 04 Jun 2025 06:00:06 PM UTC, comment #1: 

Because the case handles c being a negative integer, you might as well use gsl_pow_int(), so that in specfunc/hyperg_2F1.c this line:


int stat_e = gsl_sf_exp_mult_err_e(ln_pre_val, ln_pre_err, sg * F.val, F.err, result);


becomes this:


int stat_e = gsl_sf_exp_mult_err_e(ln_pre_val, ln_pre_err, sg * gsl_pow_int(x, -c+1) * F.val, F.err, result);



Then presumably the same should happen in the function gsl_sf_hyperg_2F1_conj_renorm_e() which appears below it, so that this line:


int stat_e = gsl_sf_exp_mult_err_e(ln_pre_val, ln_pre_err, F.val, F.err, result);


becomes this:


int stat_e = gsl_sf_exp_mult_err_e(ln_pre_val, ln_pre_err, gsl_pow_int(x, -c+1) * F.val, F.err, result);




Fermé le Lundi <fermelelundi>
Fri 11 May 2018 02:08:22 AM UTC, original submission:  

Hello,

I apologize for the poor formatting; I hope this report is simple enough that it isn't too offensive.

I was surprised today to see gsl putting out a wrong answer for the renormalized hypergeometric function. The behavior is easily observed with the following call from one of my unit tests: gsl_sf_hyperg_2F1_renorm(1, 2, -3, 0.4). This returns 2572.02, but the correct answer is 65.8436.

Upon further inspection, the cause is pretty clear: the case in question begins on line 843 of specfunc/hyperg_2F1.c, where an appeal is made to Abramowitz and Stegun (15.1.2) [see http://people.mat ... ands/page_556.htm]; this is perfectly valid, but the subsequent code omits the z^(m+1).

I would offer a patch for this directly, but I'm certain I'd do the error handling incorrectly; if not for this concern, simply changing the first argument on line 863 from "sg F.val" to "sg pow(x, -c+1) * F.val" would give the right answer.

Please let me know if there is any further assistance I can provide.

Best regards,
Charles

Charles Hussong <chussong>

 

Attached Files

This item currently has no attached files.

 

Dependencies

This item does not depend on any other items.

No items depend on this one.

 

Mail Notification Carbon-Copy List

Carbon-Copy List
  • -email is unavailable- added by fermelelundi (Posted a comment)
  • -email is unavailable- added by chussong (Submitted the item)
  •  

    History

    No changes have been made to this item

    Back to the top

    Powered by Savane 3.16-ed84.
    Corresponding source code