Integral berechnen - Grenzwert wird "ungenauer"

#1

Hey,

ich habe folgenden Algorithmus um das bestimmte Integral numerisch zu berechnen:

		double c = (a + b) / 2.0;
		double h3 = abs(b - a) / 6.0;
		double result = h3 * (f(a) + 4.0 * f(c) + f(b));
		return result;

	}
	double integrate(Function &f, double a, double b) {

		double temp = 0;
		double temp2 = 0;
		double h = 1;
		int i = 0;
		int k = 0;
		do{
            k++;
			for (double i = a; i <=b; i += h) {
				if (i + h > b) {
					temp += simpson(f, i, b);
				} else {
					temp += simpson(f, i, i + h);

				}
			}

			for (double i = a; i <=b; i += h/1.1) {
				if (i +h > b) {
					temp2 += simpson(f, i, b/1.1);
				} else {
					temp2 += simpson(f, i, i + h/1.1);

				}
			}

			h/=10;
			printf("%d ", k);
			printf("Differenz: %f " ,abs(temp-temp2));
			fflush(stdout);


		}while(abs(temp-temp2)>EPSN && k<=I_MAX);
        printf("Durchläufe: %d ", i);
		return temp;
	}```

h wird im Moment immer durch 10 geteilt. Vllt ist das etwas zu krass aber gleiche Probleme bestehen mit dem Faktor 2 etc.

Also dadurch, dass das h kleiner wird, müsste das zusammengesetzte Integral mit der Simpsonregel ja immer genauer werden. 
Allerdings passiert eher das Gegenteil.
Wenn ich damit das Integral von sin mit den Grenzen 0 und 5 ausrechnen will entwickeln sich die "Differenzen" so: 

`1 Differenz: 0.449850 2 Differenz: 0.273352 3 Differenz: 1.155551 4 Differenz: 1.603636 5 Differenz: 2.053291 6 Differenz: 2.503103 `
Danach ist h einfach schon zu klein und es dauert ewig.
Werden praktisch immer größer. Wenn man die Fläche in Rechtecke einteilt und dann addiert, dann müssten sich für so kleine "h" ja die Anzahl dieser Rechtecke ja fast nicht mehr unterscheiden und der Grenzwert gegen einen festen Wert laufen. 
Klappt leider für meinen Algorithmus nicht.

Vllt hat einer eine Idee, wie man es hinkriegt dass es genauer wird (anstatt genau das Gegenteil wie jetzt bei mir) und wie man frühzeitig abbrechen kann und nicht ewig läuft. Für kleines Abstände zwischen a und b kann man h ja vllt geschickter wählen also für a=5 und b=10.000 oder so. 

Hoffe jemand kann mir da helfen :)) 

Grüße
#2

Vielleicht ist double einfach nicht genau genug und du kommst mit [japi]BigDecimal[/japi] eher ans Gewünschte?
Oder ich habe dich einfach nur falsch verstanden bzw. keine Ahnung…

#3

Double Genauigkeit sollte für mich schon reichen.
Das komische ist allerdings:

1 Differenz: 0.449850 2 Differenz: 0.273352 3 Differenz: 1.155551 4 Differenz: 1.603636 5 Differenz: 2.053291 6 Differenz: 2.503103

Beim ersten durchlauf ist h=1. Und der 2. Wert wird ja immer Rechtecken der “Breite” h/1.1 berechnet.
Im zweiten Durchlauf ist h dann 0.1.
Im dritten 0.01 und so weiter.

h wird also immer kleiner, dementsprchend sollte die Differenz von temp und temp2 gegen null laufen.
Allerdings wird die Differenz bei mir immer größer.

Deswegen gehe ich irgendwie von einem Denkfehler in meinem Algorithmus aus.

#4

@L-ectron-X
CyborgBeta, hör auf andere Accounts zu hacken :wink:
also wenigstens doch Area C++ statt Java erkennen, wobei ich gerade auch javamäßig dazu denke


Fundamentalprinzip der Programmierung Nummer 1:
konkrete Beispiele anschauen, aufschreiben, Ergebnisse vergleichen!
die analysierte Funktion ist nichtmal genannt, wahrscheinlich aber zu kompliziert,

mit etwas einfachen anfangen, bei f(x) = x wird die Treppchenannäherung wohl immer gleich den richtigen Wert berechnen,
an sich kein gutes Beispiel, bei dir aber auch Fehler denkbar, also würde sich auch lohnen,
sonst irgendwas mit Knick in der Mitte, oder x^2 usw.

h könntest du generell von b-a abhängig machen, das war ja auch deine Frage,
(b-a) / 8 und dann immer weiter halbieren,
Faktor 10 gerade für Anfang mit Testen ungünstig schneller Anstieg

dass du nicht recht selber dran arbeitest ist etwas schade für dich, ergibt andererseits auch ein schönes Forum-Rätsel,
wollte mich gerade an Java-Nachprogrammierung dranmachen, aber dann fielen mir doch zwei Dinge auf:

  • temp/ temp2 wird nicht nach jedem Durchlauf auf 0 zurückgesetzt?!

  • der zweite Durchlauf mit /1.1 ist verbuggt:
    wenn b = 10.000 ist und h sich mit 9997, 9997.8, 9998.6 oder ähnlich langsam annähert
    wird als letzert Schritt b/1.1 genommen ~ 9091, weit vor dem letzten h, gut möglich auch weit vor a

generell ist so ein doppelter Durchlauf ärgerlich, zumal noch mit gleichen Code doppelt hingeschrieben,
vergleiche lieber mit Wert aus der vorherigen Iteration,
Faktor 10 ist dann freilich eher ungünstig, Faktor 2 dürfte reichen

b ist doch immer > a, wenn man von deinem Bug absieht?

edit:
eine entsprechende if-Prüfung mit Exception schadet im experimentellen Bau auch nicht

[quote=Tarrew;107267]h wird also immer kleiner, dementsprchend sollte die Differenz von temp und temp2 gegen null laufen.
Allerdings wird die Differenz bei mir immer größer.[/quote]
was zu antworten ist ja vielleicht angebracht auf andere Antworten,
aber die bekannten Bedenken des ersten Postings zu wiederholen muss nicht sein :wink:

#5

Au! Miiiist! Ist das peinlich… :o Das habe ich voll übersehen. Entschuldigung.

…kommt davon, wenn man Code nur mit Überschall überfliegt…

#6

Einfach ein Fehler in der Berechnung

Versuch mal das Integral von 1 von 0 bis 1, dann das Integral von x von 0 bis 1

#7

In meinem ersten Algorithmus war wohl ziemlich viel Müll.
Hab eure Hinweise mal beachtet und rausgekommen ist:

{
    double c = (a + b) / 2.0;
    double h3 = abs(b - a) / 6.0;
    double result = h3 * (f(a) + 4.0 * f(c) + f(b));
    return result;

}
double integrate(Function &f, double a, double b)
{

    double i = a;
    double h=abs((b-a)/2);
    double temp=0;
    double temp2=0;
    int k;
    do
    {
        temp2=temp;
        temp=0;
        for(i=a; i<b; i+=h)
        {
            printf("a: %f, b: %f, i: %f, h: %f 
", a,b,i,h);
            temp+= simpson(f, i, i+h);
        }
        h/=2;
        k++;
    }
    while((abs(temp2-temp))>EPSN&&k<=I_MAX);
    printf("temp1: %f, temp2: %f 
", temp, temp2);
    printf("Durchläufe: %d 
", k);
    fflush(stdout);
    if(k<I_MAX)
    {
        return temp;
    }
    return NULL;
}

Scheint nach meinen ersten Tests sehr genau zu sein :slight_smile: Werde das noch ein bisschen testen, aber denke der Code geht jetzt in Ordnung.
Danke :slight_smile: