Hi,

wahrscheinlich wiedermal einvollkommen triviales Problem für jemanden der sich auskennt mit C++:

Ich habe ein problem bei einer numerischen integration und habe daher folgendes testprogramm geschrieben um mir die funktionsweise besser klar zu machen:

#include <cstdlib>

#include <iostream>

#include "gauleg.cpp" /* gauleg ist eine funktion zur berechnung von stützstellen und gewichten für die numerische integration*/

using namespace std;

double quad(double x)

{

double z2;

z2=x*x;

return z2;

}

int main()

{

const double xu =0.0;

const double xo =1.0;

const int n_max=100;

int i;

double gew[n_max], absc[n_max];

double z1;

gauleg( xu, xo, absc, gew, n_max);

z1=0.0;

for(i=0;i<n_max;i++)

{ z1 += gew*quad(absc);

cout<<"Gewichte: "<<gew<<endl;

cout<<"Stuetzstellen: "<<absc<<endl;

cout<<"Integral von x^2: "<<z1 <<endl;

}

system("Pause");

return 0;

}

Das Problem ist daß ich für z1 enweder #QNAN oder irgendwas mit exp (294) herausbekomme. QNAN lässt ja auf eine falsche var definition schliessen. Aber was ist hier genau das problem? Liegt es daran dass ich elemente eines Arrays miteinander mult will?

Gruß,

Mesrine

Schau Dir mal die Inhalte der Arrays nach dem Gauleg an. Eventuell werden die Arrays nicht bis zum Schluß mit gültigen Werten befüllt. Vielleicht ist die Funktion auf eine bestimmte Maximalgröße beschränkt? Probiere mal testweise xu = 0.1 und xo = .9 - steht irgendwo ein ungültiger Wert liegt das Problem im Gauleg(). Das würde bedeuten, der Gauleg() füllt die Arrays vielleicht nicht mit doubles, sondern mit einem anderen Datentyp. Wie sieht die cpp dazu aus?

Danke für deine Hilfe hier ist die (Gauss-Legendre) gauleg.cpp:

void gauleg(double x1,double x2,double x[],double w[], int n)

{

int m,j,i;

double z1,z,xm,xl,pp,p3,p2,p1;

const double EPS=3.0e-7;

m=(n+1)/2;

xm=0.5*(x2+x1);

xl=0.5*(x2-x1);

for (i=1;i<=m;i++) {

z=cos(3.141592654*(i-0.25)/(n+0.5));

do {

p1=1.0;

p2=0.0;

for (j=1;j<=n;j++) {

p3=p2;

p2=p1;

p1=((2.0*j-1.0)*z*p2-(j-1.0)*p3)/j;

}

pp=n*(z*p1-p2)/(z*z-1.0);

z1=z;

z=z1-p1/pp;

} while (fabs(z-z1) > EPS);

x=xm-xl*z;

x[n+1-i]=xm+xl*z;

w=2.0*xl/((1.0-z*z)*pp*pp);

w[n+1-i]=w;

}

}

Meines Erachtens füllt gauleg die arrays schon mit doubles aber ich kann sie irgendwie nicht miteinander multiplizieren..

Mesrine

Übrigens klappt alles wenn ich n=80 nehme. Aber für n=100 kommt wieder #QNAN für den Integralwert

Hallo,

mit n_max = 100 erhälst Du Array im Bereich von 0-99 (nicht etwa bis 100), damit sind die gültigen Bereiche für gew[0..99], absc[0.99]. n_max wird als Parameter an gauleg übergeben und dort in einer Schleife verwendet:

gauleg( xu, xo, absc, gew, n_max); // n_max=100

in gauleg verwendest Du folgende Zuweisungen:

x[n+1-i]=xm+xl*z;

w[n+1-i]=w;

Für den ersten Schleifendurchlauf (also mit i=1, und n=100), erhälst Du als Index 100:

n+1-i = 100+1-1 und greifst auf x[100] bzw. w[100] und damit auf nicht reservierten Speicher zu.

Nic

Danke, jetzt is es mir klar!

Gruß,

Mesrine