fft huff

Prabhat1308 · Oct 2, 2023 · c007eb7 · c007eb7
1 parent 644c03e
commit c007eb7
Show file tree

Hide file tree

Showing 6 changed files with 350 additions and 1 deletion.
diff --git a/src/.DS_Store b/src/.DS_Store
diff --git a/src/FFT_implementation/.DS_Store b/src/FFT_implementation/.DS_Store
diff --git a/src/FFT_implementation/CONSTANT_FFT.huff b/src/FFT_implementation/CONSTANT_FFT.huff
@@ -0,0 +1,2 @@
+#define constant PI = 0x2B992DDFA23249D6
+#define constant N = 0x04
diff --git a/src/FFT_implementation/FFT.huff b/src/FFT_implementation/FFT.huff
@@ -0,0 +1,339 @@
+#include './CONSTANTS.huff'
+#include '../ComplexHUff/Trigonometry.huff'
+#include '../ComplexHUff/Complex.huff'  
+#include '../ComplexHUff/PRBMath.huff'
+#include '../ComplexHUff/Constants.huff'
+#include '../ComplexHuff/Helper.huff'
+
+#define macro PUT_THETA_T() = takes(0) returns(1) {
+    [N]
+    [PI]
+    div 
+}
+
+#define macro INIT() = takes(8) returns(0){
+    [X3] 0x00 mstore
+    0x00 0x20 mstore
+    PUT_NEG1() [X3] mul 0x40 mstore
+    0x00 0x60 mstore
+    0x00 0x80 mstore
+    [X3] 0xa0 mstore
+    0x00 0xc0 mstore
+    PUT_NEG1() [X3] mul 0xe0 mstore
+}
+
+#define macro FFT() = takes(0) returns(0) {
+   loop1:
+     [N] // [k]
+     INNER_CALC1()
+     0x01 //[1,k]
+     shr // [k>>1] 
+     dup1 // [k>>1,k>>1]
+     0x01 // [1,k>>1,k>>1]
+     lt  // [(k>>1) >1 , k>>1]
+     loop1
+     jumpi
+     pop
+}
+
+#define macro INNER_CALC1() = takes(0) returns(0) {
+    //INPUT stack => [k]
+    dup1 // [n=k,k]
+    PUT_THETA_T() // [Theta_T,n,k]
+    COS() // [cos(Theta_T),n,k]
+    dup1 // [cos(Theta_T),cos(Theta_T),n,k]
+    PUT_THETA_T() // [Theta_T,cos(Theta_T),cos(Theta_T),n,k]
+    SIN() // [sin(Theta_T),cos(Theta_T),cos(Theta_T),n,k] , [Phi_t_im,Phi_t_re,Phi_t_re_temp,..]
+    dup1 // [sin(Theta_T),sin(Theta_T),cos(Theta_T),cos(Theta_T),n,k]
+    dup1 // [sin(Theta_T),sin(Theta_T),sin(Theta_T),cos(Theta_T),cos(Theta_T),n,k]
+    mul // [sin(Theta_T)*sin(Theta_T),sin(Theta_T),cos(Theta_T),cos(Theta_T),n,k]
+    swap2 // [cos(Theta_T),sin(Theta_T)*sin(Theta_T),sin(Theta_T),cos(Theta_T),n,k]
+    dup1 // [cos(Theta_T),cos(Theta_T),sin(Theta_T)*sin(Theta_T),sin(Theta_T),cos(Theta_T),n,k]
+    mul // [cos(Theta_T)*cos(Theta_T),sin(Theta_T)*sin(Theta_T),sin(Theta_T),cos(Theta_T),n,k]
+    [X3] // [1e18,cos(Theta_T)*cos(Theta_T),sin(Theta_T)*sin(Theta_T),sin(Theta_T),cos(Theta_T),n,k]
+    swap1 // [cos(Theta_T)*cos(Theta_T),1e18,sin(Theta_T)*sin(Theta_T),sin(Theta_T),cos(Theta_T),n,k]
+    sdiv // [cos(Theta_T)*cos(Theta_T)/1e18,sin(Theta_T)*sin(Theta_T),sin(Theta_T),cos(Theta_T),n,k]
+    swap1 // [sin(Theta_T)*sin(Theta_T),cos(Theta_T)*cos(Theta_T)/1e18,sin(Theta_T),cos(Theta_T),n,k]
+    [X3] // [1e18,sin(Theta_T)*sin(Theta_T),cos(Theta_T)*cos(Theta_T)/1e18,sin(Theta_T),cos(Theta_T),n,k]
+    swap1 // [sin(Theta_T)*sin(Theta_T),1e18,cos(Theta_T)*cos(Theta_T)/1e18,sin(Theta_T),cos(Theta_T),n,k]
+    sdiv // [sin(Theta_T)*sin(Theta_T)/1e18,cos(Theta_T)*cos(Theta_T)/1e18,sin(Theta_T),cos(Theta_T),n,k]
+    swap1 // [cos(Theta_T)*cos(Theta_T)/1e18,sin(Theta_T)*sin(Theta_T)/1e18,sin(Theta_T),cos(Theta_T),n,k]
+    sub // [cos(Theta_T)*cos(Theta_T)/1e18-sin(Theta_T)*sin(Theta_T)/1e18,sin(Theta_T),cos(Theta_T),n,k] , Let cos(Theta_T)*cos(Theta_T)/1e18-sin(Theta_T)*sin(Theta_T)/1e18 = Phi_t_re
+    swap2 // [cos(Theta_T),sin(Theta_T),Phi_t_re',n,k]
+    mul // [cos(Theta_T)*sin(Theta_T),Phi_t_re',n,k]
+    0x02 // [2,cos(Theta_T)*sin(Theta_T),Phi_t_re',n,k]
+    mul // [2*cos(Theta_T)*sin(Theta_T),Phi_t_re,n,k] , Let 2*cos(Theta_T)*sin(Theta_T) = Phi_t_im
+    [X3] // [T_real,Phi_t_im,Phi_t_re,n,k]
+    0x00 // [T_im,T_real,Phi_t_im,Phi_t_re,n,k]
+    INNER_CALC2() // returns with []
+    [N] // [N]
+    LOG2() // [LOG2(N)]
+    INNER_CALC3()
+}
+
+#define macro INNER_CALC2() = takes(0) returns(0) {
+ // Input_Stack => [T_im,T_real,Phi_t_im,Phi_t_re,n,k]
+ 0x00 // [0,T_im,T_real,Phi_t_im,Phi_t_re,n,k]
+ 0x2000
+ mstore
+ LOOP:
+ INNER_CALC4() // [T_im,T_real,Phi_t_im,Phi_t_re,n,k]
+ dup2 // [T_real,T_im,T_real_temp,Phi_t_im,Phi_t_re,n,k]
+ dup5  
+ mul // [Phi_te_re*T_real,T_im,T_real_temp,Phi_t_im,Phi_t_re,n,k]
+ [X3]
+ swap1
+ sdiv // [Phi_te_re*T_real/1e18,T_im,T_real_temp,Phi_t_im,Phi_t_re,n,k]
+ dup2 // [T_im,Phi_te_re*T_real/1e18,T_im,T_real_temp,Phi_t_im,Phi_t_re,n,k]
+ dup5 // [Phi_t_im,T_im,Phi_te_re*T_real/1e18,T_im,T_real_temp,Phi_t_im,Phi_t_re,n,k]
+ mul
+ [X3]
+ swap1
+ sdiv // [Phi_t_im*T_im/1e18,Phi_te_re*T_real/1e18,T_im,T_real_temp,Phi_t_im,Phi_t_re,n,k]
+ swap1
+ sub // [Phi_te_re*T_real/1e18-Phi_t_im*T_im/1e18,T_im,T_real_temp,Phi_t_im,Phi_t_re,n,k] , Let Phi_te_re*T_real/1e18-Phi_t_im*T_im/1e18 = T_real
+ // [T_real,T_im,T_real_temp,Phi_t_im,Phi_t_re,n,k]
+ dup5 // [Phi_t_re,T_im,T_real_temp,Phi_t_im,T_real,n,k]
+ mul
+ [X3]
+ swap1
+ sdiv // [Phi_t_re*T_im/1e18,T_real_temp,Phi_t_im,Phi_t_re,T_real,n,k]
+ swap1 
+ dup3 // [Phi_t_im,T_real_temp,Phi_t_re*T_im/1e18,Phi_t_im,Phi_t_re,T_real,n,k]
+ mul
+ [X3]
+ swap1
+ sdiv
+ sub // [Phi_t_im*T_real_temp/1e18-Phi_t_re*T_im/1e18,Phi_t_im,Phi_t_re,T_real,n,k] , Let Phi_t_im*T_real_temp/1e18-Phi_t_re*T_im/1e18 = T_im
+ 0x2000
+ mload
+ 0x01
+ add 
+ dup1 // [l+1,l+1,Phi_t_im,Phi_t_re,T_real,n,k]
+ 0x2000
+ mstore // [l+1,Phi_t_im,Phi_t_re,T_real,n,k]
+ dup6
+ gt
+ jumpi
+ LOOP
+ pop pop pop pop pop pop
+}
+
+#define macro INNER_CALC3() = takes(1) returns(0) {
+    // Input Stack = [m] 
+    0x00 // [0,m]
+    Looper_calc3: // [a,m]
+    dup1 // [b=a,a,m]
+    dup1 // [b,b,a,m]
+    0xaaaaaaaa
+    and // [b&0xaaaaaaaa,b,a,m]
+    0x01
+    shr // [(b&0xaaaaaaaa)>>1,b,a,m]
+    swap1 // [b,(b&0xaaaaaaaa)>>1,a,m]
+    0x55555555 // 
+    and // [b&0x55555555,(b&0xaaaaaaaa)>>1,a,m]
+    0x01
+    shl // [(b&0x55555555)<<1,(b&0xaaaaaaaa)>>1,a,m]
+    or // [((b&0x55555555)<<1)|((b&0xaaaaaaaa)>>1),a,m]  , Let ((b&0x55555555)<<1)|((b&0xaaaaaaaa)>>1) = new_b = b
+    dup1
+    0xcccccccc
+    and
+    0x02
+    shr
+    swap1
+    0x33333333
+    and
+    0x02
+    shl 
+    or
+    dup1
+    0xf0f0f0f0
+    and
+    0x04
+    shr
+    swap1
+    0x0f0f0f0f
+    and
+    0x04
+    shl 
+    or
+    dup1
+    0xff00ff00
+    and
+    0x08
+    shr
+    swap1
+    0x00ff00ff
+    and
+    0x08
+    shl 
+    or
+    // [b,a,m]
+    dup1 // [b,b,a,m]
+    16
+    shr // [b>>16,b,a,m]
+    swap1 // [b,b>>16,a,m]
+    16
+    shl // [b<<16,b>>16,a,m]
+    or // [(b<<16)|(b>>16),a,m]
+    dup3 // [m,(b<<16)|(b>>16),a,m]
+    32 
+    sub // [(32-m),(b<<16)|(b>>16),a,m]
+    shr // [(b<<16)|(b>>16)>>(32-m),a,m]
+    COND_SOLVE() // [a,m]
+    0x01
+    add  // [a+1,m]
+    dup1
+    [N]
+    gt // [N>a+1,a+1,m]
+    Looper_calc3
+    jumpi
+}
+
+#define macro COND_SOLVE() = takes(0) returns(0){
+  // Input stack => [b,a,m]
+  dup1 // [b,b,a,m]
+  dup2 // [a,b,b,a,m]
+  lt 
+  OK
+  jumpi
+  pop
+
+  OK: // [b,a,m]
+  dup2 // 
+  dup1  // [a,a,b,a,m]
+  0x20
+  mul 
+  mload // [Re(a),a,b,a,m]
+  swap1 // [a,Re(a),b,a,m]
+  0x20
+  mul
+  0x20
+  add
+  mload // [Im(a),Re(a),b,a,m]
+  dup3
+  dup1 // [b,b,Im(a),Re(a),b,a,m]
+  0x20 
+  mul 
+  mload // [Re(b),b,Im(a),Re(a),b,a,m]
+  swap1 // [b,Re(b),Im(a),Re(a),b,a,m]
+  0x20 
+  mul
+  0x20
+  add 
+  mload // [Im(b),Re(b),Im(a),Re(b),b,a,m]
+  dup6 // [a,Im(b),Re(b),Im(a),Re(b),b,a,m]
+  0x20
+  mul
+  0x20
+  add 
+  mstore // [Re(b),Im(a),Re(b),b,a,m]
+  dup5 // [a,Re(b),Im(a),Re(b),b,a,m]
+  0x20 
+  mul
+  mstore // [Im(a),Re(b),b,a,m]
+  dup3 // [b,Im(a),Re(b),b,a,m] 
+  0x20
+  mul
+  0x20
+  add 
+  mstore // [Re(b),b,a,m]
+  dup2 
+  0x20
+  mul
+  mstore // [b,a,m]
+  pop
+
+}
+
+#define macro INNER_CALC4() = takes(0) returns(0) {
+    0x2000 // l memory location
+    mload
+    0x2200 // a memory location
+    mstore
+
+   // INPUT STACK => [T_im,T_real,Phi_t_im,Phi_t_re,n,k]
+   LOOP_last:
+    0x2200 
+    mload
+    dup7
+    add // [a+k=b,T_im,T_real,Phi_t_im,Phi_t_re,n,k]
+    dup1 // [b,b,T_im,T_real,Phi_t_im,Phi_t_re,n,k]
+    0x20
+    mul
+    mload // [Re(b),b,T_im,T_real,Phi_t_im,Phi_t_re,n,k]
+    swap1
+    0x20
+    mul
+    0x20
+    add
+    mload // [Im(b),Re(b),T_im,T_real,Phi_t_im,Phi_t_re,n,k]
+    0x2200
+    mload 
+    dup1 // [a,a,Im(b),Re(b),T_im,T_real,Phi_t_im,Phi_t_re,n,k]
+    0x20
+    mul
+    mload // [Re(a),a,Im(b),Re(b),T_im,T_real,Phi_t_im,Phi_t_re,n,k] 
+    swap1
+    0x20
+    mul
+    0x20
+    add 
+    mload // [Im(a),Re(a),Im(b),Re(b),T_im,T_real,Phi_t_im,Phi_t_re,n,k]
+    swap3 // [Re(b),Re(a),Im(b),Im(a),T_im,T_real,Phi_t_im,Phi_t_re,n,k]
+    dup4
+    dup4
+    dup4
+    dup4
+    SUB_Z() // [t_real,t_im,Re(b),Re(a),Im(b),Im(a),T_im,T_real,Phi_t_im,Phi_t_re,n,k]
+    0x2400
+    mstore
+    0x2600
+    mstore
+    ADD_Z() // [t_real',t_im',T_im,T_real,Phi_t_im,Phi_t_re,n,k]
+    0x2000
+    mload
+    0x20
+    mul
+    mstore
+    0x2000
+    mload
+    0x20
+    mul
+    0x20
+    add 
+    mload // [T_im,T_real,Phi_t_im,Phi_t_re,n,k]
+    0x2400
+    mload
+    0x2600
+    mstore // [t_im,t_real,T_im,T_real,Phi_t_im,Phi_t_re,n,k]
+    MUL_Z() // [re,im,T_im,T_real,Phi_t_im,Phi_t_re,n,k]
+    0x2200
+    mload
+    dup9
+    add 
+    0x20
+    mul 
+    mload
+    0x2200
+    mload
+    dup9
+    add
+    0x20
+    mul
+    0x20
+    add 
+    mload // [T_im,T_real,Phi_t_im,Phi_t_re,n,k]
+    0x2200 
+    mload
+    0x01 
+    add 
+    dup1
+    0x2200
+    mstore
+    [N]
+    gt
+    jumpi
+    LOOP_last
+
+}
diff --git a/.DS_Store → test/.DS_Store b/.DS_Store → test/.DS_Store
diff --git a/test/Complex.t.sol b/test/Complex.t.sol
@@ -66,7 +66,15 @@ contract ComplexTest is Test {
     }
 
     function testToPolar() public {
-        (int r, int t) = complex.toPolar(3, 4);
+        (int r, int t) = complex.toPolar(3 * scale, 4 * scale);
+        assertEq(r / scale, 5); // r = 5
+        assertEq((t * 100) / scale, 92); // T = arctan(4/3) == 0.92 rad
+    }
+
+    function testFromPolar() public {
+        (int r, int i) = complex.fromPolar(5 * scale, 92729522 * 1e10);
+        assertApproxEqAbs(r, 3 * scale, 1e15);
+        assertApproxEqAbs(i, 4 * scale, 1e15);
     }
 }