Issue with process = (1.0 + ma.EPSILON/2.0) > 1.0

[sletz]

The code process = (1.0 + ma.EPSILON/2.0) > 1.0 correctly returns 0 when compiled in -double but incorrectly returns 1 when compiled in -single.

Even if ma.EPSILON definition is dependant of -single/-double/-quad (see https://github.com/grame-cncm/faustlibraries/blob/master/maths.lib#L157), the internal computation is actually using double and has to be fixed.

[oleg-nesterov]

It is clear why -single/-double makes a difference in this case, but...

To be honest, I never really understood the purpose of ma.EPSILON. And what is it? Say,

    doubleprecision EPSILON = 2.2204460492503131e-016;

is it DBL_MIN? Or DBL_TRUE_MIN?

Nevermind. In any case I'd say it should be used with care.

process = (1.0 + ma.EPSILON/2.0) > 1.0 correctly returns 0 when compiled in -double

OK, lets suppose that "0" is correct in this case. Now consider

    import("stdfaust.lib");

    process =
            (0.0 + ma.EPSILON/2.0) > 0.0,
            (0.5 + ma.EPSILON/2.0) > 0.5,
            (1.0 + ma.EPSILON/2.0) > 1.0;

compiled in double, it outputs

    1, 1, 0

Doesn't look consistent but what can we do? and do we really care?

Oleg.

[sletz]

• see: https://en.cppreference.com/w/cpp/types/climits
• as in C:

#include <iostream>
#include <cfloat>

int main() {
    bool a = (0.0f + FLT_EPSILON / 2.0f) > 0.f;
    bool b = (0.5f + FLT_EPSILON / 2.0f) > 0.5f;
    bool c = (1.f + FLT_EPSILON / 2.0f) > 1.f;


    std::cout << "a: " << a << std::endl;
    std::cout << "b: " << b << std::endl;
    std::cout << "c: " << c << std::endl;
    
    bool d = (0.0f + DBL_EPSILON / 2.0f) > 0.f;
    bool e = (0.5f + DBL_EPSILON / 2.0f) > 0.5f;
    bool f = (1.f + DBL_EPSILON / 2.0f) > 1.f;

    std::cout << "d: " << d << std::endl;
    std::cout << "e: " << e << std::endl;
    std::cout << "f: " << f << std::endl;

    return 0;
}

[dariosanfilippo]

On Sun, 29 Oct 2023 at 22:47, oleg-nesterov @.***> wrote:

It is clear why -single/-double makes a difference in this case, but...

To be honest, I never really understood the purpose of ma.EPSILON. And what is it? Say,

doubleprecision EPSILON = 2.2204460492503131e-016;

is it DBL_MIN? Or DBL_TRUE_MIN?

I think that it is defined as the smallest quantity, let's call it EPS, for which (1.0 + EPS) > 1.0 returns true.

I'm sure that you know issues related to float comparison and numerical errors better than, so what I think it's useful for is comparisons: you can't do (.1 + .2) == .3, it will return false, so you'd normally do abs((.1 + .2) - .3) < EPS.

Nevermind. In any case I'd say it should be used with care.

process = (1.0 + ma.EPSILON/2.0) > 1.0 correctly returns 0 when compiled in -double

OK, lets suppose that "0" is correct in this case. Now consider

import("stdfaust.lib");

process =
        (0.0 + ma.EPSILON/2.0) > 0.0,
        (0.5 + ma.EPSILON/2.0) > 0.5,
        (1.0 + ma.EPSILON/2.0) > 1.0;

compiled in double, it outputs

1, 1, 0

Doesn't look consistent but what can we do? and do we really care?

I think that a common way to make EPS more useful for comparisons is to scale it down by the max abs of the two terms, but it's still not the most precise way of comparing floats. If you did that, you'd get 0, 1, 1, which is still not perfect but a but better.

I don't think that this bug is a big deal at all; I only noticed it because I was working in double and I could still hear artefacts that made me think of noise from numerical errors. So I tested that to make sure that FaustLive was taking the -double flag correctly.

Dario

Oleg.

— Reply to this email directly, view it on GitHub https://github.com/grame-cncm/faust/issues/962#issuecomment-1784233886, or unsubscribe https://github.com/notifications/unsubscribe-auth/AHG3I2BKQBSTITZHJAJWQH3YB26ENAVCNFSM6AAAAAA6U5OUUCVHI2DSMVQWIX3LMV43OSLTON2WKQ3PNVWWK3TUHMYTOOBUGIZTGOBYGY . You are receiving this because you are subscribed to this thread.Message ID: @.***>

[oleg-nesterov]

• see: https://en.cppreference.com/w/cpp/types/climits

  • as in C:

Ah. Yes, yes, I know. And I too wrote the similar program to verify.

Sorry for confusion. When I said "Doesn't look consistent" I meant "Doesn't look consistent but in fact correct".

My only point was, again, in any case ma.EPSILON should be used with care. It is not immediately clear that, again,

    process = (0.5 + ma.EPSILON/2.0) > 0.5;

and process = (1.0 + ma.EPSILON/2.0) > 1.0;

do not output the same result.

And I think we do not really care.

But I agree that the fact it depends on -single/-double make the things worse.

[oleg-nesterov]

I think that it is defined as the smallest quantity, let's call it EPS, for which (1.0 + EPS) > 1.0 returns true.

OK, quite possible, then ma.EPSILON == nextafter(1.0,2.0) - 1.0

And this means that even a more simple example can surprise the user:

    import("stdfaust.lib");

    process =
            (1.0 + ma.EPSILON) > 1.0,
            (2.0 + ma.EPSILON) > 2.0;

compiled with -double outputs

    output0[i0] = FAUSTFLOAT(1);
    output1[i0] = FAUSTFLOAT(0);

Again, this is correct in that this is how the FP math works.

But again, this means that ma.EPSILON should be used with care. And "the smallest quantity" depends on the usage, while *.lib blindly use ma.ma.EPSILON assuming it always fits the needs.

[sletz]

I guess we want ma.EPSILON to exactly behave as the FP norm says, in -single and -double mode, no more no less (and which is not currently the case).

[sletz]

Improved documentation: https://faustlibraries.grame.fr/libs/maths/#maepsilon