Fehler in Differenzieren einer Funktion



  • Wir sollen gerade einen Algorithmus zum differenzieren einer Funktion schreiben und so wie ich diesen vertehe, habe ich ihn auch mal aufgeschrieben aber da scheint wohl noch irgendwo ein Fehler zu sein. Ich habe auch testweise einmal den "normalen" Differenzierungsalgorithmus geschrieben und der funktioniert so auch (siehe weiter unten) aber diesen, den unser Prof. 4-Punkt (Approximation) genannt hat, läuft leider nicht.

    Ich bekomme immer nan raus, wenn ich bspw. die Funktionen Sinus und Tangens teste (oder auch andere). Sieht hier jemand den bzw. einen (oder mehrere) Fehler in meinem Algorithmus? Habe ich etwas weiter unten aus dem Blatt falsch verstanden? Ich hänge gerade bei der Differenzierung echt fest.

    using namespace std;
    
    #define H 1
    #define EPS 1.E-9
    
    double fourPointApprox(Function& func, int x, int h);
    double test(int n, Function& func, double x, double h);
    
    double differentiate(Function& func, double x) {
    	double h = H;
    	double lastRes = 0, result;
    
    	double dh = 0, dh_2 = 0, dh_4 = 0;
    	double f_h = 0, fh_2 = 0;
    
    	do {
    		lastRes = result;
    
    		dh = fourPointApprox(func, x, h);
    		dh_2 = fourPointApprox(func, x, 2 * h);
    		dh_4 = fourPointApprox(func, x, 4 * h);
    
    		f_h  = (((8 * dh) - dh_2) / (12 * h)) - ((pow(h, 4) / 30) * test(5, func, x, h));
    		fh_2 = ((8 * dh_2) - dh_4) / (24 * h) - (16 * (pow(h, 4) / 30) * test(5, func, x, h));
    
    		result = ((16 * f_h) - fh_2) / 15;
    		h /= 2;
    	} while(fabs(lastRes - result) > EPS);
    
    	return result;
    }
    
    double fourPointApprox(Function& func, int x, int h) {
    	int i;
    	double result = 0;
    
    	result = h * (func(x + h) - func(x - h)) / (2 * h);
    	for(i = 3; i <= 7; i += 2) {
    		result += (pow(h, 3) / (double) faculty(i)) * test(i, func, x, h);
    	}
    
    	return (2 * fabs(result));
    }
    
    double test(int n, Function& func, double x, double h) {
    	int i;
    	double result = 0;
    
    	for(i = 0; i < n; i++) {
    		result = (func(x + h) - func(x - h)) / (2 * h);
    	}
    	return result;
    }
    

    Hier ist ein Bild, von dem was wir dazu bekomemn haben aber so kannte ich die Differenzierung auch noch nicht.
    http://www.bilder-upload.eu/show.php?file=c64181-1513631337.jpg

    Wir sollen den mit der (8) gekennzeichneten Differenzierungsalgorithmus implementieren und da da ja auch die Fakultät vorkommt, habe ich mir eine Fakultätsfunktion geschrieben (faculty), da mir da keine Standard für bekannt ist. Diese habe ich aber auch getestet und zeige Sie daher hier nicht. Ist ja auch relativ leicht sowohl iterativ als auch rekursiv zu implementieren.
    Mit func( .. ) hole ich mir jeweils den Wert an der Stelle, nur um das noch einmal zu erwähnen (Klammer Operator überladen) und EPS ist die gewünschte Genauigkeit.

    Das hier ist der "normale (einfache)" Differenzierungsalgorithmus aber den sollen wir ja nicht realisieren, obwohl er auch auf dem Blatt steht aber ich wollte sichergehen, dass alles mit der Funktion etc. in Ordnung ist, da wir das ganze darum auch selber schreiben sollten und daran soll es dann wohl nicht liegen aber der Vollstängidkeit her, zeige ich ihn einfach mal.

    double differentiate(Function& func, double x) {
    	double h = H, n = 0, lastRes = 0;
    
    	do {
    		last = result;
    		result = (func(x + h) - func(x - h)) / (2 * h);
    		h /= 4;
    	} while(fabs(lastRes - result) > EPSr);
    
    	return result;
    }
    


  • Ohne den Algorithmus jetzt Inhaltlich geprüft zu haben, aber du initialisierst result nicht und weist dann den nicht initialisierten Wert dann lastRes zu.



  • Eher Adlerauge als Schlangenmensch 👍



  • Ich habe mal schnell angefangen:

    #include <iostream>
    #include <cmath>
    
    template<typename Func>
    double delta_h(Func f, double h, double x) {
        return f(x+h) - f(x-h);
    }
    
    template<typename Func>
    double f_diff_h(Func f, double h, double x) {
        return (8 * delta_h(f, h, x) - delta_h(f, 2*h, x)) / 12 / h;
    }
    
    template<typename Func>
    double f_diff(Func f, double x) {
        auto h = 0.000005;
        return (16*f_diff_h(f, h, x) - f_diff_h(f, 2*h, x))/15;
    }
    
    void test() {
        auto func = [](double x){return sin(x);};
        const auto pi = 4*atan(1);
        std::cout << f_diff(func, pi) - (-1) << " = 0\n";
        std::cout << f_diff(func, pi/2) << " = 0\n";
        std::cout << f_diff(func, pi/4) - cos(pi/4) << " = 0\n";
    }
    

    scheint zu funktionieren. Nur habe ich das h konstant gelassen und nicht geändert. Das kannst du ja gerne einbauen.

    Übrigens noch 2 andere Hinweise:

    Vermeide #define in C++, das ist hier komplett unnötig. Verwende einfach eine Konstante! (und EPS sieht so sehr nach dem gleichnamigen Format aus, da würde ich doch "epsilon" deutlich bevorzugen)

    Bemerkung 2: Übersetze richtig ins Englische, sonst wirkt es komisch. Faculty ist sowas wie die Fakultät/der Fachbereich an einer Uni in UK bzw. das akademische Personal (US). n! ist das "factorial". Aber wo kommt das denn vor?!

    Achso, und noch was:
    die Test-Funktion sollte am besten auch etwas testen. Ich habe hier mal die Ableitung von sin an 3 verschiedenen Stellen getestet (da soll also 3x0 herauskommen). Hier kannst du die Ausgabe prüfen. Wenn mal alles funktioniert, kannst du dir auch z.B. sowas wie googletest anschauen. Dann würdest du da sowas wie

    for (double x : x_test_values) 
        EXPECT_FLOAT_EQ(cos(x), f_diff(func, x));
    

    schreiben können. Vielleicht ist float_eq auch zu stark, dann eher:

    for (double x : x_test_values) 
        EXPECT_TRUE(std::abs(cos(x) - f_diff(func, x)) < epsilon);
    

    oder so ähnlich.

    Noch was:

    for(i = 0; i < n; i++) { 
            result = (func(x + h) - func(x - h)) / (2 * h); 
        }
    

    Das sieht komisch aus. Der Loop-Body hängt nicht von der Loop-Variablen ab. func soll wohl auch keine Nebeneffekte haben. Was soll also das for?


Anmelden zum Antworten