frankcox has asked for the wisdom of the Perl Monks concerning the following question:

I recently added a new test to Geo::ReadGRIB and started getting fails from CPAN Testers only on 64bit Perls where nvtype='long double' and nvsize= 12 or 16.

I found the problem is in a floating point calculation for some values of the variables. The problem code fragment is:

$out = (( ($lat - $self->La1) / $self->LaInc ) * $self->Ni)
+ (($thislong ) / $self->LoInc);

return sprintf "%d", $out;

[download]

This will work for most values but as an example of the problem, the following test will fail on 64bit uselongdouble Perls:

my $calc = sprintf "%d",(((63 - -90)/.6) * 601) + (360 /.6);
ok( $calc == 153855 ) or
    diag ("((63 - -90)/.6) * 601) + (360 /.6) = 153855 not $calc ");

[download]

This test will pass on any other CPAN Tester platform and 153855 is what I get if work this by hand.

My module will not work correctly if this code returns the wrong value. Is there any way I can modify the code to get the expected result on any version of Perl?

Another option is to detect uselongdouble in Makefile.PL and not build on the problem systems. I'd rather not do that.

Thank you for considering my problem.

Replies are listed 'Best First'.
Re: test fails on 64bit uselongdouble Perl
by ikegami (Patriarch) on Oct 29, 2009 at 19:04 UTC

    Just like some numbers are periodic in decimal (e.g. 1/3), some numbers are periodic in binary (e.g. 6/10). They cannot be represented exactly using a float since that would require infinite storage.

    You shouldn't check if two floats are equal. You should check if two floats are sufficiently identical. 

    my $calc = (((63 - -90)/.6) * 601) + (360 /.6);

ok( abs( $calc - 153855 ) < 0.000000001 )
    or diag("Expected 153855. Got $calc");

    [download]
[reply]
[d/l]

      Ahh! Good explanation of issue with .6 ikegami. I missed that

      My real problem is that I get different results with different Perl variants. I need my algorithm to work the same on all Perls.

[reply]
        Then you have to avoid all floating point numbers. But seriously, that's not a worthy goal. If your algorithms are numerically stable, the errors won't be larger than 1e-18 or so - does that really matter for your purposes?
[reply]
[d/l]

        You can't use floats if differences are a problem. Like I said, with floats, you have to be satisfied with no differences within a tolerance.

[reply]
Re: test fails on 64bit uselongdouble Perl
by moritz (Cardinal) on Oct 29, 2009 at 17:53 UTC
    I can think of several options:
    • Use explicit rounding instead of the implicit rounding that printf uses (something like int($num + 0.5))
    • Re-work your code to use integers instead of floating-point numbers
    • use rationals, for example by including bigrat
    Perl 6 - links to (nearly) everything that is Perl 6.
[reply]
[d/l]

      Thanks for the ideas moritz.

      I should have mentioned that this specific test also fails without the sprintf. (And for these values the divisions should yield whole number results.) Also, for my application there will be floating-point numbers as part of the input so it would be hard to eliminate them in this calculation

      I'll look at bigrat, thanks.

[reply]

        well, I should have tried bigrat before I replied...

        If I just 'use bigrat' I get this during make test:

        Deep recursion on subroutine "Math::BigInt::bsub"

        And then, after a very long time:

        Out of memory!

        It's possible some Math::* module will help but I don't know which one since I don't understand the underlying problem.

        One question I have is why uselongdouble gives less accurate results in this case?

        And also, shouldn't all versions of Perl return the same result for a relatively simple calculation?

[reply]
[d/l]
[select]
Re: test fails on 64bit uselongdouble Perl
by frankcox (Acolyte) on Oct 30, 2009 at 01:06 UTC

    This discussion has helped me work out this issue in my mind. I'll be reworking my algorithm and I have a better idea now about how.

    What I need to do is take a latitude/longitude pair, which may include fractions of a degree, and return an integer index into a binary data structure in a GRIB weather data file. I've been worrying about rounding to the nearest int but really I need to round to the nearest data point. The data points, projected onto a lat/long grid, are centered in a rectangle LaInc tall and LoInc wide (from my code above). These increments can also include fractions. All, I have to do is figure out which rectangle a given lat/long pair is in and return a number based on the scanning order for the type of file it is. (The above code works for west-east and south-north scanning.)

    uselongdouble really doesn't have much to do with this except that that's what whacked me up side a' the head. So, thanks for that, and thanks CPAN Testers!

[reply]

      I tested the advice of ikegami here and changed my code to round floating-point numbers rather than truncate as I was inadvertently doing. That is, I changed:

      return sprintf "%d", $out; # Wrong, truncates.

      [download]

      to:

      return sprintf "%.0f", sprintf "%.6f", $out;

      [download]

      With this change I now see consistent results on all Perl platforms including 64bit uselongdouble ones. Thanks to everyone who helped out here! I apologize for being a little slow to see the light but the end result is good.

      I still plan to rework my algorithm but now I can do that as part of a regular release rather than as part of a bug fix.

[reply]
[d/l]
[select]

Back to Seekers of Perl Wisdom