zoukankan      html  css  js  c++  java
  • 再写FFT模板

      没什么好说的,今天又考了FFT(虽然不用FFT也能过)但是确实有忘了怎么写FFT了,于是乎只有重新写一遍FFT模板练一下手了。第一部分普通FFT,第二部分数论FFT,记一下模数2^23*7*17+1

      

    #include<iostream>
    #include<cstdio>
    #include<cstring>
    #include<algorithm>
    #include<cmath>
    using namespace std;
    typedef double real;
    namespace task1
    {
            const int maxn = 110000;
            const real pi = 3.1415926535897932;
            struct Complex
            {
                    real x,y;
                    Complex(){}
                    Complex(real x,real y=0):x(x),y(y){}
            };
            Complex operator + (Complex c1,Complex c2)
            {
                    return Complex(c1.x+c2.x,c1.y+c2.y);
            }
            Complex operator - (Complex c1,Complex c2)
            {
                    return Complex(c1.x-c2.x,c1.y-c2.y);
            }
            Complex operator * (Complex c1,Complex c2)
            {
                    return Complex(c1.x*c2.x-c1.y*c2.y,c1.y*c2.x+c1.x*c2.y);
            }
            Complex operator / (Complex c1,real k)
            {
                    return Complex(c1.x/k,c1.y/k);
            }
            char s1[maxn],s2[maxn];
            Complex g1[maxn],g2[maxn],rr[maxn];
            int lst[maxn];
            void FFT(Complex g[],int l,int p)
            {
                    int x;
                    for (int i=0;i<l;i++)
                    {
                            x=0;
                            for (int j=1,k=(l>>1);j<l;j<<=1,k>>=1)
                                    if (i&j)x+=k;
                            if (x<i)
                                    swap(g[x],g[i]);
                    }
                    Complex wt,w,tmp;
                    for (int i=1;i<l;i<<=1)
                    {
                            wt=Complex(cos(pi/i*p),sin(pi/i*p));
                            for (int j=0;j<l;j+=(i<<1))
                            {
                                    Complex w=Complex(1,0);
                                    for (int k=0;k<i;k++)
                                    {
                                            tmp=g[j+k];
                                            g[j+k]=g[j+k]+g[i+j+k]*w;
                                            g[i+j+k]=tmp-g[i+j+k]*w;
                                            w=w*wt;
                                    }
                            }
                    }
            }
            void main()
            {
                    int l1,l2;
                    scanf("%s",s1);
                    scanf("%s",s2);
                    l1=(int)strlen(s1);
                    l2=(int)strlen(s2);
                    for (int i=0;i<l1;i++)
                            g1[(l1-i-1)/4]=g1[(l1-i-1)/4].x*10+s1[i]-'0';
                    for (int i=0;i<l2;i++)
                            g2[(l2-i-1)/4]=g2[(l2-i-1)/4].x*10+s2[i]-'0';
                    int l=max(l1,l2)/4+1;
                    while (l!=(l&(-l)))l-=l&(-l);
                    l<<=2;
                    FFT(g1,l,1);
                    FFT(g2,l,1);
                    for (int i=0;i<l;i++)rr[i]=g1[i]*g2[i];
                    FFT(rr,l,-1);
                    for (int i=0;i<l;i++)
                    {
                            rr[i]=rr[i]/l;
                            lst[i]=(int)(rr[i].x+0.5);
                    }
                    for (int i=0;i<l;i++)
                    {
                            lst[i+1]+=lst[i]/10000;
                            lst[i]%=10000;
                    }
                    while (l>1 && !lst[l-1])l--;
                    printf("%d",lst[l-1]);
                    for (int i=l-2;i>=0;i--)
                            printf("%04d",lst[i]);
                    printf("
    ");
            }
    }
    namespace task2
    {
            typedef long long qword;
            const int maxn = 10000;
            const qword p = 998244353;
            char s1[maxn],s2[maxn];
            qword g1[maxn],g2[maxn];;
            qword rr[maxn];
            qword pow_mod(qword x,qword y)
            {
                    qword ret=1;
                    while (y)
                    {
                            if (y&1)ret=ret*x%p;
                            x=x*x%p;
                            y>>=1;
                    }
                    return ret;
            }
            void FFT(qword g[],int l,int f)
            {
                    int x=0;
                    for (int i=1;i<l;i++)
                    {
                            x=0;
                            for (int j=1,k=(l>>1);j<l;j<<=1,k>>=1)
                                    if (i&j)x+=k;
                            if (i<x)
                                    swap(g[i],g[x]);
                    }
                    qword w,wt,tmp;
                    for (int i=1,t=1;i<l;i<<=1,t++)
                    {
                            wt=pow_mod(3,7*17*(1<<(23-t)))%p;
                            if (f==-1)
                                    //wt=-wt;
                                    wt=pow_mod(wt,p-2);
                            for (int j=0;j<l;j+=(i<<1))
                            {
                                    w=1;
                                    for (int k=0;k<i;k++)
                                    {
                                            tmp=g[j+k];
                                            g[j+k]=(g[j+k]+g[i+j+k]*w)%p;
                                            g[i+j+k]=(tmp-g[i+j+k]*w)%p;
                                            w=w*wt%p;
                                    }
                            }
                    }
            }
            void main()
            {
                    int l1,l2;
                //    for (int i=0;i<5;i++)
                //            printf("[1/%d*pi]:%d
    ",(1<<i),7*17*(1<<(23-i)));
                    scanf("%s",s1);
                    scanf("%s",s2);
                    l1=(int)strlen(s1);
                    l2=(int)strlen(s2);
                    for (int i=0;i<l1;i++)
                            g1[i]=s1[l1-i-1]-'0';
                    for (int i=0;i<l2;i++)
                            g2[i]=s2[l2-i-1]-'0';
                    int l=max(l1,l2);
                    while (l!=(l&(-l)))l-=l&(-l);
                    l<<=2;
                    FFT(g1,l,1);
                    FFT(g2,l,1);
                    for (int i=0;i<l;i++)rr[i]=g1[i]*g2[i]%p;
                    FFT(rr,l,-1);
                    qword invl=pow_mod(l,p-2);
                    for (int i=0;i<l;i++)rr[i]=(rr[i]*invl%p+p)%p;
                    for (int i=0;i<l;i++)
                    {
                            rr[i+1]+=rr[i]/10;
                            rr[i]%=10;
                    }
                    while (l>1 && ! rr[l-1])l--;
                    for (int i=l-1;i>=0;i--)
                            printf("%d",(int)rr[i]);
            }
    }
    
    int main()
    {
            freopen("input.txt","r",stdin);
            task1::main();
            task2::main();
    }
  • 相关阅读:
    Spring Boot 学习随记
    Prometheus 普罗米修斯监控
    安装VC++6.0步骤及心得
    NFS 系统搭建
    Centos 搭建邮箱系统
    搭建 RTMP 服务器
    阿里云 DTS 实践
    ELK 搭建
    Prometheus 和 Grafana 安装部署
    Centos7 Nagios 搭建
  • 原文地址:https://www.cnblogs.com/mhy12345/p/4531821.html
Copyright © 2011-2022 走看看