FFT/NTT 总结

本总结主要用于帮助个人理解,讲得不足之处,还请各位看官谅解

FFT

补充知识

\(n\)次单位复根(\(w_n\)):

使得\(z^n=1\)的一类复数,这些复数一共有\(n\)个,它们都分布在复平面的单位圆上,并且连线构成一个正\(n\)边形

点值表示:

多项式\(f(x)=\sum_{i=0}^{n-1}{a_ix^i}\)的点值表示为\(n\)个点\((x_i,y_i)\),其中\(y_i=f(x_i)\)

递归算法主要思路

由折半引理\(w_{2n}^{2k}=w_n^k\),将其代入可得:

\[f(w_{n}^{k})=\sum_{i=0}^{n-1}{a_iw_{n}^{ki}}\]

\[=\sum_{i=0}^{\frac{n}{2}-1}{a_{2i}w_n^{2ki}}+w_n^{k}\sum_{i=0}^{\frac{n}{2}-1}{a_{2i+1}w_n^{2ki}}\]

\[=\sum_{i=0}^{\frac{n}{2}-1}{a_{2i}w_{\frac{n}{2}}^{ki}}+w_n^{k}\sum_{i=0}^{\frac{n}{2}-1}{a_{2i+1}w_{\frac{n}{2}}^{ki}}\]

而且由于\(w_n^{k+\frac{n}{2}}=-w_n^k\),

所以

\[f(w_{n}^{k+\frac{n}{2}})=\sum_{i=0}^{\frac{n}{2}-1}{a_{2i}w_{\frac{n}{2}}^{ki}}+w_n^{k+\frac{n}{2}}\sum_{i=0}^{\frac{n}{2}-1}{a_{2i+1}w_{\frac{n}{2}}^{ki}}\]

\[=\sum_{i=0}^{\frac{n}{2}-1}{a_{2i}w_{\frac{n}{2}}^{ki}}-w_n^k\sum_{i=0}^{\frac{n}{2}-1}{a_{2i+1}w_{\frac{n}{2}}^{ki}}\]

逐步递归下去,最终只需要代入\(w_1^0\)的值\((1)\)并在回溯的过程中不断更新该表达式的项即可

最后每一个系数\(a_i\)都会被转换成点值表示中的\(y_i(f(w_n^{i}))\)

非递归算法主要思路

Rader:由于最后递归完毕的编号为该编号对应下标的二进制表示的倒序

故将系数预先排好序,最后再一层层向上更新即可

Rader实现的一行代码(r[i]表示i倒序后对应的编号):

for(RG int i=0;i<n;++i)r[i]=(r[i>>1]>>1)|((i&1)<<(l-1));

主要思想其实就是把\(\times2\)的操作转换成了\(>>1\),再根据情况考虑是否要在最高位\(+1\)

极大地减小了常数啊

代码

#include<bits/stdc++.h>
#include<algorithm>
#include<iostream>
#include<cstdlib>
#include<iomanip>
#include<complex>
#include<cstring>
#include<vector>
#include<cstdio>
#include<string>
#include<bitset>
#include<cmath>
#include<queue>
#include<stack>
#include<map>
#include<set>
#define isr(i) (t[t[i].fa].s[1]==i)
#define ls t[i].s[0]
#define rs t[i].s[1]
#define pb push_back
#define RG register
#define il inline
using namespace std;
const double pi=acos(-1);
const int mod=666623333;
const int N=3000010;
typedef unsigned long long ull;
typedef vector<int>VI;
typedef long long ll;
typedef double dd;
il ll read(){
    RG ll data=0,w=1;RG char ch=getchar();
    while(ch!=‘-‘&&(ch<‘0‘||ch>‘9‘))ch=getchar();
    if(ch==‘-‘)w=-1,ch=getchar();
    while(ch<=‘9‘&&ch>=‘0‘)data=data*10+ch-48,ch=getchar();
    return data*w;
}

struct Complex{double r,i;};
Complex operator +(Complex a,Complex b)
{return (Complex){a.r+b.r,a.i+b.i};}
Complex operator -(Complex a,Complex b)
{return (Complex){a.r-b.r,a.i-b.i};}
Complex operator *(Complex a,Complex b)
{return (Complex){a.r*b.r-a.i*b.i,a.i*b.r+a.r*b.i};}
typedef Complex CD;

int n,m,l,r[N];
CD A[N],B[N];
il void fft(CD *A,int n,int opt){
    for(RG int i=0;i<n;i++)if(i<r[i])swap(A[i],A[r[i]]);
    for(RG int i=2;i<=n;i<<=1){
        RG CD w=(CD){cos(2*pi/i),opt*sin(2*pi/i)};
        for(RG int j=0;j<n;j+=i){
            RG CD wn=(CD){1,0};
            for(RG int k=j;k<j+(i>>1);k++,wn=wn*w){
                RG CD x=wn*A[k+(i>>1)];
                A[k+(i>>1)]=A[k]-x;
                A[k]=A[k]+x;
            }
        }
    }
}

