diff options
| -rw-r--r-- | Physikunit.pas | 431 | ||||
| -rwxr-xr-x | Plasmapropagation | bin | 3387904 -> 3295712 bytes | |||
| -rw-r--r-- | Plasmapropagation.lpi | 6 | ||||
| -rw-r--r-- | Plasmapropagation.lpr | 8 | ||||
| -rw-r--r-- | Plasmapropagation.lps | 171 | ||||
| -rw-r--r-- | input.epost | 2 | ||||
| -rw-r--r-- | input.plap | 14 |
7 files changed, 370 insertions, 262 deletions
diff --git a/Physikunit.pas b/Physikunit.pas index 492866c..1d50c5a 100644 --- a/Physikunit.pas +++ b/Physikunit.pas @@ -8,7 +8,7 @@ unit Physikunit; interface uses - Classes, SysUtils, Math, protokollunit, matheunit, mystringlistunit, lowlevelunit, crt; + Classes, SysUtils, Math, protokollunit, matheunit, mystringlistunit, lowlevelunit; type tZeitverfahren = (zfEulerVorwaerts,zfRungeKuttaVier); @@ -24,27 +24,35 @@ type mfGamma,mfIGamma ); + tGitter = class; + { tAusgabeDatei } tAusgabeDatei = class - ableitung: boolean; - emInhalt: tEMFeldInhalt; - matInhalt: tMaterieFeldInhalt; - teilchen: longint; // wenn < 0, dann EM-Feld - name: string; - datei: file; - constructor create(feldName,prefix,suffix: string; prot: tProtokollant); + private + ableitung: boolean; + emInhalt: tEMFeldInhalt; + matInhalt: tMaterieFeldInhalt; + teilchen,nNum,tCnt: longint; // wenn teilchen < 0, dann EM-Feld + pre,suf: string; + datei: file; + pro: tProtokollant; + public + zeitAnz: longint; + sDT: extended; + constructor create(feldName,prefix,suffix: string; prot: tProtokollant; nam: string); destructor destroy; override; + function dump: string; + procedure schreibeKopf; + procedure speichereWerte(gitter: tGitter; sDX: extended); end; - tGitter = class; - { tWertePunkt } tWertePunkt = class(tObject) // repräsentiert alle Werte eines Punktes im Gitter und deren Zeitableitungen - matWerte: array of array[tMaterieFeldInhalt,boolean] of extended; + matWerte: array of array[tMaterieFeldInhalt,boolean] of double; // Materiefelder (getrennt für verschiedene Teilchenspezies) und deren Zeitableitungen - emWerte: array[tEMFeldInhalt,boolean] of extended; + emWerte: array[tEMFeldInhalt,boolean] of double; // EM-Felder und deren Zeitableitungen // A, p[xyz]? und i?Gamma haben keine sinnvolle Ableitung hier! lN,rN: tWertePunkt; @@ -53,7 +61,7 @@ type destructor destroy; override; procedure berechneGammaUndP; procedure berechneNAbleitung(iDX,pDNMax: extended); - procedure berechneAbleitungen(dT,dX,iDX: extended); + procedure berechneAbleitungen(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; @@ -64,14 +72,14 @@ type private matAnz: longint; lichters: array[boolean] of tMyStringlist; - procedure setzeRaender(iDX: extended); + procedure setzeRaender(dT,iDX,iDT: extended); public inhalt: array of tWertePunkt; par: tGitter; constructor create(groesse: longint; dichten,lichter: tMyStringList; parent: tGitter); destructor destroy; override; - procedure berechneAbleitungen(dT,dX,iDX,pDNMax: extended); + procedure berechneAbleitungen(dT,dX,iDT,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; @@ -91,11 +99,11 @@ type aktuelleFelder: longint; felders: array of tFelder; // mehrere komplette Simulationsboxen von Feldern, benötigt um Zwischenschritte für die Zeitentwicklung zu speichern - constructor create(size: longint; deltaT,deltaX,pDNMa: extended; bekannteWerte: tKnownValues; dichten, lichter: tMyStringlist; zv: tZeitverfahren; protokollant: tProtokollant); + constructor create(size: longint; deltaT,deltaX,pDNMa: extended; bekannteWerte: tKnownValues; dichten, lichter: tMyStringlist; zv: tZeitverfahren; protokollant: tProtokollant; name: string); destructor destroy; override; - procedure iteriereSchritt(dT: extended); - procedure macheAusgabe(aD: tAusgabedatei; sDX: extended); + procedure iteriereSchritt(dT,iDT: extended); function gibErhaltungsgroessen: string; + procedure berechne(was: char; teilchen: longint); end; { tSimulation } @@ -108,7 +116,7 @@ type fortschrittsAnzeige: boolean; ausgabeDateien: array of tAusgabeDatei; public - constructor create(inName: string; Protokollant: tProtokollant); + constructor create(inName: string; protokollant: tProtokollant; name: string); destructor destroy; override; function iteriereSchritt(start: extended; var zeitPhysik,zeitDatei: extended): boolean; // noch nicht zu Ende? end; @@ -125,7 +133,7 @@ const { tAusgabeDatei } -constructor tAusgabeDatei.create(feldName,prefix,suffix: string; prot: tProtokollant); +constructor tAusgabeDatei.create(feldName,prefix,suffix: string; prot: tProtokollant; nam: string); var emF: tEMFeldInhalt; maF: tMaterieFeldInhalt; @@ -133,6 +141,7 @@ var abl,gef: boolean; begin inherited create; + pro:=tProtokollant.create(prot,nam); num:=length(feldName); while (num>0) and (feldName[num] in ['0'..'9']) do dec(num); @@ -147,6 +156,9 @@ begin emInhalt:=efA; matInhalt:=mfN; feldName:=ansiUpperCase(feldName); + nNum:=0; // Header 0 wurde also noch nicht geschrieben + zeitAnz:=-1; + tCnt:=-1; gef:=false; for abl:=false to true do begin @@ -168,29 +180,108 @@ begin end; if not gef then begin - prot.schreibe('tAusgabeDatei.create kennt Feld '''+feldName+''' nicht!',true); + pro.schreibe('tAusgabeDatei.create kennt Feld '''+feldName+''' nicht!',true); + pro.destroyAll; halt(1); end; if teilchen>=0 then - name:=matFeldNamen[maF]+inttostr(teilchen) + pre:=matFeldNamen[maF]+inttostr(teilchen) else - name:=emFeldNamen[emF]; + pre:=emFeldNamen[emF]; - if abl then - name:='d'+name; - name:=prefix+name+suffix; + if ableitung then + pre:='d'+pre; - assignFile(datei,name); - rewrite(datei); + pre:=prefix+pre; + suf:=suffix; end; destructor tAusgabeDatei.destroy; begin - closeFile(datei); + if (tCnt<>zeitAnz) and ((nNum<>0) or (tCnt<>-1)) then begin + pro.schreibe('Falsche Anzahl an Zeitschritten in '+pre+'-'+inttostr(nNum-1)+suf+' geschrieben ('+inttostr(tCnt)+' statt '+inttostr(zeitAnz)+')!',true); + halt(1); + end; inherited destroy; end; +function tAusgabeDatei.dump: string; +begin + if teilchen<0 then + result:=emFeldNamen[emInhalt] + else + result:=matFeldNamen[matInhalt]+inttostr(teilchen); + if ableitung then + result:=result+''''; + result:=result+' -> '+pre+'-$i'+suf; +end; + +procedure tAusgabeDatei.schreibeKopf; +begin + if (tCnt<>zeitAnz) and ((nNum<>0) or (tCnt<>-1)) then begin + pro.schreibe('Falsche Anzahl an Zeitschritten in '+pre+'-'+inttostr(nNum-1)+suf+' geschrieben ('+inttostr(tCnt)+' statt '+inttostr(zeitAnz)+')!',true); + halt(1); + end; + tCnt:=0; + assignfile(datei,pre+'-'+inttostr(nNum)+suf); + inc(nNum); + rewrite(datei,1); + blockwrite(datei,nNum,sizeof(integer)); + blockwrite(datei,zeitAnz,sizeof(integer)); + closefile(datei); +end; + +procedure tAusgabeDatei.speichereWerte(gitter: tGitter; sDX: extended); +var i,cnt: longint; + sX,cX: extended; +begin + if (teilchen<0) and (emInhalt=efA) then + gitter.berechne('A',teilchen); + if (teilchen>=0) and (matInhalt=mfP) then + gitter.berechne('P',teilchen); + + if gitter.t>=nNum+sDT/2 then + schreibeKopf; + + cnt:=floor((length(gitter.felders[gitter.aktuelleFelder].inhalt)-1)/sDX*gitter.dX+1); + cX:=gitter.xl; + sX:=cX+(cnt-1)*sDX; + + inc(tCnt); + reset(datei,1); + seek(datei,fileSize(datei)); + blockWrite(datei,cX,sizeof(extended)); + blockWrite(datei,sX,sizeof(extended)); + blockWrite(datei,cnt,sizeof(integer)); + sX:=cX-Min(gitter.dX,sDX)/2; + if teilchen<0 then begin + for i:=0 to length(gitter.felders[gitter.aktuelleFelder].inhalt)-1 do begin + while cX>=sX do begin + dec(cnt); + blockWrite(datei,gitter.felders[gitter.aktuelleFelder].inhalt[i].emWerte[emInhalt,ableitung],sizeof(double)); + sX:=sX+sDX; + end; + cX:=cX+gitter.dX; + end; + end + else begin + for i:=0 to length(gitter.felders[gitter.aktuelleFelder].inhalt)-1 do begin + while cX>=sX do begin + dec(cnt); + blockWrite(datei,gitter.felders[gitter.aktuelleFelder].inhalt[i].matWerte[teilchen,matInhalt,ableitung],sizeof(double)); + sX:=sX+sDX; + end; + cX:=cX+gitter.dX; + end; + end; + closeFile(datei); + if cnt<>0 then begin + pro.schreibe('Falsche Anzahl an Ortsschritten geschrieben ('+inttostr(cnt)+')!',true); + halt(1); + end; +end; + { tWertePunkt } constructor tWertePunkt.create(linkerNachbar: tWertePunkt; materien: longint); @@ -275,7 +366,7 @@ begin // alpha = pDNMax * dX end; -procedure tWertePunkt.berechneAbleitungen(dT,dX,iDX: extended); // Zeitableitungen berechnen +procedure tWertePunkt.berechneAbleitungen(dX,iDX: extended); // Zeitableitungen berechnen var i: longint; begin @@ -316,7 +407,9 @@ begin end; // dA/dt = dA/dt - emWerte[efAX,true]:=emWerte[efDAXDT,false] + emWerte[efDAXDT,true]*dT; + emWerte[efAX,true]:=emWerte[efDAXDT,false];// + emWerte[efDAXDT,true]*dT; + emWerte[efAY,true]:=emWerte[efDAYDT,false];// + emWerte[efDAYDT,true]*dT; + emWerte[efAZ,true]:=emWerte[efDAZDT,false];// + emWerte[efDAZDT,true]*dT; end; procedure tWertePunkt.liKo(in1,in2: tWertePunkt; fak2: extended); // Werte werden auf (in1 + fak2*in2') gesetzt @@ -399,7 +492,7 @@ begin end; parent.kvs.rem('x'); for rechts:=false to true do begin - lichters[rechts]:=tMyStringlist.create(nil); + lichters[rechts]:=tMyStringlist.create(nil,''); lichters[rechts].text:=lichter.text; end; lichters[false].grep('^links .*'); @@ -415,9 +508,12 @@ begin inherited destroy; end; -procedure tFelder.setzeRaender(iDX: extended); +procedure tFelder.setzeRaender(dT,iDX,iDT: extended); var - emF: tEMFeldInhalt; + emF: tEMFeldInhalt; + rechts: boolean; + i: longint; + oVal: extended; begin for emF:=efAX to efAZ do begin // Vakuumrandbedingungen für das A-Feld inhalt[0].emWerte[emF,true]:= @@ -428,15 +524,18 @@ begin 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; *) + for rechts:=false to true do + for i:=0 to lichters[rechts].count-1 do begin + par.kvs.add('t',par.t-dT); + oVal:=exprToFloat(false,lichters[rechts][i],par.kvs,nil); + par.kvs.add('t',par.t); + inhalt[(length(inhalt)-1)*byte(rechts)].emWerte[efAY,true]:= + inhalt[(length(inhalt)-1)*byte(rechts)].emWerte[efAY,true] + + (exprToFloat(false,lichters[rechts][i],par.kvs,nil) - oVal)*iDT; + end; end; -procedure tFelder.berechneAbleitungen(dt,dx,iDX,pDNMax: extended); +procedure tFelder.berechneAbleitungen(dT,dX,iDT,iDX,pDNMax: extended); var i: longint; begin @@ -449,12 +548,12 @@ begin -abs(inhalt[length(inhalt)-2].matWerte[i,mfPX,false]); end; - setzeRaender(iDX); + setzeRaender(dT,iDX,iDT); 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); + inhalt[i].berechneAbleitungen(dX,iDX); end; procedure tFelder.liKo(in1,in2: tFelder; fak2: extended); // Werte werden auf (in1 + fak2*in2') gesetzt @@ -475,14 +574,14 @@ end; { tGitter } -constructor tGitter.create(size: longint; deltaT,deltaX,pDNMa: extended; bekannteWerte: tKnownValues; dichten, lichter: tMyStringlist; zv: tZeitverfahren; protokollant: tProtokollant); +constructor tGitter.create(size: longint; deltaT,deltaX,pDNMa: extended; bekannteWerte: tKnownValues; dichten, lichter: tMyStringlist; zv: tZeitverfahren; protokollant: tProtokollant; name: string); var i: longint; begin inherited create; zeitverfahren:=ZV; kvs:=bekannteWerte; - Prot:=TProtokollant.create(Protokollant,'tGitter'); + prot:=tProtokollant.create(protokollant,name); dTMaximum:=deltaT; dX:=deltaX; iDX:=1/dX; @@ -502,6 +601,7 @@ begin aktuelleFelder:=0; t:=0; + kvs.add('t',t); end; destructor tGitter.destroy; @@ -514,28 +614,30 @@ begin inherited destroy; end; -procedure tGitter.iteriereSchritt(dT: extended); +procedure tGitter.iteriereSchritt(dT,iDT: extended); var i: longint; {$IFDEF Zeitschrittueberwachung} pro: tProtokollant; j: longint; {$ENDIF} begin - felders[aktuelleFelder].berechneAbleitungen(dT,dX,iDX,pDNMax); // y' = y'(t,y(t)) + t:=t+dT; + kvs.add('t',t); + felders[aktuelleFelder].berechneAbleitungen(dT,dX,iDT,iDX,pDNMax); // y' = y'(t,y(t)) case zeitverfahren of zfEulerVorwaerts: // y(t+dt) = y(t) + y' dt felders[1-aktuelleFelder].liKo(felders[aktuelleFelder],felders[aktuelleFelder],dT); zfRungeKuttaVier: begin - 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[2].liKo(felders[aktuelleFelder],felders[aktuelleFelder],dT/2); // ya = y(t) + y' dt/2 + felders[2].berechneAbleitungen(dT/2,dX,iDT,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[3].liKo(felders[aktuelleFelder],felders[2],dT/2); // yb = y(t) + ya' dt/2 + felders[3].berechneAbleitungen(dT/2,dX,iDT,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[4].liKo(felders[aktuelleFelder],felders[3],dT); // yc = y(t) + yb' dt + felders[4].berechneAbleitungen(dT/2,dX,iDT,iDX,pDNMax); // yc' = y'(t+dt,yc) felders[1-aktuelleFelder].liKo( felders[aktuelleFelder], @@ -551,7 +653,6 @@ begin end; end{of case}; - t:=t+dT; aktuelleFelder:=1-aktuelleFelder; {$IFDEF Dichteueberwachung} with felders[aktuelleFelder] do @@ -581,58 +682,10 @@ begin {$ENDIF} end; -procedure tGitter.macheAusgabe(aD: tAusgabedatei; sDX: extended); -var i: longint; - sX,cX,tmp: extended; - emF: tEMFeldInhalt; - maF: tMaterieFeldInhalt; -begin - if (aD.teilchen<0) and (aD.emInhalt=efA) then - for i:=0 to length(felders[aktuelleFelder].inhalt)-1 do begin - tmp:=0; - 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; - 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; - 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: integer; ns: tExtendedArray; - Pro: TProtokollant; + pro: tProtokollant; begin pro:=tProtokollant.create(prot,'gibErhaltungsgroessen'); setlength(ns,felders[aktuelleFelder].matAnz); @@ -684,23 +737,53 @@ begin end; end; +procedure tGitter.berechne(was: char; teilchen: longint); +var + i: longint; + tmp: extended; + emF: tEMFeldInhalt; + maF: tMaterieFeldInhalt; +begin + if was='A' then begin + for i:=0 to length(felders[aktuelleFelder].inhalt)-1 do begin + tmp:=0; + 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; + exit; + end; + if was='P' then begin + 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[teilchen,maF,false]); + felders[aktuelleFelder].inhalt[i].matWerte[teilchen,maF,false]:=sqrt(tmp); + end; + exit; + end; + prot.schreibe('Kann '''+was+''' nicht berechnen, weil ich es nicht verstehe!',true); + halt; +end; + { tSimulation } -constructor tSimulation.create(inName: string; Protokollant: tProtokollant); +constructor tSimulation.create(inName: string; protokollant: tProtokollant; name: string); var ifile: tMyStringlist; zeitverfahren: tZeitverfahren; s,t,aSuffix,aPrefix: string; - deltaX,deltaT,breite,pDNMax: extended; - i,j: longint; + deltaX,breite,pDNMax: extended; + i: longint; kvs: tKnownValues; dichten,lichter: tMyStringlist; + pro: tProtokollant; begin inherited create; - prot:=tProtokollant.create(Protokollant,'tSimulation'); + prot:=tProtokollant.create(protokollant,name); kvs:=tKnownValues.create; - ifile:=tMyStringlist.create(prot); + ifile:=tMyStringlist.create(prot,'ifile'); ifile.loadfromfile(inName); if not ifile.unfoldMacros then begin prot.schreibe('Fehlerhafte Macros in Parameterdatei '''+inName+'''!',true); @@ -712,8 +795,9 @@ begin deltaX:=1e-2; kvs.add('λ',1); kvs.add('T',1); + kvs.add('ω',2*pi); kvs.add('dX',deltaX); - deltaT:=-1; + dT:=-1; sDT:=-1; sDX:=-1; tEnde:=100; @@ -727,18 +811,26 @@ begin setlength(ausgabeDateien,0); // Standardeinstellungen Breich 'lichtVon...' - lichter:=tMyStringlist.create(prot); + lichter:=tMyStringlist.create(prot,'lichter'); // Standardeinstellungen Breich 'teilchen...' - dichten:=tMyStringlist.create(prot); + dichten:=tMyStringlist.create(prot,'teilchen'); + + repeat + if not ifile.readln(s) then begin + prot.schreibe('Unerwartetes Dateiende in '''+inName+'''!',true); + halt(1); + end; + + if s='Dateiende' then + break; - 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; writeln(s); + end; if s='allgemeinEnde' then break; if s='runge-Kutta-4' then begin Zeitverfahren:=zfRungeKuttaVier; @@ -754,8 +846,8 @@ begin continue; end; if startetMit('zeitschritt',s) then begin - deltaT:=exprtofloat(false,s,kvs,nil); - kvs.add('dT',deltaT); + dT:=exprtofloat(false,s,kvs,nil); + kvs.add('dT',dT); continue; end; if startetMit('diffusionsterm',s) then begin @@ -812,7 +904,7 @@ begin while s<>'' do begin t:=erstesArgument(s,','); setlength(ausgabeDateien,length(ausgabeDateien)+1); - ausgabeDateien[length(ausgabeDateien)-1]:=tAusgabeDatei.create(t,aPrefix,aSuffix,prot); + ausgabeDateien[length(ausgabeDateien)-1]:=tAusgabeDatei.create(t,aPrefix,aSuffix,prot,'ausgabeDateien['+inttostr(length(ausgabeDateien)-1)+']'); end; continue; end; @@ -822,32 +914,73 @@ begin continue; end; + if startetMit('licht von ',s) then begin + if startetMit('links ',s) then begin + lichter.add('links '+s); + continue; + end; + if startetMit('rechts ',s) then begin + lichter.add('rechts '+s); + continue; + end; + prot.schreibe('Licht kann nicht von '''+erstesArgument(s)+''' kommen!',true); + halt(1); + end; + prot.schreibe('Unbekannter Befehl '''+s+''' in Parameterdatei '''+inName+'''!',true); halt(1); - end; + until false; + + ifile.free; if length(ausgabedateien)=0 then begin prot.schreibe('Du solltest irgendetwas abspeichern lassen!',true); halt(1); end; - if deltaT<0 then - deltaT:=deltaX/10; + if dT<0 then + dT:=deltaX/10; if sDT<0 then - sDT:=deltaT; + sDT:=dT; + sDT:=1/round(1/sDT); + sT:=sDT/2; 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(tEnde/sDT+1); - BlockWrite(Ausgabedateien[i].Datei,j,sizeof(longint)); + pro:=tProtokollant.create(prot,'create'); + case zeitverfahren of + zfEulerVorwaerts: + pro.schreibe('Iteration mittels Euler-Vorwärts'); + zfRungeKuttaVier: + pro.schreibe('Iteration mittels Runge-Kutta-4'); + else + pro.schreibe('Iteration mittels unbekanntem Verfahren'); + end{of case}; + pro.schreibe('dX = '+floattostr(deltaX)); + pro.schreibe('dT = '+floattostr(dT)); + pro.schreibe('breite = '+floattostr(breite)+' (= '+inttostr(round(breite/deltaX)+1)+' Punkte)'); + pro.schreibe('pDNMax = '+floattostr(pDNMax)); + pro.schreibe('bekannte Werte:'); + kvs.dump(pro,' '); + if dichten.count>0 then begin + pro.schreibe('dichten:'); + dichten.dump(pro,' '); end; - ifile.free; + if lichter.count>0 then begin + pro.schreibe('lichter:'); + lichter.dump(pro,' '); + end; + if length(ausgabeDateien)>0 then begin + pro.schreibe('ausgaben:'); + for i:=0 to length(ausgabeDateien)-1 do begin + ausgabeDateien[i].zeitAnz:=round(1/sDT); + ausgabeDateien[i].sDT:=sDT; + pro.schreibe(' '+ausgabeDateien[i].dump); + end; + end; + pro.free; - gitter:=tGitter.create(round(breite/deltaX),deltaT,deltaX,pDNMax,kvs,dichten,lichter,zeitverfahren,prot); + gitter:=tGitter.create(round(breite/deltaX)+1,dT,deltaX,pDNMax,kvs,dichten,lichter,zeitverfahren,prot,'gitter'); dichten.free; lichter.free; end; @@ -862,51 +995,33 @@ 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); + gitter.iteriereSchritt(dT,1/dT); zeitPhysik:=zeitPhysik+now; zeitDatei:=zeitDatei-now; - while gitter.t>=sT do begin + while (gitter.t>=sT) and (sT<tEnde) do begin sT:=sT+sDT; - for i:=0 to length(Ausgabedateien)-1 do - gitter.macheAusgabe(ausgabeDateien[i],sDX); + for i:=0 to length(ausgabeDateien)-1 do + ausgabeDateien[i].speichereWerte(gitter,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}; + if (floor(100*Gitter.t/tEnde) < floor(100*(Gitter.t+dT)/tEnde)) then 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; - result:=gitter.t>=tEnde; + result:=gitter.t<tEnde; end; end. diff --git a/Plasmapropagation b/Plasmapropagation Binary files differindex 98472f4..3501611 100755 --- a/Plasmapropagation +++ b/Plasmapropagation diff --git a/Plasmapropagation.lpi b/Plasmapropagation.lpi index f707ff5..28ad15a 100644 --- a/Plasmapropagation.lpi +++ b/Plasmapropagation.lpi @@ -68,7 +68,6 @@ <OverflowChecks Value="True"/> <StackChecks Value="True"/> </Checks> - <VerifyObjMethodCallValidity Value="True"/> </CodeGeneration> <Linking> <Options> @@ -78,10 +77,7 @@ </Options> </Linking> <Other> - <Verbosity> - <ShoLineNum Value="True"/> - <ShowDebugInfo Value="True"/> - </Verbosity> + <WriteFPCLogo Value="False"/> </Other> </CompilerOptions> <Debugging> diff --git a/Plasmapropagation.lpr b/Plasmapropagation.lpr index 5c914d5..5559f51 100644 --- a/Plasmapropagation.lpr +++ b/Plasmapropagation.lpr @@ -8,9 +8,7 @@ uses {$ENDIF}{$ENDIF} Classes, SysUtils, CustApp, { you can add units after this } - math, Physikunit, crt, protokollunit; - -var ErrorCode: longint; + math, Physikunit, protokollunit; type @@ -45,7 +43,7 @@ begin zeitDatei:=0; zeitPhysik:=0; - simulation:=tSimulation.create(paramstr(1),prot); + simulation:=tSimulation.create(paramstr(1),prot,'simulation'); while simulation.iteriereSchritt(start,zeitPhysik,zeitDatei) do ; @@ -66,7 +64,7 @@ begin 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; + prot.free; // stop program loop Terminate; diff --git a/Plasmapropagation.lps b/Plasmapropagation.lps index cd5e8e9..cbe3abd 100644 --- a/Plasmapropagation.lps +++ b/Plasmapropagation.lps @@ -7,74 +7,72 @@ <Unit0> <Filename Value="Plasmapropagation.lpr"/> <IsPartOfProject Value="True"/> - <TopLine Value="28"/> - <CursorPos X="39" Y="11"/> - <UsageCount Value="140"/> + <CursorPos X="34" Y="11"/> + <UsageCount Value="149"/> <Loaded Value="True"/> </Unit0> <Unit1> <Filename Value="Physikunit.pas"/> <IsPartOfProject Value="True"/> <EditorIndex Value="1"/> - <TopLine Value="699"/> - <CursorPos X="13" Y="714"/> - <UsageCount Value="81"/> + <TopLine Value="438"/> + <CursorPos Y="538"/> + <FoldState Value=" T3k60-4 piYkM0O]9j8jg0M[943jb0J5]9SkI0O3 pnBo2051g"/> + <UsageCount Value="90"/> <Bookmarks Count="1"> - <Item0 Y="525"/> + <Item0 Y="627"/> </Bookmarks> <Loaded Value="True"/> </Unit1> <Unit2> <Filename Value="../units/protokollunit.pas"/> <IsPartOfProject Value="True"/> - <EditorIndex Value="-1"/> - <TopLine Value="26"/> - <CursorPos X="103" Y="47"/> - <UsageCount Value="43"/> + <EditorIndex Value="2"/> + <TopLine Value="20"/> + <CursorPos X="22" Y="31"/> + <UsageCount Value="52"/> + <Loaded Value="True"/> </Unit2> <Unit3> <Filename Value="input.plap"/> <IsPartOfProject Value="True"/> - <IsVisibleTab Value="True"/> - <EditorIndex Value="5"/> - <CursorPos X="7" Y="19"/> - <UsageCount Value="42"/> + <EditorIndex Value="4"/> + <CursorPos Y="22"/> + <UsageCount Value="51"/> <Loaded Value="True"/> <DefaultSyntaxHighlighter Value="None"/> </Unit3> <Unit4> <Filename Value="input.epost"/> - <EditorIndex Value="-1"/> - <TopLine Value="45"/> - <CursorPos Y="69"/> - <UsageCount Value="56"/> + <IsVisibleTab Value="True"/> + <EditorIndex Value="5"/> + <CursorPos X="27" Y="21"/> + <UsageCount Value="58"/> + <Loaded Value="True"/> <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"/> + <TopLine Value="316"/> + <CursorPos Y="356"/> + <UsageCount Value="19"/> <Loaded Value="True"/> </Unit5> <Unit6> <Filename Value="../units/lowlevelunit.pas"/> - <EditorIndex Value="4"/> - <TopLine Value="373"/> - <CursorPos X="77" Y="399"/> - <UsageCount Value="10"/> - <Loaded Value="True"/> + <EditorIndex Value="-1"/> + <TopLine Value="540"/> + <CursorPos X="20" Y="573"/> + <UsageCount Value="11"/> </Unit6> <Unit7> <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"/> + <EditorIndex Value="-1"/> + <TopLine Value="367"/> + <CursorPos X="17" Y="390"/> + <FoldState Value=" T3i3075 piZjD0WQ"/> + <UsageCount Value="12"/> </Unit7> <Unit8> <Filename Value="../epost/werteunit.pas"/> @@ -107,121 +105,122 @@ <JumpHistory Count="30" HistoryIndex="29"> <Position1> <Filename Value="Physikunit.pas"/> + <Caret Line="533" Column="45" TopLine="502"/> </Position1> <Position2> <Filename Value="Physikunit.pas"/> - <Caret Line="391" Column="52" TopLine="358"/> + <Caret Line="679" Column="54" TopLine="522"/> </Position2> <Position3> - <Filename Value="Plasmapropagation.lpr"/> - <Caret Line="38" TopLine="11"/> + <Filename Value="../units/protokollunit.pas"/> + <Caret Line="28" Column="23"/> </Position3> <Position4> - <Filename Value="Plasmapropagation.lpr"/> - <Caret Line="11" Column="39" TopLine="28"/> + <Filename Value="../units/protokollunit.pas"/> + <Caret Line="29" Column="23"/> </Position4> <Position5> - <Filename Value="Physikunit.pas"/> - <Caret Line="111" Column="17" TopLine="91"/> + <Filename Value="../units/protokollunit.pas"/> + <Caret Line="76" Column="15" TopLine="58"/> </Position5> <Position6> - <Filename Value="Physikunit.pas"/> - <Caret Line="36" Column="23" TopLine="16"/> + <Filename Value="../units/protokollunit.pas"/> + <Caret Line="114" TopLine="77"/> </Position6> <Position7> - <Filename Value="Physikunit.pas"/> - <Caret Line="671" TopLine="631"/> + <Filename Value="../units/protokollunit.pas"/> + <Caret Line="29" TopLine="29"/> </Position7> <Position8> - <Filename Value="../units/mystringlistunit.pas"/> - <Caret Line="421" Column="57" TopLine="386"/> + <Filename Value="../units/protokollunit.pas"/> + <Caret Line="72" Column="23" TopLine="52"/> </Position8> <Position9> - <Filename Value="../units/mystringlistunit.pas"/> + <Filename Value="../units/protokollunit.pas"/> + <Caret Line="31" Column="24"/> </Position9> <Position10> - <Filename Value="../units/mystringlistunit.pas"/> - <Caret Line="17" Column="23"/> + <Filename Value="../units/protokollunit.pas"/> + <Caret Line="114" Column="25" TopLine="80"/> </Position10> <Position11> - <Filename Value="../units/mystringlistunit.pas"/> - <Caret Line="27" Column="23"/> + <Filename Value="Plasmapropagation.lpr"/> + <Caret Line="35" Column="63"/> </Position11> <Position12> - <Filename Value="../units/mystringlistunit.pas"/> - <Caret Line="28" Column="23"/> + <Filename Value="Plasmapropagation.lpr"/> + <Caret Line="26" Column="11" TopLine="15"/> </Position12> <Position13> - <Filename Value="../units/mystringlistunit.pas"/> - <Caret Line="83" Column="10" TopLine="44"/> + <Filename Value="Physikunit.pas"/> + <Caret Line="509" TopLine="369"/> </Position13> <Position14> - <Filename Value="../units/mystringlistunit.pas"/> - <Caret Line="81" TopLine="42"/> + <Filename Value="Physikunit.pas"/> + <Caret Line="11" Column="84"/> </Position14> <Position15> - <Filename Value="../units/mystringlistunit.pas"/> - <Caret Line="27" Column="54" TopLine="16"/> + <Filename Value="Physikunit.pas"/> + <Caret Line="1010" Column="69" TopLine="990"/> </Position15> <Position16> <Filename Value="Physikunit.pas"/> - <Caret Line="732" Column="56" TopLine="702"/> + <Caret Line="1017" TopLine="764"/> </Position16> <Position17> <Filename Value="Physikunit.pas"/> - <Caret Line="402" Column="52" TopLine="382"/> + <Caret Line="937" Column="40" TopLine="915"/> </Position17> <Position18> - <Filename Value="../units/mystringlistunit.pas"/> - <Caret Line="264" Column="39" TopLine="256"/> + <Filename Value="Physikunit.pas"/> + <Caret Line="262" Column="50" TopLine="242"/> </Position18> <Position19> - <Filename Value="../units/mystringlistunit.pas"/> - <Caret Line="25" Column="6" TopLine="8"/> + <Filename Value="Physikunit.pas"/> + <Caret Line="55" Column="53" TopLine="23"/> </Position19> <Position20> - <Filename Value="../units/mystringlistunit.pas"/> - <Caret Line="83" Column="3" TopLine="46"/> + <Filename Value="Physikunit.pas"/> </Position20> <Position21> - <Filename Value="../units/mystringlistunit.pas"/> - <Caret Line="276" Column="48" TopLine="256"/> + <Filename Value="Physikunit.pas"/> + <Caret Line="15" Column="57"/> </Position21> <Position22> - <Filename Value="../units/mystringlistunit.pas"/> - <Caret Line="340" Column="30" TopLine="320"/> + <Filename Value="Physikunit.pas"/> + <Caret Line="410" Column="48" TopLine="382"/> </Position22> <Position23> - <Filename Value="../units/mystringlistunit.pas"/> - <Caret Line="27" Column="54" TopLine="17"/> + <Filename Value="Physikunit.pas"/> + <Caret Line="528" Column="149" TopLine="415"/> </Position23> <Position24> <Filename Value="Physikunit.pas"/> - <Caret Line="402" Column="47" TopLine="382"/> + <Caret Line="511"/> </Position24> <Position25> <Filename Value="Physikunit.pas"/> - <Caret Line="728" Column="38" TopLine="709"/> + <Caret Line="529" Column="35" TopLine="473"/> </Position25> <Position26> - <Filename Value="../units/matheunit.pas"/> - <Caret Line="102" Column="33" TopLine="25"/> + <Filename Value="Physikunit.pas"/> + <Caret Line="511" Column="19" TopLine="475"/> </Position26> <Position27> - <Filename Value="../units/matheunit.pas"/> - <Caret Line="45" Column="21" TopLine="25"/> + <Filename Value="Physikunit.pas"/> + <Caret Line="75" Column="31" TopLine="61"/> </Position27> <Position28> - <Filename Value="../units/matheunit.pas"/> - <Caret Line="347" Column="21" TopLine="314"/> + <Filename Value="Physikunit.pas"/> + <Caret Line="534" Column="61" TopLine="472"/> </Position28> <Position29> - <Filename Value="../units/matheunit.pas"/> - <Caret Line="443" Column="34" TopLine="425"/> + <Filename Value="Physikunit.pas"/> + <Caret Line="64" Column="35" TopLine="47"/> </Position29> <Position30> - <Filename Value="../units/lowlevelunit.pas"/> - <Caret Line="54" Column="18" TopLine="46"/> + <Filename Value="Physikunit.pas"/> + <Caret Line="369" Column="43" TopLine="344"/> </Position30> </JumpHistory> </ProjectSession> diff --git a/input.epost b/input.epost index 8f9dc9e..68a9f8f 100644 --- a/input.epost +++ b/input.epost @@ -18,7 +18,7 @@ Ende Threadanzahl: 10 # !Schleife: $Feld: dDPHIDX dDAYDT dAY dAX -!Schleife: $Feld: AX +!Schleife: $Feld: AY DAYDT Daten einlesen Genauigkeit: extended @@ -3,9 +3,9 @@ allgemein runge-Kutta-4 ortsschritt 10^-2 * λ - zeitschritt 10^-3 * T - zeit 100 * T - breite 10 * λ + zeitschritt 10^-2 * T + zeit 20 * T + breite 5 * λ mit Fortschrittsanzeige allgemeinEnde @@ -15,12 +15,12 @@ ausgaben # felder AX,AY,dAYDT,N,dPhiDX,dPsiDX,Gamma ausgabenEnde -!setze $tFwhm: 5 T -!setze $tMitte: 20 T +!setze $tFwhm: (2.5 * T) +!setze $tMitte: (1 * T) -lichtVonLinks 0.5 * 2^(-2*((t-$tMitte)/$tFwhm)^2) * Sin(ω*t) +licht von links 0.5 * 2^(-2*((t-$tMitte)/$tFwhm)^2) * sin(ω*t) -dateiEnde +Dateiende teilchen1 ladung -q |