URL: <https://savannah.gnu.org/bugs/?67728>
Summary: gsl_sf_psi_n_e yields domain error
Group: GNU Scientific Library
Submitter: fermelelundi
Submitted: Sun 23 Nov 2025 05:21:01 PM UTC
Category: Accuracy problem
Severity: 3 - Normal
Priority: 5 - Normal
Operating System:
Status: None
Assigned to: None
Open/Closed: Open
Discussion Lock: Any
Release: 2.8
_______________________________________________________
Follow-up Comments:
-------------------------------------------------------
Date: Sun 23 Nov 2025 05:21:01 PM UTC By: Fermé le Lundi <fermelelundi>
As per gsl_sf_psi_n_e yields domain error
(https://lists.gnu.org/archive/html/bug-gsl/2025-11/msg00000.html), attached
is a drop-in replacement for the function gsl_sf_psi_n_e in file
specfunc/psi.c
The given test case now passes.
Limitations:
- Not recommended for use of values of n >= 5; in fact, sign is incorrect for
certain sections of x for these domain values of n.
- Close to singularities from the left (-d.1, -d.01, -d.001 etc) accuracy is
particularly poor.
Test cases (specfunc/test_sf.c):
// x < 0
TEST_SF(s, gsl_sf_psi_n_e, (0, -0.5, &r), 2.0 - M_EULER - log(4), TEST_TOL0,
GSL_SUCCESS);
TEST_SF(s, gsl_sf_psi_n_e, (0, -1.5, &r), 8/3.0 - M_EULER - log(4),
TEST_TOL0, GSL_SUCCESS);
TEST_SF(s, gsl_sf_psi_n_e, (1, -0.5, &r), 4.0 + (M_PISQR / 2.0), TEST_TOL0,
GSL_SUCCESS);
TEST_SF(s, gsl_sf_psi_n_e, (1, -1.5, &r), 40/9.0 + M_PISQR / 2, TEST_TOL0,
GSL_SUCCESS);
TEST_SF(s, gsl_sf_psi_n_e, (2, -0.5, &r), 16 - 14*ZETA3, TEST_TOL0,
GSL_SUCCESS);
TEST_SF(s, gsl_sf_psi_n_e, (2, -1.5, &r), (448/27.0) - 14*ZETA3, TEST_TOL2,
GSL_SUCCESS);
TEST_SF(s, gsl_sf_psi_n_e, (2, -2.5, &r), (56432/3375.0) - 14*ZETA3,
TEST_TOL2, GSL_SUCCESS);
TEST_SF(s, gsl_sf_psi_n_e, (2, -4503599224717311 / 536870912.0, &r),
124.02510962421644344751, TEST_TOL0, GSL_SUCCESS);
TEST_SF(s, gsl_sf_psi_n_e, (3, -1.5, &r), 194.594276219187622, TEST_TOL0,
GSL_SUCCESS);
TEST_SF(s, gsl_sf_psi_n_e, (4, -1.5, &r), -0.3137559995067250, TEST_TOL4,
GSL_SUCCESS);
// TEST_SF(s, gsl_sf_psi_n_e, (5, -1.5, &r), -1.5380336530580968e04, 1000,
GSL_SUCCESS); wrong sign
// x < 0 and near pole -1
TEST_SF(s, gsl_sf_psi_n_e, (0, -0.9, &r), -9.312643829299965684057,
TEST_TOL0, GSL_SUCCESS);
TEST_SF(s, gsl_sf_psi_n_e, (0, -0.99, &r), -99.55078444776766439647,
TEST_TOL3, GSL_SUCCESS);
TEST_SF(s, gsl_sf_psi_n_e, (0, -0.999, &r), -999.57457093080929947,
TEST_TOL4, GSL_SUCCESS);
TEST_SF(s, gsl_sf_psi_n_e, (0, -1.1, &r), 10.154163959143857699, TEST_TOL2,
GSL_SUCCESS);
TEST_SF(s, gsl_sf_psi_n_e, (0, -1.01, &r), 100.39631270584539495, TEST_TOL3,
GSL_SUCCESS);
TEST_SF(s, gsl_sf_psi_n_e, (0, -1.001, &r), 1000.42013919689235357,
TEST_TOL6, GSL_SUCCESS);
TEST_SF(s, gsl_sf_psi_n_e, (1, -0.9, &r), 102.66786705202732671845,
TEST_TOL1, GSL_SUCCESS);
TEST_SF(s, gsl_sf_psi_n_e, (1, -0.99, &r), 10002.64151757892030103,
TEST_TOL1, GSL_SUCCESS);
TEST_SF(s, gsl_sf_psi_n_e, (1, -0.999, &r), 1.00000264453619987398e06,
TEST_TOL1, GSL_SUCCESS);
TEST_SF(s, gsl_sf_psi_n_e, (1, -1.1, &r), 102.74898624046893905366,
TEST_TOL3, GSL_SUCCESS);
TEST_SF(s, gsl_sf_psi_n_e, (1, -1.01, &r), 10002.64960015047874540,
TEST_TOL5, GSL_SUCCESS);
TEST_SF(s, gsl_sf_psi_n_e, (1, -1.001, &r), 1.00000264534442778204e06,
1000*TEST_TOL6, GSL_SUCCESS);
TEST_SF(s, gsl_sf_psi_n_e, (2, -0.9, &r), -1999.117973153378299866,
TEST_TOL1, GSL_SUCCESS);
TEST_SF(s, gsl_sf_psi_n_e, (2, -0.99, &r), -2.00000027917837282778e06,
TEST_TOL1, GSL_SUCCESS);
TEST_SF(s, gsl_sf_psi_n_e, (2, -0.999, &r), -2.0000000003916202697e09,
TEST_TOL1, GSL_SUCCESS);
TEST_SF(s, gsl_sf_psi_n_e, (2, -1.1, &r), 1998.3008325580084355940,
TEST_TOL5, GSL_SUCCESS);
TEST_SF(s, gsl_sf_psi_n_e, (2, -1.01, &r), 1.999999470862128567269e06,
1e-05, GSL_SUCCESS);
// TEST_SF(s, gsl_sf_psi_n_e, (2, -1.001, &r), 1.99999999958339177079e09,
TEST_TOL6, GSL_SUCCESS);
TEST_SF(s, gsl_sf_psi_n_e, (3, -0.9, &r), 60013.657824206819061857,
TEST_TOL1, GSL_SUCCESS);
TEST_SF(s, gsl_sf_psi_n_e, (3, -0.99, &r), 6.000000124971840069947e08,
TEST_TOL1, GSL_SUCCESS);
TEST_SF(s, gsl_sf_psi_n_e, (3, -0.999, &r), 6.00000000001249317417e12,
TEST_TOL1, GSL_SUCCESS);
// TEST_SF(s, gsl_sf_psi_n_e, (3, -1.1, &r), 6.0000001251491133342506e08,
TEST_TOL6, GSL_SUCCESS);
// TEST_SF(s, gsl_sf_psi_n_e, (3, -1.01, &r), 6.000000125149113334250e08,
TEST_TOL6, GSL_SUCCESS);
// TEST_SF(s, gsl_sf_psi_n_e, (3, -1.001, &r), 6.00000000001249494671e12,
TEST_TOL6, GSL_SUCCESS);
_______________________________________________________
File Attachments:
Name: psi_dropin.c Size: 1.8KiB
<https://file.savannah.gnu.org/file/psi_dropin.c?file_id=57842>
AGPL NOTICE
These attachments are served by Savane. You can download the corresponding
source code of Savane at
https://savannah.gnu.org/source/savane-ed84fe80348165e3d4f98c86300132395c6fd7e1.tar.gz
_______________________________________________________
Reply to this item at:
<https://savannah.gnu.org/bugs/?67728>
_______________________________________________
Message sent via Savannah
https://savannah.gnu.org/
signature.asc
Description: PGP signature
