HDU

来源:互联网 发布:javascript长度函数 编辑:程序博客网 时间:2024/06/03 20:27

题目链接点这里

过程看这里http://blog.csdn.net/u013368721/article/details/53001532

,,然而,,这位大佬的程序也gg了,,可能改过数据了,,不过思路是这样的,,

然而我与标程对拍,,答案总是差0.00几。。反正,,我写的fft自带低精度,大内存,高常数。。,,以后再看吧

#include<iostream>#include<cstdio>#include<math.h>#include<algorithm>#include<map>#include<set>#include<bitset>#include<stack>#include<queue>#include<string.h>#include<cstring>#include<vector>#include<time.h>#include<stdlib.h>using namespace std;#define INF 0x3f3f3f3f#define INFLL 0x3f3f3f3f3f3f3f3f#define FIN freopen("input.txt","r",stdin)#define mem(x,y) memset(x,y,sizeof(x))typedef unsigned long long ULL;typedef long long LL;#define fuck(x) cout<<"x"<<endl;#define lson l,m,rt<<1#define rson m+1,r,rt<<1|1typedef pair<pair<int,int>,int> PIII;typedef pair<int,int> PII;const double eps=1e-100;const int MX=1111111;int P;int G;double a[MX];//fftconst double pi = acos(-1.0);int len,res[MX],mx;//开大4倍struct Complex{    double r,i;    Complex(double r=0,double i=0):r(r),i(i) {};    Complex operator+(const Complex &rhs)    {        return Complex(r + rhs.r,i + rhs.i);    }    Complex operator-(const Complex &rhs)    {        return Complex(r - rhs.r,i - rhs.i);    }    Complex operator*(const Complex &rhs)    {        return Complex(r*rhs.r - i*rhs.i,i*rhs.r + r*rhs.i);    }} va[MX],vb[MX];void rader(Complex F[],int len)   //len = 2^M,reverse F[i] with  F[j] j为i二进制反转{    int j = len >> 1;    for(int i = 1; i < len - 1; ++i)    {        if(i < j) swap(F[i],F[j]);  // reverse        int k = len>>1;        while(j>=k)        {            j -= k;            k >>= 1;        }        if(j < k) j += k;    }}void FFT(Complex F[],int len,int t){    rader(F,len);    for(int h=2; h<=len; h<<=1)    {        Complex wn(cos(-t*2*pi/h),sin(-t*2*pi/h));        for(int j=0; j<len; j+=h)        {            Complex E(1,0); //旋转因子            for(int k=j; k<j+h/2; ++k)            {                Complex u = F[k];                Complex v = E*F[k+h/2];                F[k] = u+v;                F[k+h/2] = u-v;                E=E*wn;            }        }    }    if(t==-1)   //IDFT        for(int i=0; i<len; ++i)            F[i].r/=len;}/*inline void FFT(Complex *a,int n,int r){    rader(a,len);    for (int i=1; i<n; i<<=1)    {        Complex wn(cos(pi/i),r*sin(pi/i));        for (int j=0; j<n; j+=(i<<1))        {            Complex w(1,0);            for (int k=0; k<i; k++,w=w*wn)            {                Complex x=a[j+k],y=w*a[j+k+i];                a[j+k]=x+y;                a[j+k+i]=x-y;            }        }    }    if (r==-1) for (int i=0; i<n; i++) a[i].r/=n;}*/void Conv(Complex a[],Complex b[],int len)   //求卷积{    FFT(a,len,1);    FFT(b,len,1);    for(int i=0; i<len; ++i) a[i] = a[i]*b[i];    FFT(a,len,-1);}int change[MX];double ans[MX];double ss(int ret){    return pow(2,sin(2*pi*ret/P)*sin(2*pi*ret/P)*sin(2*pi*ret/P));}void gao(){    mx=P-2+P-2;    len=1;    while(len<mx)len<<=1;//mx为卷积后最大下标    int ret=1;    va[P-2-0].r=a[1];    vb[0].r=ss(ret);    change[0]=1;    for(int i=1; i<P-1; i++)    {        ret=(ret*G)%P;        change[i]=ret;        va[P-2-i].r=a[ret];        va[P-2-i].i=0;        vb[i].r=ss(ret);        vb[i].i=0;    }    for(int i=P-1; i<len; i++)va[i].i=va[i].r=vb[i].i=vb[i].r=0;    Conv(va,vb,len);    ans[1]=va[P-2].r;    ans[0]=a[P-1];    for(int i=0; i<P-2; i++)    {        ans[change[i+1]]=va[i].r+va[i+P-1].r;        ans[0]+=a[i+1];    }    for(int i=0; i<=P-1; i++)    {        printf("%.3f ",ans[i]+a[0]);    }    puts("");}int main(){    FIN;    freopen("output1.txt","w",stdout);    while(~scanf("%d",&P))    {        if(P==13)G=2;        else if(P==103)G=5;        else G=2;        for(int i=0; i<P; i++)scanf("%lf",&a[i]);        gao();    }    return 0;}


原创粉丝点击