CODEVS3123 a*b problem plus (FFT)

 1 type xh=record
 2         x,y:double;
 3         end;
 4       arr=array[0..1000008] of xh;
 5 var n,m:longint;
 6     s1,s2:ansistring;
 7     a,b,g,w:arr;
 8     ch:char;
 9 operator -(a,b:xh) c:xh;
10 begin
11     c.x:=a.x-b.x;
12     c.y:=a.y-b.y;
13 end;
14 operator +(a,b:xh) c:xh;
15 begin
16     c.x:=a.x+b.x;
17     c.y:=a.y+b.y;
18 end;
19 operator *(a,b:xh) c:xh;
20 begin
21     c.x:=a.x*b.x-a.y*b.y;
22     c.y:=a.x*b.y+a.y*b.x;
23 end;
24 procedure dft(var a:arr;s,t:longint);//a待处理数组 s初始位置 1<<t长度
25 var i,p:longint;
26     cnt:xh;
27 begin
28     if n>>t=1 then exit;
29     dft(a,s,t+1); dft(a,s+1<<t,t+1);
30     for i:=0 to n>>t>>1-1 do
31     begin
32         p:=i<<t<<1+s;
33         cnt:=w[i<<t]*a[p+1<<t];
34         g[i]:=a[p]+cnt;
35         g[i+n>>t>>1]:=a[p]-cnt;
36     end;
37     for i:=0 to n>>t-1 do a[s+i<<t]:=g[i];
38 end;
39 procedure clr(var a:arr);
40 begin
41     fillchar(a,sizeof(a),0);
42 end;
43 procedure FFT;
44 var i,len:longint;
45     lx:longint;
46 begin
47     n:=1;
48     while n<m<<1 do n:=n<<1;
49     for i:=0 to n-1 do w[i].x:=cos(pi*2*i/n);
50     for i:=0 to n-1 do w[i].y:=sin(pi*2*i/n);
51     dft(a,0,0); dft(b,0,0);
52     for i:=0 to n-1 do a[i]:=a[i]*b[i];
53     for i:=0 to n-1 do w[i].y:=-w[i].y;  //!!!!
54     dft(a,0,0);
55     for i:=m<<1-1 downto 0 do a[i].x:=a[i].x/n;//!!!!
56     for i:=0 to m<<1-1 do
57     begin
58         lx:=round(a[i].x);
59         a[i+1].x:=a[i+1].x+lx div 10;
60         a[i].x:=lx mod 10;
61     end;
62         len:=m<<1-1;
63         while a[len].x<1e-12 do dec(len);
64     for i:=len downto 0 do write(a[i].x:0:0);
65     writeln;
66 end;
67 procedure main;
68 var i:longint;
69 begin
70     read(ch);
71     while ord(ch)>=48 do begin s1:=s1+ch; read(ch); end;
72     readln(s2);
73     clr(a); clr(b);
74     for i:=length(s1) downto 1 do a[length(s1)-i].x:=ord(s1[i])-48;
75     for i:=length(s2) downto 1 do b[length(s2)-i].x:=ord(s2[i])-48;
76     if length(s1)>length(s2) then m:=length(s1) else m:=length(s2);
77     FFT;
78 end;
79 begin
80     main;
81 end.

 

posted @ 2015-03-11 19:55  rpSebastian  阅读(202)  评论(0编辑  收藏  举报