New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
Enhancements for mathmodule #45981
Comments
The patch adds several small enhancements to the math module and pyport.h.
Together with http://bugs.python.org/issue1635 it implements most of PEP-754. |
Cool! If there's a move to add functions to the math module, there are |
Minor typo. Should be IEEE: |
Mark Dickinson wrote:
I've added the inverse hyperbolic, log1p and expm1. copysign is too low Do you want some more functions? Feel free to provide a patch :) Christian |
The patch fixes some compatibility problems on Windows and implements |
Unfortunately, implementing asinh, acosh, atanh correctly is not as simple as just using the formulas asinh(x) = log(x + sqrt(1.+(x*x))) has significant error (i.e. much more than just a few ulps) when x is small (so that the argument of So either some serious work is required here, or code should be borrowed in an appropriately legal Mark |
i suggest abandoning any attempts at implementing math ourselves. I |
Guido van Rossum wrote:
How do you like a compromise? I rip out the new functions for the math Mark is right. My naive functions are a numerical nightmare for small We could reuse the functions from the BSD libm. They are used in Christian |
The trunk_math_sign_inf_nan patch contains just three new method |
One nit: you added a blank line to the end of test_math.py. One bigger issue: the sign() function doesn't seem to work properly for >>> inf = 1e1000
>>> nan = inf/inf
>>> mnan = -nan
>>> math.sign(nan)
-1
>>> math.sign(mnan)
1 IOW a positive nan is considered negative and vice versa. (This is I'm also curious why math.sign(0.0) returns 1 -- this is going to cause This is also missing doc patches. |
Guido van Rossum wrote:
*grr* stupid white space check hook
If I recall the definition correctly NaNs don't have a sign. The content
math.sign(0.0) == 1 and math.sign(-0.0) == -1 is the main purpose of the Because 0.0 == -0.0 it's not possible the distinguish a positive from a Christian |
No, you edited a line that didn't need editing. :-)
Perhaps you recall wrong, as negating the nan returns one with the
Hm, OK, but then passing a zero of some other type (e.g. int) should Perhaps it would be better to have a function math.isneg() that I suggest that you check in the isinf() and isnan() functions and |
Guido van Rossum wrote:
I'm fine with a isneg() function but I wouldn't "fix" it for NaN. It has "Note: Since NaN is "not a number," there really isn't a "+" or "-" sign Christian |
Why not implement copysign? It's a standard, familiar function with well- |
Good idea. Since you seem to like providing patches, can you create |
Guido van Rossum wrote:
Don't forget it's copysign() on Unix but _copysign() on Windows. #if defined(MS_WINDOWS) || defined(HAVE_COPYSIGN)
#ifdef MS_WINDOWS
FUNC2(copysign, _copysign,
#else
FUNC2(copysign, copysign,
#endif
"doc");
#endif should work on all systems. Christian |
Well, the Python API in the math module should always be called copysign(). :-) And what to do if neither is present? Are there any systems without it? |
Guido van Rossum wrote:
takes care of it. It's a macro to define a function which accepts two On Windows _copysign is always defined and on other systems we can use Christian |
OK, just check it in then, but do add docs! On Jan 3, 2008 2:03 PM, Christian Heimes <report@bugs.python.org> wrote:
|
The functionality of what was called (and I agree confusingly so) http://www.opengroup.org/onlinepubs/000095399/functions/signbit.html The 754 standard says NaNs have sign bits, but assigns no meaning to |
Is there still interest in implementing the inverse hyperbolic trig |
George: I'm certainly still interested in having asinh, acosh and atanh in math---I'm The right thing to do would be to wrap the libm functions when they exist, but provide There's a fallback version of asinh already provided in the patch to issue bpo-1381, which For acosh and atanh, you'd have to find the right way to deal with domain errors (e.g. |
I'm +1 in adding fallbacks for important functions like copysign, asinh, But why write our own version if we can reuse existing implementations?
In bpo-1381 I suggested two new files Python/pymath.c and Include/pymath.h. |
The patch implements replacements for copysign, log1p (taken from Mark's copysign is now defined on Windows, too. I also found a bug in configure.in (my fault). The C99 function is Christian |
Excellent! I'll take a look. |
One question: is there a policy on what should happen for singularities and domain errors? If For domain errors (e.g. sqrt(-1), log(-1), acosh(0)), the obvious two options are:
For singularities (e.g. log(0), atanh(1)), the options are basically the same:
I suspect there are use-cases for both types of behaviour here. Of course, the *right* thing to do, in some sense, would be to have a thread-local floating-point A few months ago I definitely would have said that an exception should be raised in both cases, as Tim, if you're listening: any words of wisdom? Should I ask this on python-dev? |
whoops. OverflowError should be something else in the previous post; of course, Currently, on Linux I get:
On OS X I get:
|
Mark Dickinson wrote:
Windows raises the same exceptions as Linux. |
Okay: for now, I guess we just follow the pattern that already exists on I think the OS X sqrt(-1) behaviour is a bug. |
Just a quick addition here regarding the singularities to these |
No: IEEE-754r and the C99 standard both say clearly that atanh(1) should be infinity The general idea is that it's meaningful to set atanh(1) = infinity because that's |
I misunderstood the rationale for the function returning infinite at |
George: I think my last post was a bit rude. I apologize if it came across that way. Mathematical rigor and IEEE-754 recommendations aren't necessarily in conflict here, At any rate, I agree with you that log(0) and atanh(1) should raise Python exceptions, at |
Mark, these are the "spirit of 754" rules:
And these are what Python has historically /tried/ to do (but not always For #1, raise OverflowError. For #2, return a zero (with the appropriate sign if that happens by For #3 and #4, raise ValueError. It may have made sense to raise You're right that you can't make everyone happy, so may as well stick |
Thanks, Tim! Dare I suggest extending these rules to encompass things like sqrt(NaN), log(inf), etc.,
So e.g. cos(infinity) should give a ValueError, while log(infinity) and exp(infinity) |
Mark Dickinson wrote:
The matter should be discussed in a proper PEP and targeted for Python I still don't like the idea of math.atanh(1) == inf. Why? See for yourself: 18.714973875118524
>>> math.atanh(.99999999999999999)
Traceback (most recent call last):
File "<stdin>", line 1, in <module>
ValueError: math domain error (Linux) Christian |
The mail interface screwed up the message: >>> math.atanh(.9999999999999999)
18.714973875118524
>>> math.atanh(.99999999999999999)
Traceback (most recent call last):
File "<stdin>", line 1, in <module>
ValueError: math domain error |
I've added your complex patch and its tests to my patch. The tests are AssertionError: atanh0022:atanh(-1.000000) returned -18.36840028483855, On Linux the tests are passing just fine. |
It's most probably the fault of log1p(): AssertionError: atanh0022:atanh(-0.99999999999999989) returned |
Christian: I'm definitely not proposing atanh(1) = inf: it should raise The new thing here is that since you've made inf and nan more accessible and I think I'm missing the point of your math.atanh(.999...) example. >>> x = .99999999999999999
>>> x == 1.0
True The atanh0022 result is definitely a bug: it looks like either asinh or log1p |
Christian: would it be possible for you to tell me what the argument of the On OS X I get: Input to log1p is 3.2451855365842669075e+32 It's hard to imagine that there's anything wrong with log1p here, since all |
Okay: here's an attempted guess at what's happening on Windows: Near the end of c_atanh, there's a line that reads: r.real = log1p(4.*z.real/((1-z.real)*(1-z.real) + ay*ay))/4.; Somehow, when z.real is 0.99999999999999989 and ay is 0, the argument to log1p is ending up one-quarter of the size that it should be. I Christian: if you have any time to play with this, could you try replacing this line with something like: double temp = 1-z.real
printf("z.real is %.19e\n", z.real);
r.real = log1p(4.*z.real/(temp*temp + ay*ay))/4.; and see if either the problem goes away or if you can confirm that temp is coming out to be 2.2204460492503131e-16 rather than 1.1102230246251565e-16. |
Sorry: those lines should have been: double temp = 1-z.real;
printf("temp is %.19e\n", temp);
r.real = log1p(4.*z.real/(temp*temp + ay*ay))/4.; |
Mark Dickinson wrote:
You got me wrong (and I didn't explain it properly). All complex Christian |
Sorry: I should have read more carefully. So math.atanh works on Linux It's still rather puzzling though. I still suspect that it's the |
Mark Dickinson wrote:
It uses t = 0.5 * log1p((x + x) / (1.0 - x)) for t > 0.5. I presume the Christian |
I disabled the tests for atanh0022 and atanh0023. Other changes: |
The problem with atanh is that it should be using absx throughout, rather than x. So in "if (absx < 0.5)" branch and the following branch, replace all occurrences of x with One other comment: asinh shouldn't directly set errno for a NaN. It should do the same as This makes asinh(float("nan")) return a nan, which makes it consistent with acosh and I think asinh should also return x+x for x an infinity. This again should make it |
Also, for the C-level routines, atanh(1.0) and atanh(-1.0) should definitely return Note that this is not an argument about what Python should do: Python will still raise a |
One more comment: the Py_IS_NAN test in acosh is likely never reached, since the |
Also: -1.0 shouldn't be returned at this level to indicate an error; these are pure C functions Just so I'm doing something a little more constructive than yelling criticisms from the sideline,
|
Mark Dickinson wrote:
I tried your patch. It fixes the problem with atanh0022 and 0023 test The libm of uclibc *does* set errno for signaling NaN. The relevant code w_acosh.c:
double acosh(double x) /* wrapper acosh */
{
#ifdef _IEEE_LIBM
return __ieee754_acosh(x);
#else
double z;
z = __ieee754_acosh(x);
if(_LIB_VERSION == _IEEE_ || isnan(x)) return z;
if(x<1.0) {
return __kernel_standard(x,x,29); /* acosh(x<1) */
} else
return z;
#endif
}
kernel_standard:
case 129:
/* acosh(x<1) */
exc.type = DOMAIN;
exc.name = type < 100 ? "acosh" : "acoshf";
exc.retval = zero/zero;
if (_LIB_VERSION == _POSIX_)
errno = EDOM;
else if (!matherr(&exc)) {
if (_LIB_VERSION == _SVID_) {
(void) WRITE2("acosh: DOMAIN error\n", 20);
}
errno = EDOM;
}
break; |
Hmmm. For atanh(1): raising OverflowError is actually consistent with what currently happens. For acosh(0): you're right, and I'm wrong---this should definitely return a NaN and set I do still think that asinh(nan), atanh(nan) and acosh(nan) should return nan and not I guess I don't really care about asinh(+/-inf), etc: an infinite return value will be |
Mark Dickinson wrote:
But on Linux atanh(1) sets errno = EDOM and raises a ValueError. Either
I'm not an expert on IEEE754 and NaN but I guess "Signaling NaN" also
Yes, it makes perfectly sense and glibc's libm already does the right I inserted the four lines in mathmodule.c:math_1() and something similar
in math_2() to fix Windows.
#ifndef __GNUC__ /* Windows et al */
if (Py_IS_NAN(x))
return PyFloat_FromDouble(x+x);
#endif
Are you fine with OverflowError? I'm going to create a branch for your work. It's getting tedious sending Christian |
The patch was part of the merge of Mark's and my trunk-math branch. It |
Note: these values reflect the state of the issue at the time it was migrated and might not reflect the current state.
Show more details
GitHub fields:
bugs.python.org fields:
The text was updated successfully, but these errors were encountered: