exp() Funktion selber implementieren



  • Ich habe mal umgesetzt, dass es nicht immer neu berechnet wird. Das hat super geholfen:

    long double ehoch(double x){
        if(x == 0)
            return 1.;
        bool invers = false;
        if(x < 0){
            invers = true;
            x *= -1;
        }
        long double erg = 0;
        long long counter = 1;
        double epsilon = 0.00001;
        long double numerator = 1;
        long double denominator = 1;
        while(epsilon < (numerator/denominator)){
            erg += numerator/denominator;
            numerator *=  x;
            denominator *= counter;
            counter++;
        }
        return invers ? 1/erg : erg;
    }
    

    Für negative Exponenten lasse ich einfach Exp(|x|) berechnen und gebe dann die Inverse zurück. Also 1/exp(|x|).
    Jetzt funktioniert es mit der gewünschten Genauigkeit.

    Edit: Ich nutze jetzt noch long double. Dann funktioniert das auch für große Exponenten.



  • Ich frage mich aber weiterhin wieso das komplette Neuberechnen von pow() und fakultaet() zur nummerischen Auslöschung führt und diese Lösung nicht. Letztlich werden bei der Neuberechnung doch einfach alle Schritte "wiederholt" und dann das nächste multipliziert und zurückgegeben. In meinem jetzigen Fall spare ich mich doch nur die Wiederholungen und multipliziere direkt das neue dran. Wo ist denn nun der Trick versteckt?


  • Mod

    Mit dem "Trick" kommt keines deiner Zwischenergebnisse jemals in ungünstige Wertebereiche. Wenn du pow durch die Fakultät teilst, dann sind zwischendurch der Rückgabewert von pow und der Rückgabewert der Fakultät beide sehr groß. Das ist numerisch sehr ungünstig. Auch wenn das Ergebnis am Ende klein sein mag, ist es dann unnötig ungenau. Bei der Fakultät sprengst du sogar ganz schnell den Wertebereich sogar von long long und bekommst dann ein komplett falsches Ergebnis.

    Als Faustregel willst du immer mit möglichst gleich großen und allgemein mit nicht zu großen Zahlen arbeiten (außer bei Subtraktion).

    Damit ergibt sich auch direkt ein weiterer Tipp: Du addierst hier lauter kleine Zahlen zu einer immer größer werdenden Summe. Das ist auch ungünstig, weil diese Addition dadurch mit jedem Schritt numerisch ungenauer wird. Idealerweise möchte man mit dem Addieren der kleinen Zahlen anfangen und dann zu den größeren Zahlen kommen, so dass die laufende Summe und der neue Faktor stets ungefähr gleich groß sind.

    Das ist hier leider nicht so ganz einfach, aber als einfache und effektive Maßnahme betrachte mal den ersten Faktor: Der ist immer 1. Das ist ungünstig, denn 1 plus etwas kleines ist numerisch viel ungenauer als 0 plus etwas kleines. Daher: Lass den ersten Faktor weg und berechne erst exp(x) -1. Beim return am Ende addierst du dann die fehlende 1. So hast du gerade bei kleinen Exponenten durch einfaches Umstellen der Rechnung deutlich an Genauigkeit gewonnen.



  • Hier noch meine Lösung: Ideone-Code
    Bei Reduzierung von epsilon kommen auf 10 Stellen exakt die Werte von exp(x) heraus:

    0: 1 <-> 1
    1: 2.718281828 <-> 2.718281828
    2: 7.389056099 <-> 7.389056099
    3: 20.08553692 <-> 20.08553692
    4: 54.59815003 <-> 54.59815003
    5: 148.4131591 <-> 148.4131591
    6: 403.4287935 <-> 403.4287935
    7: 1096.633158 <-> 1096.633158
    8: 2980.957987 <-> 2980.957987
    9: 8103.083928 <-> 8103.083928
    

  • Mod

    @Th69 sagte in exp() Funktion selber implementieren:

    Hier noch meine Lösung: Ideone-Code
    Bei Reduzierung von epsilon kommen auf 10 Stellen exakt die Werte von exp(x) heraus:

    0: 1 <-> 1
    1: 2.718281828 <-> 2.718281828
    2: 7.389056099 <-> 7.389056099
    3: 20.08553692 <-> 20.08553692
    4: 54.59815003 <-> 54.59815003
    5: 148.4131591 <-> 148.4131591
    6: 403.4287935 <-> 403.4287935
    7: 1096.633158 <-> 1096.633158
    8: 2980.957987 <-> 2980.957987
    9: 8103.083928 <-> 8103.083928
    

    Und mit meiner Verbesserung bekommen wir fast noch eine Stelle mehr geschenkt, ohne das epsilon zu verkleinern zu müssen. Einzig exp(1) kommt durch eine ungünstige Rundung zu einem anderen Ergebnis:
    https://ideone.com/ZRvXHM

    0: 1 <-> 1
    1: 2.7182818284 <-> 2.7182818285
    2: 7.3890560989 <-> 7.3890560989
    3: 20.085536923 <-> 20.085536923
    4: 54.598150033 <-> 54.598150033
    5: 148.4131591 <-> 148.4131591
    6: 403.42879349 <-> 403.42879349
    7: 1096.6331584 <-> 1096.6331584
    8: 2980.957987 <-> 2980.957987
    9: 8103.0839276 <-> 8103.0839276
    


  • @Asse sagte in exp() Funktion selber implementieren:

    Edit: fakultaet habe ich auch selber implementiert das funktioniert meines Erachtens aber einwandfrei:
    Aus Verzweiflung arbeite ich überall schon mit den größten mir bekannten Datentypen.

    long long fakultaet(long long x){  // long long für x ist ganz sicher übertrieben (char würde reicht)
        if(x < 0){                     // ok ein unsigned
            //Fehler
            return -1;
        }
     ....
    
    

    Das Problem liegt auch daran, dass schon ab 22! der Wertebereich vom long long (+/- 63-Bit) überschritten wird.
    unsigned long long macht es auch nicht besser.



  • @SeppJ: Ich habe mal beide Varianten auf 16 Stellen erweitert (sowie epsilon auf 1e-16 gesetzt):

        0: 1 <-> 1
        1: 2.718281828459045 <-> 2.718281828459045
        2: 7.38905609893065 <-> 7.38905609893065
        3: 20.08553692318767 <-> 20.08553692318767
        4: 54.59815003314424 <-> 54.59815003314424
        5: 148.4131591025766 <-> 148.4131591025766
        6: 403.4287934927351 <-> 403.4287934927351
        7: 1096.633158428459 <-> 1096.633158428459
        8: 2980.957987041728 <-> 2980.957987041728
        9: 8103.083927575384 <-> 8103.083927575384
    

    Auch ohne das Ändern von epsilon waren die beiden Ausgaben schon identisch (nur eben noch in den letzten 4-5 Stellen abweichend von exp(x)).



  • Ich danke euch für die Hilfe!!

    Als finale Version nehme ich folgenden Code.

    long double ehoch(double x){
        if(x == 0)
            return 1.;
        bool invers = false;
        if(x < 0){
            invers = true;
            x *= -1;
        }
        long double erg = 0;
        long long counter = 2;
        long double epsilon = 1e-16;
        long double next = x;
        while(epsilon < abs(next)){
            erg += next;
            next *= x / counter;
            counter++;
        }
        return invers ? 1/(erg + 1) : (erg + 1);
    }
    

    Der Funktioniert dann auch für negative Exponenten und ist auch für große Exponenten recht genau.



  • @SeppJ Ok das verstehe ich. Super Tipps, danke!

    Also wäre es nummerisch besser, wenn man erst die Summanden "von klein nach groß" addiert?
    Die Reihe, die wir hier verwenden liefert für die ersten paar Indizes noch relativ kleiner Werte, wächst dann stark an und wird ab einen gewissen Index immer kleiner und geht gegen Null. Müsste es dann helfen die ganze Rechnung "rückwärts" zu machen? Z.B in einer schleife den Summanden next aufzubauen ohne an erg zu addieren und dann ab dem Eintauchindex wo expsilon > abs(next) ist rückwärts anzufangen auf erg zu addieren und dann statt next *= x / counter -> next *= i / x wobei i Index der Schleife vom Eintauchindex bis 0 ist?!

    Ich hoffe die Idee ist verständlich.

    Edit: Eine andere Idee wäre es die einzelnen Summanden in einem Array oder Vector zu speichern, dann der größe nach zu sortieren und dann nacheinander addieren?!


  • Mod

    Hier ist es nicht so ganz wichtig, darauf zu achten, weil mit dem epsilon eine gewisse zu erreichende Genauigkeit vorgegeben ist und genau so lange gerechnet wird, bis diese erreicht ist. Ob die Stellen danach mit gleichem Aufwand auch noch hätten besser sein können, ist dann weniger interessant. Meistens hat man es eher umgekehrt, dass man die Zahl der Reihenglieder vorgibt, und dann eben hofft, dass das Ergebnis so genau wie möglich ist. Dann ist es wichtiger, dass man mit dem vorgegebenen Rechenaufwand ien möglichst gutes Ergebnis zu erzielen.



  • @SeppJ Ok, danke für deine Hilfe!!

    Mit meinem oben genannten Ergebnis bin ich auf jeden Fall zufrieden.


  • Gesperrt

    Falls das hier jemanden interessiert: So haben das die alten Taschenrechner gemacht:
    https://de.wikipedia.org/wiki/CORDIC


Anmelden zum Antworten