{This program computes the product of two binary numbers
in multiple precision by two methods.
The first method computes the FFT of each sequence, multiply
the FFTs component-wise, then computes the inverse FFT,
resulting in the convolution.
The second computes the convolution by the naive method.
In either method, the convolution is converted to the binary product 
by scanning from the lower positions in O(n) time.
}
program ex(input,output);
const pai = 3.141592;
type complex =
  record
    re : real;
    im : real
  end;
 t = array[0..100000] of complex;
var i,j,m,n,c1,c2:integer;
    x,y,z,u,v,s,w:t;
    a, b,c: array[0..100000] of integer;
procedure add(var z, x, y : complex);
begin z.re:=x.re+y.re; z.im:=x.im+y.im end;
procedure sub(var z, x, y : complex);
begin z.re:=x.re-y.re; z.im:=x.im-y.im end;
procedure mul(var z, x, y : complex);
begin z.re:=x.re*y.re-x.im*y.im; z.im:=x.re*y.im+x.im*y.re end;

procedure fft(var u:t; x:t; n:integer; sign:integer);
{ sign = 1 : FFT, sign = -1 : inverse FFT }
var
    i,l,n1,k1,k,d,e,i0,i1,i2,jL,jR,kL,kR,kk,c1,c2 : integer;
    z0,z1 : complex;
    y : t;
begin
{c1:=clock;}
  for i:=0 to n-1 do w[i].re:=cos(i*2*pai/n);
  for i:=0 to n-1 do w[i].im:=sign*sin(i*2*pai/n);
  k:=1; m:=0;
  while k<n do begin k:=2*k; m:=m+1 end; { m = log n }
  n1:=n div 2; k1:=n1; d:=1;
  for l:=1 to m do
  begin
    jL:=-1; e:=n div (2*d);
    for i0:=0 to e-1 do begin
      i2:=i0*(2*d);
      for i1:=0 to d-1 do begin
        jL:=jL+1; jR:=jL+n1;
        kL:=i1+i2; kR:=kL+d;
        z0:=x[jL]; mul(z1,x[jR],w[k1*i1]);
        add(y[kL],z0,z1); sub(y[kR],z0,z1);
      end;
    end;
    d:=d*2; k1:=k1 div 2;
    for i:=0 to n-1 do x[i]:=y[i]
  end;
  for i:=0 to n-1 do begin u[i].re:=x[i].re; u[i].im:=x[i].im end;
  if sign<0 then 
    for i:=0 to n-1 do begin u[i].re:=u[i].re/n; u[i].im:=u[i].im/n end
{c2:=clock; writeln('time is ',c2-c1)}
end;
begin
  readln(n);
  for i:=0 to n-1 do begin {read(x[i].re);}x[i].re:=1; x[i].im:=0; 
                           x[i+n].re:=0; x[i+n].im:=0 end;
  for i:=0 to n-1 do begin {read(y[i].re);}y[i].re:=1; y[i].im:=0; 
                           y[i+n].re:=0; y[i+n].im:=0 end;
 { c1:=clock;}
  fft(u,x,2*n,1);
  fft(v,y,2*n,1);
  for i:=0 to 2*n-1 do mul(s[i],u[i],v[i]);
  fft(z,s,2*n,-1);
  for i:=0 to 2*n-1 do c[i]:=trunc(z[i].re+0.5);
  for i:=0 to 2*n-2 do begin
    c[i+1]:=c[i+1]+c[i] div 2;
    c[i]:=c[i] mod 2
  end;
{  c2:=clock;}
  write('n = ', n:5,'  ');
  writeln('time is ',c2-c1);
  {for i:=0 to 9 do writeln(z[i].re);}
  for i:=2*n-1 downto 0 do write(c[i]:4);
  writeln;
{  c1:=clock;}
  for i:=0 to n-1 do a[i]:=trunc(x[i].re);
  for i:=0 to n-1 do b[i]:=trunc(y[i].re);
  for i:=0 to 2*n-1 do c[i]:=0;
  for i:=0 to n-1 do
    for j:=0 to n-1 do
      c[i+j]:=c[i+j] + a[i]*b[j];
  for i:=0 to 2*n-2 do begin
    c[i+1]:=c[i+1]+c[i] div 2;
    c[i]:=c[i] mod 2
  end;
{  c2:=clock;}
  writeln('time is ',c2-c1);
  for i:=2*n-1 downto 0 do write(c[i]:4);
  writeln;
end.

{          FFT                        School method
n =   256  time is        283         time is   100
n =   512  time is        533         time is   417
n =  1024  time is        733         time is   783
n =  2048  time is       1433         time is  6916
n =  4096  time is       2550         time is 27867

sample session
n=4

        x     0 0 0 0 1 1 1 1       decimal 15
        y     0 0 0 0 1 1 1 1       decimal 15
    -------------------------------
                      1 1 1 1
                    1 1 1 1
                  1 1 1 1 
                1 1 1 1
    -------------------------------
              0 1 2 3 4 3 2 1    convolution
                          0 1
       carry            1 
                        0 0 1
       carry        1 
                      0 0 0 1
       carry      1
                    0 0 0 0 1
       carry     1
                  1 0 0 0 0 1
       carry    1
                1 1 0 0 0 0 1
       carry  1 
            0 1 1 1 0 0 0 0 1       decimal 225
}