From: Jarkko Hietaniemi Date: Sat, 16 Mar 2002 19:57:03 +0000 (+0000) Subject: Thwart IRIX long doubles and sloppy pow() (in Perl, **); X-Git-Url: http://git.shadowcat.co.uk/gitweb/gitweb.cgi?a=commitdiff_plain;h=58d76dfdf30230046d41beb912124c3381ef2bc8;p=p5sagit%2Fp5-mst-13.2.git Thwart IRIX long doubles and sloppy pow() (in Perl, **); from Nicholas Clark. Fixes lib/integer.t and t/op/pack.t which assume that 2**someinteger is accurate. p4raw-id: //depot/perl@15267 --- diff --git a/MANIFEST b/MANIFEST index 3416aa4..6bcb8ba 100644 --- a/MANIFEST +++ b/MANIFEST @@ -2354,6 +2354,7 @@ t/op/override.t See if operator overriding works t/op/pack.t See if pack and unpack work t/op/pat.t See if esoteric patterns work t/op/pos.t See if pos works +t/op/pow.t See if ** works t/op/push.t See if push and pop work t/op/pwent.t See if getpw*() functions work t/op/qq.t See if qq works diff --git a/pp.c b/pp.c index ead07f0..15bf351 100644 --- a/pp.c +++ b/pp.c @@ -877,10 +877,98 @@ PP(pp_postdec) PP(pp_pow) { dSP; dATARGET; tryAMAGICbin(pow,opASSIGN); +#ifdef PERL_PRESERVE_IVUV + /* ** is implemented with pow. pow is floating point. Perl programmers + write 2 ** 31 and expect it to be 2147483648 + pow never made any guarantee to deliver a result to 53 (or whatever) + bits of accuracy. Which is unfortunate, as perl programmers expect it + to, and on some platforms (eg Irix with long doubles) it doesn't in + a very visible case. (2 ** 31, which a regression test uses) + So we'll implement power-of-2 ** +ve integer with multiplies, to avoid + these problems. */ { - dPOPTOPnnrl; - SETn( Perl_pow( left, right) ); - RETURN; + SvIV_please(TOPm1s); + if (SvIOK(TOPm1s)) { + bool baseuok = SvUOK(TOPm1s); + UV baseuv; + + if (baseuok) { + baseuv = SvUVX(TOPm1s); + } else { + IV iv = SvIVX(TOPm1s); + if (iv >= 0) { + baseuv = iv; + baseuok = TRUE; /* effectively it's a UV now */ + } else { + baseuv = -iv; /* abs, baseuok == false records sign */ + } + } + SvIV_please(TOPs); + if (SvIOK(TOPs)) { + UV power; + + if (SvUOK(TOPs)) { + power = SvUVX(TOPs); + } else { + IV iv = SvIVX(TOPs); + if (iv >= 0) { + power = iv; + } else { + goto float_it; /* Can't do negative powers this way. */ + } + } + /* now we have integer ** positive integer. + foo & (foo - 1) is zero only for a power of 2. */ + if (!(baseuv & (baseuv - 1))) { + /* We are raising power-of-2 to postive integer. + The logic here will work for any base (even non-integer + bases) but it can be less accurate than + pow (base,power) or exp (power * log (base)) when the + intermediate values start to spill out of the mantissa. + With powers of 2 we know this can't happen. + And powers of 2 are the favourite thing for perl + programmers to notice ** not doing what they mean. */ + NV result = 1.0; + NV base = baseuok ? baseuv : -(NV)baseuv; + int n = 0; + + /* The logic is this. + x ** n === x ** m1 * x ** m2 where n = m1 + m2 + so as 42 is 32 + 8 + 2 + x ** 42 can be written as + x ** 32 * x ** 8 * x ** 2 + I can calculate x ** 2, x ** 4, x ** 8 etc trivially: + x ** 2n is x ** n * x ** n + So I loop round, squaring x each time + (x, x ** 2, x ** 4, x ** 8) and multiply the result + by the x-value whenever that bit is set in the power. + To finish as soon as possible I zero bits in the power + when I've done them, so that power becomes zero when + I clear the last bit (no more to do), and the loop + terminates. */ + for (; power; base *= base, n++) { + /* Do I look like I trust gcc with long longs here? + Do I hell. */ + UV bit = (UV)1 << (UV)n; + if (power & bit) { + result *= base; + /* Only bother to clear the bit if it is set. */ + power &= ~bit; + } + } + SP--; + SETn( result ); + RETURN; + } + } + } + } + float_it: +#endif + { + dPOPTOPnnrl; + SETn( Perl_pow( left, right) ); + RETURN; } } diff --git a/t/op/pow.t b/t/op/pow.t new file mode 100644 index 0000000..2e1d29f --- /dev/null +++ b/t/op/pow.t @@ -0,0 +1,46 @@ +#!./perl -w +# Now they'll be wanting biff! and zap! tests too. + +BEGIN { + chdir 't' if -d 't'; + @INC = '../lib'; + require './test.pl'; +} + +# This calcualtion ought to be within 0.001 of the right answer. +my $bits_in_uv = int (0.001 + log (~0+1) / log 2); + +# 3**30 < 2**48, don't trust things outside that range on a Cray +# Likewise other 3 should not overflow 48 bits if I did my sums right. +my @pow = ([3,30,1e-14], [4,32,0], [5,20,1e-14], [2.5, 10,,1e-14], [-2, 69,0]); +my $tests; +$tests += $_->[1] foreach @pow; + +plan tests => 1 + $bits_in_uv + $tests; + +# Ought to be 32, 64, 36 or something like that. + +my $remainder = $bits_in_uv & 3; + +cmp_ok ($remainder, '==', 0, 'Sanity check bits in UV calculation') + or printf "# ~0 is %d (0x%d) which gives $bits_in_uv bits\n", ~0, ~0; + +# These are a lot of brute force tests to see how accurate $m ** $n is. +# Unfortunately rather a lot of perl programs expect 2 ** $n to be integer +# perfect, forgetting that it's a call to floating point pow() which never +# claims to deliver perfection. +foreach my $n (0..$bits_in_uv - 1) { + my $exp = 2 ** $n; + my $int = 1 << $n; + cmp_ok ($exp, '==', $int, "2 ** $n vs 1 << $n"); +} + +foreach my $pow (@pow) { + my ($base, $max, $range) = @$pow; + my $fp = 1; + foreach my $n (0..$max-1) { + my $exp = $base ** $n; + within ($exp, $fp, $range, "$base ** $n [$exp] vs $base * $base * ..."); + $fp *= $base; + } +} diff --git a/t/test.pl b/t/test.pl index 351963a..91daf1a 100644 --- a/t/test.pl +++ b/t/test.pl @@ -148,6 +148,70 @@ sub isnt { _ok($pass, _where(), $name, @mess); } +sub cmp_ok { + my($got, $type, $expected, $name, @mess) = @_; + + my $pass; + { + local $^W = 0; + local($@,$!); # don't interfere with $@ + # eval() sometimes resets $! + $pass = eval "\$got $type \$expected"; + } + unless ($pass) { + # It seems Irix long doubles can have 2147483648 and 2147483648 + # that stringify to the same thing but are acutally numerically + # different. Display the numbers if $type isn't a string operator, + # and the numbers are stringwise the same. + # (all string operators have alphabetic names, so tr/a-z// is true) + # This will also show numbers for some uneeded cases, but will + # definately be helpful for things such as == and <= that fail + if ($got eq $expected and $type !~ tr/a-z//) { + unshift @mess, "# $got - $expected = " . ($got - $expected) . "\n"; + } + unshift(@mess, "# got "._q($got)."\n", + "# expected $type "._q($expected)."\n"); + } + _ok($pass, _where(), $name, @mess); +} + +# Check that $got is within $range of $expected +# if $range is 0, then check it's exact +# else if $expected is 0, then $range is an absolute value +# otherwise $range is a fractional error. +# Here $range must be numeric, >= 0 +# Non numeric ranges might be a useful future extension. (eg %) +sub within { + my ($got, $expected, $range, $name, @mess) = @_; + my $pass; + if (!defined $got or !defined $expected or !defined $range) { + # This is a fail, but doesn't need extra diagnostics + } elsif ($got !~ tr/0-9// or $expected !~ tr/0-9// or $range !~ tr/0-9//) { + # This is a fail + unshift @mess, "# got, expected and range must be numeric\n"; + } elsif ($range < 0) { + # This is also a fail + unshift @mess, "# range must not be negative\n"; + } elsif ($range == 0) { + # Within 0 is == + $pass = $got == $expected; + } elsif ($expected == 0) { + # If expected is 0, treat range as absolute + $pass = ($got <= $range) && ($got >= - $range); + } else { + my $diff = $got - $expected; + $pass = abs ($diff / $expected) < $range; + } + unless ($pass) { + if ($got eq $expected) { + unshift @mess, "# $got - $expected = " . ($got - $expected) . "\n"; + } + unshift@mess, "# got "._q($got)."\n", + "# expected "._q($expected)." (within "._q($range).")\n"; + } + _ok($pass, _where(), $name, @mess); +} + # Note: this isn't quite as fancy as Test::More::like(). sub like { my ($got, $expected, $name, @mess) = @_;