int main()
{
    n=read();m=read();
    for(RG int i=0;i<=n;i++)A[i].r=read();
    for(RG int i=0;i<=m;i++)B[i].r=read();
    for(m+=n,n=1;n<=m;n<<=1)l++;
    for(RG int i=0;i<n;i++)r[i]=(r[i>>1]>>1)|((i&1)<<(l-1));

    fft(A,n,1);fft(B,n,1);
    for(RG int i=0;i<n;i++)A[i]=A[i]*B[i];
    fft(A,n,-1);
    for(RG int i=0;i<=m;i++)printf("%d ",(int)(A[i].r/n+0.5));

    return 0;
}

NTT

补充知识

原根:

(个人理解)对于一个给定的质数\(p\),它的原根\(g\)指使得\(g^k(k=0,1,2,...,p-2)\)在\(mod\) \(p\)意义下两两不同的最小值

由于在FFT时我们需要求出\(w_2,w_4,w_8,...w_{\frac{n}{2}},w_n\),因此我们希望\(p-1\)是一个\(r\times2^k\)的形式

正好,\(998244353=119*2^{23}+1\),在\(n\leq20\)的情况下够我们用的了

主要思路

---

可以知道在FFT下的单位复根\(w_i\)就是在NTT下的\(g^{\frac{p-1}{i}}(i=2^k)\)

其他过程和FFT大致相同,其正变换的旋转向量为\(g\),逆变换的旋转向量设为\(g^{-1}\)即可

代码

与上面FFT代码不同的地方已在代码中标出

#include<bits/stdc++.h>
#include<algorithm>
#include<iostream>
#include<cstdlib>
#include<iomanip>
#include<cstring>
#include<complex>
#include<vector>
#include<cstdio>
#include<string>
#include<bitset>
#include<cmath>
#include<queue>
#include<stack>
#include<map>
#include<set>
#define mp make_pair
#define pb push_back
#define RG register
#define il inline
using namespace std;
typedef unsigned long long ull;
typedef vector<int>VI;
typedef long long ll;
typedef double dd;
const int mod=998244353;//模数
const int rev=332748118;//原根的逆元
const int N=3000010;
const dd eps=1e-10;
const int g=3;//原根
il ll read(){
    RG ll data=0,w=1;RG char ch=getchar();
    while(ch!=‘-‘&&(ch<‘0‘||ch>‘9‘))ch=getchar();
    if(ch==‘-‘)w=-1,ch=getchar();
    while(ch<=‘9‘&&ch>=‘0‘)data=data*10+ch-48,ch=getchar();
    return data*w;
}

il void file(){
    //freopen(".in","r",stdin);
    //freopen(".out","w",stdout);
}

int n,m,l,r[N],A[N],B[N];
il ll poww(ll a,ll b){
    RG ll ret=1;
    for(a%=mod;b;b>>=1,a=a*a%mod)
        if(b&1)ret=ret*a%mod;
    return ret;
}

il void NTT(int *A,int n,bool opt){
    for(RG int i=0;i<n;i++)if(i<r[i])swap(A[i],A[r[i]]);
    for(RG int i=2;i<=n;i<<=1){
        RG int w=poww((opt?g:rev),(mod-1)/i);//这里变成了原根的次幂

        for(RG int j=0;j<n;j+=i){
            RG int wn=1;

            for(RG int k=j;k<j+(i>>1);k++,wn=1ll*w*wn%mod){
                RG int x=1ll*wn*A[k+(i>>1)]%mod;

                A[k+(i>>1)]=(A[k]-x+mod)%mod;
                A[k]=(A[k]+x)%mod;
            }
        }
    }
}

int main()
{
    n=read();m=read();
    for(RG int i=0;i<=n;i++)A[i]=read();
    for(RG int i=0;i<=m;i++)B[i]=read();
    for(m+=n,n=1;n<=m;n<<=1)l++;
    for(RG int i=0;i<n;i++)r[i]=(r[i>>1]>>1)|((i&1)<<(l-1));

    NTT(A,n,1);NTT(B,n,1);
    for(RG int i=0;i<n;i++)A[i]=1ll*A[i]*B[i]%mod;
    NTT(A,n,0);

    for(RG int i=0,rv=poww(n,mod-2);i<=m;i++)
        printf("%lld ",1ll*A[i]*rv%mod);
    //这里除n的时候记得变成乘逆元
    puts("");

    return 0;
}

原文地址：https://www.cnblogs.com/cjfdf/p/8428344.html

时间： 2024-10-18 13:53:07

FFT/NTT 总结的相关文章

FFT/NTT做题方法与调试技巧(+提高码题效率的一些想法)

(其实本文应该写成了"结合FFT讨论的调试技巧+码题方法",语文不好一写文章就偏题QAQ) 有意见欢迎提出,有遗漏欢迎补充! FFT(快速傅里叶变换)/NTT(数论变换)是卷积运算常见而实用的优化 但是FFT/NTT的处理过程并不像暴力运算(差不多是多项式乘法)那样能够直观地反映卷积结果的实时变化. 因此在使用FFT时将会或多或少地加大调试的难度. 如果调试程序时直接跟踪变量,每步手算结果比对,不仅会耽误大量时间,而且效果可能并不理想. 直接肉眼查错效率可能也不太高. 但也正由于FFT

