[Untitled]

Hallo zusammen. Ich besuche momentan eine Vorlesung namens Numerik I und muss dazu auch Sachen in C bzw C++ programmieren. Leider haben die nicht so wirklich gute Scripte bzw Tutorien dazu angeboten, dafür die Aufgaben aber recht schwer gestaltet ;-)

Nun bräuchte ich bei der aktuellen ein wenig Hilfe

Hier das Blatt
http://www.uni-due.de/ingmath/downlo...s09/prog05.pdf

Kann das irgend jemand von eich, oder kann mir gute Tipps dazu geben?

[Hisoka999]

1. Romberg-Integration verstehen: http://de.wikipedia.org/wiki/Romberg-Integration
2. Diese Vorgehensweise in Code umsetzen.
==> nicht viel mehr als ne Schleife und ein paar zerlegte Formeln.

[Untitled]

1. Romberg-Integration verstehen: http://de.wikipedia.org/wiki/Romberg-Integration
2. Diese Vorgehensweise in Code umsetzen.
==> nicht viel mehr als ne Schleife und ein paar zerlegte Formeln.

Das ist der Punkt. Das ganze umsetzen ;-) Ich habe in C leider eher rudimentäre Kenntnisse und beiße mir gerade ein bisschen die Zähne daran aus.

[Xorlosch]

In der englischen Wikipedia gibt es eine C-Implementierung der Romber-Integration, vielleicht hilft dir die.

[Thoronador]

Das ist der Punkt. Das ganze umsetzen ;-) Ich habe in C leider eher rudimentäre Kenntnisse und beiße mir gerade ein bisschen die Zähne daran aus.

Wie weit reichen die Kenntnisse denn bzw. wo hängt es?

Im Prinzip gibt es in C/C++ zwei (bzw. drei) für diese Art von Aufgabe relevante Datentypen: int und float (sowie double). Der Typ int kann ganze Zahlen aufnehmen, der Typ float ist für Gleitkommazahlen (sprich reelle Zahlen bzw. eine Teilmenge davon). double ist im Grunde erst einmal das gleiche wie float, nur dass man mit double eine größere Genauigkeit erhält, weil hier mehr Speicher für eine Zahl verwendet wird. Gleitkommazahlen sind nämlich aufgrund der Art, wie sie gespeichert werden, nicht immer genau sondern enthalten meist einen gewissen Rundungsfehler. Der Wert 2,5 ließe sich noch genau darstellen, der Wert 2,4 z.B. nicht mehr. Das führt dazu, dass bei Rechenoperationen mit Fließkommazahlen kleinere (und manchmal auch größere) Ungenauigkeiten auftreten und ein Algorithmus, der algebraisch betrachtet 2,0 als Wert liefern sollte, am Ende einen Wert wie 1,999997 oder 2,000003 liefert. (Und das sind noch die angenehmeren Fälle.)
Daher ist es bei Gleitkommazahlen meistens empfehlenswert, den Datentyp double zu nehmen. (Eine Ausnahme könnten Matrizen oder ähnliche Objekte sein, wo man sehr viele Zahlen hat und float den Vorzug gibt, um Speicherplatz zu sparen.)

Weiterhin brauchst du hier (und bei vielen anderen numerischen Aufgaben/ Verfahren auch) die Zählschleife/ for-Schleife. Grundlegender Aufbau derselbigen ist folgender:

int i; /* Zählvariable */

for (i=0; i<50; i=i+1)
{
  Anweisung;
  Noch_eine_Anweisung;
}

Oder allgemeiner:

for (Startanweisung; Bedingung; Schleifenanweisung)
{
  Anweisungsblock;
}

Dabei ist Startanweisung eine Anweisung, die einmalig ganz zu Anfang der Schleife ausgeführt wird. Normalerweise wird das für die Initialisierung der Zählvariable genutzt, damit diese einen definierten Anfangswert hat. Mit Bedingung ist im Prinzip ein Ausdruck gemeint, der sich entweder zu wahr oder falsch auswerten lässt. Beispiele wären 1<5 (wahr), 2.4>=variable oder variable==andere_variable (möglicherweise wahr). Man könnte da auch noch Funktionsausrufe reinbauen oder was es sonst noch gibt. Solange diese Bedingung wahr ist, wird die Schleife ausgeführt, d.h. der Anweisungsblock wird immer wieder abgearbeitet. Geprüft wird die Bedingung immer zu Beginn eines Schleifendurchlaufes. Es wäre also auch denkbar, dass der Anweisungsblock nie ausgeführt wird, weil die Bedingung schon ganz am Anfang nicht erfüllt wird. Oft ist die Bedingung derart, dass sie prüft, ob die Zählvariable einen gewissen Wert erreicht hat bzw. über- oder unterschreitet. Und schließlich ist Schleifenanweisung eine Anweisung, welche je einmal nach Abarbeitung des Anweisungsblocks ausgeführt wird. Meistens benutzt man dies, um die Zählvariable um eins (oder einen beliebigen anderen Wert) zu erhöhen.
Im obigen Beispiel würden also die beiden Zeilen im Schleifenrumpf 50 mal hintereinander ausgeführt. Im Falle der Rombergintegration könnte man so eine Schleife beispielsweise nutzen, um auf eine Variable mehrmals verschiedene Werte draufzuaddieren.

