FFT

雪藏到NOIP2018结束

JU_force

#include <bits/stdc++.h>
const int oo=2147483000;
#define sc(x) scanf("%d",&x)
const int N=2001000;
#define PI 3.1415926535897932384626
using namespace std;
typedef complex <double> cd;
int rev[N];
void get_rev (int bit)
{
    for (int i=0; i<(1<<bit); i++)
        rev[i]=(rev[i>>1]>>1)|((i&1)<<(bit-1));
}
void FFT (cd *a,int n,int kind)
{
    for (int i=0; i<n; i++)
        if (i<rev[i]) swap (a[i],a[rev[i]]);
    for (int i=1; i<n; i<<=1)
    {
        cd wn=exp(cd(0,kind*PI/i));
        for (int j=0; j<n; j+=i<<1)
        {
            cd wnk(1,0);
            for (int k=j; k<j+i; k++)
            {
                cd x=a[k],y=wnk*a[k+i];
                a[k]=x+y;
                a[k+i]=x-y;
                wnk*=wn;
            }
        }
    }
    if (kind==-1) for (int i=0; i<n; i++) a[i]/=n;
}
cd a[N],b[N];
char s1[N],s2[N];
int ans[N];
int main()
{
    //freopen (".in","r",stdin);
    //freopen (".out","w",stdout);
    int sbsbsbsbs; sc(sbsbsbsbs); 
    scanf ("%s%s",s1,s2);
    int l1=strlen (s1),l2=strlen (s2);
    int s=2,bit=1;
    for (bit=1; (1<<bit)<l1+l2-1; bit++) s<<=1;
    for (int i=0; i<l1; i++) a[i]=(double)(s1[l1-i-1]-'0');
    for (int i=0; i<l2; i++) b[i]=(double)(s2[l2-i-1]-'0');
    get_rev (bit);
    FFT (a,s,1);
    FFT (b,s,1);
    for (int i=0; i<s; i++) a[i]*=b[i];
    FFT (a,s,-1);
    for (int i=0; i<s; i++)
    {
        ans[i]+=(int)(a[i].real()+0.5);
        ans[i+1]+=ans[i]/10;
        ans[i]%=10;
    }
    int len=l1+l2;
    while (!ans[len]&&len>=0) len--;
    if (len==-1) printf ("0");
    for (int i=len; i>=0; i--) printf ("%d",ans[i]);
    return 0;
}

LZH_force

 

#include <iostream>
#include <fstream>
#include <cstring>
#include <string>
#include <cstdlib>
#include <cstdio>
#include <algorithm>
#include <cmath>
#include <bitset>
#include <ctime>
#include <map>
#include <queue>
#include <complex>
#define pi M_PI
using namespace std;

struct Node
{
    complex <double> a[400005];
    int len;
}q,g,s1,s2,A,ans1,ans2,a;
int n;
double x;
int done(int x)
{
    int g=1;
    for (;g<x;g<<=1);
    return g;
}
int rev(int k,int v)
{
    int g=0;
    for (int i=1;i<v;i<<=1)
    {
        g<<=1;
        g+=k&1;
        k>>=1;
    }
    return g;
}
int fft(int op)
{
    complex <double > wm,w,t,u;
    A.len=a.len;
    for (int i=0;i<a.len;i++)
        A.a[rev(i,a.len)]=a.a[i];
    int n=a.len;
    int xd=(int)(log(n)/log(2));
    for (int i=1;i<=xd;i++)
    {
        int m=1 << i;
        wm=complex <double> (cos(2*pi*op/m),sin(2*pi*op/m));
        for (int k=0;k<n;k+=m)
        {
            w=1.0;
            for (int j=0;j<m/2;j++)
            {
                t=w*A.a[k+j+m/2];
                u=A.a[k+j];
                A.a[k+j]=u+t;
                A.a[k+j+m/2]=u-t;
                w*=wm;
            }
        }
    }
    return 0;
}
int main()
{
    freopen("force.in","r",stdin);
    freopen("force.out","w",stdout);
    scanf("%d",&n);
    for (int i=0;i<n;i++) scanf("%lf",&x),q.a[i]=x;
    q.len=g.len=done(n)*2;
    for (int i=1;i<n;i++) g.a[i]=1.0/i/i;
    a=q;
    fft(1);
    s1=A;
    a=g;
    fft(1);
    s2=A;
    g=s1;
    for (int i=0;i<g.len;i++) g.a[i]*=s2.a[i];
    a=g;
    fft(-1);
    s1=A;
    ans1=s1;
    for (int i=0;i<ans1.len;i++) ans1.a[i]/=ans1.len;
    for (int i=0;i<n/2;i++) swap(q.a[i],q.a[n-i-1]);
    for (int i=0;i<ans1.len;i++) g.a[i]=0;
    for (int i=1;i<n;i++) g.a[i]=-1.0/i/i;
    a=q;
    fft(1);
    s1=A;
    a=g;
    fft(1);
    s2=A;
    g=s1;
    for (int i=0;i<g.len;i++) g.a[i]*=s2.a[i];
    a=g;
    fft(-1);
    s1=A;
    ans2=s1;
    for (int i=0;i<ans2.len;i++) ans2.a[i]/=ans2.len;
    for (int i=0;i<n/2;i++) swap(ans2.a[i],ans2.a[n-i-1]);
    for (int i=0;i<n;i++) 
       printf("%lf\n",(ans1.a[i]+ans2.a[i]).real());
    return 0;
}