4. problem mit function

problem mit function

  • K
     

    //temporaere Impulse 
        temp_p[i][0] = pold[i][0] + delta/2 * F1_p[i][0]; 
        temp_p[i][1] = pold[i][1] + delta/2 * F1_p[i][1]; 
        temp_p[i][2] = pold[i][2] + delta/2 * F1_p[i][2]; 

        temp_omz = omzold + delta/2 * F1_omz; 
//p(t+delta/2) 

/*x-Kompon.*/    F2_p[i][0] = temp_F[i][0] - gam * temp_p[i][1] 
                 - zet(temp_F,temp_p) * temp_p[i][0] 
                 - nustir/temp * Fdotp(temp_F,temp_p) * temp_omz * temp_p[i][1];

    ich denke ich weiß jetzt, wo der hund begraben liegt: du du benutzt uninitialisierte arrayfelder.

    du denkst warscheinlich, du kannst zet(temp_F,temp_p) berechnen, weil die arrays oben bereits initialisiert worden sind. das ist falsch! du initialisierst tatsächlich nur ein feld von temp_p! die felder temp_p[2][] werden erst in den nächsten schleifendurchläufen berechnet. und temp_p[0][] berechnest du überhaupt nicht!

    tut mir leid, aber der ganze code ist murks 😃 kannst du weitgehend neu schreiben.

    0
  • *

    😡 oh ja indem ich folgendes tue

    for (i=start; i<particles;i++)
		{

//r(t+delta/2)
		F2_r[i][0] = (pold[i][0]+ delta/2 * F1_p[i][0]) + gam * (rold[i][1]+ delta/2 * F1_r[i][1]) ;
		F2_r[i][1] = pold[i][1] + delta/2 * F1_p[i][1];
		F2_r[i][2] = pold[i][2] + delta/2 * F1_p[i][2];

		//temporaere Teilchenkraft 
		temp_F[i][0] = F[i][0]+ delta/2 * tempKraft[i][0];
		temp_F[i][1] = F[i][1]+ delta/2 * tempKraft[i][1];
		temp_F[i][2] = F[i][2]+ delta/2 * tempKraft[i][2];
		//temporaere Impulse
		temp_p[i][0] = pold[i][0] + delta/2 * F1_p[i][0]; 
		temp_p[i][1] = pold[i][1] + delta/2 * F1_p[i][1];
		temp_p[i][2] = pold[i][2] + delta/2 * F1_p[i][2];

		temp_omz = omzold + delta/2 * F1_omz;
//p(t+delta/2)

		F2_p[i][0] = temp_F[i][0] - gam * temp_p[i][1] 
			     - zet(temp_F,temp_p) * temp_p[i][0] 
			     - nustir/temp * Fdotp(temp_F,temp_p) * temp_omz * temp_p[i][1];
		F2_p[i][1] = temp_F[i][1]			- zet(temp_F,temp_p) * temp_p[i][1] 
							+ nustir/temp * Fdotp(temp_F,temp_p) * temp_omz * temp_p[i][0];
		F2_p[i][2] = temp_F[i][2]			- zet(temp_F,temp_p) * temp_p[i][2] ;

//omz(t+delta/2)
		F2_omz = omz_dt(temp_F,temp_p);

		}

    ... mache ich tatsächlich einen riesen Fehler ... die Größen zet, Fdotp und omz_dt hängen natürlich von allen temp_F und temp_p ab. Ich benötige zuerst das gesamte Array von temp_F und temp_p um F2_p berechnen zu können. Alles mit 0 zu initialiesieren hat deshalb geholfen, weil die 0 eine gutartige Zahl ist, die die Beiträge, die ich falsch berechne einfach löscht. Ohne Initialisierung mit 0 stand in diesen Arrays zu Beginn nur Datenmüll, was das sofortige Auftreten von nan erklärt.

    Das ist alles ein Copy and paste Fehler denn bei der Berechnung von F1_p müsste es wie oben gut gehen, da ich dort nur rold, pold und F benötige, die als vollständig gefüllte Arrays vorliegen.

    Aber dass der gesamte Code murks ist würd ich jetzt nicht behaupten. Meine Änderung sollte die folgende sein:

    //Kraft am Ort F1_r ausrechnen	
	Gesamtkraft(F1_r,tempKraft); //tempKraft liegt als vollständiger Satz vor :)

	//Schleife über alle Teilchen	
	for (i=start; i<particles;i++)
		{

//alle F1_p, und F1_r, sowie pold liegen vollständig vor

		F2_r[i][0] = (pold[i][0]+ delta/2 * F1_p[i][0]) + gam * (rold[i][1]+ delta/2 * F1_r[i][1]) ;
		F2_r[i][1] = pold[i][1] + delta/2 * F1_p[i][1];
		F2_r[i][2] = pold[i][2] + delta/2 * F1_p[i][2];

//auch alle tempKraft habe ich		
		//temporaere Teilchenkraft 
		temp_F[i][0] = F[i][0]+ delta/2 * tempKraft[i][0];
		temp_F[i][1] = F[i][1]+ delta/2 * tempKraft[i][1];
		temp_F[i][2] = F[i][2]+ delta/2 * tempKraft[i][2];
		//temporaere Impulse
		temp_p[i][0] = pold[i][0] + delta/2 * F1_p[i][0]; 
		temp_p[i][1] = pold[i][1] + delta/2 * F1_p[i][1];
		temp_p[i][2] = pold[i][2] + delta/2 * F1_p[i][2];
//F1_omz gibt es aus dem vorherigen Schritt, ist auch kein Array		
		temp_omz = omzold + delta/2 * F1_omz;

		} //hier endet die Schleife über alle Teilchen

//hier beginnt wieder Schleife über alle Teilchen	
// ich besitze temp_F und temp_p vollständig, da schleife oben wieder geschlossen
//ich arbeite also mit den physikalisch richtigen Werten
	for (i=start; i<particles;i++)
		{

		F2_p[i][0] = temp_F[i][0] - gam * temp_p[i][1] 
			     - zet(temp_F,temp_p) * temp_p[i][0] 
			     - nustir/temp * Fdotp(temp_F,temp_p) * temp_omz * temp_p[i][1];
		F2_p[i][1] = temp_F[i][1]			- zet(temp_F,temp_p) * temp_p[i][1] 
							+ nustir/temp * Fdotp(temp_F,temp_p) * temp_omz * temp_p[i][0];
		F2_p[i][2] = temp_F[i][2]			- zet(temp_F,temp_p) * temp_p[i][2] ;

//omz(t+delta/2)
		F2_omz = omz_dt(temp_F,temp_p);

		}

    Danke für den Tip ... ich wär da so schnell nicht drauf gekommen ... mal schauen ob's dann geht ...

    jesus (schämt sich ein wenig ... )

    edit: temp_p[0] und so brauche ich eigentlich nicht, ich setze auf das Teilchen 0 meinen Koordinatenursprung. In meinem mitbewegte Koordinatensystem dürfen keine Kräfte auf Teilchen 0 wirken (sonst wär das Koordinatensystem nicht fest auf 0 ) deshalb durchlaufe ich alle Schleifen mit dem Startwert 1 (#define start 1). Ich hab das so gewählt, weil ich später auch mal die Dynamik des Systems gegenüber einem festen Beobachtungsort sehen wollte, dann muss ich nur start auf 0 setzen und schon wird alles für Teilchen 0 mitberechnet. solnage ich das nciht tue bleiben alle Parameter für Teilchen 0 wie sie anfangs initialisiert wurden

    0
Anmelden zum Antworten