关于fft的一点总结

好吧，其实我并没有深入运用fft，只会优化卷积

听说fft经常和生成函数结合在一起………………oi真是迅猛发展，我真是与时代脱节了……

关于fft的学习推荐直接去看算法导论，写得非常清楚

主要弄懂n次单位根的相关性质定理（消去定理，折半定理）即可，当然也可以直接背代码……

bzoj2179

模板题，fft可以优化高精度乘法

顺便说一句，pascal可以定义operator，但跑得慢

这题我跑了10s……

  1 uses math;
  2 type point=record
  3        x,y:double;
  4      end;
  5      arr=array[0..150010] of point;
  6 
  7 var a,b,c:arr;
  8     d,r:array[0..150010] of longint;
  9     i,n,m,l:longint;
 10     s:ansistring;
 11 
 12 operator +(a,b:point)c:point;
 13   begin
 14     c.x:=a.x+b.x;
 15     c.y:=a.y+b.y;
 16   end;
 17 
 18 operator -(a,b:point)c:point;
 19   begin
 20     c.x:=a.x-b.x;
 21     c.y:=a.y-b.y;
 22   end;
 23 
 24 operator *(a,b:point)c:point;
 25   begin
 26     c.x:=a.x*b.x-a.y*b.y;
 27     c.y:=a.x*b.y+a.y*b.x;
 28   end;
 29 
 30 procedure fft(var a:arr; ty:longint);
 31   var i,j,s,k:longint;
 32       w,p,u,v:point;
 33   begin
 34     for i:=0 to n-1 do
 35       if i<r[i] then
 36       begin
 37         w:=a[r[i]];
 38         a[r[i]]:=a[i];
 39         a[i]:=w;
 40       end;
 41     i:=1;
 42     while i<n do
 43     begin
 44       p.x:=cos(pi/i);
 45       p.y:=ty*sin(pi/i);
 46       s:=i shl 1;
 47       j:=0;
 48       while j<n do
 49       begin
 50         w.x:=1; w.y:=0;
 51         k:=0;
 52         while k<i do
 53         begin
 54           u:=a[j+k];
 55           v:=w*a[j+k+i];
 56           a[j+k]:=u+v;
 57           a[j+k+i]:=u-v;
 58           w:=w*p;
 59           inc(k);
 60         end;
 61         inc(j,s);
 62       end;
 63       i:=i shl 1;
 64     end;
 65   end;
 66 
 67 begin
 68   readln(n);      
 69   readln(s);
 70   for i:=0 to n-1 do
 71     a[n-1-i].x:=ord(s[i+1])-48;
 72   readln(s);
 73   for i:=0 to n-1 do
 74     b[n-1-i].x:=ord(s[i+1])-48;
 75   dec(n); 
 76   m:=2*n;
 77   n:=1;
 78   l:=0;
 79   while n<=m do
 80   begin
 81     n:=n shl 1;
 82     inc(l);
 83   end;
 84   for i:=0 to n-1 do
 85     r[i]:=(r[i shr 1] shr 1) or ((i and 1) shl (l-1));
 86   fft(a,1); fft(b,1);
 87   for i:=0 to n-1 do
 88     c[i]:=a[i]*b[i];
 89   fft(c,-1);
 90   i:=0;
 91   while i<=m do
 92   begin
 93     d[i]:=d[i]+trunc(round(c[i].x/n));
 94     if d[i]>=10 then
 95     begin
 96       d[i+1]:=d[i+1]+d[i] div 10;
 97       d[i]:=d[i] mod 10;
 98     end;
 99     if d[m+1]>0 then inc(m);
100     inc(i);
101   end;
102   for i:=m downto 0 do
103     write(d[i]);
104   writeln;
105 end.

View Code