From 1edf88972a1eb55c50b33bc37e707077f93dae5e Mon Sep 17 00:00:00 2001 From: Erich Eckner Date: Thu, 12 Nov 2015 11:03:27 +0100 Subject: Berechnung der Ableitung durch FFT korrigiert --- Physikunit.pas | 90 +++++++++++++++++++++++++++++++++++++++++++++------------- 1 file changed, 71 insertions(+), 19 deletions(-) diff --git a/Physikunit.pas b/Physikunit.pas index 5beac01..bf76b6f 100644 --- a/Physikunit.pas +++ b/Physikunit.pas @@ -659,15 +659,67 @@ begin fillchar(ffts,sizeof(ffts),#0); // Vorsicht, die Planung der FFTs überschreibt den Input, auch wenn die eigentlichen FFTs das nicht tun! 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_exhaustive]); - 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_destroy_input,fftw_exhaustive]); - - 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 div 2 + 1,[fftw_preserve_input,fftw_exhaustive]); - ffts[i,1,true]:= // Planung der Rücktransformationen über p - fftw_plan_many_dft_c2r(1,@_aP,aX,fftTmp,nil,1,aP div 2 + 1,impulsRaumGradient[i,1],nil,1,aP,[fftw_destroy_input,fftw_exhaustive]); + ffts[i,0,false]:= // Planung der Hintransformationen über x: + fftw_plan_many_dft_r2c( + 1, // jeweils 1-d + @_aX, // jeweils _aX Werte + aP, // aP Mal + impulsraum[i,false], // von dort + nil, // kein Padding + aP, // Abstand benachbarter Werte + 1, // Verschiebung der Werte zwischen den Transformationen + fftTmp, // nach dort + nil, // kein Padding + aP, // Abstand benachbarter Werte + 1, // Verschiebung der Werte zwischen den Transformationen + [fftw_preserve_input,fftw_exhaustive] // input nicht überschreiben, erschöpfende Suche nach FFT + ); + ffts[i,0,true]:= // Planung der Rücktransformationen über x: + fftw_plan_many_dft_c2r( + 1, // jeweils 1-d + @_aX, // jeweils _aX Werte + aP, // aP Mal + fftTmp, // von hier + nil, // kein Padding + aP, // Abstand benachbarter Werte + 1, // Verschiebung der Werte zwischen den Transformationen + impulsRaumGradient[i,0], // nach dort + nil, // kein Padding + aP, // Abstand benachbarter Werte + 1, // Verschiebung der Werte zwischen den Transformationen + [fftw_destroy_input,fftw_exhaustive] // input überschreiben, erschöpfende Suche nach FFT + ); + + ffts[i,1,false]:= // Planung der Hintransformationen über p: + fftw_plan_many_dft_r2c( + 1, // jeweils 1-d + @_aP, // jeweils _aP Werte + aX, // aX Mal + impulsraum[i,false], // von dort + nil, // kein Padding + 1, // Abstand benachbarter Werte + aP, // Verschiebung der Werte zwischen den Transformationen + fftTmp, // nach dort + nil, // kein Padding + 1, // Abstand benachbarter Werte + aP div 2 + 1, // Verschiebung der Werte zwischen den Transformationen + [fftw_preserve_input,fftw_exhaustive] // input nicht überschreiben, erschöpfende Suche nach FFT + ); + ffts[i,1,true]:= // Planung der Rücktransformationen über p: + fftw_plan_many_dft_c2r( + 1, // jeweils 1-d + @_aP, // jeweils _aP Werte + aX, // aX Mal + fftTmp, // von dort + nil, // kein Padding + 1, // Abstand benachbarter Werte + aP div 2 + 1, // Verschiebung der Werte zwischen den Transformationen + impulsRaumGradient[i,1], // nach dort + nil, // kein Padding + 1, // Abstand benachbarter Werte + aP, // Verschiebung der Werte zwischen den Transformationen + [fftw_destroy_input,fftw_exhaustive] // input überschreiben, erschöpfende Suche nach FFT + ); end; {$IFDEF exzessiveArrayBereichsTests} @@ -1161,11 +1213,11 @@ begin pruefeArrayEnden('Vorwärts-FFT ('+inttostr(i)+',0) schreibt beim Ausführen über das Ende hinaus!'); {$ENDIF} - for k:=0 to aP-1 do - for j:=0 to aX div 2 do begin // *-i*k*2*pi; - tmp:=(fftTmp+j+k*(aX div 2 + 1))^.re; - (fftTmp+j+k*(aX div 2 + 1))^.re:= 2*pi/dX/aX * j * (fftTmp+j+k*(aX div 2 + 1))^.im/aX * byte(3*j