Ach ja, dann gibt es in C/C++ die Funktion pow() von engl. power zum Berechnen von Potenzen, welche zwei Argumente (Gleitkommazahlen) entgegennimmt, Basis und Exponent. Der Aufruf pow(a,b) würde die Potenz a^b zurückliefern. Damit die Funktion verfügbar ist, muss man aber erst den zugehörigen Header einbinden, was über das Hinzufügen der Zeile

#include <cmath>

am Anfang des Codes passiert. Ist vielleicht ganz hilfreich für die Schrittweitenberechnung.

[Untitled]

Danke für das Grundtutorium ;-)
Ich habe in der Schule (ist mittlerweile auch 5 Jahre her) ein bisschen java gemacht; ist ja alles recht ähnlich.

Wie krieg ich das denn mit den Integralen hin? Wie kann ich die in C beschreiben?

[Thoronador]

Ah, gut, dann kennst du vermutlich auch if-Verzweigungen aus Java, das ist in C/C++ ähnlich:

if (Bedingung)
{
  Anweisungsblock1;
}
else
{
  Anweisungsblock2;
}

Prinzip sollte aus Java bekannt sein: wenn Bedingung ausgewertet wahr ergibt, dann wird der Anweisungsblock1 ausgeführt, andernfalls der Anweisungsblock2. Den Teil ab dem Schlüsselwort else kann man auch weglassen, wenn man etwas nur machen möchte, falls die Bedingung zutrifft und im anderen Fall nichts geschehen soll.

Was könnte man noch brauchen? Funktionen vermutlich. Eine (zugegeben etwas einfache/ sinnlose) Funktion zum Addieren zweier Gleitkommawerte mit doppelter Genauigkeit sähe zum Beispiel so aus:

double Add(double a, double b)
{
  return a+b;
}

Allgemeiner:

Rückgabetyp Funktionsname(Parameterliste)
{
  Anweisungen;
}

Dabei ist Rückgabetyp der Datentyp des Wertes, welchen die Funktion zurückgeben soll. Funktionsname ist selbsterklärend, und Parameterliste ist eine (möglicherweise auch leere) Liste von Parametern, welche an die Funktion übergeben werden. Parameter werden dabei vom Prinzip her wie Variablendeklarationen angegeben, also Typ gefolgt vom Namen des Parameters.
In den Anweisungen im Funktionsrumpf kann alles mögliche vorkommen, auch wieder Verzweigungen oder Schleifen. Wichtig ist aber, dass mindestens einmal irgendwo eine Anweisung mit dem Schlüsselwort return erscheint. return weißt die Funktion an, den Wert zurückzugeben, der direkt auf return folgt, also im obigen Beispiel die Summe der Parameter a und b. Man kann in einer Funktion auch mehrere return-Anweisungen haben, allerdings wird die Funktion beendet, sobald die erste return-Anweisung erreicht wird.


Jetzt zu der Frage, wie man das mit den Integralen macht: gar nicht. Es gibt in C oder C++ keine vordefinierte Funktion, der man einfach zwei Werte a und b sowie eine Funktion übergibt (ja, auch Funktionen können Parameter einer Funktion sein), und die dann das bestimmte Integral dieser Funktion von a bis b zurückliefert. Also muss man sich selbst etwas einfallen lassen, und hier kommt in deinem Fall die Rombergintegration ins Spiel.

Ein recht einfaches Beispiel unter Nutzung der Trapezregel wäre folgendes:

double Integration(double a, double b, double (*f)(double))
{
  double h;
  h = b-a;
  return 0.5*h*(f(a)+f(b));
}

Der Wert wird natürlich erstmal nicht sehr genau sein, aber das soll hier keine Rolle spielen. Hat man irgendwo im Code noch eine Funktion definiert, die man damit "integrieren" möchte, nennen wir sie FunktionZwei, dann könnte das in C++ so aussehen:

#include <iostream>

double Integration(double a, double b, double (*f)(double))
{
  double h;
  h = b-a;
  return 0.5*h*(f(a)+f(b));
}

double FunktionZwei(double x)
{
  return 1.0f/(1.0f+x*x);
}

