URL: <https://savannah.gnu.org/bugs/?66993>
Summary: Bug: pow_int issues
Group: GNU Scientific Library
Submitter: fermelelundi
Submitted: Sat 05 Apr 2025 09:04:13 PM UTC
Category: Accuracy problem
Severity: 3 - Normal
Operating System:
Status: None
Assigned to: None
Open/Closed: Open
Discussion Lock: Any
Release: 2.8
_______________________________________________________
Follow-up Comments:
-------------------------------------------------------
Date: Sat 05 Apr 2025 09:04:13 PM UTC By: Fermé le Lundi <fermelelundi>
There are multiple functions in GSL which support raising a double to an int:
1- gsl_pow_int.h: has inline functions such as gsl_pow_2(const double x) {
return x*x; }
2- sys/pow_int.c: functions gsl_pow_int(double x, int n) and
gsl_pow_uint(double x, int n)
3- specfunc/pow_int.c: function gsl_sf_pow_int(const double x, const int n)
To reduce function duplication we could rationalise a few, ie start to mark as
'deprecated'; suggest: gsl_sf_pow_int(x, n).
There are a couple of issues with these functions:
1- Source code documentation in specfunc/pow_int.c states "Does not check for
overflow/underflow.", but this should read: "Does not check for underflow, but
does check for overflow when x = 0.0 and n negative."
2- All functions report the same values, and some are incorrect:
(a) Error compounding: (+0.5)^(-9) = 512.0000000000010232. Unclear how to
eliminate, as this looks like the result of machine imprecision; perhaps a
note in the source code and/or documentation that results are correct up to
the tenth decimal for -9 <= n <= 9?
(b) Raising to power of 1: (-0.3)^1 = -0.3000000000000002. We could benefit
from handling raising to 1 as a special case.
(c) Raising zero from below to a negative int: (-0.0)^(-1) =
-7205759403792794.0000000000000000; where it should report "Division by zero"
error, ie GSL_EZERODIV.
(d) Raising zero to a negative int: sys/gsl_pow_int(0.0, n) reports +inf for n
<= -1, but this should be GSL_EZERODIV.
Note: complex functions gsl_complex_pow (gsl_complex a, gsl_complex b) and
gsl_complex_pow_real (gsl_complex a, double b) probably warrant a review of
their own.
_______________________________________________________
Reply to this item at:
<https://savannah.gnu.org/bugs/?66993>
_______________________________________________
Message sent via Savannah
https://savannah.gnu.org/
signature.asc
Description: PGP signature