FFT NTT 错误总结(持续更新)

FFT NTT错误总结 1 处理\(r\)数组时忘记赋值 r[i] = (r[i >> 1] >> 1) | ((i & 1) << (l - 1)); 2 负数重载运算符 point operator * (point a,point b){ return point(a.x * b.x - a.y * b.y,a.x * b.y + a.y * b.x); } 3 欧拉公式记不清楚 point Wn = point(cos(Pi / mid),type *

FFT/NTT基础题总结

在学各种数各种反演之前把以前做的$FFT$/$NTT$的题整理一遍 还请数论$dalao$口下留情 T1快速傅立叶之二 题目中要求求出 $c_k=\sum\limits_{i=k}^{n-1}a_i*b_{i-k}$ 首先可以把$a$翻转, $c_k=\sum\limits_{i=k}^{n-1}a_{n-1-i}*b_{i-k}$ $c_k=\sum\limits_{i=0}^{n-k-1}a_{n-k-1-i}*b_{i}$ 也就是说对新的$a$,$b$数组做一遍$FFT$得到的便是$c$数

多项式$fft$,$ntt$,$fwt$初步

最近在学多项式和生成函数. 上课听$lnc$大神讲还是$mengbier$. 作为多项式的前置芝士,$fft,ntt$等是必学的. 在此记录一些关于$fft,ntt,fwt$的知识及例题... FFT: 应用在处理$\sum _{i+j=k} f[i]*g[j]$的卷积上. 看网上大佬的博客,基本入了门吧. 自己的关于原理的一些见解: 多项式有系数表示和点值表示,两种表示方法可以相互转化. FFT可以在$O(n*longn)$内解决多项式乘法. 具体是,点值表示的方法应用在多项式上是$O(n)$

HDU-4609(FFT/NTT)

HDU-4609(FFT/NTT) 题意: 给出n个木棒,现从中不重复地选出3根来,求能拼出三角形的概率. 计算合法概率容易出现重复,所以建议计算不合法方案数 枚举选出的最大边是哪条,然后考虑剩下两条边之和小于等于它 两条边之和为\(x\)的方案数可以\(FFT/NTT\)得到,是一个简单的构造 即\(f(x)=\sum x^{length_i}\),求出\(f(x)^2\),就能得到和的方案数,但是会重复,包括自己和自己算,一对算两次 处理一下前缀和即可 #include<bits/stdc+

HDU-5730(CDQ+FFT/NTT)

HDU-5730(CDQ+FFT/NTT) 题意:将长度为\(n\)的序列分成若干段,每段\([l,r]\)的权值为\(a_{r-l+1}\),一种分法的权值为所有段的乘积,求所有可能的分法的权值和 根据题意可以得到简单\(dp\) \(dp_0=1,dp_i=\sum_0^{i-1}dp_j \cdot a_{i-j}\) 可以看到是一个\(i-j\)形式的作差卷积 但是直接卷积我们无法保证先求出了\(dp_j\),所有可以用\(CDQ\)分治优化,复杂度\(n\log^2 n\) 具体实现见

[拉格朗日反演][FFT][NTT][多项式大全]详解

1.多项式的两种表示法 1.系数表示法 我们最常用的多项式表示法就是系数表示法,一个次数界为\(n\)的多项式\(S(x)\)可以用一个向量\(s=(s_0,s_1,s_2,\cdots,s_n-1)\)系数表示如下:\[S(x)=\sum_{k=0}^{n-1}s_kx^k\] 系数表示法很适合做加法,可以在\(O(n)\)的时间复杂度内完成,表达式为:\[S(x)=A(x)+B(x)=\sum_{k=0}^{n-1}(a_k+b_k)x^k\] 当中\[s_k=a_k+b_k\] 但是,系数

多项式FFT/NTT模板(含乘法/逆元/log/exp/求导/积分/快速幂)

自己整理出来的模板 存在的问题: 1.多项式求逆常数过大(尤其是浮点数FFT) 2.log只支持f[0]=1的情况,exp只支持f[0]=0的情况 有待进一步修改和完善 FFT: 1 #include<bits/stdc++.h> 2 using namespace std; 3 typedef long long ll; 4 typedef double db; 5 const db pi=acos(-1); 6 const int N=4e5+10,M=1e6+10,mod=9982443

FFT&amp;NTT模板

博主菜鸡,只会背模板,根本解释不来嘤嘤嘤~ FFT: #include<bits/stdc++.h> #define IL inline #define LF double using namespace std; const int N=5e6+5; const LF Pi=acos(-1.0); struct com{ LF x,y; com(LF xx=0,LF yy=0){x=xx,y=yy;} com operator+(const com &a) const{return