FFT与NTT模版存档

这是一个悲伤的故事。。

打算从今天开始重学FFT,然后顺便学NTT
因为很多基本操作都不记得了。。
一开始是想在学校教了复数之后再学的。。
但是似乎并没有这个时间了。。
于是就去自己看了高中的复数内容。。
然后发现和高中的微积分类容一样,都是没什么用的东西。。
和我想要的根本不是一个东西,好吧。。
然后去苦战了一波算法导论
发现我看不懂。。
投降。。
于是被下几个结论,就和当时第二次学的时候一样草草收尾了。。
这是第三次了,都没有学会TAT
于是根据结论,我还是只会递归版本

Update4.18

递归版的实现虽然复杂度是一样的,但是常数实在是太大了太大了太大了。实在受不了,于是搞了一些迭代的版本
一般来说,速度可能是3,4,5,6倍等等

迭代版的

NTT

void ntt(LL *a,LL n,LL op)
{
    for (LL u=0;u<n;u++) bin[u]=(bin[u>>1]>>1)|((u&1)*(n>>1));
    for (LL u=0;u<n;u++) if (u<bin[u]) swap(a[u],a[bin[u]]);
    for (LL u=1;u<n;u<<=1)
    {
        LL wn=pow(op==1?g:gi,(MOD-1)/(u<<1)),w,t;
        for (LL i=0;i<n;i=i+(u<<1))
        {
            w=1;
            for (LL k=0;k<u;k++)
            {
                t=w*a[u+i+k]%MOD;
                a[u+i+k]=(a[i+k]-t+MOD)%MOD;
                a[i+k]=(a[i+k]+t)%MOD;
                w=w*wn%MOD;
            }
        }
    }
    if(op==-1)
    {
        LL Inv=pow(n,MOD-2);
        for(LL i=0;i<n;i++) a[i]=a[i]*Inv%MOD;
    }
}

FFT

void fft (qq *a,int n,int o)
{
    for (int u=0;u<n;u++) bin[u]=(bin[u>>1]>>1)|((u&1)*(n>>1));
    for (int u=1;u<n;u++)
        if (u<bin[u])
            swap(a[u],a[bin[u]]);
    for (int u=1;u<n;u<<=1)
    {
        qq wn=qq(cos(2*pi/(u<<1)),o*sin(2*pi/(u<<1)));
        for (int i=0;i<n;i+=(u<<1))
        {
            qq w=qq(1,0);
            for (int j=0;j<u;j++)
            {
                qq t=w*a[u+i+j];
                a[u+i+j]=a[i+j]-t;
                a[i+j]=a[i+j]+t;
                w=w*wn;
            }
        }
    }
}

递归版的

FFT

#include<cstdio>
#include<iostream>
#include<cstring>
#include<cmath>
#include<complex>
using namespace std;
const double pi=M_PI;
const int N=1000005;
int n,m;
complex<double> a[N],b[N];
void fft (complex<double> *a,int n,int o)
{
    if (n==1) return ;
    int k=(n>>1);
    complex<double> w=1,wn(cos(2*pi/n),o*sin(2*pi/n)),a0[k],a1[k];
    for (int u=0;u<k;u++)
    {
        int i=u*2;
        a0[u]=a[i];
        a1[u]=a[i+1];
    }
    fft(a0,k,o);fft(a1,k,o);
    for (int u=0;u<k;u++)
    {
        a[u]=a0[u]+w*a1[u];
        a[u+k]=a0[u]-w*a1[u];
        w=w*wn;
    }
}
int main()
{
    scanf("%d%d",&n,&m);
    for (int u=0;u<=n;u++) scanf("%lf",&a[u]);
    for (int u=0;u<=m;u++) scanf("%lf",&b[u]);
    m=m+n;n=1;while (n<=m) n<<=1;
    fft(a,n,1);fft(b,n,1);
    for (int u=0;u<=n;u++) a[u]*=b[u];
    fft(a,n,-1);
    for (int u=0;u<=m;u++)
        printf("%d ",(int)(a[u].real()/n+0.5));
    return 0;
}

NTT

就是用原根代替单位复数根
使用范围一点都不广。。
什么?你问我原根怎么找?
暴力找啊!一般都不大的
基本上和上面这个代码是一样的
CODE:

#include<cstdio>
#include<iostream>
#include<cstring>
#include<cmath>
#include<complex>
using namespace std;
typedef long long LL;
const LL MOD=998244353,g=3,gi=332748118;
const LL N=1000005;
LL n,m;
LL a[N],b[N];
LL pow (LL x,LL y)
{
    if (y==1) return x;
    LL lalal=pow(x,y>>1);
    lalal=lalal*lalal%MOD;
    if (y&1) lalal=lalal*x%MOD;
    return lalal;
}
void ntt (LL *a,LL n,LL o)
{
    if (n==1) return ;
    LL k=(n>>1);
    LL w=1,wn=pow(o==1?g:gi,(MOD-1)/n),a0[k],a1[k];
    for (LL u=0;u<k;u++)
    {
        LL i=u*2;
        a0[u]=a[i];
        a1[u]=a[i+1];
    }
    ntt(a0,k,o);ntt(a1,k,o);
    for (LL u=0;u<k;u++)
    {
        a[u]=a0[u]+w*a1[u]%MOD;
        a[u]=(a[u]%MOD+MOD)%MOD;
        a[u+k]=a0[u]-w*a1[u];
        a[u+k]=(a[u+k]%MOD+MOD)%MOD;
        w=w*wn%MOD;
    }
}
int main()
{
    scanf("%I64d%I64d",&n,&m);
    for (LL u=0;u<=n;u++) scanf("%I64d",&a[u]);
    for (LL u=0;u<=m;u++) scanf("%I64d",&b[u]);
    m=m+n;n=1;while (n<=m) n<<=1;
    ntt(a,n,1);ntt(b,n,1);
    for (LL u=0;u<=n;u++) a[u]*=b[u];
    ntt(a,n,-1);
    LL inv=pow(n,MOD-2);
    for (LL u=0;u<=m;u++)   printf("%I64d ",a[u]*inv%MOD);
    return 0;
}
  • 1
    点赞
  • 0
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值