[Fredrik Johansson] >>> I'd rather like to see a well implemented math.nthroot. 64**(1/3.0) >>> gives 3.9999999999999996, and this error could be avoided.
[Steven D'Aprano] >> >>> math.exp(math.log(64)/3.0) >> 4.0 >> >> Success!!! [Tom Anderson] > Eeeeeeenteresting. I have no idea why this works. Given that math.log is > always going to be approximate for numbers which aren't rational powers of > e (which, since e is transcendental, is all rational numbers, and > therefore all python floats, isn't it?), i'd expect to get the same > roundoff errors here as with exponentiation. Is it just that the errors > are sufficiently smaller that it looks exact? Writing exp(log(x)*y) rather than x**y is in _general_ a terrible idea, but in the example it happens to avoid the most important rounding error entirely: 1./3. is less than one-third, so 64**(1./3.) is less than 64 to the one-third. Dividing by 3 instead of multiplying by 1./3. is where the advantage comes from here: >>> 1./3. # less than a third 0.33333333333333331 >>> 64**(1./3.) # also too small 3.9999999999999996 >>> exp(log(64)/3) # happens to be on the nose 4.0 If we feed the same roundoff error into the exp+log method in computing 1./3., we get a worse result than pow got: >>> exp(log(64) * (1./3.)) # worse than pow's 3.9999999999999991 None of this generalizes usefully -- these are example-driven curiousities. For example, let's try 2000 exact cubes, and count how often "the right" answer is delivered: c1 = c2 = 0 for i in range(1, 2001): p = i**3 r1 = p ** (1./3.) r2 = exp(log(p)/3) c1 += r1 == i c2 += r2 == i print c1, c2 On my box that prints 3 284 so "a wrong answer" is overwhelmingly more common either way. Fredrik is right that if you want a library routine that can guarantee to compute exact n'th roots whenever possible, it needs to be written for that purpose. ... > YES! This is something that winds me up no end; as far as i can tell, > there is no clean programmatic way to make an inf or a NaN; All Python behavior in the presence of infinities, NaNs, and signed zeroes is a platform-dependent accident, mostly inherited from that all C89 behavior in the presence of infinities, NaNs, and signed zeroes is a platform-dependent crapshoot. > in code i write which cares about such things, i have to start: > > inf = 1e300 ** 1e300 > nan = inf - inf That would be much more portable (== would do what you intended by accident on many more platforms) if you used multiplication instead of exponentiation in the first line. ... > And then god forbid i should actually want to test if a number is NaN, > since, bizarrely, (x == nan) is true for every x; instead, i have to > write: > > def isnan(x): > return (x == 0.0) and (x == 1.0) The result of that is a platform-dependent accident too. Python 2.4 (but not eariler than that) works hard to deliver _exactly_ the same accident as the platform C compiler delivers, and at least NaN comparisons work "as intended" (by IEEE 754) in 2.4 under gcc and MSVC 7.1 (because those C implementations treat NaN comparisons as intended by IEEE 754; note that MSVC 6.0 did not): >>> inf = 1e300 * 1e300 >>> nan == nan >>> nan = inf - inf >>> nan == 1.0 False >>> nan < 1.0 False >>> nan > 1.0 False >>> nan == nan False >>> nan < nan False >>> nan > nan False >>> nan != nan True So at the Python level you can do "x != x" to see whether x is a NaN in 2.4+(assuming that works in the C with which Python was compiled; it does under gcc and MSVC 7.1). > The IEEE spec actually says that (x == nan) should be *false* for every x, > including nan. I'm not sure if this is more or less stupid than what > python does! Python did nothing "on purpose" here before Python 2.4. > And while i'm ranting, how come these expressions aren't the same: > > 1e300 * 1e300 > 1e300 ** 2 Because all Python behavior in the presence of infinities, NaNs and signed zeroes is a platform-dependent accident. > And finally, does Guido know something about arithmetic that i don't, Probably yes, but that's not really what you meant to ask <wink>. > or is this expression: > > -1.0 ** 0.5 > > Evaluated wrongly? Read the manual for the precedence rules. -x**y groups as -(x**y). -1.0 is the correct answer. If you intended (-x)**y, then you need to insert parentheses to force that order. -- http://mail.python.org/mailman/listinfo/python-list