boost::float_next mittels Lookup-Table?



  • Hallo,

    in C++ habe ich ja (zum Beispiel) mittels boost::math::float_next die Möglichkeit mir sehr einfach den nächsten repräsentierbaren Wert eines floats zu holen. Genau diese Funktionalität brauche ich jetzt in einer Anwendung die ich schreibe, dort habe ich aber keinen Zugriff auf diese Funktion.

    Hintergrund: Ich hab ein (closed-source) Program welches in C++ geschrieben wurde. Dieses Program bietet über eine symbolische Schaltbild-Sprache die Möglichkeit an, Abläufe zu automatisieren. Die Möglichkeiten dieser "Programmiersprache" sind ziemlich eingeschränkt.

    Unter anderem gibt es den Datentyp "float" in dieser Umgebung, welcher direkt einem C++ float entspricht. Ich brauche jetzt in dieser Umgebung eine Funktion wie float_next , kann aber keine C/C++ Funktionen einbinden. Daher würde ich jetzt versuchen mir diese Funktion selbst zu schreiben.

    Meine Idee wäre, eine Art Lookup Tabelle zu implementieren, in die je nach Größe der aktuellen Zahl, die minimale Differenz zur nächsten Zahl festgelegt ist. Also vom Prinzip her:

    float lookup(float f)
    {
        if (f >= 2 && f <= 4)
            return 0.0000002;
        else if (...)
            // usw für alle relevanten Bereiche
    }
    
    float next_float(float f)
    {
        float min_diff = lookup(f);
        return f + min_diff; // Nächste float Zahl
    }
    

    Jetzt meine Frage, ist das so praktikabel? Oder gibt es viel zu viele Genauigkeitsbereiche für eine Lookup Table? Und gibt es irgendwo eine Tabelle wo man vielleicht sehen kann, welche minimale Differenz in welchem Bereich für float in C++ gültig ist? Oder gibts eine noch viel bessere Methode?



  • Also es gibt maximal 2^32, ist also praktisch schon eine Moeglichkeit. Wenn auch nicht besonders elegant.



  • TGGC schrieb:

    Also es gibt maximal 2^32, ist also praktisch schon eine Moeglichkeit. Wenn auch nicht besonders elegant.

    Naja, das wär schon übel die alle manuell zu schreiben.

    Es müssten ja auch wesentlich weniger sein, weil die minimale Differenz ja immer für bestimmte feste Bereiche gilt und diese Bereiche ja exponentiell wachsen. Also etwa so:

    Bereich       Minimale Differenz
    --------------------------------
    [1, 2]        1.1920929e-7
    [2, 4]        2.3841858e-7
    [4, 8]        4.7683716e-7
    [8, 16]       9.5367432e-7
    etc.
    

    Also wären es ja nur 2^8 == 256 Bereiche die abgedeckt werden müssten, ist das korrekt?

    Der Teil ist leider Performance-kritisch, deswegen suche ich nach einer möglichst schnellen Variante.



  • Deine Bereiche kannst du ja sehr simpel mit boost::math::float_next ueberpruefen, indem du die 2^32 Moeglichkeiten durchgehst und entsprechend zusammenfasst.



  • Die Boost-Funktion sieht mir eigentlich nicht so kompliziert aus, als dass ich einen Lookup-Table gerechtfertigt halte:
    https://github.com/boostorg/math/blob/master/include/boost/math/special_functions/next.hpp#L108
    Ich glaube man hat sich schneller in die nötigen Details zu Floating-Point-Zahlen eingearbeitet und das selbst implementiert, als dass man sich eine LUT-Lösung zusammenbastelt,
    welche durch ihre Speicherzugriffe wahrscheinlich auch noch nicht besonders effizient sein wird.

    Mein Tip: float_next selbst implementieren da wahrscheinlich

    - schneller fertig
    - effizienter
    - man weiss nachher einiges mehr über Floating-Point-Zahlen (kann nie schaden!) 😃

    Finnegan


  • Mod

    Die nächste repräsentierbare Fließkommazahl ist, zumindest mit IEEE 754, trivial zu holen: Dort sind Floats nämlich schon lexikographisch sortiert, man braucht nur mit einem integralen Typ zu aliasen:

    static_assert( sizeof(float) == sizeof(std::uint32_t) );
    ++reinterpret_cast<std::uint32_t&>(f);
    

    Alternativ kann man einfach std::nexttowards verwenden.



  • Arcoth schrieb:

    Die nächste repräsentierbare Fließkommazahl ist, zumindest mit IEEE 754, trivial zu holen: Dort sind Floats nämlich schon lexikographisch sortiert, man braucht nur mit einem integralen Typ zu aliasen:

    static_assert( sizeof(float) == sizeof(std::uint32_t) );
    ++reinterpret_cast<std::uint32_t&>(f);
    

    Alternativ kann man einfach std::nexttowards verwenden.

    Das Problem ist, ich muss mit der Funktionalität der Schaltbild-Sprache auskommen. Diese erlaubt zwar das definieren von eigenen Blöcken (== Funktionen) und hat Konstrukte wie if , for , while oder goto sowie Variablen und ein paar simple Funktionen, ist aber ansonsten sehr eingeschränkt. Etwas vergleichbares zu reinterpret_cast gibt es z.B. nicht, weswegen das leider schonmal wegfällt.

    Finnegan schrieb:

    Mein Tip: float_next selbst implementieren

    Hm, wenn das nur mit den oben genannten Mitteln auskommt wäre das natürlich schön... muss mich mal durch den Source-Code arbeiten ob ich das so nachbauen kann...

    Auf jeden Fall verwendet das ja schonmal ldexp , was ja mit doubles arbeitet... ich habe aber nur float zur Verfügung, macht das was? Oder kann ich einfach mit einer float Variante von ldexp arbeiten ohne dass das Ergebnis beeinträchtigt wird?



  • Arcoth schrieb:

    Die nächste repräsentierbare Fließkommazahl ist, zumindest mit IEEE 754, trivial zu holen: Dort sind Floats nämlich schon lexikographisch sortiert, man braucht nur mit einem integralen Typ zu aliasen:

    static_assert( sizeof(float) == sizeof(std::uint32_t) );
    ++reinterpret_cast<std::uint32_t&>(f);
    

    Das ruft undefiniertes Verhalten hervor, weil es die "strikte Aliasingregel" verletzt. Genau denselben Bug findet man auch in dem viel zitiertem Quake3-Invers-Wurzel-Trick, der dann beim GCC ab -O2 nicht mehr das tut, was er soll.

    Außerdem erhöhst Du da den Betrag, womit Du Dich bei einer negativen Zahl Richtung Richting Minus Unendlich bewegst und aus Unendlich (positiv wie negativ) ein NAN machst. Da könnte man eine domain_error Ausnahme schmeißen. Den Fall +/- 0 müsste man auch noch abfangen. IMHO, sollte float_next(-0.0f) == float_next(0.0f) sein. std::numeric_limits<float>::is_iec559 könnte man noch mit in die static_assert ion einbauen, um sicherzustellen, dass man es tatsächlich mit IEEE754 floats zu tun hat.


  • Mod

    Das ruft undefiniertes Verhalten hervor, weil es die "strikte Aliasingregel" verletzt.

    Das ist klar, anders ist diese effiziente Methode aber (ohne Lib-Funktionen wie nexttoward oder builtins etc.) nicht zu haben. Natürlich -fno-strict-aliasing mitgeben.

    Außerdem erhöhst Du da den Betrag, womit Du Dich bei einer negativen Zahl Richtung Richting Minus Unendlich bewegst und aus Unendlich (positiv wie negativ) ein NAN machst. Da könnte man eine domain_error Ausnahme schmeißen. Den Fall +/- 0 müsste man auch noch abfangen. IMHO, sollte float_next(-0.0f) == float_next(0.0f) sein. std::numeric_limits<float>::is_iec559 könnte man noch mit in die static_assert ion einbauen, um sicherzustellen, dass man es tatsächlich mit IEEE754 floats zu tun hat.

    Ja, obiges war nur ein grobes PoC. Kann man alles noch verfeinern.


  • Mod

    Hier ist eine etwas ausgereiftere Skizze (die annimmt, das wenigstens isinf und isnan gegeben sind):

    #include <climits>
    #include <cmath>
    #include <cstdint>
    
    #include <limits>
    #include <stdexcept>
    #include <type_traits>
    
    template <typename T> struct identity {using type=T;};
    
    template <typename Float>
    Float next( Float f ) {
    	static_assert( std::numeric_limits<Float>::is_iec559, "" );
    	using IntAlias = typename std::conditional_t<sizeof(Float)*CHAR_BIT == 32, identity<std::uint32_t>,
    	                   std::conditional_t<sizeof(Float)*CHAR_BIT == 64, identity<std::uint64_t>,
    #ifdef __SIZEOF_INT128__
    	                     std::conditional_t<sizeof(Float)*CHAR_BIT == 128, identity<unsigned __int128>, std::enable_if<false>>
    #else
    	                     std::enable_if<false>
    #endif // defined
    	                  >>::type;
    
    	if (std::isinf(f) || std::isnan(f))
    		throw std::domain_error("");
    
    	auto& alias = reinterpret_cast<IntAlias&>(f);
    	static const auto signbit_mask = IntAlias(1) << (sizeof(alias)*CHAR_BIT - 1);
    	if (alias & signbit_mask) {
    		if (f != 0) {
    			--alias;
    			if (f == 0) // If negative zero has emerged…
    				f = 0;  // correct to positive one
    		}
    		else { // f is negative zero
    			alias &= ~signbit_mask;
    			++alias;
    		}
    	}
    	else
    		++alias;
    
    	return f;
    }
    

    Ich hab versucht irgendwo Attribute oder Pragmas reinzuhauen, die strict aliasing unterbinden; Jedoch scheinen diese nicht zu funktionieren. Daher ist das oben genannte Flag hier Pflicht.

    Edit: sizeof(alias) * CHAR_BIT-1 ist nicht die schlauste Verteilung von Leerraum…



  • Arcoth schrieb:

    Das ruft undefiniertes Verhalten hervor, weil es die "strikte Aliasingregel" verletzt.

    Das ist klar, anders ist diese effiziente Methode aber (ohne Lib-Funktionen wie nexttoward oder builtins etc.) nicht zu haben.

    🙂
    Doch, ist sie:

    #include <limits>
    #include <cstring>
    #include <cstdint>
    
    std::uint32_t reinterpret_float_as_uint(float x) {
        static_assert(sizeof(float) == sizeof(std::uint32_t),
                      "float/uint32_t incompatible");
        std::uint32_t result;
        std::memcpy(&result, &x, sizeof result);
        return result;
    }
    
    float reinterpret_uint_as_float(std::uint32_t x) {
        static_assert(sizeof(float) == sizeof(std::uint32_t),
                      "float/uint32_t incompatible");
        float result;
        std::memcpy(&result, &x, sizeof result);
        return result;
    }
    

    Im Kompilat wird bei eingeschalteter Optimierung auch nichts mehr von memcpy auftauchen. Der Compiler versteht, was da passiert und macht genau das, was man will ohne Overhead. Habe ich gerade per GCC mit "-O2" auf für die x86_64-Architektur ausprobiert und mir den Assemblercode angesehen.



  • Arcoth schrieb:

    Das ruft undefiniertes Verhalten hervor, weil es die "strikte Aliasingregel" verletzt.

    Das ist klar, anders ist diese effiziente Methode aber (ohne Lib-Funktionen wie nexttoward oder builtins etc.) nicht zu haben. Natürlich -fno-strict-aliasing mitgeben.

    Unsinn. Man kann das sehr wohl so machen, dass es funktioniert:

    std::uint32_t i;
    std::memcpy(&i, &f, sizeof(i));
    ++i;
    std::memcpy(&f, &i, sizeof(i));
    

    Man muss einfach nur hinschreiben, was man vom Compiler will. In C++. Nicht in "ich bin so super schlau"-Pseudo-Assembly-Style mit wilden Casts.


  • Mod

    TyRoXx schrieb:

    Arcoth schrieb:

    Das ruft undefiniertes Verhalten hervor, weil es die "strikte Aliasingregel" verletzt.

    Das ist klar, anders ist diese effiziente Methode aber (ohne Lib-Funktionen wie nexttoward oder builtins etc.) nicht zu haben. Natürlich -fno-strict-aliasing mitgeben.

    Unsinn. Man kann das sehr wohl so machen, dass es funktioniert:

    std::uint32_t i;
    std::memcpy(&i, &f, sizeof(i));
    ++i;
    std::memcpy(&f, &i, sizeof(i));
    

    Man muss einfach nur hinschreiben, was man vom Compiler will. In C++. Nicht in "ich bin so super schlau"-Pseudo-Assembly-Style mit wilden Casts.

    Dein Code erzeugt genauso UB wie meiner.



  • Leider habe ich diese ganzen Casts/Funktionen nicht zur Verfügung. Ich brauche eine Lösung ohne reinterpret_cast und ohne Sachen wie memcpy . Ich hab nur eine Art static_cast .



  • happystudent schrieb:

    Leider habe ich diese ganzen Casts/Funktionen nicht zur Verfügung. Ich brauche eine Lösung ohne reinterpret_cast und ohne Sachen wie memcpy . Ich hab nur eine Art static_cast .

    Was genau hast du denn zur Verfügung?

    Arcoth schrieb:

    Dein Code erzeugt genauso UB wie meiner.

    In welcher Zeile und warum?



  • TyRoXx schrieb:

    Was genau hast du denn zur Verfügung?

    Ich habe drei Datentypen: int, float (beide 32 bit) und eine Art string (der aber eigentlich nur zur Text-Ausgabe dient). Für diese 3 Typen habe ich Konvertierungsfunktionen in alle Richtungen (für float -> int zum Beispiel entspricht das ja dann etwa einem static_cast). Die Typen kann ich auch als Array deklarieren.

    An Operatoren für diese Datentypen habe ich +-*/, zum Vergleichen ==, !=, <, >, >=, <=. Sonst nix (also auch keine bit-shift operatoren o.ä.).

    Darüber hinaus kann ich selbst erstellen: Funktionen (inklusive Rekursion), Schleifen, Verzweigungen, Variablen (von den oben genannten Typen). Und das wars dann im Großen und Ganzen auch schon. Sonst gibts nur noch ein paar Mathematische Grundfunktionen (Sinus, Cosinus, Pow, etc.).

    Das sind quasi die Vorgaben. Und mit diesen Mitteln will ich jetzt eine (möglichst effiziente) Version von float_next erstellen.


  • Mod

    Arcoth schrieb:

    Dein Code erzeugt genauso UB wie meiner.

    In welcher Zeile und warum?

    Ab ++i . Laut [expr.pre.incr] und [expr.ass]/7 ist dieser Ausdruck zu i = i + 1; (bis auf die doppelte Auswertung von i ) äquivalent. Es wird also eine l-t-r Konvertierung auf i durchgeführt.

    Jedoch ist dort kein Objekt vom Typ uint32_t mehr existent; [basic.life]/(1.4):

    The lifetime of an object of type T ends when:
    (1.3) — if T is a class type with a non-trivial destructor (12.4), the destructor call starts, or
    (1.4) — the storage which the object occupies is reused or released.

    Du hast per memcpy den Speicher reused, daher wurde die lifetime des ursprünglichen Objektes beendet, und demnach ist die l-t-r Konvertierung von i - welches auf das ursprüngliche Objekt verweist - undefiniert:

    [basic.life]/6 schrieb:

    Similarly, […] after the lifetime of an object has ended and before the storage which the object occupied is reused or released, any glvalue that refers to the original object may be used but only in limited ways. […] The program has undefined behavior if:
    (6.1) — an lvalue-to-rvalue conversion (4.1) is applied to such a glvalue,

    Die Regel in [basic.types]/2 greift hier ebenfalls nicht.

    Die aktuellen Regeln bezüglich lifetime sind zwar ziemlich schrott, aber obiges ist eben der Status Quo.



  • Ich muss gerade feststellen, dass der Standard ganz schön schwammig bei so einem wichtigen Thema ist. Ich weiß nur, dass man seit C die Bytes von PODs manipulieren darf. Das Ergebnis ist eben implementierungsabhängig, aber nicht undefiniert.
    Was mich gerade stört, ist, dass nach deiner Argumentation selbst das undefiniert wäre:

    float f = 2;
    float i;
    std::memcpy(&i, &f, sizeof(i));
    

    Das ist offensichtlich nicht der Fall.


  • Mod

    TyRoXx schrieb:

    I
    Was mich gerade stört, ist, dass nach deiner Argumentation selbst das undefiniert wäre:

    float f = 2;
    float i;
    std::memcpy(&i, &f, sizeof(i));
    

    Das ist offensichtlich nicht der Fall.

    Nein. Hier greift [basic.types]/3:

    For any trivially copyable type T , if two pointers to T point to distinct T objects obj1 and obj2 , where neither obj1 nor obj2 is a base-class subobject, **if the underlying bytes (1.7) making up obj1 are copied into obj2 , obj2 shall subsequently hold the same value as obj1 **.



  • Arcoth schrieb:

    Nein. Hier greift [basic.types]/3:

    For any trivially copyable type T , if two pointers to T point to distinct T objects obj1 and obj2 , where neither obj1 nor obj2 is a base-class subobject, **if the underlying bytes (1.7) making up obj1 are copied into obj2 , obj2 shall subsequently hold the same value as obj1 **.

    Die Regel besagt nur, dass der "value" gleich ist. Das bedeutet, dass kein versteckter Zustand existieren darf, den memcpy nicht mitnehmen kann. Ich kann daraus jetzt nicht erkennen, warum das bei ungleichen Typen sofort undefiniertes Verhalten ist.

    Man kann das Problem auf das hier reduzieren, richtig?

    std::uint32_t t = 2;
    *reinterpret_cast<unsigned char *>(&t) += 1;
    

    Ich dachte immer, dass so etwas erlaubt ist und keineswegs die Lebenszeit von t betrifft.


Log in to reply