/* Hauptfunktion des Programmes */
int main()
{
    std::cout << "Integration von 0 bis 1: " << Integration(0, 1, &FunktionZwei);
    return 0;
}

Die Aufgabe könnte nun also sein, die Funktion Integration() so abzuändern, dass sie anstelle der Trapezregel das Verfahren von Romberg nutzt, um das bestimmte Integral zu berechnen.

[Untitled]

Kann mir jemand sagen, wo der Fehler im Programm ist? Bei Aufruf stürzzt es ab
Es ist ein C Programm

#include <stdio.h>
#include <stdlib.h>
#include <math.h>


double betrag(double x)
{
  if (x>=0) 
  { return x;}
  else
  { return (-x);}
}

double wurzel(double x)
{
    x = pow(x, 0.5);
	return (x);
}

void Romberg_Integration(double a, double b, int imax, double (*f)(double), double T[imax][imax],double wert)
{
       int i,j;
       double h[imax];
       double x[imax];
       int n[imax];
       double summe;
       double F[imax][imax];
       
       /*definiere Schrittweitenfolge*/
       h[0] = b-a;
       for(i=1; i<imax; i++)
       {
          h[i] = pow(2,(-i));
          h[i] = h[i] * h[0];
       }
       
       /*definiere Anzahl der Stützstellen*/
       for(i=0; i<imax; i++)
       {
          n[i] = pow(2,i) + 1;
       }
       
       /*definiere Stützstellen*/
       x[0]= a;
       for(i=1; i<n[i]; i++)
       {
          x[i] = a + i*h[i];
          x[n[i]] = b;
       }
       
       for(i=1; i<n[i]-1; i++)
       {
          summe += (*f) (x[i]);
       }
       
       /*definiere erste Spalte der Matrix(zusammengesetzte Trapezregel)*/
       for(i=0; i<imax; i++)
       {
          T[i][0] = (*f) (a) + (*f) (b) + 2*summe;
          T[i][0] *= (h[i]/2);
       }
       
       /*definiere übrigen Näherungen des Schemas*/
       for(i=1; i<imax; i++)
       {
                for(j=1; j<=i; j++)       
                {
                   T[i][j] = pow(4,j) * T[i][j-1] - T[i-1][j-1] ;
                   T[i][j] = T[i][j] / (pow(4,j)-1) ;
                }
       }
       
       /*Extrapolationsschema ausgeben*/
       printf("Die Nährungen für das Integral von a nach b von f(x)dx lauten: \n");
       for(i=0; i<imax; i++)
       {
          for(j=0; j<=i; j++)
          {
              printf(" %lf ", T[i][j], i,j);
          }
          printf("\n");
        }
        printf("\n");
        
        /*absolute Fehler der Näherungen*/
        for(i=0; i<imax; i++)
        {
           for(j=0; j<=i; j++)
           {
              F[i][j] = T[i][j] - wert;
              F[i][j] = betrag(F[i][j]);
           }
        }
     
        printf("Die absoluten Fehler der Näherungen betragen: \n");
        for(i=0; i<imax; i++)
        {
            for(j=0; j<=i; j++)
            {
               printf(" %lf ", F[i][j], i,j);
            }
            printf("\n");
        }
        printf("\n");
       
}

double funktion1(double x)
{
       x = pow(x,1/8);
       return x;
}

double funktion2(double x)
{
       x = 1 + x*x;
       x = 1/x;
       
       return x;
}

/*double funktion3(double x)
{
       double const pi = 3.141592654;
       
       x = ( sin(x) / (wurzel(2) );
       x = exp (x);
       x = x * ( 1/(2*pi));
       
       return x;
}*/
       
int main()
{
    double const pi = 3.141592654;
    
    int imax = 5;
    double a1 = 0;
    double b1 = 1;
    double wert1 = 8/9;
    
    double a2 = -1;
    double b2 = 1;
    double wert2 = pi/2;
    
    double a3 = 0;
    double b3 = 2*pi;
    double wert3 = 1.1289609;
    
    double T[imax][imax];
    
    Romberg_Integration(a1,b1,5,funktion1,T,wert1);
    Romberg_Integration(a2,b2,5,funktion2,T,wert2);
   // Romberg_Integration(a3,b3,5,funktion3,T,wert3);
    
    system("PAUSE");
    return 0;
    
}

[Thoronador]

Es hat zwar (vermutlich) nichts mit dem Absturz zu tun, aber in der Funktion Romberg_Integration() hast du die Variable summe nicht initialisiert. Nach der Deklarationszeile am Anfang ist die nächste Zeile, in der summe erscheint, diese:
summe += (*f) (x[i]);
Die Addition klappt zwar, aber da du nirgends einen Anfangswert für summe definiert hast, könnte da am Anfang auch theoretisch jeder beliebige Wert drinstehen, was dann das Ergebnis verfälschen würde. Daher sollte man noch irgendwo davor die Zeile
summe = 0.0;
einfügen.

