Berechnung der Ableitung durch FFT korrigiert

1 files 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<aX);
-        (fftTmp+j+k*(aX div 2 + 1))^.im:=- 2*pi/dX/aX * j * tmp/aX                             * byte(3*j<aX);
+    for k:=0 to aX div 2 do          // Iteration über die Werte einer Transformation
+      for j:=0 to aP-1 do begin      // Iteration über die Transformationen
+        tmp:=(fftTmp+j+k*aP)^.re;    //    *-i*k*2*pi;
+        (fftTmp+j+k*aP)^.re:=  2*pi/dX/aX * k * (fftTmp+j+k*aP)^.im/aX * byte(3*k<aX);
+        (fftTmp+j+k*aP)^.im:=- 2*pi/dX/aX * k * tmp/aX                 * byte(3*k<aX);
       end;
 
     {$IFDEF exzessiveArrayBereichsTests}
@@ -1198,11 +1250,11 @@ begin
     pruefeArrayEnden('Vorwärts-FFT ('+inttostr(i)+',1) schreibt beim Ausführen über das Ende hinaus!');
     {$ENDIF}
 
-    for k:=0 to aP div 2 do
-      for j:=0 to aX-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 * byte(4*k < aP);
-        (fftTmp+j+k*aX)^.im:=- 2*pi/dP/aP * k * tmp/aP                 * byte(4*k < aP);
+    for k:=0 to aX-1 do                         // Iteration über die Transformationen
+      for j:=0 to aP div 2 do begin             // Iteration über die Werte einer Transformation
+        tmp:=(fftTmp+j+k*(aP div 2 + 1))^.re;   //    *-i*k*2*pi;
+        (fftTmp+j+k*(aP div 2 + 1))^.re:=  2*pi/dP/aP * j * (fftTmp+j+k*(aP div 2 + 1))^.im/aP * byte(3*j < aP);
+        (fftTmp+j+k*(aP div 2 + 1))^.im:=- 2*pi/dP/aP * j * tmp/aP                             * byte(3*j < aP);
       end;
 
     {$IFDEF exzessiveArrayBereichsTests}