in reply to what did I just see..?

No one seems to have actually told you how to fix this. Basically if you are subtracting two floating point numbers (which you are) and then rounding (which is what print does - downward) then the upper bound of the error in the result will be twice the machines epsilon. So to fix this your code need to be this:

use Machine::Epsilon; for (my $x=0.8; $x > 0.1; $x -= 0.01) { print "".($x+(2*machine_epsilon()))."\n"; }

Replies are listed 'Best First'.
Re^2: what did I just see..?
by pryrt (Abbot) on Mar 22, 2021 at 17:37 UTC
    Sorry, sectokia, but I don't think that adding twice the machine_epsilon is the cure-all you really think it is.

    Per the documentation, machine_epsilon is the maximum relative error ("Machine::Epsilon - The maximum relative error while rounding a floating point number"), but you are using it in an absolute fashion, which means that you are not properly scaling the error. The error possible is different for 1_000_000 vs 1 vs 0.000_001... it's even different for 2 vs 1 vs 1/2. So that won't help you "correct" the cumulative floating point errors. Since you are using values from 0.8 down to 0.1, you will be in the ranges [0.5,1.0), [0.25,0.50), [0.125,0.25), and [0.0625,0.125), which actually have four different ULP sizes.

    Besides that, the cumulative errors will keep getting bigger, until your x+2*epsilon is no longer enough to increase it. Starting at 0.8 and decreasing by 0.01 each loop, by 70 iterations, the stringification of the "real" number and the stringification of your x+2*epsilon will not be matching each other.

    The example in the spoiler shows both of those issues in more detail (where "diff" represents the real number, start - n*step, "x" represents the multiple individual subtractions, and "sectokia" represents the x+2epsilon).

    #!perl use 5.012; # strict, // use warnings; use Machine::Epsilon; use Data::IEEE754::Tools qw/ulp/; use Test::More; sub dprint { my ($name, $value) = @_; diag sprintf '%-16s => %24.17e', + $name, $value } my $eps = machine_epsilon(); # note the docs say that it's relative +, _not_ absolute my $ulp1 = ulp(1); # ULP = https://en.wikipedia.org/wiki/ +Unit_in_the_last_place is $eps, $ulp1, 'machine epsilon is same as ULP(1)'; dprint eps => $eps; dprint 'ulp(1)' => $ulp1; isnt $eps, ulp(2), 'machine epsilon is not same as ULP(2)'; dprint 'ulp(2)'=> ulp(2); is $eps*2, ulp(2), 'machine epsilon scaled by 2 is actually the same a +s ULP(2)'; dprint '2*eps', 2*$eps; isnt $eps*2.1, ulp(2.1), 'machine epsilon is _not_ scaled by 2.1 to gi +ve the same as ULP(2)'; dprint '2.1*eps', 2.1*$eps; dprint 'ulp(2.1)' => ulp(2.1); isnt $eps, ulp(1/10), 'machine epsilon is not the same as ULP(1/10)'; dprint 'ulp(1/10)', ulp(1/10); diag 'using machine epsilon does not give the results you seem to thin +k'; my $start = 1e5; my $x = $start; $x -= (1/10) for 1..1e6; # one million subtractions of 0.1 my $diff = $start - 1e6*(1/10); # subtract one million times the fl +oating point 0.1 isnt $x, $diff, '1 million single steps is not the same as removing on +e million times the step-value'; dprint 'x' => $x; dprint 'diff' => $diff; my $sectokia = $x + 2*$eps; # this is the fix that sectokia ap +plied isnt $sectokia, $diff; dprint 'sectokia' => $sectokia; diag "stringify sectokia = ".($sectokia)."\n"; diag "in case you don't believe mine is equivalent:"; $start = 0.8; $_ = 0; $x = $start; $sectokia = $x + 2*$eps; # this is adding the two epsilons... $diff = $start - $_ * 0.01; diag sprintf "loop %4d:\n\tx = %24.17e = '%s',\n\tsecto +kia = %24.17e = '%s',\n\tdiff = %24.17e = '%s'", $_, $x, $x, $sec +tokia, $sectokia, $diff, $diff; for (1..1000) { $x -= 0.01; # same as your $x -= 0.01 in the for l +oop $sectokia = $x + 2*$eps; # this is adding the two epsilons... $diff = $start - $_ * 0.01; if( $_ < 10 or 0 == $_ % 100 ) { diag sprintf "loop %4d:\n\tx = %24.17e = '%s',\n\tsecto +kia = %24.17e = '%s',\n\tdiff = %24.17e = '%s'", $_, $x, $x, $sec +tokia, $sectokia, $diff, $diff; } if( "$sectokia" ne "$diff" ) { # if the stringified diag sprintf "loop %4d:\n\tx = %24.17e = '%s',\n\tsecto +kia = %24.17e = '%s',\n\tdiff = %24.17e = '%s'", $_, $x, $x, $sec +tokia, $sectokia, $diff, $diff; diag "mismatch ends loop at $_\n\n"; last; } } done_testing;

    edit: add my output from above code

    ok 1 - machine epsilon is same as ULP(1) # eps => 2.22044604925031308e-16 # ulp(1) => 2.22044604925031308e-16 ok 2 - machine epsilon is not same as ULP(2) # ulp(2) => 4.44089209850062616e-16 ok 3 - machine epsilon scaled by 2 is actually the same as ULP(2) # 2*eps => 4.44089209850062616e-16 ok 4 - machine epsilon is _not_ scaled by 2.1 to give the same as ULP( +2) # 2.1*eps => 4.66293670342565767e-16 # ulp(2.1) => 4.44089209850062616e-16 ok 5 - machine epsilon is not the same as ULP(1/10) # ulp(1/10) => 1.38777878078144568e-17 # using machine epsilon does not give the results you seem to think ok 6 - 1 million single steps is not the same as removing one million +times the step-value # x => -1.33288113454699264e-06 # diff => 0.00000000000000000e+00 ok 7 # sectokia => -1.33288113410290343e-06 # stringify sectokia = -1.3328811341029e-06 # in case you don't believe mine is equivalent: # loop 0: # x = 8.00000000000000044e-01 = '0.8', # sectokia = 8.00000000000000488e-01 = '0.8', # diff = 8.00000000000000044e-01 = '0.8' # loop 1: # x = 7.90000000000000036e-01 = '0.79', # sectokia = 7.90000000000000480e-01 = '0.79', # diff = 7.90000000000000036e-01 = '0.79' # loop 2: # x = 7.80000000000000027e-01 = '0.78', # sectokia = 7.80000000000000471e-01 = '0.78', # diff = 7.80000000000000027e-01 = '0.78' # loop 3: # x = 7.70000000000000018e-01 = '0.77', # sectokia = 7.70000000000000462e-01 = '0.77', # diff = 7.70000000000000018e-01 = '0.77' # loop 4: # x = 7.60000000000000009e-01 = '0.76', # sectokia = 7.60000000000000453e-01 = '0.76', # diff = 7.60000000000000009e-01 = '0.76' # loop 5: # x = 7.50000000000000000e-01 = '0.75', # sectokia = 7.50000000000000444e-01 = '0.75', # diff = 7.50000000000000000e-01 = '0.75' # loop 6: # x = 7.39999999999999991e-01 = '0.74', # sectokia = 7.40000000000000435e-01 = '0.74', # diff = 7.39999999999999991e-01 = '0.74' # loop 7: # x = 7.29999999999999982e-01 = '0.73', # sectokia = 7.30000000000000426e-01 = '0.73', # diff = 7.29999999999999982e-01 = '0.73' # loop 8: # x = 7.19999999999999973e-01 = '0.72', # sectokia = 7.20000000000000417e-01 = '0.72', # diff = 7.20000000000000084e-01 = '0.72' # loop 9: # x = 7.09999999999999964e-01 = '0.71', # sectokia = 7.10000000000000409e-01 = '0.71', # diff = 7.10000000000000075e-01 = '0.71' # loop 70: # x = 9.99999999999994643e-02 = '0.0999999999999995', # sectokia = 9.99999999999999084e-02 = '0.0999999999999999', # diff = 9.99999999999999778e-02 = '0.1' # mismatch ends loop at 70 # 1..7
    /edit

    Further, your x+2*epsilon assumes that your subtraction step size is slightly bigger than the real value; that is true for a step of 0.1 (1.00000000000000006e-01) or 0.01 (1.00000000000000002e-02), but for a step of 0.03 (2.99999999999999989e-02), the step size is smaller, so now your adjustment takes the string in the wrong direction. (Not shown in the spoiler code.)

      For the range 0 to 1, then epsilon will always be greater than the error (as it would only scale smaller), but of course you are correct in that it should be scaled both up and down for a normalized solution.

      You are also right about needing to apply it to each subtraction operation.

      I don't agree with the bit around it being in the wrong direction if the step happens to be just under the desired ideal value. Print rounds down. If the float is +/- epsilon from ideal, then adding an epsiol brings it into range of 0 to +2 epsilon from ideal, which will round down to ideal. It doesn't matter if you started +ve or -ve from ideal.

        For the range 0 to 1, then epsilon will always be greater than the error

        So while we have the example of 0.8 down to 0.1, the difference between the epsilon and the ULP won't be huge. But who knows whether the OP will eventually go beyond that range, and get upset when numbers near 2e-16 start displaying as near 4e-16. That's one of the reasons I was cautioning against applying epsilon in this case, because it might later be generalized to a condition where it's not even close to the real absolute error.

        Print rounds down

        that's not true, as syphilis showed. Here's another example.

        C:\usr\local\share>perl -le "printf qq(%.17e => %s\n), $_, $_ for 2.99 +999999999999989e-02, 2.99999999999999503e-02, 2.99999999999999468e-02 +" 2.99999999999999989e-02 => 0.03 2.99999999999999503e-02 => 0.03 2.99999999999999468e-02 => 0.0299999999999999

        Print rounds to nearest to the precision that print "likes" rounding to.

        If the float is +/- epsilon from ideal, then adding an epsiol brings it into range of 0 to +2 epsilon from ideal, which will round down to ideal.

        Again, no.

        C:\usr\local\share>perl -MMachine::Epsilon -le "for (2.999999999999999 +89e-02, 2.99999999999999503e-02, 2.99999999999999468e-02) { printf qq +(%.17e => %-30s x+2eps = %s\n), $_, $_, $_+2*machine_epsilon}" 2.99999999999999989e-02 => 0.03 x+2eps = 0.0 +300000000000004 2.99999999999999503e-02 => 0.03 x+2eps = 0.0 +300000000000004 2.99999999999999468e-02 => 0.0299999999999999 x+2eps = 0.0 +300000000000004

        Notice that x+2epsilon prints a number bigger than 0.03, whereas the value of x is slightly less than 0.03 (by 1 ULP).

        Print rounds down

        As a generalization this is not true.
        Even when it is true for some value $x, it will be untrue for -$x.
        However, there are times when print() rounds up for positive values.

        Consider the following (perl-5.32.0, nvtype is "double"):
        use strict; use warnings; my $d = 0.4444444444444446; printf "%a\n", $d; print "$d\n"; print "print() has rounded upwards\n" if "$d" > $d; __END__ Outputs: 0x1.c71c71c71c71fp-2 0.444444444444445 print() has rounded upwards
        And how do I calculate the value of this "epsilon" that is being mentioned ?

        Cheers,
        Rob
Re^2: what did I just see..?
by LanX (Saint) on Mar 22, 2021 at 02:12 UTC
    > No one seems to have actually told you how to fix this.

    I did, I said calculate in cents if you want that precision in the end

    Even when using one single division to a float in the final output, it'll print correctly.

    DB<3> for ($x=80; $x > 10; $x -= 1) { say $x/100; } 0.8 0.79 0.78 0.77 0.76 0.75 0.74 0.73 0.72 0.71 0.7 ... yadda 0.21 0.2 0.19 0.18 0.17 0.16 0.15 0.14 0.13 0.12 0.11 DB<4>

    Cheers Rolf
    (addicted to the Perl Programming Language :)
    Wikisyntax for the Monastery

      Thats not a solution, its sticking your head in the sand. You'll end up coming unstuck at certain values. If are really are going to use cents, use it, and do integer division.
        > You'll end up coming unstuck at certain values

        Show me!

        Cheers Rolf
        (addicted to the Perl Programming Language :)
        Wikisyntax for the Monastery