Weiter geht's mit den Zeilen zur Ausgabe der Näherung und des Fehlers. Dort steht u.a.:

printf(" %lf ", T[i][j], i,j);

Wie man sieht, willst du hier drei Werte ausgeben lassen, aber der Formatstring ist nur für einen Wert ausgelegt. Das müsste man noch ändern. Möglicherweise ist das auch der Punkt, wo das Programm sagt "Nein, danke!" und sich verabschiedet.



Der eigentliche Fehler mit Absturz als Folge liegt vermutlich in diesem Teil von Romberg_Integration:

       /*definiere Anzahl der Stützstellen*/
       for(i=0; i<imax; i++)
       {
          n[i] = pow(2,i) + 1;
       }
       
       /*definiere Stützstellen*/
       x[0]= a;
       for(i=1; i<n[i]; i++)
       {
          x[i] = a + i*h[i];
          x[n[i]] = b;
       }

In der ersten Schleife weißt du n[i] um eins erhöhte Potenzen von 2 zu, die bei i=4 den Wert 2^4+1=17 erreichen können. In der nächsten Schleife nutzt du diese Werte als Arrayindex für x, der Indexbereich von x läuft aber nur von 0 bis 4, wenn ich das richtig gesehen habe. Damit gibt es dort einen Speicherzugriffsfehler und das Programm verabschiedet sich.

[Untitled]

Der eigentliche Fehler mit Absturz als Folge liegt vermutlich in diesem Teil von Romberg_Integration:

       /*definiere Anzahl der Stützstellen*/
       for(i=0; i<imax; i++)
       {
          n[i] = pow(2,i) + 1;
       }
       
       /*definiere Stützstellen*/
       x[0]= a;
       for(i=1; i<n[i]; i++)
       {
          x[i] = a + i*h[i];
          x[n[i]] = b;
       }

In der ersten Schleife weißt du n[i] um eins erhöhte Potenzen von 2 zu, die bei i=4 den Wert 2^4+1=17 erreichen können. In der nächsten Schleife nutzt du diese Werte als Arrayindex für x, der Indexbereich von x läuft aber nur von 0 bis 4, wenn ich das richtig gesehen habe. Damit gibt es dort einen Speicherzugriffsfehler und das Programm verabschiedet sich.

Ah, ich verstehe was du meinst.

Wie müsste ich es denn umändern?

[Thoronador]

Du müsstest verhindern, dass auf einen Index zugegriffen wird, der die Größe des Arrays überschreitet. Im Wesentlichen gibt es da zwei Methoden: entweder man legt den Array von Anfang an mit einer Größe an, die dem maximalen Index entspricht, oder aber man schränkt den Index auf die Größe des Arrays.

[Untitled]

Gut, das Programm stürzt schonmal nicht mehr ab! Sehr gut, danke ;-)

Jetzt muss ich nur noch die Rechenfehler bzw Programmfehler finden, denn es kommen nicht die richtigen Ergebnisse raus

Kann ich die Funktionen eigentlich so definieren? Ich meine die Integrale. Wie die eigentlich aussehen sollen, sieht man ja auf dem verlinkten Arbeitsblatt..

[Thoronador]

Du meinst funktion1(), funktion2() sowie die auskommentierte funktion3() in deinem Code? Ja, kann man, zumindest solange der Compiler das ohne Murren akzeptiert. Aber das war wohl nicht die Frage.
Im Prinzip kann man das so schreiben, aber sofern man die Klammern richtig setzt, könnte man die Funktionen auch in weniger Zeilen schreiben.

double funktion1(double x)
{
  return pow(x, 0.125);
}

double funktion2(double x)
{
 return (1.0f/(1.0f+x*x));
}

double funktion3(double x)
{
  double const pi = 3.141592654;
  return (1.0f/(2*pi) * exp(sin(x)/sqrt(2)));
}

Wie man sieht, habe ich in der dritten Funktion deine Funktion wurzel() durch sqrt() ersetzt. Die Funktion sqrt() ist Bestandteil von math.h und liefert die Quadratwurzel der angegebenen Zahl zurück. Im Englischen "square root", daher kommt der Funktionsname. Wenn man Funktion 3 noch ein wenig beschleunigen will, kann man die Wurzel von 2 auch ausrechnen und ähnlich wie Pi als Konstante definieren. Damit verhindert man, dass der Wert bei jedem Aufruf der Funktion erneut berechnet werden muss.