From b2b8170095b90d3944096e76077d0fa9bef1187b Mon Sep 17 00:00:00 2001 From: Erich Eckner Date: Tue, 29 Sep 2015 16:41:27 +0200 Subject: extended -> double scheint viele Probleme zu lösen :-/ MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit --- Physikunit.pas | 212 +++++++++++++++++++++++++------------------------- Plasmapropagation.lpr | 43 +--------- Plasmapropagation.lps | 128 ++++++++++++++++-------------- 3 files changed, 176 insertions(+), 207 deletions(-) diff --git a/Physikunit.pas b/Physikunit.pas index b49675b..f95dc32 100644 --- a/Physikunit.pas +++ b/Physikunit.pas @@ -13,11 +13,11 @@ unit Physikunit; interface uses - Classes, SysUtils, Math, protokollunit, matheunit, mystringlistunit, lowlevelunit, baseUnix, fftw_l; + Classes, SysUtils, Math, protokollunit, matheunit, mystringlistunit, lowlevelunit, baseUnix, fftw; type tZeitverfahren = (zfEulerVorwaerts,zfRungeKuttaDreiAchtel,zfRungeKuttaVier,zfRungeKuttaZehn,zfRungeKuttaZwoelf,zfRungeKuttaVierzehn); - tVerteilungsfunktion = function(x: extended): extended; + tVerteilungsfunktion = function(x: double): double; tEMFeldGroesze = ( efAX,efAY, efDAXDT,efDAYDT, @@ -29,7 +29,7 @@ type tMaterieFeldGroesze = (mfPY); tMaterieSpeicherGroesze = (msPX,msPXSqr,msPY,msVX,msVY,msN,msPXRipple); tTeilchenSpeziesGroeszen = (tsgMasse,tsgIQdrMasse,tsgLadung); - tEMQuellen = array[tEMQuellGroesze] of extended; + tEMQuellen = array[tEMQuellGroesze] of double; const erstesEMFMitAbleitung = efAX; @@ -60,36 +60,36 @@ type procedure schreibeBuffer(minBufLen: longint); public zeitAnz: longint; - sDT: extended; + sDT: double; constructor create(feldName,prefix,suffix: string; prot: tProtokollant; nam: string); destructor destroy; override; function dump: string; procedure schreibeKopf; - procedure speichereWerte(gitter: tGitter; sT, sDX: extended); + procedure speichereWerte(gitter: tGitter; sT, sDX: double); end; { tTeilchenSpezies } tTeilchenSpezies = class private - ortsGrenzen: array of extended; + ortsGrenzen: array of double; dichteStuecke: array of string; public - nMin,nMax,breite: extended; // minimale und maximale Massendichte in nc; Breite der Impulsverteilung - eigenschaften: array[tTeilchenSpeziesGroeszen] of extended; // q in e; m in me - diffusion: array[0..1] of extended; // Diffusionsterm in x- und p-Richtung (= maximales d(ln n)/d[xp] * iD[XP]) + nMin,nMax,breite: double; // minimale und maximale Massendichte in nc; Breite der Impulsverteilung + eigenschaften: array[tTeilchenSpeziesGroeszen] of double; // q in e; m in me + diffusion: array[0..1] of double; // Diffusionsterm in x- und p-Richtung (= maximales d(ln n)/d[xp] * iD[XP]) constructor create; overload; constructor create(original: tTeilchenSpezies); overload; destructor destroy; override; procedure dump(prot: tProtokollant; prefix: string); - function gibDichte(x: extended; kvs: tKnownValues): extended; // Massendichte + function gibDichte(x: double; kvs: tKnownValues): double; // Massendichte procedure liesDichteFunktion(rest: string; ifile: tMyStringList; kvs: tKnownValues; teilchen: array of tTeilchenSpezies; prot: tProtokollant); end; { TFelder } (* - procedure berechneAbleitungen(pX: extended); inline; - function nichtnegativieren: extended; inline; + procedure berechneAbleitungen(pX: double); inline; + function nichtnegativieren: double; inline; *) tFelder = class(tObject) // repräsentiert eine Simulationsbox von Feldern inklusive deren Ableitungen @@ -99,35 +99,35 @@ type teilchen: array of tTeilchenSpezies; lichters: array[boolean] of tMyStringlist; - massen: array of extended; // Gesamtmassen je Teilchenspezies zur Überwachung der Erhaltungsgrößen - gesamtDefizit: extended; - impulsRaumGradient: array of array[0..1] of pExtended; // Ableitung des Impulsraumes nach x und px - iMGamma: array of pExtended; // 1/m/gamma (= v/p) + massen: array of double; // Gesamtmassen je Teilchenspezies zur Überwachung der Erhaltungsgrößen + gesamtDefizit: double; + impulsRaumGradient: array of array[0..1] of pDouble; // Ableitung des Impulsraumes nach x und px + iMGamma: array of pDouble; // 1/m/gamma (= v/p) _aX,_aP: longint; - _dX,_dP,iDX: extended; - ffts: array of array[0..1,boolean] of fftw_plan_extended; // fft-Algorithmen - fftTmp: Pcomplex_extended; // Zwischenspeicher für ffts - procedure initialisiereDichte(ort,tlc: longint; breite,n: extended); inline; + _dX,_dP,iDX: double; + ffts: array of array[0..1,boolean] of fftw_plan_double; // fft-Algorithmen + fftTmp: Pcomplex_double; // Zwischenspeicher für ffts + procedure initialisiereDichte(ort,tlc: longint; breite,n: double); inline; public - emFelder: array[tEMFeldGroesze,boolean] of pExtended; // EM-Felder und deren Zeitableitungen + emFelder: array[tEMFeldGroesze,boolean] of pDouble; // EM-Felder und deren Zeitableitungen // dPhiDX und B[xyz] haben keine sinnvolle Ableitung hier! - emQuellen: array[tEMQuellGroesze] of pExtended; // Quellen für die EM-Felder - matFelder: array of array[tMaterieFeldGroesze] of pExtended; // px-unabhängige Felder, nach Teilchenspezies sortiert - impulsraum: array of array[boolean] of pExtended; // px-abhängige Felder, nach Teilchenspezies sortiert, sowie deren Ableitung + emQuellen: array[tEMQuellGroesze] of pDouble; // Quellen für die EM-Felder + matFelder: array of array[tMaterieFeldGroesze] of pDouble; // px-unabhängige Felder, nach Teilchenspezies sortiert + impulsraum: array of array[boolean] of pDouble; // px-abhängige Felder, nach Teilchenspezies sortiert, sowie deren Ableitung gitter: tGitter; - constructor create(anzX,anzP: longint; deltaX,deltaP: extended; _teilchen: array of tTeilchenSpezies; lichter: tMyStringList; parent: tGitter); + constructor create(anzX,anzP: longint; deltaX,deltaP: double; _teilchen: array of tTeilchenSpezies; lichter: tMyStringList; parent: tGitter); destructor destroy; override; procedure berechneAbleitungen; procedure setzeNull; inline; procedure dumpErhaltungsgroeszen(prot: tProtokollant); property aX: longint read _aX; property aP: longint read _aP; - property dX: extended read _dX; - property dP: extended read _dP; - function impulsIntegral(ort,tlc: longint; emQ: tEMQuellGroesze): extended; overload; inline; - function impulsIntegral(ort,tlc: longint; maF: tMaterieSpeicherGroesze): extended; overload; inline; + property dX: double read _dX; + property dP: double read _dP; + function impulsIntegral(ort,tlc: longint; emQ: tEMQuellGroesze): double; overload; inline; + function impulsIntegral(ort,tlc: longint; maF: tMaterieSpeicherGroesze): double; overload; inline; {$DEFINE LiKoInterface} {$INCLUDE linearkombinationen.inc} {$UNDEF LiKoInterface} @@ -141,7 +141,7 @@ type prot: tProtokollant; zeitverfahren: tZeitverfahren; abbruch: boolean; - dX,iDX,xl,t: extended; + dX,iDX,xl,t: double; kvs: tKnownValues; procedure abbrechen; public @@ -149,9 +149,9 @@ type aktuelleFelder: longint; felders: array of tFelder; // mehrere komplette Simulationsboxen von Feldern, benötigt um Zwischenschritte für die Zeitentwicklung zu speichern - constructor create(derBesitzer: tSimulation; aX,aP: longint; deltaX,deltaP: extended; bekannteWerte: tKnownValues; teilchen: array of tTeilchenSpezies; lichter: tMyStringlist; zv: tZeitverfahren; protokollant: tProtokollant; name: string); + constructor create(derBesitzer: tSimulation; aX,aP: longint; deltaX,deltaP: double; bekannteWerte: tKnownValues; teilchen: array of tTeilchenSpezies; lichter: tMyStringlist; zv: tZeitverfahren; protokollant: tProtokollant; name: string); destructor destroy; override; - procedure iteriereSchritt(dT: extended); + procedure iteriereSchritt(dT: double); procedure dumpErhaltungsgroeszen; end; @@ -161,14 +161,14 @@ type private prot: tProtokollant; gitter: tGitter; - dT,tEnde,sT,sDT,sDX: extended; + dT,tEnde,sT,sDT,sDX: double; fortschrittsAnzeige: boolean; ausgabeDateien: array of tAusgabeDatei; public gotSighup: boolean; constructor create(inName: string; protokollant: tProtokollant; name: string); destructor destroy; override; - function iteriereSchritt(start: extended; var zeitPhysik,zeitDatei: extended): boolean; // noch nicht zu Ende? + function iteriereSchritt(start: double; var zeitPhysik,zeitDatei: double): boolean; // noch nicht zu Ende? end; procedure SignalCapture(signal : longint); cdecl; @@ -349,9 +349,9 @@ begin closefile(datei); end; -procedure tAusgabeDatei.speichereWerte(gitter: tGitter; sT, sDX: extended); +procedure tAusgabeDatei.speichereWerte(gitter: tGitter; sT, sDX: double); var i,cnt: longint; - sX,cX,val: extended; + sX,cX,val: double; begin if (teilchen>=length(gitter.felders[gitter.aktuelleFelder].teilchen)) then begin pro.schreibe('Teilchen '+inttostr(teilchen+1)+' gibt es nicht, da kann ich auch nichts speichern!',true); @@ -367,19 +367,19 @@ begin sX:=cX+(cnt-1)*sDX; inc(tCnt); - if bufLen-bufPos < sizeof(integer) + (2+cnt)*sizeof(extended) then - schreibeBuffer(sizeof(integer) + (2+cnt)*sizeof(extended)); + if bufLen-bufPos < sizeof(integer) + (2+cnt)*sizeof(double) then + schreibeBuffer(sizeof(integer) + (2+cnt)*sizeof(double)); - if bufLen-bufPos < sizeof(integer) + (2+cnt)*sizeof(extended) then begin - pro.schreibe('Schreibbuffer ist zu klein ('+inttostr(bufLen)+' statt mindestens '+inttostr(sizeof(integer) + (2+cnt)*sizeof(extended))+' Bytes)!'); + if bufLen-bufPos < sizeof(integer) + (2+cnt)*sizeof(double) then begin + pro.schreibe('Schreibbuffer ist zu klein ('+inttostr(bufLen)+' statt mindestens '+inttostr(sizeof(integer) + (2+cnt)*sizeof(double))+' Bytes)!'); pro.destroyall; halt(1); end; - move(cX,(buf+bufPos)^,sizeof(extended)); - inc(bufPos,sizeof(extended)); - move(sX,(buf+bufPos)^,sizeof(extended)); - inc(bufPos,sizeof(extended)); + move(cX,(buf+bufPos)^,sizeof(double)); + inc(bufPos,sizeof(double)); + move(sX,(buf+bufPos)^,sizeof(double)); + inc(bufPos,sizeof(double)); move(cnt,(buf+bufPos)^,sizeof(integer)); inc(bufPos,sizeof(integer)); @@ -389,8 +389,8 @@ begin for i:=0 to gitter.felders[gitter.aktuelleFelder].aX-1 do begin while cX>=sX do begin dec(cnt); - move((gitter.felders[gitter.aktuelleFelder].emFelder[emFeld,ableitung]+i)^,(buf+bufPos)^,sizeof(extended)); - inc(bufPos,sizeof(extended)); + move((gitter.felders[gitter.aktuelleFelder].emFelder[emFeld,ableitung]+i)^,(buf+bufPos)^,sizeof(double)); + inc(bufPos,sizeof(double)); sX:=sX+sDX; end; cX:=cX+gitter.dX; @@ -400,8 +400,8 @@ begin while cX>=sX do begin dec(cnt); val:=gitter.felders[gitter.aktuelleFelder].impulsIntegral(i,teilchen,emQuelle); - move(val,(buf+bufPos)^,sizeof(extended)); - inc(bufPos,sizeof(extended)); + move(val,(buf+bufPos)^,sizeof(double)); + inc(bufPos,sizeof(double)); sX:=sX+sDX; end; cX:=cX+gitter.dX; @@ -411,8 +411,8 @@ begin while cX>=sX do begin dec(cnt); val:=gitter.felders[gitter.aktuelleFelder].impulsIntegral(i,teilchen,matFeld); - move(val,(buf+bufPos)^,sizeof(extended)); - inc(bufPos,sizeof(extended)); + move(val,(buf+bufPos)^,sizeof(double)); + inc(bufPos,sizeof(double)); sX:=sX+sDX; end; cX:=cX+gitter.dX; @@ -484,7 +484,7 @@ begin end; end; -function tTeilchenSpezies.gibDichte(x: extended; kvs: tKnownValues): extended; +function tTeilchenSpezies.gibDichte(x: double; kvs: tKnownValues): double; var i: longint; begin @@ -547,7 +547,7 @@ begin end; (* -function tRaumPunkt.nichtnegativieren: extended; // Dichten nicht negativ machen +function tRaumPunkt.nichtnegativieren: double; // Dichten nicht negativ machen var i: longint; begin @@ -559,11 +559,11 @@ end; *) // tFelder ********************************************************************* -constructor tFelder.create(anzX,anzP: longint; deltaX,deltaP: extended; _teilchen: array of tTeilchenSpezies; lichter: tMyStringList; parent: tGitter); +constructor tFelder.create(anzX,anzP: longint; deltaX,deltaP: double; _teilchen: array of tTeilchenSpezies; lichter: tMyStringList; parent: tGitter); var i,j: longint; rechts,abl: boolean; - dens: extended; + dens: double; emF: tEMFeldGroesze; emQ: tEMQuellGroesze; maF: tMaterieFeldGroesze; @@ -585,39 +585,39 @@ begin setlength(impulsRaumGradient,length(teilchen)); for i:=0 to length(impulsRaumGradient)-1 do for j:=0 to length(impulsRaumGradient[i])-1 do begin - fftw_getmem(impulsRaumGradient[i,j],aP*aX*sizeof(extended)); - fillchar(impulsRaumGradient[i,j]^,aP*aX*sizeof(extended),#0); + fftw_getmem(impulsRaumGradient[i,j],aP*aX*sizeof(double)); + fillchar(impulsRaumGradient[i,j]^,aP*aX*sizeof(double),#0); end; fillchar(iMGamma,sizeof(iMGamma),#0); setlength(iMGamma,length(teilchen)); for i:=0 to length(iMGamma)-1 do begin - fftw_getmem(iMGamma[i],aP*aX*sizeof(extended)); - fillchar(iMGamma[i]^,aP*aX*sizeof(extended),#0); + fftw_getmem(iMGamma[i],aP*aX*sizeof(double)); + fillchar(iMGamma[i]^,aP*aX*sizeof(double),#0); end; for emF:=low(tEMFeldGroesze) to high(tEMFeldGroesze) do begin emFelder[emF,true]:=nil; for abl:=false to not(emF in [efEX,efBZ]) do begin - fftw_getmem(emFelder[emF,abl],aX*sizeof(extended)); - fillchar(emFelder[emF,abl]^,aX*sizeof(extended),#0); + fftw_getmem(emFelder[emF,abl],aX*sizeof(double)); + fillchar(emFelder[emF,abl]^,aX*sizeof(double),#0); end; end; for emQ:=low(tEMQuellGroesze) to high(tEMQuellGroesze) do begin - fftw_getmem(emQuellen[emQ],aX*sizeof(extended)); - fillchar(emQuellen[emQ]^,aX*sizeof(extended),#0); + fftw_getmem(emQuellen[emQ],aX*sizeof(double)); + fillchar(emQuellen[emQ]^,aX*sizeof(double),#0); end; fillchar(matFelder,sizeof(matFelder),#0); setlength(matFelder,length(teilchen)); for i:=0 to length(matFelder)-1 do for maF:=low(tMaterieFeldGroesze) to high(tMaterieFeldGroesze) do begin - fftw_getmem(matFelder[i,maF],aX*sizeof(extended)); - fillchar(matFelder[i,maF]^,aX*sizeof(extended),#0); + fftw_getmem(matFelder[i,maF],aX*sizeof(double)); + fillchar(matFelder[i,maF]^,aX*sizeof(double),#0); end; fillchar(impulsraum,sizeof(impulsraum),#0); setlength(impulsraum,length(teilchen)); for i:=0 to length(impulsraum)-1 do for abl:=false to true do begin - fftw_getmem(impulsraum[i,abl],aX*aP*sizeof(extended)); - fillchar(impulsraum[i,abl]^,aX*aP*sizeof(extended),#0); + fftw_getmem(impulsraum[i,abl],aX*aP*sizeof(double)); + fillchar(impulsraum[i,abl]^,aX*aP*sizeof(double),#0); end; fillchar(massen,sizeof(massen),#0); @@ -639,21 +639,22 @@ begin lichters[true].grep('^rechts .*'); lichters[true].subst('^rechts *',''); - fftw_getmem(fftTmp,(aP+1)*(aX+1)*sizeof(complex_extended)); + fftw_getmem(fftTmp,(aP+1)*(aX+1)*sizeof(complex_double)); + fillchar(fftTmp^,(aP+1)*(aX+1)*sizeof(complex_double),#0); fillchar(ffts,sizeof(ffts),#0); - setlength(ffts,length(teilchen)); write('.'); + setlength(ffts,length(teilchen)); for i:=0 to length(ffts)-1 do begin -// ffts[i,0,false]:= // Planung der Hintransformationen über x -// fftw_plan_many_dft_r2c(1,@_aX,aP,impulsraum[i,false],nil,aP,1,fftTmp,nil,aP,1,[fftw_preserve_input,fftw_measure]); write('1'); + ffts[i,0,false]:= // Planung der Hintransformationen über x + fftw_plan_many_dft_r2c(1,@_aX,aP,impulsraum[i,false],nil,aP,1,fftTmp,nil,aP,1,[fftw_preserve_input,fftw_measure]); ffts[i,0,true]:= // Planung der Rücktransformationen über x - fftw_plan_many_dft_c2r(1,@_aX,aP,fftTmp,nil,aP,1,impulsRaumGradient[i,0],nil,aP,1,[fftw_measure]); write('2'); + fftw_plan_many_dft_c2r(1,@_aX,aP,fftTmp,nil,aP,1,impulsRaumGradient[i,0],nil,aP,1,[fftw_measure]); -// ffts[i,1,false]:= // Planung der Hintransformationen über p -// fftw_plan_many_dft_r2c(1,@_aP,aX,impulsraum[i,false],nil,1,aP,fftTmp,nil,1,aP,[fftw_preserve_input,fftw_measure]); write('3'); -// ffts[i,1,true]:= // Planung der Rücktransformationen über p -// fftw_plan_many_dft_c2r(1,@_aP,aX,fftTmp,nil,1,aP,impulsRaumGradient[i,1],nil,1,aP,[fftw_measure]); write('4'); - end; write('>'); + ffts[i,1,false]:= // Planung der Hintransformationen über p + fftw_plan_many_dft_r2c(1,@_aP,aX,impulsraum[i,false],nil,1,aP,fftTmp,nil,1,aP,[fftw_preserve_input,fftw_measure]); + ffts[i,1,true]:= // Planung der Rücktransformationen über p + fftw_plan_many_dft_c2r(1,@_aP,aX,fftTmp,nil,1,aP,impulsRaumGradient[i,1],nil,1,aP,[fftw_measure]); + end; (* fftw_plan_many_dft_r2c(int rank, const int *n, int howmany, double *in, const int *inembed, @@ -713,10 +714,10 @@ begin inherited destroy; end; -procedure tFelder.initialisiereDichte(ort,tlc: longint; breite,n: extended); +procedure tFelder.initialisiereDichte(ort,tlc: longint; breite,n: double); var i: longint; - ges: extended; + ges: double; begin ges:=0; for i:=0 to aP-1 do begin @@ -734,7 +735,7 @@ var emQ: tEMQuellGroesze; emF: tEMFeldGroesze; rechts: boolean; - tmp,err:extended; + tmp: double; begin // PX ggf. korrigieren @@ -824,35 +825,32 @@ begin // dA/dt = dA/dt - E for i:=1 to aX-2 do (emFelder[efAX,true]+i)^:=(emFelder[efDAXDT,false]+i)^ - (emFelder[efEX,false]+i)^; - move(emFelder[efDAYDT,false]^,emFelder[efAY,true]^,aX*sizeof(extended)); + move(emFelder[efDAYDT,false]^,emFelder[efAY,true]^,aX*sizeof(double)); // Gradienten der Phasenraumbesetzungsdichte berechnen for i:=0 to length(ffts)-1 do begin -// fftw_execute(ffts[i,0,false]); - // *i*k*2*pi; -// fftw_execute(ffts[i,0,true]); - write('A'); + fftw_execute(ffts[i,0,false]); + for j:=0 to aX-1 do + for k:=0 to aP-1 do begin // *i*k*2*pi; + tmp:=(fftTmp+j+k*aX)^.re; + (fftTmp+j+k*aX)^.re:=- 2*pi/dX/aX * j * (fftTmp+j+k*aX)^.im/aX; + (fftTmp+j+k*aX)^.im:= 2*pi/dX/aX * j * tmp/aX; + end; + + fftw_execute(ffts[i,0,true]); fftw_execute(ffts[i,1,false]); - write('B'); - // *i*k*2*pi; + + for j:=0 to aX-1 do + for k:=0 to aP-1 do begin // *i*k*2*pi; + tmp:=(fftTmp+j+k*aX)^.re; + (fftTmp+j+k*aX)^.re:=- 2*pi/dP/aP * k * (fftTmp+j+k*aX)^.im/aP; + (fftTmp+j+k*aX)^.im:= 2*pi/dP/aP * k * tmp/aP; + end; + fftw_execute(ffts[i,1,true]); - write('C'); end; - tmp:=0; - err:=0; - for i:=0 to length(impulsraum)-1 do - for j:=0 to aP*aX-1 do begin - tmp:=tmp+abs((impulsraum[i,false]+j)^)+abs((impulsRaumGradient[i,1]+j)^); - err:=err+abs((impulsraum[i,false]+j)^-(impulsRaumGradient[i,1]+j)^); -// tmp:=tmp+abs((impulsraum[i,false]+j)^)+abs((impulsRaumGradient[i,0]+j)^)+abs((impulsRaumGradient[i,1]+j)^); -// err:=err+abs((impulsraum[i,false]+j)^-(impulsRaumGradient[i,0]+j)^) -// +abs((impulsraum[i,false]+j)^-(impulsRaumGradient[i,1]+j)^); - end; - writeln(tmp,' ',err); - halt; - // Zeitableitung der Phasenraumbesetzungsdichte berechnen for i:=0 to length(teilchen)-1 do @@ -902,7 +900,7 @@ end; procedure tFelder.dumpErhaltungsgroeszen(prot: tProtokollant); var i,j: longint; - dens: extended; + dens: double; begin dens:=0; for i:=0 to length(massen)-1 do @@ -917,7 +915,7 @@ begin end; end; -function tFelder.impulsIntegral(ort,tlc: longint; emQ: tEMQuellGroesze): extended; +function tFelder.impulsIntegral(ort,tlc: longint; emQ: tEMQuellGroesze): double; begin if tlc<0 then begin result:=(emQuellen[emQ]+ort)^; // das ist leicht :-) @@ -934,7 +932,7 @@ begin result:=result * teilchen[tlc].eigenschaften[tsgLadung]; end; -function tFelder.impulsIntegral(ort,tlc: longint; maF: tMaterieSpeicherGroesze): extended; +function tFelder.impulsIntegral(ort,tlc: longint; maF: tMaterieSpeicherGroesze): double; var i: longint; begin @@ -999,7 +997,7 @@ end; // tGitter ********************************************************************* -constructor tGitter.create(derBesitzer: tSimulation; aX,aP: longint; deltaX,deltaP: extended; bekannteWerte: tKnownValues; teilchen: array of tTeilchenSpezies; lichter: tMyStringlist; zv: tZeitverfahren; protokollant: tProtokollant; name: string); +constructor tGitter.create(derBesitzer: tSimulation; aX,aP: longint; deltaX,deltaP: double; bekannteWerte: tKnownValues; teilchen: array of tTeilchenSpezies; lichter: tMyStringlist; zv: tZeitverfahren; protokollant: tProtokollant; name: string); var i: longint; begin @@ -1061,7 +1059,7 @@ begin felders[i].setzeNull; end; -procedure tGitter.iteriereSchritt(dT: extended); +procedure tGitter.iteriereSchritt(dT: double); {$IFDEF negativeDichteueberwachung} var i,j: longint; pro: tProtokollant; @@ -1157,7 +1155,7 @@ var ifile: tMyStringlist; zeitverfahren: tZeitverfahren; s,t,aSuffix,aPrefix: string; - dX,breite,dP,pMax: extended; + dX,breite,dP,pMax: double; i: longint; kvs: tKnownValues; teilchen: array of tTeilchenSpezies; @@ -1520,7 +1518,7 @@ begin inherited destroy; end; -function tSimulation.iteriereSchritt(start: extended; var zeitPhysik,zeitDatei: extended): boolean; // noch nicht zu Ende? +function tSimulation.iteriereSchritt(start: double; var zeitPhysik,zeitDatei: double): boolean; // noch nicht zu Ende? var i: longint; begin diff --git a/Plasmapropagation.lpr b/Plasmapropagation.lpr index d3b8f0a..e1346c5 100644 --- a/Plasmapropagation.lpr +++ b/Plasmapropagation.lpr @@ -8,7 +8,7 @@ uses {$ENDIF}{$ENDIF} Classes, SysUtils, CustApp, { you can add units after this } - math, Physikunit, protokollunit, fftw_l; + math, Physikunit, protokollunit, lowlevelunit; type @@ -27,49 +27,10 @@ type procedure TPlasmapropagation.DoRun; var simulation: tSimulation; - start,zeitPhysik,zeitDatei: extended; + start,zeitPhysik,zeitDatei: double; prot: tProtokollant; s,t,u: string; - - xlen,ylen,i,j: longint; - ein: Pcomplex_extended; - aus: Pextended; - plan: fftw_plan_extended; begin - - xlen:=4; - ylen:=4; - - fftw_getmem(ein,xlen*ylen*sizeof(complex_extended)); - fftw_getmem(aus,xlen*ylen*sizeof(extended)); - - plan:=fftw_plan_many_dft_c2r(1,@xlen,ylen,ein,nil,ylen,1,aus,nil,ylen,1,[fftw_measure]); - - for i:=0 to xlen-1 do - for j:=0 to ylen-1 do begin - (ein+j+i*ylen)^.re:=j; - (ein+j+i*ylen)^.im:=i; - end; - for i:=0 to xlen-1 do begin - for j:=0 to ylen-1 do - write((ein+j+i*ylen)^.re,' + ',(ein+j+i*ylen)^.im,' *I ; '); - writeln; - end; - - fftw_execute(plan); - - for i:=0 to xlen-1 do begin - for j:=0 to ylen-1 do - write((aus+j+i*ylen)^,' '); - writeln; - end; - - fftw_destroy_plan(plan); - fftw_freemem(ein); - fftw_freemem(aus); - - halt; - prot:=tProtokollant.create('error'); if paramcount<>1 then begin diff --git a/Plasmapropagation.lps b/Plasmapropagation.lps index eee3207..d9cdf52 100644 --- a/Plasmapropagation.lps +++ b/Plasmapropagation.lps @@ -7,19 +7,19 @@ - - - + + + - - - + + + @@ -28,14 +28,14 @@ - + - + - + @@ -44,18 +44,18 @@ - + - + - + @@ -64,19 +64,19 @@ - + - + - - - + + + @@ -126,146 +126,156 @@ - - - - - + + + + - - + + + - - + + + - + - + - + - + - + - + - + - + - + - + - + - + - + - + - + - + - + - + - - + - + - - + - + - + - - + + - - + + - + - + - - + + + + + + + + + + + + + -- cgit v1.2.3-70-g09d2