Hi,
[For conclusive proposal, see below.]
>> 1) Internal FPU precision on x86.
>
> Do you have any test cases show the error?
Sure. Consider the following C code:
#include <stdio.h>
int main (int argc, char **argv) {
volatile double v = 1000000.0;
printf ("%.35f\n%.35f\n", 0.002877, 2877.0 / v);
return 0;
}
According to the C standard, both values should give
0.00287699999..something (the closest double representation of 0.002877).
On most x86 platforms without SSE, the above code will give:
0.00287699999999999978320119886632256
0.00287700000000000021688206786052433
The first is the closest double representation of 0.002877, the second
is the closest double-extended representation of 0.002877 trucated back
to double.
Current PHP (5.2.6) gives me (on my Linux platform):
php -r 'printf ("%.35f\n%.35f\n", 0.002877, 2877.0 / 1000000.0);'
0.00287700000000000021688206786052433
0.00287700000000000021688206786052433
This is due to the fact that PHP has an own implementation of strtod() -
which uses double-extended instead of double for the entire internal
calculation and just returns a double.
PHP 5.2.6 on Windows on the other hand gives me the correct double
representation (because MSVC/Windows always sets the FPU into double
precision mode anyway).
PHP on x86_64 or PHP on x86 Macs (where the compiler defaults to SSE) or
PHP on PPC or Sparc etc. also have the correct double representation.
Thus, FP semantics of PHP are highly platform and (!) compiler
dependent. Worse, since PHP only exposes the double data type, the
internal precision of x87 FPUs can never actually be profited from in
PHP itself - consider $a * $b / $c. This will first calculate $a * $b
(with double-extended precision), truncate that to double and store it
in a temp var. Then it will calculate temp var / $c (double-extended)
and truncate that. So the extra precision is never actually used within
PHP and programmers can never actually take advantage of it (unlike in C
where the variables remain in FPU registers by default [2]).
On the other hand, having the above situation where the division of two
exactly representable (!) FP numbers does not yield the closest FP
representation of the actual result is problematic in my eyes.
[1] http://www.wrcad.com/linux_numerics.txt
[2] And even that is problematic in C, since for example before function
calls, FPU registers are stored back to memory and loaded again after
the function call returns, so one can't even be sure when extended
precision is used and when not.
>> 4) round (1.255, 2) should give 1.26 but gives 1.25. The FUZZ stuff
>> tries to resolve this issue (but not the other three) but I've
>> come to the conclusion that the FUZZ is actually the wrong solution
>> for the problem.
>
> The current implementation does return 1.26
Errr, yes, of course. Wrong grammar. I meant »gave 1.25« on certain
platforms. Your solution works for the cases discussed here and is
certainly much better than the previous code.
But it doesn't always work: Take 32769.255 for example. The FUZZ version
gives 32769.25 where 32769.26 would be expected behaviour (previous PHP
versions were of course never better).
> The proposed solution seems fairly complex and wide-reaching to me, I am
> also concerned with the overheads it introduces as that was a problem
> with many in-C implementation libc folks have tried and rejected.
The libc folks didn't have that many problems. C's round() function
doesn't accept a precision parameter, thus *always* rounds to integer.
This means, if you have a look at my previous points:
1) Internal FPU precision
2) Spec problem (rounding mode)
3) Dividing/multiplying by >= 10^23 inexact
4) round(1.255, 2)
1, 3 and 4 immediately disappear because they ONLY applie if you
actually have a precision parameter (FPU precision is irrelevant for
rounding to integer - if the original number with integer and fraction
part could be represented, the integer part alone can also be
represented - and no division / multiplication is necessary for integer
rounding, so 3 is irrelevant. Since there's no precision parameter,
rounding 1.255 to two places with C's round is irrelevant anyway).
2 also disappears since the C 99 standard defines round() to behave like
arithmetic rounding (round-to-nearest, round-half-up).
On the other hand, PHP has this precision parameter. We can't just drop
it. The best thing would probably have been to never add it but we are
about 8 years or so late for that. So we have to live with it.
My proposal adresses all the four problems. 1 is addressed by the FPU
precision macros (my two patches for ZE2 alone), 3 and 4 are also
addressed (and my implementation correctly handles the above example I
provided) - for implementation details see my RFC.
2 is addressed by the mere fact that I add a third optional parameter to
round() that specifies the rounding mode. Thus the user of the round
function can decide if he/she wants arithmetic or banker's rounding (and
for completeness I've added two others).
As for performance, I've already done some tests. With compiler
optimization, my implemetation is about 14% slower than the previous
implementation for the VAST majority of the use cases. This is certainly
a speed setback but on my laptop (Pentium M 2 GHz) I can still achieve
nearly 1.4 million rounding operations per second with PHP (previously
neraly 1.6 million). I don't really see a problem with this slowdown.
For special cases where somebody wants to round to 23 or 24 digits
precision (i.e. forcing to-string-and-back conversion in the proposed
algorithm), my algorithm is actually about 5 times slower. But I
honestly believe that somebody who writes round ($num, 30) is by far
more interested in results as exact as possible than execution speed.
----------------------------- snip ------------------------------------
Anyway, could we agree that we use your implementation for PHP 5.2.7
(it's better than the current situation and although it isn't perfect it
doesn't change the entire logic in just a bugfix-only revision) and have
my proposal used for PHP >= 5.3?
Regards,
Christian
--
PHP Internals - PHP Runtime Development Mailing List
To unsubscribe, visit: http://www.php.net/unsub.php
