Beefy Boxes and Bandwidth Generously Provided by pair Networks
Problems? Is your data what you think it is?
 
PerlMonks  

Re^2: Portability of floor(log(N))

by toolic (Bishop)
on Mar 13, 2018 at 14:46 UTC ( #1210827=note: print w/replies, xml ) Need Help??


in reply to Re: Portability of floor(log(N))
in thread Portability of floor(log(N))

Thank you for the suggestions.

I agree with you that, before I pester the tester, I should try to have a pretty good plan. Using diag for all the intermediate values is a great idea. I was unaware of %a. Now all I have to do is play with it until I have a grasp at what the hex fp format means. There's a ton I don't understand.

I like the approach of releasing a new debug version of the module, too. However, in this case, it may not be practical since it could be months before I get any feedback (it did take 4 months to land on this one machine). I'd feel like an astronaut in one of those sci-fi movies transmitting a message back to earth and waiting 6 months for a response :)

This is all more involved than I was expecting (and hoping). Creating portable software is hard! Based on your reply and that of syphilis as well, I'm struggling to decide whether my code has a fundamental flaw or whether this one machine is just badly broken in some way. Is my code bound to fail on other machines given enough variety in input values? Since the input is the set of all real numbers, it's impossible to exhaustively test it. I should try to test more thoroughly (but intelligently).

Or, is there a way to change my code to avoid this error? Can I pre-process the value, say with something as simple as sprintf? Or, do I need to use one of syphilis's suggestions, like Data::Float? Or, can I post-process the value? And, if I do make this type of change, will it break backward-compatibility in a significant way? Or, will it be minor enough to justify the change? Lot's of questions for me to mull over.

Replies are listed 'Best First'.
Re^3: Portability of floor(log(N))
by pryrt (Monsignor) on Mar 13, 2018 at 15:57 UTC

    I guess I'd say the next steps depend on whether it bothers you (or your end users) to get 1000e-36 instead of 1e-33. If it does, post-processing is probably your best bet -- but I wouldn't use sprintf. I'd look at the number of digits before the "e" (or the decimal point), and if it's four digits, I'd divide by 1000 (and re-round, if necessary) and adjust the exponent. syphilis's suggestion of generating the 1e-33 using Data::Float's hex_float() would work to "fix" your test -- make the test pass; but if you or your end user would consider it wrong to ever see 1000e-36, then forcing a "correct" 1e-33 won't show the actual error.

    I was able to make it show 1000e-36 by running

    C:\WINDOWS\system32>perl -MNumber::FormatEng=format_pref -le "print fo +rmat_pref($_) for 0.9999999994e-33, 0.9999999995e-33, 0.9999999996e-3 +3" 999.999999e-36 999.999999e-36 1000e-36
    . Your code $num = 0 + sprintf '%.6f', ($num / $mult); only allows 6 digits after the decimal, and you don't do any more checks after that, so numbers that round from 999.999999x to 1000.000000 are going to show up wrong. In your library, you need to fix that. In your test suite, you need to make sure you test the edge cases of what your code does with that (for example, the %.6f introduces edge cases) -- you don't have to test every representable floating value in binary32, binary64, and binary128, just the edge cases that your specific code introduces. I'd probably do a simple test after that line that if(abs($num)>=1000) { ... $num /= 1000, $e++ } or similar. <code>
Re^3: Portability of floor(log(N))
by syphilis (Bishop) on Mar 13, 2018 at 22:50 UTC
    Thank you for the suggestions.

    I think they were completely irrelevant to the problem you're facing.
    Silly me didn't think to actually check what the correct calculation of floor( log($num) / log(1000) ) should yield under 'long double' (64-bit) precision.

    When I do check, I find that the correct answer is -12.
    The error is in your test script ... either that, or in Math::MPFR:
    use strict; use warnings; use Math::MPFR qw(:mpfr); # Set precision for extended # precision 'long double' Rmpfr_set_default_prec(64); my $num = Math::MPFR->new('1e-33'); my $thou = Math::MPFR->new('1000'); # Convert $num & $thou to their # respective log() Rmpfr_log($num, $num, MPFR_RNDN); Rmpfr_log($thou, $thou, MPFR_RNDN); print $num / $thou; # prints -1.10000000000000000009e1
    This indicates that when the floating point precision is 64 bits, the correct floor() result is -12.
    And, IIUC, the test script should accept '1000e-36' as correct when NV precision is 64 bits (as is the case with those 2 FAIL reports).

    Interestingly, with full quad precision of 113 bits, the calculation reverts to returning -11.

    Cheers,
    Rob
      And, IIUC, the test script should accept '1000e-36' as correct when NV precision is 64 bits

      Alternatively (and since you're already using POSIX) you should get the result you expected, if you replace the log() calls with POSIX::log10(). That way both numerator and divisor should be exact (for these particular inputs), as should also be the result of the division ... irrespective of the precision of the NV:
      $ perl -MPOSIX -le 'print POSIX::floor(POSIX::log10(1e-33) / POSIX::lo +g10(1000));' -11
      Cheers,
      Rob
        Using POSIX::log10() looks like the simplest solution. I will attempt to pester the tester with a modified version of FormatEng.pm.

        Do you think I can simplify it by replacing POSIX::log10(1000) with 3, or will it be better to keep your code as-is?

        my $e = floor( POSIX::log10($num) / 3 );

Log In?
Username:
Password:

What's my password?
Create A New User
Node Status?
node history
Node Type: note [id://1210827]
help
Chatterbox?
and the web crawler heard nothing...

How do I use this? | Other CB clients
Other Users?
Others examining the Monastery: (5)
As of 2020-10-29 03:01 GMT
Sections?
Information?
Find Nodes?
Leftovers?
    Voting Booth?
    My favourite web site is:












    Results (267 votes). Check out past polls.

    Notices?