标签:
好吧,其实我并没有深入运用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.
bzoj2194
模板题,倒过来就变成卷积了
不用operator就跑的快多了
1 type point=record 2 x,y:double; 3 end; 4 arr=array[0..300010] of point; 5 6 var a,b,c:arr; 7 r:array[0..300010] of longint; 8 h,i,n,l,m:longint; 9 10 procedure fft(var a:arr; ty:longint); 11 var i,j,k,s:longint; 12 u,v,p,w:point; 13 begin 14 for i:=0 to n-1 do 15 if i<r[i] then 16 begin 17 w:=a[r[i]]; 18 a[r[i]]:=a[i]; 19 a[i]:=w; 20 end; 21 22 i:=1; 23 while i<n do 24 begin 25 p.x:=cos(pi/i); 26 p.y:=ty*sin(pi/i); 27 s:=i shl 1; 28 j:=0; 29 while j<n do 30 begin 31 w.x:=1; 32 w.y:=0; 33 for k:=0 to i-1 do 34 begin 35 u:=a[j+k]; 36 v.x:=w.x*a[j+k+i].x-w.y*a[j+k+i].y; 37 v.y:=w.x*a[j+k+i].y+w.y*a[j+k+i].x; 38 a[j+k].x:=u.x+v.x; 39 a[j+k].y:=u.y+v.y; 40 a[j+k+i].x:=u.x-v.x; 41 a[j+k+i].y:=u.y-v.y; 42 u:=w; 43 w.x:=u.x*p.x-u.y*p.y; 44 w.y:=u.x*p.y+u.y*p.x; 45 end; 46 inc(j,s); 47 end; 48 i:=i shl 1; 49 end; 50 end; 51 52 53 begin 54 readln(n); 55 h:=n; 56 for i:=0 to n-1 do 57 readln(a[i].x,b[n-1-i].x); 58 m:=2*n-2; 59 n:=1; 60 while n<=m do 61 begin 62 n:=n*2; 63 inc(l); 64 end; 65 for i:=0 to n-1 do 66 r[i]:=(r[i shr 1] shr 1) or ((i and 1) shl (l-1)); 67 fft(a,1); 68 fft(b,1); 69 for i:=0 to n-1 do 70 begin 71 c[i].x:=a[i].x*b[i].x-a[i].y*b[i].y; 72 c[i].y:=a[i].x*b[i].y+a[i].y*b[i].x; 73 end; 74 fft(c,-1); 75 for i:=h-1 to m do 76 writeln(round(c[i].x/n)); 77 78 end.
bzoj3527
给一下题面吧,, 令 Ei=Fi/qi,求Ei
带进去稍微弄弄就是卷积啊……
1 type point=record 2 x,y:extended; 3 end; 4 arr=array[0..300010] of point; 5 6 var a,b,c:arr; 7 ans,d:array[0..300010] of extended; 8 r:array[0..300010] of longint; 9 i,n,m,l:longint; 10 11 procedure fft(var a:arr; ty:longint); 12 var i,j,k,s:longint; 13 w,p,u,v:point; 14 begin 15 for i:=0 to n-1 do 16 if i<r[i] then 17 begin 18 w:=a[r[i]]; 19 a[r[i]]:=a[i]; 20 a[i]:=w; 21 end; 22 23 i:=1; 24 while i<n do 25 begin 26 p.x:=cos(pi/i); 27 p.y:=ty*sin(pi/i); 28 s:=i shl 1; 29 j:=0; 30 while j<n do 31 begin 32 w.x:=1; 33 w.y:=0; 34 for k:=0 to i-1 do 35 begin 36 u:=a[j+k]; 37 v.x:=w.x*a[j+k+i].x-w.y*a[j+k+i].y; 38 v.y:=w.x*a[j+k+i].y+w.y*a[j+k+i].x; 39 a[j+k].x:=u.x+v.x; 40 a[j+k].y:=u.y+v.y; 41 a[j+k+i].x:=u.x-v.x; 42 a[j+k+i].y:=u.y-v.y; 43 u:=w; 44 w.x:=u.x*p.x-u.y*p.y; 45 w.y:=u.x*p.y+u.y*p.x; 46 end; 47 inc(j,s); 48 end; 49 i:=i shl 1; 50 end; 51 if ty=-1 then 52 begin 53 for i:=0 to n-1 do 54 a[i].x:=a[i].x/extended(n); 55 end; 56 end; 57 58 begin 59 readln(n); 60 m:=n; 61 n:=1; 62 while n<=2*(m-1) do 63 begin 64 n:=n shl 1; 65 inc(l); 66 end; 67 for i:=0 to n-1 do 68 r[i]:=(r[i shr 1] shr 1) or ((i and 1) shl (l-1)); 69 for i:=0 to m-1 do 70 readln(d[i]); 71 for i:=0 to m-1 do 72 begin 73 a[i].x:=d[i]; 74 if i<>0 then b[i].x:=1/i/i; 75 end; 76 fft(a,1); fft(b,1); 77 for i:=0 to n-1 do 78 begin 79 c[i].x:=a[i].x*b[i].x-a[i].y*b[i].y; 80 c[i].y:=a[i].x*b[i].y+a[i].y*b[i].x; 81 end; 82 fft(c,-1); 83 for i:=0 to m-1 do 84 ans[i]:=c[i].x; 85 86 for i:=0 to n-1 do 87 begin 88 a[i].x:=0; 89 a[i].y:=0; 90 if i<m then a[i].x:=d[m-1-i]; 91 end; 92 fft(a,1); 93 for i:=0 to n-1 do 94 begin 95 c[i].x:=a[i].x*b[i].x-a[i].y*b[i].y; 96 c[i].y:=a[i].x*b[i].y+a[i].y*b[i].x; 97 end; 98 fft(c,-1); 99 for i:=0 to m-1 do 100 begin 101 ans[i]:=ans[i]-c[m-i-1].x; 102 writeln(ans[i]:0:3); 103 end; 104 end.
bzoj3160
好题,首先合法方案=不连续间隔相等的回文-连续的回文
连续回文的方案显然可以manacher搞出来
连续的呢?考虑穷举中心点,我们只要快算计算中间点两边的相等间隔字符相等的对数s
则此中心点贡献的方案数即为2^s-1
设中心点为k,则s=sigama([a(i)=a(k-i)]) 这里k代表的是manacher翻倍后的下标,i为原串的下标
考虑字符只有a,b两种,我们只要分别统计即可
当统计a时,把a标为1,b标为0,则[a(i)=a(k-i)]=a(i)*a(k-i) 卷积就出现了,fft搞搞即可
注意inline可以加快operator的速度,但要慎用,有时候会出错
1 uses math; 2 const mo=1000000007; 3 type point=record 4 x,y:double; 5 end; 6 arr=array[0..300010] of point; 7 8 var a,b:arr; 9 p,d,r:array[0..300010] of longint; 10 c:array[0..300010] of char; 11 s:ansistring; 12 x,n,l,len,i,ans:longint; 13 14 operator +(a,b:point)c:point;inline; 15 begin 16 c.x:=a.x+b.x; 17 c.y:=a.y+b.y; 18 end; 19 20 operator -(a,b:point)c:point;inline; 21 begin 22 c.x:=a.x-b.x; 23 c.y:=a.y-b.y; 24 end; 25 26 operator *(a,b:point)c:point;inline; 27 begin 28 c.x:=a.x*b.x-a.y*b.y; 29 c.y:=a.x*b.y+a.y*b.x; 30 end; 31 32 function min(a,b:longint):longint; 33 begin 34 if a>b then exit(b) else exit(a); 35 end; 36 37 function manacher:longint; 38 var w,t,i,right,k:longint; 39 begin 40 c[0]:=‘$‘; 41 t:=0; 42 for i:=1 to len do 43 begin 44 inc(t); c[t]:=‘#‘; 45 inc(t); c[t]:=s[i]; 46 end; 47 c[t+1]:=‘#‘; 48 right:=0; k:=0; 49 w:=0; 50 for i:=1 to t do 51 begin 52 if right>i then p[i]:=min(p[2*k-i],right-i) 53 else p[i]:=1; 54 while c[i+p[i]]=c[i-p[i]] do inc(p[i]); 55 w:=(w+p[i] div 2) mod mo; 56 if p[i]+i>right then 57 begin 58 right:=p[i]+i; 59 k:=i; 60 end; 61 end; 62 exit(w); 63 end; 64 65 procedure fft(var a:arr; ty:longint); 66 var i,j,k,s:longint; 67 w,p,u,v:point; 68 begin 69 for i:=0 to n-1 do 70 if i<r[i] then 71 begin 72 w:=a[r[i]]; 73 a[r[i]]:=a[i]; 74 a[i]:=w; 75 end; 76 i:=1; 77 while i<n do 78 begin 79 p.x:=cos(pi/i); 80 p.y:=ty*sin(pi/i); 81 s:=i shl 1; 82 j:=0; 83 while j<n do 84 begin 85 w.x:=1; w.y:=0; 86 for k:=0 to i-1 do 87 begin 88 u:=a[j+k]; 89 v:=w*a[j+k+i]; 90 a[j+k]:=u+v; 91 a[j+k+i]:=u-v; 92 w:=w*p; 93 end; 94 inc(j,s); 95 end; 96 i:=i shl 1; 97 end; 98 99 end; 100 begin 101 readln(s); 102 len:=length(s); 103 d[0]:=1; 104 for i:=1 to len do 105 d[i]:=d[i-1]*2 mod mo; 106 n:=1; 107 l:=0; 108 while n<=(len-1) shl 1 do 109 begin 110 n:=n*2; 111 inc(l); 112 end; 113 for i:=0 to n-1 do 114 r[i]:=(r[i shr 1] shr 1) or ((i and 1) shl (l-1)); 115 for i:=1 to len do 116 if s[i]=‘a‘ then a[i-1].x:=1; 117 fft(a,1); 118 for i:=0 to n-1 do 119 begin 120 b[i]:=a[i]*a[i]; 121 a[i].x:=0; a[i].y:=0; 122 end; 123 for i:=1 to len do 124 if s[i]=‘b‘ then a[i-1].x:=1; 125 fft(a,1); 126 for i:=0 to n-1 do 127 b[i]:=b[i]+a[i]*a[i]; 128 fft(b,-1); 129 for i:=0 to n-1 do 130 begin 131 x:=trunc(round(b[i].x/n)); 132 ans:=(ans+d[(x+1) shr 1]-1) mod mo; 133 end; 134 writeln((ans-manacher+mo) mod mo); 135 end.
标签:
原文地址:http://www.cnblogs.com/phile/p/4665374.html