summaryrefslogtreecommitdiff
diff options
context:
space:
mode:
Diffstat
-rw-r--r--Physikunit.pas1016
-rwxr-xr-xPlasmapropagationbin3344888 -> 3387904 bytes
-rw-r--r--Plasmapropagation.lpr73
-rw-r--r--Plasmapropagation.lps186
-rw-r--r--error521
-rw-r--r--input.plap16
6 files changed, 750 insertions, 1062 deletions
diff --git a/Physikunit.pas b/Physikunit.pas
index 71e6a79..492866c 100644
--- a/Physikunit.pas
+++ b/Physikunit.pas
@@ -2,75 +2,99 @@ unit Physikunit;
{$mode objfpc}{$H+}
-{ $DEFINE Zeitschrittueberwachung}
-{ $DEFINE Dichteueberwachung}
+{$DEFINE Zeitschrittueberwachung}
+{$DEFINE Dichteueberwachung}
interface
uses
- Classes, SysUtils, Math, protokollunit, matheunit, mystringlistunit, lowlevelunit;
+ Classes, SysUtils, Math, protokollunit, matheunit, mystringlistunit, lowlevelunit, crt;
type
- TZeitverfahren = (zfEulerVorwaerts,zfRungeKuttaVier);
- TVerteilungsfunktion = function(x: extended): extended;
- TFeldInhalt = (
- fiA,fiAX,fiAY,fiAZ,
- fidAXDT,fidAYDT,fidAZDT,
- fiN,fidPhiDX,fidPsiDX,
- fiP,fiPX,fiPY,fiPZ,
- fiGamma,fiiGamma);
-
- TAusgabeDatei = record
+ tZeitverfahren = (zfEulerVorwaerts,zfRungeKuttaVier);
+ tVerteilungsfunktion = function(x: extended): extended;
+ tEMFeldInhalt = (
+ efA,efAX,efAY,efAZ,
+ efDAXDT,efDAYDT,efDAZDT,
+ efDPhiDX
+ );
+ tMaterieFeldInhalt = (
+ mfN,mfDPsiDX,
+ mfP,mfPX,mfPY,mfPZ,
+ mfGamma,mfIGamma
+ );
+
+ { tAusgabeDatei }
+
+ tAusgabeDatei = class
ableitung: boolean;
- inhalt: tFeldInhalt;
+ emInhalt: tEMFeldInhalt;
+ matInhalt: tMaterieFeldInhalt;
+ teilchen: longint; // wenn < 0, dann EM-Feld
name: string;
datei: file;
+ constructor create(feldName,prefix,suffix: string; prot: tProtokollant);
+ destructor destroy; override;
end;
- tFelder = class;
+ tGitter = class;
+
+ { tWertePunkt }
- { tWerteGitter }
+ tWertePunkt = class(tObject) // repräsentiert alle Werte eines Punktes im Gitter und deren Zeitableitungen
+ matWerte: array of array[tMaterieFeldInhalt,boolean] of extended;
+ // Materiefelder (getrennt für verschiedene Teilchenspezies) und deren Zeitableitungen
+ emWerte: array[tEMFeldInhalt,boolean] of extended;
+ // EM-Felder und deren Zeitableitungen
+ // A, p[xyz]? und i?Gamma haben keine sinnvolle Ableitung hier!
+ lN,rN: tWertePunkt;
- tWerteGitter = class(tObject) // repräsentiert ein Gitter von Werten und deren Zeitableitungen
- werte: array[boolean] of array of extended;
- par: tFelder;
+ constructor create(linkerNachbar: tWertePunkt; materien: longint);
+ destructor destroy; override;
+ procedure berechneGammaUndP;
+ procedure berechneNAbleitung(iDX,pDNMax: extended);
+ procedure berechneAbleitungen(dT,dX,iDX: extended);
+ procedure liKo(in1,in2: tWertePunkt; fak2: extended); overload; // Werte werden auf (in1 + fak2*in2') gesetzt
+ procedure liKo(in1,in2,in3,in4,in5: tWertePunkt; fak2,fak3,fak4,fak5: extended); overload; // Werte werden auf (in1 + \sum_i faki * ini') gesetzt
end;
{ TFelder }
tFelder = class(tObject) // repräsentiert eine Simulationsbox von Feldern inklusive deren Ableitungen
+ private
+ matAnz: longint;
+ lichters: array[boolean] of tMyStringlist;
+ procedure setzeRaender(iDX: extended);
public
- materieFelder: array of tWerteGitter; // Materiefelder (getrennt für verschiedene Teilchenspezies) und deren Zeitableitungen
- elektroMagnetischeFelder: tWerteGitter; // EM-Felder und deren Zeitableitungen
- // A, p[xyz]? und i?Gamma haben keine sinnvolle Ableitung hier!
- dX,iDX,x,pDNMax: extended;
+ inhalt: array of tWertePunkt;
+ par: tGitter;
- constructor create(deltaX,pDNMa: extended);
- procedure berechneGammaUndP;
- procedure berechneNAbleitung;
- procedure berechneAbleitungen(dT: extended);
+ constructor create(groesse: longint; dichten,lichter: tMyStringList; parent: tGitter);
+ destructor destroy; override;
+ procedure berechneAbleitungen(dT,dX,iDX,pDNMax: extended);
+ procedure liKo(in1,in2: tFelder; fak2: extended); overload; // Werte werden auf (in1 + fak2*in2') gesetzt
+ procedure liKo(in1,in2,in3,in4,in5: tFelder; fak2,fak3,fak4,fak5: extended); overload; // Werte werden auf (in1 + \sum_i faki * ini') gesetzt
end;
{ tGitter }
tGitter = class(tObject)
private
- groesse: longint;
- prot: tProtokollant;
- procedure berechneAbleitungen(Felder: longint; dt: extended); overload;
- procedure berechneAbleitungen(Felder: longint; dt: extended; out dTMax: extended); overload;
- procedure setzeRaender(Felder: longint);
+ prot: tProtokollant;
+ zeitverfahren: tZeitverfahren;
+ pDNMaxDynamisch: boolean;
+ dX,iDX,pDNMax,dTMaximum,xl,t: extended;
+ kvs: tKnownValues;
+ procedure diffusionsTermAnpassen(pro: tProtokollant);
public
+
aktuelleFelder: longint;
felders: array of tFelder; // mehrere komplette Simulationsboxen von Feldern, benötigt um Zwischenschritte für die Zeitentwicklung zu speichern
- dTMaximum,dX,iDX,xl,xr,t,pDNMax: extended;
- pDNMaxDynamisch: boolean;
- Zeitverfahren: TZeitverfahren;
- Al,Ar: TVerteilungsfunktion;
- constructor create(size: longint; deltaT,deltaX,pDNMa: extended; dichten, lichter: tMyStringlist; zv: tZeitverfahren; protokollant: tProtokollant);
+
+ constructor create(size: longint; deltaT,deltaX,pDNMa: extended; bekannteWerte: tKnownValues; dichten, lichter: tMyStringlist; zv: tZeitverfahren; protokollant: tProtokollant);
destructor destroy; override;
- procedure iteriereSchritt(var dT: extended);
- procedure macheAusgabe(AD: TAusgabedatei; sDX: extended);
+ procedure iteriereSchritt(dT: extended);
+ procedure macheAusgabe(aD: tAusgabedatei; sDX: extended);
function gibErhaltungsgroessen: string;
end;
@@ -78,59 +102,157 @@ type
tSimulation = class(tObject)
private
- prot: tProtokollant;
- gitter: tGitter;
- kvs: tKnownValues;
+ prot: tProtokollant;
+ gitter: tGitter;
+ dT,tEnde,sT,sDT,sDX: extended;
+ fortschrittsAnzeige: boolean;
+ ausgabeDateien: array of tAusgabeDatei;
public
constructor create(inName: string; Protokollant: tProtokollant);
destructor destroy; override;
+ function iteriereSchritt(start: extended; var zeitPhysik,zeitDatei: extended): boolean; // noch nicht zu Ende?
end;
+implementation
+
const
- OptionNamen: array[TFeldInhalt] of string = (
- 'A','AX','AY','AZ','dAXDT','dAYDT','dAZDT','N','dPhiDX','dPsiDX','P','PX','PY','PZ','Gamma','iGamma');
- AusgabeNamen: array[TFeldInhalt] of string = (
- 'fiA','fiAX','fiAY','fiAZ','fidAXDT','fidAYDT','fidAZDT','fiN','fidPhiDX','fidPsiDX','fiP','fiPX','fiPY','fiPZ','fiGamma','fiiGamma');
- AusgabeHilfe = '--d?A[XYZ]? --d?dA[XYZ]dt --d?N --d?dPhiDX --d?dPsiDX --P[XYZ]? --Gamma --iGamma';
+ emFeldNamen: array[tEMFeldInhalt] of string = (
+ 'A','AX','AY','AZ','DAXDT','DAYDT','DAZDT','DPHIDX'
+ );
+ matFeldNamen: array[tMaterieFeldInhalt] of string = (
+ 'N','DPSIDX','P','PX','PY','PZ','GAMMA','IGAMMA'
+ );
-implementation
+{ tAusgabeDatei }
+
+constructor tAusgabeDatei.create(feldName,prefix,suffix: string; prot: tProtokollant);
+var
+ emF: tEMFeldInhalt;
+ maF: tMaterieFeldInhalt;
+ num: longint;
+ abl,gef: boolean;
+begin
+ inherited create;
+ num:=length(feldName);
+ while (num>0) and (feldName[num] in ['0'..'9']) do
+ dec(num);
+ inc(num);
+ if num<=length(feldName) then begin
+ teilchen:=strtoint(copy(feldName,num,length(feldName)));
+ delete(feldName,num,length(feldName));
+ end
+ else
+ teilchen:=0;
+
+ emInhalt:=efA;
+ matInhalt:=mfN;
+ feldName:=ansiUpperCase(feldName);
+
+ gef:=false;
+ for abl:=false to true do begin
+ if not gef then
+ for emF:=low(tEMFeldInhalt) to high(tEMFeldInhalt) do
+ if copy('D',1,byte(abl))+emFeldNamen[emF] = feldName then begin
+ emInhalt:=emF;
+ teilchen:=-1;
+ gef:=true;
+ break;
+ end;
+ if not gef then
+ for maF:=low(tMaterieFeldInhalt) to high(tMaterieFeldInhalt) do
+ if copy('D',1,byte(abl))+matFeldNamen[maF] = feldName then begin
+ matInhalt:=maF;
+ gef:=true;
+ break;
+ end;
+ end;
+
+ if not gef then begin
+ prot.schreibe('tAusgabeDatei.create kennt Feld '''+feldName+''' nicht!',true);
+ halt(1);
+ end;
+
+ if teilchen>=0 then
+ name:=matFeldNamen[maF]+inttostr(teilchen)
+ else
+ name:=emFeldNamen[emF];
-{ TFeld }
+ if abl then
+ name:='d'+name;
+ name:=prefix+name+suffix;
-constructor TFeld.create(deltaX,pDNMa: extended);
-var FI: TFeldinhalt;
- Abl: Boolean;
+ assignFile(datei,name);
+ rewrite(datei);
+end;
+
+destructor tAusgabeDatei.destroy;
+begin
+ closeFile(datei);
+ inherited destroy;
+end;
+
+{ tWertePunkt }
+
+constructor tWertePunkt.create(linkerNachbar: tWertePunkt; materien: longint);
+var
+ emF: tEMFeldInhalt;
+ maF: tMaterieFeldInhalt;
+ abl: boolean;
+ i: longint;
begin
inherited create;
- lN:=nil;
+ fillchar(matWerte,sizeof(matWerte),#0);
+ setlength(matWerte,materien);
+ for i:=0 to length(matWerte)-1 do begin
+ for maF:=low(tMaterieFeldInhalt) to high(tMaterieFeldInhalt) do
+ for abl:=false to true do
+ matWerte[i,maF,abl]:=0;
+ matWerte[i,mfGamma,false]:=0;
+ matWerte[i,mfIGamma,false]:=0;
+ end;
+ for emF:=low(tEMFeldInhalt) to high(tEMFeldInhalt) do
+ for abl:=false to true do
+ emWerte[emF,abl]:=0;
+ lN:=linkerNachbar;
rN:=nil;
- for Abl:=False to True do
- for FI:=low(TFeldinhalt) to high(TFeldinhalt) do
- Groessen[Abl,FI]:=0;
- dX:=deltaX;
- iDX:=1/dX;
- x:=0;
- pDNMax:=pDNMa;
- Groessen[false,fiN]:=1;
- Groessen[false,fiGamma]:=1;
- Groessen[false,fiiGamma]:=1;
+ if assigned(lN) then
+ lN.rN:=self;
+end;
+
+destructor tWertePunkt.destroy;
+begin
+ setlength(matWerte,0);
+ if assigned(lN) then
+ lN.rN:=nil;
+ if assigned(rN) then
+ rN.lN:=nil;
+ inherited destroy;
end;
-procedure TFeld.berechneGammaUndP; // aus aktuellem dPsi/dX und A
-var b: byte;
- tmp: extended;
+procedure tWertePunkt.berechneGammaUndP; // aus aktuellem dPsi/dX und A
+var
+ i: longint;
begin
- tmp:=0;
- for b:=0 to 2 do begin
- Groessen[false,TFeldinhalt(Byte(fiPX)+b)]:=
- Groessen[false,TFeldinhalt(Byte(fiAX)+b)] + Groessen[false,fidPsiDX]*Byte(b=0);
- tmp:=tmp+sqr(Groessen[false,TFeldinhalt(Byte(fiPX)+b)]);
+ for i:=0 to length(matWerte)-1 do begin
+ matWerte[i,mfPX,false]:= emWerte[efAX,false] + matWerte[i,mfDPsiDX,false];
+ matWerte[i,mfPY,false]:= emWerte[efAY,false];
+ matWerte[i,mfPZ,false]:= emWerte[efAZ,false];
+
+ matWerte[i,mfGamma,false]:=
+ sqrt(
+ 1 +
+ sqr(matWerte[i,mfPX,false]) +
+ sqr(matWerte[i,mfPY,false]) +
+ sqr(matWerte[i,mfPZ,false])
+ );
+ matWerte[i,mfIGamma,false]:=
+ 1/matWerte[i,mfGamma,false];
end;
- Groessen[false,fiGamma]:=sqrt(1+tmp);
- Groessen[false,fiiGamma]:=1/Groessen[false,fiGamma];
end;
-procedure TFeld.berechneNAbleitung; // Zeitableitung von n berechnen
+procedure tWertePunkt.berechneNAbleitung(iDX,pDNMax: extended); // Zeitableitung von n berechnen
+var
+ i: longint;
begin
(* // dn/dt = -d(n/Gamma p_x)/dx
Groessen[true,fiN]:=
@@ -141,104 +263,226 @@ begin
)*iDX/6;
// es wird über die beiden nächsten linken bzw. rechten Nachbarn gemittelt *)
- // dn/dt = -d(n/Gamma p_x)/dx + (p_max*dx) * Laplace n
- Groessen[true,fiN]:=
- ( rN^.Groessen[false,fiN]*rN^.Groessen[false,fiiGamma]*(pDNMax - rN^.Groessen[false,fiPX])
- - lN^.Groessen[false,fiN]*lN^.Groessen[false,fiiGamma]*(-pDNMax - lN^.Groessen[false,fiPX])
- - Groessen[false,fiN]*Groessen[false,fiiGamma]*2*pDNMax)*iDX/2
- // Der zweite Summand entspricht (in etwa) einer thermischen Verschmierung
- // und soll die numerische Stabilität bis zu p * dn/dx / n von 2*pDNMax gewährleisten.
- // Es handelt sich um einen Diffusionsterm dn/dt = alpha * Laplace n mit
- // alpha = pDNMax * dX
+ for i:=0 to length(matWerte)-1 do
+ // dn/dt = -d(n/Gamma p_x)/dx + (p_max*dx) * Laplace n
+ matWerte[i,mfN,true]:=
+ ( rN.matWerte[i,mfN,false] * rN.matWerte[i,mfIGamma,false] * (pDNMax - rN.matWerte[i,mfPX,false])
+ - ln.matWerte[i,mfN,false] * ln.matWerte[i,mfIGamma,false] * (-pDNMax - rN.matWerte[i,mfPX,false])
+ - matWerte[i,mfN,false] * matWerte[i,mfIGamma,false] * 2 * pDNMax) * iDX/2;
+ // Der zweite Summand entspricht (in etwa) einer thermischen Verschmierung
+ // und soll die numerische Stabilität bis zu p * dn/dx / n von 2*pDNMax gewährleisten.
+ // Es handelt sich um einen Diffusionsterm dn/dt = alpha * Laplace n mit
+ // alpha = pDNMax * dX
end;
-procedure TFeld.berechneAbleitungen(dT: extended); // Zeitableitungen berechnen
-var b: byte;
+procedure tWertePunkt.berechneAbleitungen(dT,dX,iDX: extended); // Zeitableitungen berechnen
+var
+ i: longint;
begin
// d2Psi/dxdt = dPhi/dx - dGamma/dx
- Groessen[true,fidPsiDX]:=
- Groessen[false,fidPhiDX]
- - (rN^.Groessen[false,fiGamma] - lN^.Groessen[false,fiGamma])*iDX/2;
+ for i:=0 to length(matWerte)-1 do
+ matWerte[i,mfDPsiDX,true]:=
+ emWerte[efDPhiDX,false]
+ - (rN.matWerte[i,mfGamma,false] - lN.matWerte[i,mfGamma,false]) * iDX/2;
+
// d3Phi/dx2dt = dn/dt ... hier muss integriert werden:
// d2Phi/dxdt = d2Phi/dxdt(x- deltax) + <dn/dt> * deltax
- Groessen[true,fidPhiDX]:=
- lN^.Groessen[true,fidPhiDX]
- + (lN^.Groessen[true,fiN] + Groessen[true,fiN])*dX/2;
- // d2A/dt2 = - n/gamma p + Laplace(A) ...
- for b:=0 to 2 do
- Groessen[true,TFeldinhalt(Byte(fidAXDT)+b)]:=
- - (Groessen[false,fiN]*Groessen[false,fiiGamma]*Groessen[false,TFeldinhalt(Byte(fiPX)+b)])
- + (rN^.Groessen[false,TFeldinhalt(Byte(fiAX)+b)]
- - 2*Groessen[false,TFeldinhalt(Byte(fiAX)+b)]
- + lN^.Groessen[false,TFeldinhalt(Byte(fiAX)+b)])*sqr(iDX);
- Groessen[true,fidAXDT]:=Groessen[true,fidAXDT] - Groessen[false,fidPhiDX]; // ... - Nabla(phi)
+ emWerte[efDPhiDX,true]:=
+ lN.emWerte[efDPhiDX,true];
+ for i:=0 to length(matWerte)-1 do
+ emWerte[efDPhiDX,true]:=
+ emWerte[efDPhiDX,true]
+ + (lN.matWerte[i,mfN,true] + matWerte[i,mfN,true])*dX/2;
+
+ // d2A/dt2 = Laplace(A) - Nabla(phi) ...
+ emWerte[efDAXDT,true]:=
+ ( rN.emWerte[efAX,false] - 2*emWerte[efAX,false] + lN.emWerte[efAX,false] )*sqr(iDX)
+ - emWerte[efDPhiDX,false];
+ emWerte[efDAYDT,true]:=
+ ( rN.emWerte[efAY,false] - 2*emWerte[efAY,false] + lN.emWerte[efAY,false] )*sqr(iDX);
+ emWerte[efDAZDT,true]:=
+ ( rN.emWerte[efAZ,false] - 2*emWerte[efAZ,false] + lN.emWerte[efAZ,false] )*sqr(iDX);
+ // ... - n/gamma p
+ for i:=0 to length(matWerte)-1 do begin
+ emWerte[efDAXDT,true]:=
+ emWerte[efDAXDT,true]
+ - (matWerte[i,mfN,false] * matWerte[i,mfIGamma,false] * matWerte[i,mfPX,false]);
+ emWerte[efDAYDT,true]:=
+ emWerte[efDAYDT,true]
+ - (matWerte[i,mfN,false] * matWerte[i,mfIGamma,false] * matWerte[i,mfPY,false]);
+ emWerte[efDAZDT,true]:=
+ emWerte[efDAZDT,true]
+ - (matWerte[i,mfN,false] * matWerte[i,mfIGamma,false] * matWerte[i,mfPZ,false]);
+ end;
+
// dA/dt = dA/dt
- for b:=0 to 2 do
- Groessen[true,TFeldInhalt(Byte(fiAX)+b)]:=
- Groessen[false,TFeldInhalt(Byte(fidAXDT)+b)] +
- Groessen[true,TFeldInhalt(Byte(fidAXDT)+b)]*dT;
+ emWerte[efAX,true]:=emWerte[efDAXDT,false] + emWerte[efDAXDT,true]*dT;
end;
-{ TGitter }
+procedure tWertePunkt.liKo(in1,in2: tWertePunkt; fak2: extended); // Werte werden auf (in1 + fak2*in2') gesetzt
+var
+ emF: tEMFeldInhalt;
+ maF: tMaterieFeldInhalt;
+ i: longint;
+begin
+(* tEMFeldInhalt = (
+ efA,efAX,efAY,efAZ,
+ efDAXDT,efDAYDT,efDAZDT,
+ efDPhiDX
+ ); *)
+ for emF:=efAX to efDPhiDX do // alles außer efA, welchen Ableitung ja nicht berechnet wurde
+ emWerte[emF,false]:= in1.emWerte[emF,false] + fak2 * in2.emWerte[emF,true];
+
+(* tMaterieFeldInhalt = (
+ mfN,mfDPsiDX,
+ mfP,mfPX,mfPY,mfPZ,
+ mfGamma,mfIGamma
+ ); *)
+ for i:=0 to length(matWerte)-1 do // siehe oben
+ for maF:=mfN to mfDPsiDX do
+ matWerte[i,maF,false]:= in1.matWerte[i,maF,false] + fak2 * in2.matWerte[i,maF,true];
+end;
-procedure TGitter.berechneAbleitungen(Felder: longint; dT: extended; out dTMax: extended);
-var i: longint;
+procedure tWertePunkt.liKo(in1,in2,in3,in4,in5: tWertePunkt; fak2,fak3,fak4,fak5: extended); // Werte werden auf (in1 + \sum_i faki * ini') gesetzt
+var
+ emF: tEMFeldInhalt;
+ maF: tMaterieFeldInhalt;
+ i: longint;
begin
- berechneAbleitungen(Felder,dT);
- dTMax:=dTMaximum;
- for i:=0 to groesse+3 do
- if Felders[Felder,i].Groessen[true,fiN]<0 then begin
- dTMax:=min(dTMax,Abs(Felders[Felder,i].Groessen[false,fiN]/Felders[Felder,i].Groessen[true,fiN]));
- end;
+(* tEMFeldInhalt = (
+ efA,efAX,efAY,efAZ,
+ efDAXDT,efDAYDT,efDAZDT,
+ efDPhiDX
+ ); *)
+ for emF:=efAX to efDPhiDX do // alles außer efA, welchen Ableitung ja nicht berechnet wurde
+ emWerte[emF,false]:=
+ in1.emWerte[emF,false]
+ + fak2 * in2.emWerte[emF,true]
+ + fak3 * in3.emWerte[emF,true]
+ + fak4 * in4.emWerte[emF,true]
+ + fak5 * in5.emWerte[emF,true];
+
+(* tMaterieFeldInhalt = (
+ mfN,mfDPsiDX,
+ mfP,mfPX,mfPY,mfPZ,
+ mfGamma,mfIGamma
+ ); *)
+ for i:=0 to length(matWerte)-1 do // siehe oben
+ for maF:=mfN to mfDPsiDX do
+ matWerte[i,maF,false]:=
+ in1.matWerte[i,maF,false]
+ + fak2 * in2.matWerte[i,maF,true]
+ + fak3 * in3.matWerte[i,maF,true]
+ + fak4 * in4.matWerte[i,maF,true]
+ + fak5 * in5.matWerte[i,maF,true];
end;
-procedure TGitter.berechneAbleitungen(Felder: longint; dT: extended);
-var i: longint;
+{ tFelder }
+
+constructor tFelder.create(groesse: longint; dichten,lichter: tMyStringList; parent: tGitter);
+var
+ i,j: longint;
+ rechts: boolean;
begin
- for i:=0 to groesse+3 do begin
-(* Felders[Felder,i].Groessen[false,fiN]:=
- max(Felders[Felder,i].Groessen[false,fiN],0); // n >= 0 *)
- Felders[Felder,i].berechneGammaUndP;
+ inherited create;
+ par:=parent;
+ matAnz:=dichten.count;
+ fillchar(inhalt,sizeof(inhalt),#0);
+ setlength(inhalt,groesse+4); // zwei Felder links und rechts extra für Randbedingungen
+ inhalt[0]:=tWertePunkt.create(nil,matAnz);
+ for i:=1 to length(inhalt)-1 do
+ inhalt[i]:=tWertePunkt.create(inhalt[i-1],matAnz);
+ for i:=0 to length(inhalt)-1 do begin
+ parent.kvs.add('x',par.xl+(i-2)/groesse*par.dX);
+ for j:=0 to dichten.count-1 do
+ inhalt[i].matWerte[j,mfN,false]:=exprToFloat(false,dichten[j],parent.kvs,nil);
+ end;
+ parent.kvs.rem('x');
+ for rechts:=false to true do begin
+ lichters[rechts]:=tMyStringlist.create(nil);
+ lichters[rechts].text:=lichter.text;
end;
- Felders[Felder,1].Groessen[false,fiPX]:= // Teilchen werden reflektiert
- Abs(Felders[Felder,1].Groessen[false,fiPX]);
- Felders[Felder,groesse+2].Groessen[false,fiPX]:=
- -Abs(Felders[Felder,groesse+2].Groessen[false,fiPX]);
- setzeRaender(Felder);
- for i:=1 to groesse+2 do
- Felders[Felder,i].berechneNAbleitung;
- for i:=1 to groesse+2 do
- Felders[Felder,i].berechneAbleitungen(dT);
+ lichters[false].grep('^links .*');
+ lichters[false].subst('^links *','');
+ lichters[true].grep('^rechts .*');
+ lichters[true].subst('^rechts *','');
end;
-procedure TGitter.setzeRaender(Felder: longint);
-var FI: TFeldInhalt;
+destructor tFelder.destroy;
begin
- for FI:=fiAX to fiAZ do begin // Vakuumrandbedingungen für das A-Feld
- Felders[Felder,0].Groessen[true,FI]:=
- (Felders[Felder,1].Groessen[false,FI] -
- Felders[Felder,0].Groessen[false,FI])*iDX;
- Felders[Felder,groesse+3].Groessen[true,FI]:=
- (Felders[Felder,groesse+2].Groessen[false,FI] -
- Felders[Felder,groesse+3].Groessen[false,FI])*iDX;
- end; // (ein bisschen wird noch reflektiert, vmtl. durch numerische Fehler)
- Felders[Felder,0].Groessen[true,fiAy]:=
- Felders[Felder,0].Groessen[true,fiAy]
+ lichters[false].free;
+ lichters[true].free;
+ inherited destroy;
+end;
+
+procedure tFelder.setzeRaender(iDX: extended);
+var
+ emF: tEMFeldInhalt;
+begin
+ for emF:=efAX to efAZ do begin // Vakuumrandbedingungen für das A-Feld
+ inhalt[0].emWerte[emF,true]:=
+ (inhalt[1].emWerte[emF,false] -
+ inhalt[0].emWerte[emF,false])*iDX;
+ inhalt[length(inhalt)-1].emWerte[emF,true]:=
+ (inhalt[length(inhalt)-2].emWerte[emF,false] -
+ inhalt[length(inhalt)-1].emWerte[emF,false])*iDX;
+ end; // (ein bisschen wird noch reflektiert, vmtl. durch numerische Fehler bzw. nicht beachtete numerische Dispersion)
+
+ (* inhalt[0].emWerte[efAY,true]:=
+ inhalt[0].emWerte[fiAY,true]
+ (Al(t+dTMaximum)-Al(t))/dTMaximum;
Felders[Felder,groesse+3].Groessen[true,fiAy]:=
Felders[Felder,groesse+3].Groessen[true,fiAy]
- + (Ar(t+dTMaximum)-Ar(t))/dTMaximum;
+ + (Ar(t+dTMaximum)-Ar(t))/dTMaximum; *)
end;
-constructor TGitter.create(size: longint; deltaT,deltaX,pDNMa: extended; dichten, lichter: tMyStringlist; zv: tZeitverfahren; protokollant: tProtokollant);
-var i,j: longint;
+procedure tFelder.berechneAbleitungen(dt,dx,iDX,pDNMax: extended);
+var
+ i: longint;
+begin
+ for i:=0 to length(inhalt)-1 do
+ inhalt[i].berechneGammaUndP;
+ for i:=0 to matAnz-1 do begin // Teilchen werden reflektiert
+ inhalt[1].matWerte[i,mfPX,false]:=
+ abs(inhalt[1].matWerte[i,mfPX,false]);
+ inhalt[length(inhalt)-2].matWerte[i,mfPX,false]:=
+ -abs(inhalt[length(inhalt)-2].matWerte[i,mfPX,false]);
+ end;
+
+ setzeRaender(iDX);
+
+ for i:=1 to length(inhalt)-2 do
+ inhalt[i].berechneNAbleitung(iDX,pDNMax);
+ for i:=1 to length(inhalt)-2 do
+ inhalt[i].berechneAbleitungen(dT,dX,iDX);
+end;
+
+procedure tFelder.liKo(in1,in2: tFelder; fak2: extended); // Werte werden auf (in1 + fak2*in2') gesetzt
+var
+ i: longint;
+begin
+ for i:=0 to length(inhalt)-1 do
+ inhalt[i].liKo(in1.inhalt[i],in2.inhalt[i],fak2);
+end;
+
+procedure tFelder.liKo(in1,in2,in3,in4,in5: tFelder; fak2,fak3,fak4,fak5: extended); // Werte werden auf (in1 + \sum_i faki * ini') gesetzt
+var
+ i: longint;
+begin
+ for i:=0 to length(inhalt)-1 do
+ inhalt[i].liKo(in1.inhalt[i],in2.inhalt[i],in3.inhalt[i],in4.inhalt[i],in5.inhalt[i],fak2,fak3,fak4,fak5);
+end;
+
+{ tGitter }
+
+constructor tGitter.create(size: longint; deltaT,deltaX,pDNMa: extended; bekannteWerte: tKnownValues; dichten, lichter: tMyStringlist; zv: tZeitverfahren; protokollant: tProtokollant);
+var
+ i: longint;
begin
inherited create;
- Ar:=@nullfunktion;
- Al:=@nullfunktion;
- Zeitverfahren:=ZV;
- groesse:=size;
- Prot:=TProtokollant.create(Protokollant,'TGitter');
+ zeitverfahren:=ZV;
+ kvs:=bekannteWerte;
+ Prot:=TProtokollant.create(Protokollant,'tGitter');
dTMaximum:=deltaT;
dX:=deltaX;
iDX:=1/dX;
@@ -252,252 +496,192 @@ begin
zfRungeKuttaVier: Setlength(Felders,5);
end{of Case};
xl:=dX/2;
- for i:=0 to length(Felders)-1 do begin
- setlength(Felders[i],size+4); // auf jeder Seite je zwei Felder Rand extra
- for j:=0 to length(felders[i])-1 do begin
- felders[i,j]:=TFeld.create(dX,pDNMax);
- felders[i,j].x:=xl+j*dX;
- knownValues.add('x',felders[i,j].x);
- felders[i,j].Groessen[false,fiN]:=n0(felders[i,j].x);
- felders[i,j].Groessen[false,fidPsiDX]:=0.0001*Cos(j/(length(felders[i])-1)*pi*15);
- end;
- felders[i,0].lN:=@felders[i,0];
- for j:=1 to length(felders[i])-1 do begin
- felders[i,j].lN:=@felders[i,j-1];
- felders[i,j-1].rN:=@felders[i,j];
- end;
- felders[i,length(felders[i])-1].rN:=@felders[i,length(felders[i])-1];
- end;
- xr:=xl+dX*(length(Felders[0])-1);
+
+ for i:=0 to length(felders)-1 do
+ felders[i]:=tFelder.create(size,dichten,lichter,self);
+
aktuelleFelder:=0;
t:=0;
end;
-destructor TGitter.destroy;
-var i,j: longint;
+destructor tGitter.destroy;
+var
+ i: longint;
begin
- for i:=0 to length(Felders)-1 do begin
- for j:=0 to length(Felders[i])-1 do
- Felders[i,j].free;
- setlength(Felders[i],0);
- end;
- setlength(Felders,0);
+ for i:=0 to length(Felders)-1 do
+ felders[i].free;
+ setlength(felders,0);
inherited destroy;
end;
-procedure TGitter.iteriereSchritt(var dT: extended);
+procedure tGitter.iteriereSchritt(dT: extended);
var i: longint;
- FI: TFeldInhalt;
- dTMax: extended;
{$IFDEF Zeitschrittueberwachung}
- Pro: TProtokollant;
+ pro: tProtokollant;
+ j: longint;
{$ENDIF}
-const zeitObergrenze = 1/4;
- zeitSollwert = 1/8;
- zeitUntergrenze = 1/16;
begin
- berechneAbleitungen(aktuelleFelder,dT,dTMax); // y' = y'(t,y(t))
- {$IFDEF Zeitschrittueberwachung}
- Pro:=TProtokollant.create(Prot,'iteriereSchritt');
- {$ENDIF}
- while dT>=dTMax*zeitObergrenze do begin
- dT:=dTMax*zeitSollwert;
- berechneAbleitungen(aktuelleFelder,dT,dTMax); // y' = y'(t,y(t))
- {$IFDEF Zeitschrittueberwachung}
- Pro.schreibe('+ neues dT: '+floattostr(dT)+' (t='+floattostr(t)+')');
- {$ENDIF}
- end;
- {$IFDEF Zeitschrittueberwachung}
- if dT<1e-30 then begin
- for i:=1 to groesse+2 do
- if (Felders[aktuelleFelder,i].Groessen[true,fiN]<0) and
- (Abs(Felders[aktuelleFelder,i].Groessen[false,fiN]/
- Felders[aktuelleFelder,i].Groessen[true,fiN])<1e-15) then
- Prot.schreibe(
- floattostr(t)+' '+
- inttostr(i)+' '+
- floattostr(Felders[aktuelleFelder,i-1].Groessen[false,fidPsiDX]*
- Felders[aktuelleFelder,i-1].Groessen[false,fiiGamma])+' '+
- floattostr(Felders[aktuelleFelder,i-1].Groessen[false,fiN])+' '+
- floattostr(Felders[aktuelleFelder,i-1].Groessen[true,fiN])+' '+
- floattostr(Felders[aktuelleFelder,i].Groessen[false,fidPsiDX]*
- Felders[aktuelleFelder,i].Groessen[false,fiiGamma])+' '+
- floattostr(Felders[aktuelleFelder,i].Groessen[false,fiN])+' '+
- floattostr(Felders[aktuelleFelder,i].Groessen[true,fiN])+' '+
- floattostr(Felders[aktuelleFelder,i+1].Groessen[false,fidPsiDX]*
- Felders[aktuelleFelder,i+1].Groessen[false,fiiGamma])+' '+
- floattostr(Felders[aktuelleFelder,i+1].Groessen[false,fiN])+' '+
- floattostr(Felders[aktuelleFelder,i+1].Groessen[true,fiN]),
- true);
- Pro.destroyall;
- halt(1);
- end;
- {$ENDIF}
+ felders[aktuelleFelder].berechneAbleitungen(dT,dX,iDX,pDNMax); // y' = y'(t,y(t))
- case Zeitverfahren of
- zfEulerVorwaerts:
- begin // y(t+dt) = y(t) + y' dt
- for i:=0 to groesse+3 do
- for FI:=low(TFeldinhalt) to fidPsiDX do
- Felders[1-aktuelleFelder,i].Groessen[false,FI]:=
- Felders[aktuelleFelder,i].Groessen[false,FI]
- + Felders[aktuelleFelder,i].Groessen[true,FI]*dT;
- end;
+ case zeitverfahren of
+ zfEulerVorwaerts: // y(t+dt) = y(t) + y' dt
+ felders[1-aktuelleFelder].liKo(felders[aktuelleFelder],felders[aktuelleFelder],dT);
zfRungeKuttaVier:
begin
- for i:=0 to groesse+3 do // ya = y(t) + y' dt/2
- for FI:=low(TFeldinhalt) to fidPsiDX do
- Felders[2,i].Groessen[false,FI]:=
- Felders[aktuelleFelder,i].Groessen[false,FI]
- + Felders[aktuelleFelder,i].Groessen[true,FI]*dT/2;
- berechneAbleitungen(2,dT/2,dTMax); // ya' = y'(t+dt/2,ya)
- if dT/2>dTMax*zeitObergrenze then begin
- {$IFDEF Zeitschrittueberwachung}
- Pro.schreibe('1. '+floattostr(dT/2)+'>'+floattostr(dTMax*zeitObergrenze)+' (t='+floattostr(t)+')',true);
- Pro.free;
- {$ENDIF}
- dT:=dTMax*zeitSollwert;
- exit;
- end;
- for i:=0 to groesse+3 do // yb = y(t) + ya' dt/2
- for FI:=low(TFeldinhalt) to fidPsiDX do
- Felders[3,i].Groessen[false,FI]:=
- Felders[aktuelleFelder,i].Groessen[false,FI]
- + Felders[2,i].Groessen[true,FI]*dT/2;
- berechneAbleitungen(3,dT/2,dTMax); // yb' = y'(t+dt/2,yb)
- if dT/2>dTMax*zeitObergrenze then begin
- {$IFDEF Zeitschrittueberwachung}
- Pro.schreibe('2. '+floattostr(dT/2)+'>'+floattostr(dTMax*zeitObergrenze)+' (t='+floattostr(t)+')',true);
- Pro.free;
- {$ENDIF}
- dT:=dTMax*zeitSollwert;
- exit;
- end;
- for i:=0 to groesse+3 do // yc = y(t) + yb' dt
- for FI:=low(TFeldinhalt) to fidPsiDX do
- Felders[4,i].Groessen[false,FI]:=
- Felders[aktuelleFelder,i].Groessen[false,FI]
- + Felders[3,i].Groessen[true,FI]*dT;
- berechneAbleitungen(4,dT/2,dTMax); // yc' = y'(t+dt,yc)
- if dT/2>dTMax*zeitObergrenze then begin
- {$IFDEF Zeitschrittueberwachung}
- Pro.schreibe('3. '+floattostr(dT/2)+'>'+floattostr(dTMax*zeitObergrenze)+' (t='+floattostr(t)+')',true);
- Pro.free;
- {$ENDIF}
- dT:=dTMax*zeitSollwert;
- exit;
- end;
- for i:=0 to groesse+3 do // y(t+dt) = y(t) + (y' + 2(ya' + yb') + yc') dt/6
- for FI:=low(TFeldinhalt) to fidPsiDX do
- Felders[1-aktuelleFelder,i].Groessen[false,FI]:=
- Felders[aktuelleFelder,i].Groessen[false,FI]
- + (Felders[aktuelleFelder,i].Groessen[true,FI] +
- 2*(Felders[2,i].Groessen[true,FI] +
- Felders[3,i].Groessen[true,FI]) +
- Felders[4,i].Groessen[true,FI])*dT/6;
+ felders[2].liKo(felders[aktuelleFelder],felders[aktuelleFelder],dT/2); // ya = y(t) + y' dt/2
+ felders[2].berechneAbleitungen(dT/2,dX,iDX,pDNMax); // ya' = y'(t+dt/2,ya)
+
+ felders[3].liKo(felders[aktuelleFelder],felders[2],dT/2); // yb = y(t) + ya' dt/2
+ felders[3].berechneAbleitungen(dT/2,dX,iDX,pDNMax); // yb' = y'(t+dt/2,yb)
+
+ felders[4].liKo(felders[aktuelleFelder],felders[3],dT); // yc = y(t) + yb' dt
+ felders[4].berechneAbleitungen(dT/2,dX,iDX,pDNMax); // yc' = y'(t+dt,yc)
+
+ felders[1-aktuelleFelder].liKo(
+ felders[aktuelleFelder],
+ felders[aktuelleFelder],
+ felders[2],
+ felders[3],
+ felders[4],
+ dT/6,
+ dT/3,
+ dT/3,
+ dT/6
+ ); // y(t+dt) = y(t) + (y' + 2(ya' + yb') + yc') dt/6
end;
end{of case};
t:=t+dT;
aktuelleFelder:=1-aktuelleFelder;
- if dT<dTMax*zeitUntergrenze then begin
- dT:=dTMax*zeitSollwert;
- {$IFDEF Zeitschrittueberwachung}
- Pro.schreibe('- neues dT: '+floattostr(dT)+' (t='+floattostr(t)+')');
- {$ENDIF}
- end;
- {$IFDEF Zeitschrittueberwachung}
- for i:=0 to groesse+3 do
- if Felders[aktuelleFelder,i].Groessen[false,fiN]<0 then begin
- Pro.schreibe(
- 'n<=0 bei '+
- floattostr(t)+' '+
- inttostr(i)+' '+
- floattostr(Felders[aktuelleFelder,i-1].Groessen[false,fidPsiDX]*
- Felders[aktuelleFelder,i-1].Groessen[false,fiiGamma])+' '+
- floattostr(Felders[aktuelleFelder,i-1].Groessen[false,fiN])+' '+
- floattostr(Felders[aktuelleFelder,i-1].Groessen[true,fiN])+' '+
- floattostr(Felders[aktuelleFelder,i].Groessen[false,fidPsiDX]*
- Felders[aktuelleFelder,i].Groessen[false,fiiGamma])+' '+
- floattostr(Felders[aktuelleFelder,i].Groessen[false,fiN])+' '+
- floattostr(Felders[aktuelleFelder,i].Groessen[true,fiN])+' '+
- floattostr(Felders[aktuelleFelder,i+1].Groessen[false,fidPsiDX]*
- Felders[aktuelleFelder,i+1].Groessen[false,fiiGamma])+' '+
- floattostr(Felders[aktuelleFelder,i+1].Groessen[false,fiN])+' '+
- floattostr(Felders[aktuelleFelder,i+1].Groessen[true,fiN]),true);
- Pro.destroyall;
- halt(1);
- end;
+ {$IFDEF Dichteueberwachung}
+ with felders[aktuelleFelder] do
+ for i:=0 to length(inhalt)-1 do
+ for j:=0 to matAnz-1 do
+ if inhalt[i].matWerte[j,mfN,false]<0 then begin
+ pro:=tProtokollant.create(prot,'iteriereSchritt');
+ pro.schreibe(
+ 'n<=0 bei '+
+ floattostr(t)+' '+
+ inttostr(i)+' '+
+ floattostr(inhalt[i-1].matWerte[j,mfDPsiDX,false]*
+ inhalt[i-1].matWerte[j,mfIGamma,false])+' '+
+ floattostr(inhalt[i-1].matWerte[j,mfN,false])+' '+
+ floattostr(inhalt[i-1].matWerte[j,mfN,true])+' '+
+ floattostr(inhalt[i].matWerte[j,mfDPsiDX,false]*
+ inhalt[i].matWerte[j,mfIGamma,false])+' '+
+ floattostr(inhalt[i].matWerte[j,mfN,false])+' '+
+ floattostr(inhalt[i].matWerte[j,mfN,true])+' '+
+ floattostr(inhalt[i+1].matWerte[j,mfDPsiDX,false]*
+ inhalt[i+1].matWerte[j,mfIGamma,false])+' '+
+ floattostr(inhalt[i+1].matWerte[j,mfN,false])+' '+
+ floattostr(inhalt[i+1].matWerte[j,mfN,true]),true);
+ pro.destroyall;
+ halt(1);
+ end;
{$ENDIF}
end;
-procedure TGitter.macheAusgabe(AD: TAusgabedatei; sDX: extended);
-var i: longint;
- tmp,sX: extended;
- b: byte;
+procedure tGitter.macheAusgabe(aD: tAusgabedatei; sDX: extended);
+var i: longint;
+ sX,cX,tmp: extended;
+ emF: tEMFeldInhalt;
+ maF: tMaterieFeldInhalt;
begin
- if AD.Inhalt in [fiA,fiP] then
- for i:=0 to length(Felders[aktuelleFelder])-1 do begin
+ if (aD.teilchen<0) and (aD.emInhalt=efA) then
+ for i:=0 to length(felders[aktuelleFelder].inhalt)-1 do begin
tmp:=0;
- for b:=1 to 3 do
- tmp:=tmp + sqr(Felders[aktuelleFelder,i].Groessen[AD.Ableitung,TFeldinhalt(Byte(AD.Inhalt)+b)]);
- Felders[aktuelleFelder,i].Groessen[AD.Ableitung,AD.Inhalt]:=sqrt(tmp);
+ for emF:=efAX to efAZ do
+ tmp:=tmp+sqr(felders[aktuelleFelder].inhalt[i].emWerte[emF,false]);
+ felders[aktuelleFelder].inhalt[i].emWerte[efA,false]:=sqrt(tmp);
end;
- i:=floor((xr-xl)/sDX+1);
- tmp:=xl+(i-1)*sDX;
- BlockWrite(AD.Datei,xl,sizeof(extended));
- BlockWrite(AD.Datei,tmp,sizeof(extended));
- BlockWrite(AD.Datei,i,sizeof(longint));
+ if (aD.teilchen>=0) and (aD.matInhalt=mfP) then
+ for i:=0 to length(felders[aktuelleFelder].inhalt)-1 do begin
+ tmp:=0;
+ for maF:=mfPX to mfPZ do
+ tmp:=tmp+sqr(felders[aktuelleFelder].inhalt[i].matWerte[aD.teilchen,maF,false]);
+ felders[aktuelleFelder].inhalt[i].matWerte[aD.teilchen,maF,false]:=sqrt(tmp);
+ end;
+
+ i:=floor((length(felders[aktuelleFelder].inhalt)-1)/sDX*dX+1);
+ sX:=xl+(i-1)*sDX;
+ blockWrite(aD.datei,xl,sizeof(extended));
+ blockWrite(aD.datei,sX,sizeof(extended));
+ blockWrite(aD.datei,i,sizeof(longint));
sX:=xl-Min(dX,sDX)/2;
- for i:=0 to length(Felders[aktuelleFelder])-1 do
- while Felders[aktuelleFelder,i].x>=sX do begin
- BlockWrite(AD.Datei,Felders[aktuelleFelder,i].Groessen[AD.Ableitung,AD.Inhalt],sizeof(extended));
- sX:=sX+sDX;
+ cX:=xl;
+ if aD.teilchen<0 then begin
+ for i:=0 to length(felders[aktuelleFelder].inhalt)-1 do begin
+ while cX>=sX do begin
+ blockWrite(aD.datei,felders[aktuelleFelder].inhalt[i].emWerte[aD.emInhalt,aD.ableitung],sizeof(extended));
+ sX:=sX+sDX;
+ end;
+ cX:=cX+dX;
end;
+ end
+ else begin
+ for i:=0 to length(felders[aktuelleFelder].inhalt)-1 do begin
+ while cX>=sX do begin
+ blockWrite(aD.datei,felders[aktuelleFelder].inhalt[i].matWerte[aD.teilchen,aD.matInhalt,aD.ableitung],sizeof(extended));
+ sX:=sX+sDX;
+ end;
+ cX:=cX+dX;
+ end;
+ end;
end;
-function TGitter.gibErhaltungsgroessen: string;
-var i,j,k: integer;
- n: extended;
- Pro: TProtokollant;
+function tGitter.gibErhaltungsgroessen: string;
+var i,j: integer;
+ ns: tExtendedArray;
+ Pro: TProtokollant;
begin
- Pro:=TProtokollant.create(Prot,'gibErhaltungsgroessen');
- n:=0;
- for i:=0 to groesse+3 do
- n:=n+Felders[aktuelleFelder,i].Groessen[false,fiN];
- for i:=1 to groesse+2 do
- if abs(Felders[aktuelleFelder,i].Groessen[false,fiPx]*
- (Felders[aktuelleFelder,i+1].Groessen[false,fiN]-Felders[aktuelleFelder,i-1].Groessen[false,fiN])/
- Max(Felders[aktuelleFelder,i+1].Groessen[false,fiN]+Felders[aktuelleFelder,i-1].Groessen[false,fiN],1e-10))
- > pDNMax then begin
- if pDNMaxDynamisch then begin
- pDNMax:=
- 2*
- abs(Felders[aktuelleFelder,i].Groessen[false,fiPx]*
- (Felders[aktuelleFelder,i+1].Groessen[false,fiN]-Felders[aktuelleFelder,i-1].Groessen[false,fiN])/
- Max(Felders[aktuelleFelder,i+1].Groessen[false,fiN]+Felders[aktuelleFelder,i-1].Groessen[false,fiN],1e-10));
- for j:=0 to length(felders)-1 do
- for k:=0 to length(felders[j])-1 do
- felders[j,k].pDNMax:=pDNMax;
- Pro.schreibe('Neues maximales px * dn/dx / n: '+floattostr(pDNMax)+' (t='+floattostr(t)+')');
- end
- else
- if pDNMax>0 then begin
- Pro.schreibe('Warnung: maximaler Impuls * dn/dx / n von '+floattostr(pDNMax)+' wurde überschritten (t='+floattostr(t)+
- 'T), die numerische Stabilität ist nicht mehr gewährleistet!',true);
- Pro.schreibe(' Lösung: größeren Diffusionsterm wählen (-D)',true);
- Pro.schreibe(' außerdem empfohlen: Ortsauflösung in gleichem Maße verbessern (-x)',true);
- pDNMax:=-1;
- end;
+ pro:=tProtokollant.create(prot,'gibErhaltungsgroessen');
+ setlength(ns,felders[aktuelleFelder].matAnz);
+ result:='';
+ for i:=0 to length(ns)-1 do
+ ns[i]:=0;
+ for i:=0 to length(felders[aktuelleFelder].inhalt)-1 do
+ for j:=0 to length(ns)-1 do
+ ns[j]:=ns[j] + felders[aktuelleFelder].inhalt[i].matWerte[j,mfN,false];
+ for i:=0 to length(ns)-1 do begin
+ result:=result + ' n['+inttostr(i)+']='+floattostr(ns[i]);
+ {$IFDEF Dichteueberwachung}
+ if ns[i]>1000 then begin
+ errorCode:=2;
+ pro.schreibe(' n > 1000, es scheinen sehr viele neue Teilchen entstanden zu sein. Die Simulation wird abgebrochen. (t='+floattostr(t)+')');
end;
- result:='n='+floattostr(n);
- {$IFDEF Dichteueberwachung}
- if n>1000 then begin
- errorCode:=2;
- Pro.schreibe(' n > 1000, es scheinen sehr viele neue Teilchen entstanden zu sein. Die Simulation wird abgebrochen. (t='+floattostr(t)+')');
+ {$ENDIF}
end;
- {$ENDIF}
- Pro.free;
+ delete(result,1,1);
+ pro.free;
+end;
+
+procedure tGitter.diffusionsTermAnpassen(pro: tProtokollant);
+var
+ i,j: longint;
+begin
+ for i:=1 to length(felders[aktuelleFelder].inhalt)-2 do
+ for j:=0 to felders[aktuelleFelder].matAnz-1 do
+ if abs(felders[aktuelleFelder].inhalt[i].matWerte[j,mfPx,false]*
+ (felders[aktuelleFelder].inhalt[i+1].matWerte[j,mfN,false] - felders[aktuelleFelder].inhalt[i-1].matWerte[j,mfN,false])/
+ max(felders[aktuelleFelder].inhalt[i+1].matWerte[j,mfN,false] + felders[aktuelleFelder].inhalt[i-1].matWerte[j,mfN,false],1e-10))
+ > pDNMax then begin
+ if pDNMaxDynamisch then begin
+ pDNMax:=
+ 2*
+ abs(felders[aktuelleFelder].inhalt[i].matWerte[j,mfPX,false]*
+ (felders[aktuelleFelder].inhalt[i+1].matWerte[j,mfN,false] - felders[aktuelleFelder].inhalt[i-1].matWerte[j,mfN,false])/
+ Max(felders[aktuelleFelder].inhalt[i+1].matWerte[j,mfN,false] + felders[aktuelleFelder].inhalt[i-1].matWerte[j,mfN,false],1e-10));
+ pro.schreibe('Neues maximales px * dn/dx / n: '+floattostr(pDNMax)+' (t='+floattostr(t)+')');
+ end
+ else
+ if pDNMax>0 then begin
+ pro.schreibe('Warnung: maximaler Impuls * dn/dx / n von '+floattostr(pDNMax)+' wurde überschritten (t='+floattostr(t)+
+ 'T), die numerische Stabilität ist nicht mehr gewährleistet!',true);
+ pro.schreibe(' Lösung: größeren Diffusionsterm wählen (-D)',true);
+ pro.schreibe(' außerdem empfohlen: Ortsauflösung in gleichem Maße verbessern (-x)',true);
+ pDNMax:=-1;
+ end;
+ end;
end;
{ tSimulation }
@@ -507,11 +691,10 @@ var
ifile: tMyStringlist;
zeitverfahren: tZeitverfahren;
s,t,aSuffix,aPrefix: string;
- deltaX,deltaT,endzeit,breite,sDT,pDNMax: extended;
- ausgabeDateien: array of tAusgabeDatei;
- abl: boolean;
- fi: tFeldinhalt;
+ deltaX,deltaT,breite,pDNMax: extended;
i,j: longint;
+ kvs: tKnownValues;
+ dichten,lichter: tMyStringlist;
begin
inherited create;
prot:=tProtokollant.create(Protokollant,'tSimulation');
@@ -527,25 +710,35 @@ begin
// Standardeinstellungen Bereich 'allgemein'
zeitverfahren:=zfRungeKuttaVier;
deltaX:=1e-2;
- kvs.add('λ',1/deltaX);
+ kvs.add('λ',1);
+ kvs.add('T',1);
+ kvs.add('dX',deltaX);
deltaT:=-1;
sDT:=-1;
- endzeit:=100;
+ sDX:=-1;
+ tEnde:=100;
breite:=10.0;
pDNMax:=-1;
+ fortschrittsAnzeige:=false;
// Standardeinstellungen Bereich 'ausgaben'
aPrefix:=extractfilepath(paramstr(0));
aSuffix:='.dat';
setlength(ausgabeDateien,0);
- while ifile.readln(s) do begin
+ // Standardeinstellungen Breich 'lichtVon...'
+ lichter:=tMyStringlist.create(prot);
+
+ // Standardeinstellungen Breich 'teilchen...'
+ dichten:=tMyStringlist.create(prot);
+
+ while ifile.readln(s) do begin writeln(s);
if s='allgemein' then begin
repeat
if not ifile.readln(s) then begin
prot.schreibe('Unerwartetes Dateiende in Parameterdatei '''+inName+''' im Bereich allgemein!',true);
halt(1);
- end;
+ end; writeln(s);
if s='allgemeinEnde' then break;
if s='runge-Kutta-4' then begin
Zeitverfahren:=zfRungeKuttaVier;
@@ -557,7 +750,7 @@ begin
end;
if startetMit('ortsschritt',s) then begin
deltaX:=exprtofloat(false,s,kvs,nil);
- kvs.add('λ',1/deltaX);
+ kvs.add('dX',deltaX);
continue;
end;
if startetMit('zeitschritt',s) then begin
@@ -570,13 +763,21 @@ begin
continue;
end;
if startetMit('zeit',s) then begin
- endzeit:=exprtofloat(false,s,kvs,nil);
+ tEnde:=exprtofloat(false,s,kvs,nil);
continue;
end;
if startetMit('breite',s) then begin
breite:=exprtofloat(false,s,kvs,nil);
continue;
end;
+ if s='mit Fortschrittsanzeige' then begin
+ fortschrittsAnzeige:=true;
+ continue;
+ end;
+ if s='ohne Fortschrittsanzeige' then begin
+ fortschrittsAnzeige:=false;
+ continue;
+ end;
prot.schreibe('Unbekannter Befehl '''+s+''' in Parameterdatei '''+inName+''' im Bereich allgemein!',true);
halt(1);
until false;
@@ -602,20 +803,16 @@ begin
sDT:=exprtofloat(false,s,kvs,nil);
continue;
end;
+ if startetMit('ortsschritt',s) then begin
+ sDX:=exprtofloat(false,s,kvs,nil);
+ continue;
+ end;
if startetMit('felder',s) then begin
s:=s+',';
while s<>'' do begin
t:=erstesArgument(s,',');
- for abl:=false to true do
- for fI:=low(tFeldInhalt) to high(tFeldInhalt) do
- if not (abl and (fI in [fiGamma,fiP])) then
- if t=copy('d',1,byte(abl))+optionNamen[fI] then begin
- setlength(ausgabedateien,length(ausgabedateien)+1);
- ausgabedateien[length(ausgabedateien)-1].inhalt:=fI;
- ausgabedateien[length(ausgabedateien)-1].ableitung:=abl;
- ausgabedateien[length(ausgabedateien)-1].name:=aPrefix+t+aSuffix;
- assignFile(ausgabedateien[length(ausgabedateien)-1].datei,ausgabedateien[length(ausgabedateien)-1].name);
- end;
+ setlength(ausgabeDateien,length(ausgabeDateien)+1);
+ ausgabeDateien[length(ausgabeDateien)-1]:=tAusgabeDatei.create(t,aPrefix,aSuffix,prot);
end;
continue;
end;
@@ -638,26 +835,79 @@ begin
deltaT:=deltaX/10;
if sDT<0 then
sDT:=deltaT;
+ if sDX<0 then
+ sDX:=deltaX;
for i:=0 to length(ausgabedateien)-1 do begin
rewrite(ausgabedateien[i].datei,1);
j:=0;
BlockWrite(Ausgabedateien[i].Datei,j,sizeof(longint));
- j:=floor(Endzeit/sDT+1);
+ j:=floor(tEnde/sDT+1);
BlockWrite(Ausgabedateien[i].Datei,j,sizeof(longint));
end;
ifile.free;
- gitter:=tGitter.create(Breite,deltaT,deltaX,pDNMax,@Anfangsdichte,Zeitverfahren,Prot);
+ gitter:=tGitter.create(round(breite/deltaX),deltaT,deltaX,pDNMax,kvs,dichten,lichter,zeitverfahren,prot);
+ dichten.free;
+ lichter.free;
end;
destructor tSimulation.destroy;
begin
gitter.free;
- kvs.free;
prot.free;
inherited destroy;
end;
+function tSimulation.iteriereSchritt(start: extended; var zeitPhysik,zeitDatei: extended): boolean; // noch nicht zu Ende?
+var
+ i: longint;
+ c: char;
+begin
+ result:=false;
+
+ zeitPhysik:=zeitPhysik-now;
+ if errorCode<2 then
+ gitter.iteriereSchritt(dT);
+ zeitPhysik:=zeitPhysik+now;
+ zeitDatei:=zeitDatei-now;
+ while gitter.t>=sT do begin
+ sT:=sT+sDT;
+ for i:=0 to length(Ausgabedateien)-1 do
+ gitter.macheAusgabe(ausgabeDateien[i],sDX);
+ end;
+ zeitDatei:=zeitDatei+now;
+
+ gitter.iteriereSchritt(dT);
+
+ if fortschrittsAnzeige then begin
+ if (floor(100*Gitter.t/tEnde) < floor(100*(Gitter.t+dT)/tEnde)) or keypressed then begin
+ if keypressed then c:=readkey
+ else c:=#0;
+ case c of
+ #27,'q': begin
+ errorCode:=3;
+ exit;
+ end;
+ ' ': begin
+ writeln(' ... Pause (beliebige Taste drücken um fortzufahren) ...');
+ readkey;
+ writeln(' ... weiter geht''s ...');
+ end;
+ else begin
+ Prot.schreibe(inttostr(round(100*Gitter.t/tEnde))+'% (t='+floattostr(Gitter.t)+'T)',true);
+ Prot.schreibe(timetostr(now-start)+' ('+floattostr(zeitPhysik/max(1e-11,zeitPhysik+zeitDatei))+')',true);
+ Prot.schreibe('ETA: '+timetostr((now-start)*(tEnde-Gitter.t)/max(Gitter.t,dT)),true);
+ Prot.schreibe('aktueller Zeitschritt: '+floattostr(dT)+'T',true);
+ if errorCode<2 then
+ Prot.schreibe(gitter.gibErhaltungsgroessen);
+ end;
+ end{of case};
+ end;
+ end;
+
+ result:=gitter.t>=tEnde;
+end;
+
end.
diff --git a/Plasmapropagation b/Plasmapropagation
index 013329a..98472f4 100755
--- a/Plasmapropagation
+++ b/Plasmapropagation
Binary files differ
diff --git a/Plasmapropagation.lpr b/Plasmapropagation.lpr
index b861d2e..5c914d5 100644
--- a/Plasmapropagation.lpr
+++ b/Plasmapropagation.lpr
@@ -28,13 +28,10 @@ type
procedure TPlasmapropagation.DoRun;
var
- Breite,i,j: longint;
simulation: tSimulation;
- start,zeitDatei,zeitPhysik: extended;
- FI: TFeldInhalt;
- Abl: Boolean;
- Prot: TProtokollant;
- c: char;
+ start,zeitPhysik,zeitDatei: extended;
+ prot: tProtokollant;
+ s,t,u: string;
begin
prot:=tProtokollant.create('error');
@@ -48,71 +45,29 @@ begin
zeitDatei:=0;
zeitPhysik:=0;
- sT:=-Min(deltaT,sDT)/2;
- simulation:=tSimulation.create(paramstr(1));
-// Gitter.Al:=@InputFeld;
- while Gitter.t<Endzeit do begin
- if HasOption('F','Fortschrittsanzeige') then begin
- if errorCode<2 then
- s:=Gitter.gibErhaltungsgroessen;
- if (floor(100*Gitter.t/Endzeit) < floor(100*(Gitter.t+deltaT)/Endzeit)) or keypressed then begin
- if keypressed then c:=readkey
- else c:=#0;
- case c of
- #27,'q': begin
- errorCode:=3;
- break;
- end;
- ' ': begin
- writeln(' ... Pause (beliebige Taste drücken um fortzufahren) ...');
- readkey;
- writeln(' ... weiter geht''s ...');
- end;
- else begin
- Prot.schreibe(inttostr(round(100*Gitter.t/Endzeit))+'% (t='+floattostr(Gitter.t)+'T)',true);
- Prot.schreibe(timetostr(now-start)+' ('+floattostr(zeitPhysik/max(1e-11,zeitPhysik+zeitDatei))+')',true);
- Prot.schreibe('ETA: '+timetostr((now-start)*(Endzeit-Gitter.t)/max(Gitter.t,deltaT)),true);
- Prot.schreibe('aktueller Zeitschritt: '+floattostr(deltaT)+'T',true);
- Prot.schreibe(s);
- end;
- end{of case};
- end;
- end;
-
- zeitPhysik:=zeitPhysik-now;
- if errorCode<2 then
- Gitter.iteriereSchritt(deltaT);
- zeitPhysik:=zeitPhysik+now;
- zeitDatei:=zeitDatei-now;
- while Gitter.t>=sT do begin
- sT:=sT+sDT;
- for j:=0 to length(Ausgabedateien)-1 do
- Gitter.macheAusgabe(Ausgabedateien[j],sDX);
- end;
- zeitDatei:=zeitDatei+now;
- end;
+ simulation:=tSimulation.create(paramstr(1),prot);
+
+ while simulation.iteriereSchritt(start,zeitPhysik,zeitDatei) do ;
case errorCode of
- 2: Prot.schreibe('Simulation wurde auf halbem Wege abgebrochen wegen Überlaufs.',true);
- 3: Prot.schreibe('Simulation wurde auf halbem Wege abgebrochen wegen Benutzereingriffs.',true);
+ 2: prot.schreibe('Simulation wurde auf halbem Wege abgebrochen wegen Überlaufs.',true);
+ 3: prot.schreibe('Simulation wurde auf halbem Wege abgebrochen wegen Benutzereingriffs.',true);
end{of case};
- for i:=0 to length(Ausgabedateien)-1 do
- CloseFile(Ausgabedateien[i].Datei);
+ simulation.free;
+ prot.schreibe('fertig!',true);
- Prot.schreibe('fertig!',true);
s:=timetostr(now-start);
t:=timetostr(zeitDatei);
u:=timetostr(zeitPhysik);
while length(s)<max(length(t),length(u)) do s:=' '+s;
while length(t)<max(length(s),length(u)) do t:=' '+t;
while length(u)<max(length(s),length(t)) do u:=' '+u;
- Prot.schreibe('Das hat '+s+' gedauert,',true);
- Prot.schreibe(' davon '+t+' für Dateizugriffe',true);
- Prot.schreibe('und nur '+u+' für die eigentliche Physik!',true);
+ prot.schreibe('Das hat '+s+' gedauert,',true);
+ prot.schreibe(' davon '+t+' für Dateizugriffe',true);
+ prot.schreibe('und nur '+u+' für die eigentliche Physik!',true);
+ prot.Free;
- Gitter.Free;
- Prot.Free;
// stop program loop
Terminate;
if errorCode=1 then errorCode:=0;
diff --git a/Plasmapropagation.lps b/Plasmapropagation.lps
index f9ed3f9..cd5e8e9 100644
--- a/Plasmapropagation.lps
+++ b/Plasmapropagation.lps
@@ -7,21 +7,20 @@
<Unit0>
<Filename Value="Plasmapropagation.lpr"/>
<IsPartOfProject Value="True"/>
- <TopLine Value="101"/>
- <CursorPos Y="133"/>
- <UsageCount Value="129"/>
+ <TopLine Value="28"/>
+ <CursorPos X="39" Y="11"/>
+ <UsageCount Value="140"/>
<Loaded Value="True"/>
</Unit0>
<Unit1>
<Filename Value="Physikunit.pas"/>
<IsPartOfProject Value="True"/>
- <IsVisibleTab Value="True"/>
<EditorIndex Value="1"/>
- <TopLine Value="136"/>
- <CursorPos X="41" Y="252"/>
- <UsageCount Value="70"/>
+ <TopLine Value="699"/>
+ <CursorPos X="13" Y="714"/>
+ <UsageCount Value="81"/>
<Bookmarks Count="1">
- <Item0 Y="301"/>
+ <Item0 Y="525"/>
</Bookmarks>
<Loaded Value="True"/>
</Unit1>
@@ -31,195 +30,198 @@
<EditorIndex Value="-1"/>
<TopLine Value="26"/>
<CursorPos X="103" Y="47"/>
- <UsageCount Value="32"/>
+ <UsageCount Value="43"/>
</Unit2>
<Unit3>
<Filename Value="input.plap"/>
<IsPartOfProject Value="True"/>
- <EditorIndex Value="3"/>
- <CursorPos Y="19"/>
- <UsageCount Value="31"/>
+ <IsVisibleTab Value="True"/>
+ <EditorIndex Value="5"/>
+ <CursorPos X="7" Y="19"/>
+ <UsageCount Value="42"/>
<Loaded Value="True"/>
<DefaultSyntaxHighlighter Value="None"/>
</Unit3>
<Unit4>
- <Filename Value="../epost/epostunit.pas"/>
- <TopLine Value="575"/>
- <FoldState Value=" T3jg0C2 pk4kM0{A1O"/>
- <UsageCount Value="1"/>
- </Unit4>
- <Unit5>
<Filename Value="input.epost"/>
<EditorIndex Value="-1"/>
<TopLine Value="45"/>
<CursorPos Y="69"/>
- <UsageCount Value="57"/>
+ <UsageCount Value="56"/>
<DefaultSyntaxHighlighter Value="None"/>
+ </Unit4>
+ <Unit5>
+ <Filename Value="../units/matheunit.pas"/>
+ <EditorIndex Value="3"/>
+ <TopLine Value="75"/>
+ <CursorPos Y="347"/>
+ <FoldState Value=" T3i804A pibjM034]9aj60U6 plFmG04^"/>
+ <UsageCount Value="18"/>
+ <Loaded Value="True"/>
</Unit5>
<Unit6>
- <Filename Value="../units/matheunit.pas"/>
- <EditorIndex Value="2"/>
- <TopLine Value="82"/>
- <CursorPos X="33" Y="102"/>
- <UsageCount Value="13"/>
+ <Filename Value="../units/lowlevelunit.pas"/>
+ <EditorIndex Value="4"/>
+ <TopLine Value="373"/>
+ <CursorPos X="77" Y="399"/>
+ <UsageCount Value="10"/>
<Loaded Value="True"/>
</Unit6>
<Unit7>
- <Filename Value="../units/lowlevelunit.pas"/>
- <EditorIndex Value="-1"/>
- <TopLine Value="28"/>
- <CursorPos Y="37"/>
+ <Filename Value="../units/mystringlistunit.pas"/>
+ <EditorIndex Value="2"/>
+ <TopLine Value="37"/>
+ <CursorPos Y="223"/>
+ <FoldState Value=" T3i2076 pjFjc095 plcmG0H1b"/>
<UsageCount Value="11"/>
+ <Loaded Value="True"/>
</Unit7>
<Unit8>
- <Filename Value="../units/mystringlistunit.pas"/>
- <EditorIndex Value="-1"/>
- <TopLine Value="371"/>
- <CursorPos Y="399"/>
- <FoldState Value=" T3i2077 pjIk009411j"/>
- <UsageCount Value="11"/>
- </Unit8>
- <Unit9>
<Filename Value="../epost/werteunit.pas"/>
<EditorIndex Value="-1"/>
<TopLine Value="950"/>
<CursorPos X="30" Y="1054"/>
- <UsageCount Value="10"/>
- </Unit9>
- <Unit10>
+ <UsageCount Value="9"/>
+ </Unit8>
+ <Unit9>
<Filename Value="../epost/typenunit.pas"/>
<EditorIndex Value="-1"/>
<TopLine Value="347"/>
<CursorPos X="62" Y="358"/>
- <UsageCount Value="10"/>
- </Unit10>
- <Unit11>
+ <UsageCount Value="9"/>
+ </Unit9>
+ <Unit10>
<Filename Value="../units/systemunit.pas"/>
<EditorIndex Value="-1"/>
<CursorPos X="3" Y="79"/>
- <UsageCount Value="10"/>
+ <UsageCount Value="9"/>
+ </Unit10>
+ <Unit11>
+ <Filename Value="/usr/lib/fpc/src/rtl/inc/objpash.inc"/>
+ <EditorIndex Value="-1"/>
+ <TopLine Value="232"/>
+ <CursorPos X="23" Y="192"/>
+ <UsageCount Value="9"/>
</Unit11>
</Units>
<JumpHistory Count="30" HistoryIndex="29">
<Position1>
<Filename Value="Physikunit.pas"/>
- <Caret Line="498" Column="37" TopLine="482"/>
</Position1>
<Position2>
- <Filename Value="Plasmapropagation.lpr"/>
- <Caret Line="92" TopLine="80"/>
+ <Filename Value="Physikunit.pas"/>
+ <Caret Line="391" Column="52" TopLine="358"/>
</Position2>
<Position3>
<Filename Value="Plasmapropagation.lpr"/>
- <Caret Line="82" TopLine="62"/>
+ <Caret Line="38" TopLine="11"/>
</Position3>
<Position4>
<Filename Value="Plasmapropagation.lpr"/>
- <Caret Line="113" Column="22" TopLine="53"/>
+ <Caret Line="11" Column="39" TopLine="28"/>
</Position4>
<Position5>
<Filename Value="Physikunit.pas"/>
- <Caret Line="626" Column="60" TopLine="598"/>
+ <Caret Line="111" Column="17" TopLine="91"/>
</Position5>
<Position6>
<Filename Value="Physikunit.pas"/>
- <Caret Line="601" Column="124" TopLine="582"/>
+ <Caret Line="36" Column="23" TopLine="16"/>
</Position6>
<Position7>
<Filename Value="Physikunit.pas"/>
- <Caret Line="600" Column="62" TopLine="508"/>
+ <Caret Line="671" TopLine="631"/>
</Position7>
<Position8>
- <Filename Value="Physikunit.pas"/>
- <Caret Line="27" Column="6"/>
+ <Filename Value="../units/mystringlistunit.pas"/>
+ <Caret Line="421" Column="57" TopLine="386"/>
</Position8>
<Position9>
- <Filename Value="Physikunit.pas"/>
- <Caret Line="504" Column="45" TopLine="486"/>
+ <Filename Value="../units/mystringlistunit.pas"/>
</Position9>
<Position10>
- <Filename Value="Plasmapropagation.lpr"/>
- <Caret Line="113" TopLine="93"/>
+ <Filename Value="../units/mystringlistunit.pas"/>
+ <Caret Line="17" Column="23"/>
</Position10>
<Position11>
- <Filename Value="Physikunit.pas"/>
- <Caret Line="628" Column="4" TopLine="609"/>
+ <Filename Value="../units/mystringlistunit.pas"/>
+ <Caret Line="27" Column="23"/>
</Position11>
<Position12>
- <Filename Value="Physikunit.pas"/>
- <Caret Line="626" Column="76" TopLine="609"/>
+ <Filename Value="../units/mystringlistunit.pas"/>
+ <Caret Line="28" Column="23"/>
</Position12>
<Position13>
- <Filename Value="Plasmapropagation.lpr"/>
- <Caret Line="109" TopLine="96"/>
+ <Filename Value="../units/mystringlistunit.pas"/>
+ <Caret Line="83" Column="10" TopLine="44"/>
</Position13>
<Position14>
- <Filename Value="Plasmapropagation.lpr"/>
- <Caret Line="64" TopLine="53"/>
+ <Filename Value="../units/mystringlistunit.pas"/>
+ <Caret Line="81" TopLine="42"/>
</Position14>
<Position15>
- <Filename Value="Physikunit.pas"/>
- <Caret Line="637" Column="11" TopLine="608"/>
+ <Filename Value="../units/mystringlistunit.pas"/>
+ <Caret Line="27" Column="54" TopLine="16"/>
</Position15>
<Position16>
<Filename Value="Physikunit.pas"/>
- <Caret Line="504" Column="44" TopLine="483"/>
+ <Caret Line="732" Column="56" TopLine="702"/>
</Position16>
<Position17>
<Filename Value="Physikunit.pas"/>
- <Caret Line="637" Column="55" TopLine="610"/>
+ <Caret Line="402" Column="52" TopLine="382"/>
</Position17>
<Position18>
- <Filename Value="Physikunit.pas"/>
- <Caret Line="561" Column="20" TopLine="511"/>
+ <Filename Value="../units/mystringlistunit.pas"/>
+ <Caret Line="264" Column="39" TopLine="256"/>
</Position18>
<Position19>
- <Filename Value="input.plap"/>
- <Caret Line="19" Column="19" TopLine="10"/>
+ <Filename Value="../units/mystringlistunit.pas"/>
+ <Caret Line="25" Column="6" TopLine="8"/>
</Position19>
<Position20>
- <Filename Value="Physikunit.pas"/>
- <Caret Line="642" Column="55" TopLine="522"/>
+ <Filename Value="../units/mystringlistunit.pas"/>
+ <Caret Line="83" Column="3" TopLine="46"/>
</Position20>
<Position21>
- <Filename Value="Physikunit.pas"/>
- <Caret Line="523" Column="30"/>
+ <Filename Value="../units/mystringlistunit.pas"/>
+ <Caret Line="276" Column="48" TopLine="256"/>
</Position21>
<Position22>
- <Filename Value="Physikunit.pas"/>
- <Caret Line="642" Column="55" TopLine="615"/>
+ <Filename Value="../units/mystringlistunit.pas"/>
+ <Caret Line="340" Column="30" TopLine="320"/>
</Position22>
<Position23>
- <Filename Value="Physikunit.pas"/>
- <Caret Line="530" TopLine="502"/>
+ <Filename Value="../units/mystringlistunit.pas"/>
+ <Caret Line="27" Column="54" TopLine="17"/>
</Position23>
<Position24>
<Filename Value="Physikunit.pas"/>
- <Caret Line="641" Column="55" TopLine="548"/>
+ <Caret Line="402" Column="47" TopLine="382"/>
</Position24>
<Position25>
<Filename Value="Physikunit.pas"/>
- <Caret Line="61" Column="24" TopLine="40"/>
+ <Caret Line="728" Column="38" TopLine="709"/>
</Position25>
<Position26>
- <Filename Value="Physikunit.pas"/>
- <Caret Line="224" Column="21" TopLine="204"/>
+ <Filename Value="../units/matheunit.pas"/>
+ <Caret Line="102" Column="33" TopLine="25"/>
</Position26>
<Position27>
- <Filename Value="Physikunit.pas"/>
- <Caret Line="50" Column="15" TopLine="21"/>
+ <Filename Value="../units/matheunit.pas"/>
+ <Caret Line="45" Column="21" TopLine="25"/>
</Position27>
<Position28>
- <Filename Value="Physikunit.pas"/>
- <Caret Line="224" Column="154" TopLine="204"/>
+ <Filename Value="../units/matheunit.pas"/>
+ <Caret Line="347" Column="21" TopLine="314"/>
</Position28>
<Position29>
<Filename Value="../units/matheunit.pas"/>
- <Caret Line="106" Column="34" TopLine="81"/>
+ <Caret Line="443" Column="34" TopLine="425"/>
</Position29>
<Position30>
- <Filename Value="Physikunit.pas"/>
- <Caret Line="228" Column="22" TopLine="219"/>
+ <Filename Value="../units/lowlevelunit.pas"/>
+ <Caret Line="54" Column="18" TopLine="46"/>
</Position30>
</JumpHistory>
</ProjectSession>
diff --git a/error b/error
deleted file mode 100644
index 166a407..0000000
--- a/error
+++ /dev/null
@@ -1,521 +0,0 @@
- Simulationsbox
- 1500 Felder breit (1.5 lambda)
- 1 T lang (etwa 10000 Schritte)
- 0.0001 Zeitschritt und 0.001 Ortsschritt
- Mache Ausgabe fiAX.
- Mache Ausgabe fiAY.
- Mache Ausgabe fidAYDT.
- Mache Ausgabe fiN.
- Mache Ausgabe fidPhiDX.
- Mache Ausgabe fidPsiDX.
- Mache Ausgabe fiGamma.
- Maximales px * dn/dx / n wird automatisch ermittelt.
- Iterationsverfahren: Runge-Kutta
-TGitter.gibErhaltungsgroessen Neues maximales px * dn/dx / n: 4.90867701035798E-13 (t=0.000025)
-TGitter.gibErhaltungsgroessen Neues maximales px * dn/dx / n: 1.9210083367696E-6 (t=0.000025)
-TGitter.gibErhaltungsgroessen Neues maximales px * dn/dx / n: 6.95677493143965E-6 (t=0.000025)
-TGitter.gibErhaltungsgroessen Neues maximales px * dn/dx / n: 0.00001391665818 (t=0.000025)
- 1% (t=0.0099875T)
- 4s (0.93572594126268)
- ETA: 6min 12s
- aktueller Zeitschritt: 0.0000125T
- n=1100.1504001997
- 2% (t=0.0199875T)
- 7s (0.93845056986265)
- ETA: 6min 4s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040039965
- 3% (t=0.0299875T)
- 11s (0.9396443307214)
- ETA: 5min 59s
- aktueller Zeitschritt: 0.0000125T
- n=1100.1504005996
- 4% (t=0.04T)
- 15s (0.93995505725415)
- ETA: 5min 55s
- aktueller Zeitschritt: 0.0000125T
- n=1100.1504007998
- 5% (t=0.05T)
- 18s (0.94021909059327)
- ETA: 5min 51s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040099976
- 6% (t=0.06T)
- 22s (0.94028101686857)
- ETA: 5min 47s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040119971
- 7% (t=0.07T)
- 26s (0.94036230559588)
- ETA: 5min 43s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040139966
- 8% (t=0.08T)
- 30s (0.94045357808056)
- ETA: 5min 39s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040159961
- 9% (t=0.09T)
- 33s (0.94041907955508)
- ETA: 5min 36s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040179956
- 10% (t=0.1T)
- 37s (0.93998029343796)
- ETA: 5min 32s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040199952
- 11% (t=0.11T)
- 41s (0.9391274466381)
- ETA: 5min 29s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040219947
- 12% (t=0.12T)
- 44s (0.93928366492497)
- ETA: 5min 25s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040239942
- 13% (t=0.13T)
- 48s (0.93936119390127)
- ETA: 5min 21s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040259938
- 14% (t=0.14T)
- 52s (0.93947196769724)
- ETA: 5min 17s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040279933
- 15% (t=0.1499875T)
- 55s (0.93947754328833)
- ETA: 5min 14s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040299903
- 16% (t=0.1599875T)
- 59s (0.93955884945007)
- ETA: 5min 10s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040319899
- 17% (t=0.1699875T)
- 1min 3s (0.93965263560525)
- ETA: 5min 6s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040339894
- 18% (t=0.1799875T)
- 1min 6s (0.9397448828285)
- ETA: 5min 2s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040359889
- 19% (t=0.1899875T)
- 1min 10s (0.93915432836246)
- ETA: 4min 59s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040379885
- 20% (t=0.1999875T)
- 1min 14s (0.93900073347529)
- ETA: 4min 55s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040399879
- 21% (t=0.2099875T)
- 1min 17s (0.93909659887595)
- ETA: 4min 51s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040419873
- 22% (t=0.2199875T)
- 1min 21s (0.93917337012531)
- ETA: 4min 48s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040439866
- 23% (t=0.2299875T)
- 1min 25s (0.93758208132779)
- ETA: 4min 44s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040459857
- 24% (t=0.2399875T)
- 1min 29s (0.93768875699546)
- ETA: 4min 41s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040479846
- 25% (t=0.2499875T)
- 1min 32s (0.9378121327816)
- ETA: 4min 37s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040499831
- 26% (t=0.2599875T)
- 1min 36s (0.93788794049542)
- ETA: 4min 33s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040519814
- 27% (t=0.2699875T)
- 1min 39s (0.93470907399499)
- ETA: 4min 29s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040539794
- 28% (t=0.2799875T)
- 1min 43s (0.93392779228268)
- ETA: 4min 24s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040559772
- 29% (t=0.2899875T)
- 1min 46s (0.9341559070417)
- ETA: 4min 21s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040579748
- 30% (t=0.2999875T)
- 1min 50s (0.93429738634912)
- ETA: 4min 17s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040599723
- 31% (t=0.3099875T)
- 1min 54s (0.93446650346284)
- ETA: 4min 13s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040619699
- 32% (t=0.3199875T)
- 1min 57s (0.93465115233278)
- ETA: 4min 10s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040639676
- 33% (t=0.3299875T)
- 2min 1s (0.93480925657284)
- ETA: 4min 6s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040659655
- 34% (t=0.3399875T)
- 2min 5s (0.93496938818957)
- ETA: 4min 2s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040679637
- 35% (t=0.3499875T)
- 2min 8s (0.93499199359356)
- ETA: 3min 59s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040699623
- 36% (t=0.3599875T)
- 2min 12s (0.93513974146828)
- ETA: 3min 55s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040719612
- 37% (t=0.3699875T)
- 2min 16s (0.93462522712332)
- ETA: 3min 51s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040739603
- 38% (t=0.3799875T)
- 2min 20s (0.93477500648)
- ETA: 3min 48s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040759598
- 39% (t=0.3899875T)
- 2min 23s (0.93418072020113)
- ETA: 3min 43s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040779595
- 40% (t=0.3999875T)
- 2min 26s (0.93434368089496)
- ETA: 3min 39s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040799592
- 41% (t=0.4099875T)
- 2min 29s (0.93438719449992)
- ETA: 3min 35s
- aktueller Zeitschritt: 0.0000125T
- n=1100.1504081959
- 42% (t=0.4199875T)
- 2min 33s (0.93448808312355)
- ETA: 3min 31s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040839587
- 43% (t=0.4299875T)
- 2min 37s (0.93459425789043)
- ETA: 3min 28s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040859582
- 44% (t=0.4399875T)
- 2min 40s (0.93472965099255)
- ETA: 3min 24s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040879575
- 45% (t=0.4499875T)
- 2min 44s (0.93483912389575)
- ETA: 3min 20s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040899565
- 46% (t=0.4599875T)
- 2min 48s (0.93495962549421)
- ETA: 3min 17s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040919551
- 47% (t=0.4699875T)
- 2min 51s (0.93495055947613)
- ETA: 3min 13s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040939535
- 48% (t=0.4799875T)
- 2min 55s (0.93509540698992)
- ETA: 3min 10s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040959517
- 49% (t=0.4899875T)
- 2min 59s (0.93522625662536)
- ETA: 3min 6s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040979497
- 50% (t=0.4999875T)
- 3min 2s (0.93534449434932)
- ETA: 3min 2s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15040999477
- 51% (t=0.5099875T)
- 3min 6s (0.93508263291742)
- ETA: 2min 59s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15041019456
- 52% (t=0.5199875T)
- 3min 10s (0.93509288526884)
- ETA: 2min 55s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15041039437
- 53% (t=0.5299875T)
- 3min 13s (0.93488834098026)
- ETA: 2min 52s
- aktueller Zeitschritt: 0.0000125T
- n=1100.1504105942
- 54% (t=0.5399875T)
- 3min 17s (0.93499188128487)
- ETA: 2min 48s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15041079406
- 55% (t=0.5499875T)
- 3min 21s (0.93480126139187)
- ETA: 2min 44s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15041099395
- 56% (t=0.5599875T)
- 3min 25s (0.93490491597985)
- ETA: 2min 41s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15041119387
- 57% (t=0.5699875T)
- 3min 28s (0.93500136747388)
- ETA: 2min 37s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15041139383
- 58% (t=0.5799875T)
- 3min 32s (0.93465203973783)
- ETA: 2min 34s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15041159381
- 59% (t=0.5899875T)
- 3min 36s (0.93472936944875)
- ETA: 2min 30s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15041179381
- 60% (t=0.6T)
- 3min 39s (0.9348087023095)
- ETA: 2min 26s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15041199407
- 61% (t=0.61T)
- 3min 43s (0.93491070711469)
- ETA: 2min 23s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15041219408
- 62% (t=0.62T)
- 3min 46s (0.93462517792379)
- ETA: 2min 19s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15041239408
- 63% (t=0.63T)
- 3min 50s (0.93465224645404)
- ETA: 2min 15s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15041259406
- 64% (t=0.64T)
- 3min 53s (0.93472351649562)
- ETA: 2min 11s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15041279402
- 65% (t=0.65T)
- 3min 56s (0.9346948061847)
- ETA: 2min 7s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15041299395
- 66% (t=0.66T)
- 3min 59s (0.93474317500594)
- ETA: 2min 3s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15041319385
- 67% (t=0.67T)
- 4min 2s (0.93478852332337)
- ETA: 1min 59s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15041339372
- 68% (t=0.68T)
- 4min 5s (0.93477508176243)
- ETA: 1min 55s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15041359357
- 69% (t=0.69T)
- 4min 9s (0.93486059208118)
- ETA: 1min 52s
- aktueller Zeitschritt: 0.0000125T
- n=1100.1504137934
- 70% (t=0.7T)
- 4min 13s (0.93484459753493)
- ETA: 1min 48s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15041399322
- 71% (t=0.71T)
- 4min 16s (0.93492848469232)
- ETA: 1min 45s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15041419304
- 72% (t=0.72T)
- 4min 20s (0.93471191856157)
- ETA: 1min 41s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15041439287
- 73% (t=0.73T)
- 4min 24s (0.9348012495454)
- ETA: 1min 38s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15041459273
- 74% (t=0.74T)
- 4min 27s (0.9348810582785)
- ETA: 1min 34s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15041479261
- 75% (t=0.75T)
- 4min 31s (0.93466154687072)
- ETA: 1min 30s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15041499252
- 76% (t=0.76T)
- 4min 35s (0.93474800649981)
- ETA: 1min 27s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15041519246
- 77% (t=0.77T)
- 4min 39s (0.93464517725866)
- ETA: 1min 23s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15041539244
- 78% (t=0.78T)
- 4min 42s (0.93471589587541)
- ETA: 1min 20s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15041559244
- 79% (t=0.79T)
- 4min 46s (0.93479319538545)
- ETA: 1min 16s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15041579245
- 80% (t=0.8T)
- 4min 50s (0.93487046513783)
- ETA: 1min 12s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15041599248
- 81% (t=0.81T)
- 4min 53s (0.93494272921598)
- ETA: 1min 9s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15041619251
- 82% (t=0.82T)
- 4min 57s (0.93502453059895)
- ETA: 1min 5s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15041639253
- 83% (t=0.83T)
- 5min 1s (0.93509745330369)
- ETA: 1min 2s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15041659253
- 84% (t=0.84T)
- 5min 4s (0.93517156164012)
- ETA: 58s
- aktueller Zeitschritt: 0.0000125T
- n=1100.1504167925
- 85% (t=0.85T)
- 5min 8s (0.9352390216625)
- ETA: 54s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15041699245
- 86% (t=0.86T)
- 5min 12s (0.9352517564134)
- ETA: 51s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15041719236
- 87% (t=0.87T)
- 5min 16s (0.9351394704682)
- ETA: 47s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15041739224
- 88% (t=0.88T)
- 5min 19s (0.93520735055799)
- ETA: 44s
- aktueller Zeitschritt: 0.0000125T
- n=1100.1504175921
- 89% (t=0.89T)
- 5min 23s (0.93512915953451)
- ETA: 40s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15041779194
- 90% (t=0.9T)
- 5min 27s (0.93519826581303)
- ETA: 36s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15041799177
- 91% (t=0.91T)
- 5min 31s (0.93483492969017)
- ETA: 33s
- aktueller Zeitschritt: 0.0000125T
- n=1100.1504181916
- 92% (t=0.92T)
- 5min 34s (0.93490573860459)
- ETA: 29s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15041839144
- 93% (t=0.93T)
- 5min 38s (0.93496816016958)
- ETA: 25s
- aktueller Zeitschritt: 0.0000125T
- n=1100.1504185913
- 94% (t=0.94T)
- 5min 42s (0.9348235620777)
- ETA: 22s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15041879118
- 95% (t=0.95T)
- 5min 45s (0.93488922414053)
- ETA: 18s
- aktueller Zeitschritt: 0.0000125T
- n=1100.1504189911
- 96% (t=0.96T)
- 5min 49s (0.93492794655297)
- ETA: 15s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15041919104
- 97% (t=0.97T)
- 5min 53s (0.93495188615352)
- ETA: 11s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15041939102
- 98% (t=0.98T)
- 5min 56s (0.93481669414105)
- ETA: 7s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15041959102
- 99% (t=0.99T)
- 6min 0s (0.93488083192356)
- ETA: 4s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15041979103
- 100% (t=1T)
- 6min 4s (0.93494126481698)
- ETA: 0s
- aktueller Zeitschritt: 0.0000125T
- n=1100.15041999106
- fertig!
- Das hat 6min 4s gedauert,
- davon 23s für Dateizugriffe
- und nur 5min 29s für die eigentliche Physik!
diff --git a/input.plap b/input.plap
index 2e5ecc4..d2debd1 100644
--- a/input.plap
+++ b/input.plap
@@ -2,19 +2,21 @@
allgemein
runge-Kutta-4
- ortsschritt 10^-2 λ
- zeitschritt 10^-3 T
- zeit 100 T
- breite 10 λ
+ ortsschritt 10^-2 * λ
+ zeitschritt 10^-3 * T
+ zeit 100 * T
+ breite 10 * λ
+ mit Fortschrittsanzeige
allgemeinEnde
ausgaben
suffix _test.dat
- felder AX,AY,dAYDT,N,dPhiDX,dPsiDX,Gamma
+ felder AX,AY,dAYDT,dPhiDX
+# felder AX,AY,dAYDT,N,dPhiDX,dPsiDX,Gamma
ausgabenEnde
-!set $tFwhm: 5 T
-!set $tMitte: 20 T
+!setze $tFwhm: 5 T
+!setze $tMitte: 20 T
lichtVonLinks 0.5 * 2^(-2*((t-$tMitte)/$tFwhm)^2) * Sin(ω*t)
generated by cgit v1.2.3-70-g09d2 (git 2.47.0) at 2024-11-21 16:19:11 +0000