FFT
struct node{
long double x,y;
node(){}
node(long double xx,long double yy){x=xx; y=yy; return ;}
node operator +(const node &a)const{return node(x+a.x,y+a.y);}
node operator -(const node &a)const{return node(x-a.x,y-a.y);}
node operator *(const node &a)const{return node(x*a.x-y*a.y,x*a.y+y*a.x);}
inline node conj(){return node(x,-y);}
}A[N],B[N];
inline void FFT(node *f,int lim,int opt){
for(reg int i=0;i<lim;++i) if(i<r[i]) swap(f[i],f[r[i]]);
for(reg int p=2;p<=lim;p<<=1){
int len=p>>1; node tmp(cos(pi/len),opt*sin(pi/len));
for(reg int k=0;k<lim;k+=p){
node buf(1,0); for(reg int l=k;l<k+len;++l){
node tt=buf*f[l+len];
f[len+l]=f[l]-tt; f[l]=f[l]+tt; buf=buf*tmp;
}
}
} return ;
}
NTT
inline void NTT(int *f,int lim,int fl){
for(reg int i=0;i<lim;++i) if(i<r[i]) swap(f[i],f[r[i]]);
for(reg int p=2;p<=lim;p<<=1){
int len=p>>1,tmp=ksm(G,(mod-1)/p); if(fl==-1) tmp=ksm(tmp,mod-2);
for(reg int k=0;k<lim;k+=p){
int buf=1; for(reg int l=k;l<k+len;++l){
int tt=mul(buf,f[l+len]); f[l+len]=del(f[l],tt); f[l]=add(f[l],tt);
buf=mul(buf,tmp);
}
}
} if(fl==-1) for(reg int i=0;i<lim;++i) f[i]=mul(f[i],ksm(lim,mod-2)); return ;
}
所以有了这个就可以去找点基础计数来做
分治FFT
还是 \(cdq\) 的那一套
每次对于当前区间处理当前左边和 \(g\) 的卷积,把答案累加到右侧
按照左中右的顺序进行处理
有些细节:开始的时候把序列长度补成 \(2^x\) 的形式,给右边添加的时候的起点需要思考
这样就会做快速求第二类斯特林数的一列的形式
inline void NTT(int *f,int lim,int fl){
for(reg int i=0;i<lim;++i) if(i<R[i]) swap(f[i],f[R[i]]);
for(reg int p=2;p<=lim;p<<=1){
int len=p>>1,tmp=ksm(3,(mod-1)/p); if(fl==-1) tmp=ksm(tmp,mod-2);
for(reg int k=0;k<lim;k+=p){
int buf=1; for(reg int l=k;l<k+len;++l){
int tt=mul(buf,f[l+len]); f[l+len]=del(f[l],tt); f[l]=add(f[l],tt);
buf=mul(buf,tmp);
}
}
} if(fl==-1) for(reg int i=0;i<lim;++i) f[i]=mul(f[i],ksm(lim,mod-2)); return ;
}
inline void solve(int l,int r){
if(l+1==r) return ; int mid=l+((r-l)>>1); solve(l,mid);
int len=r-l; for(reg int i=0;i<len;++i) R[i]=R[i>>1]>>1|((i&1)?(len>>1):0);
for(reg int i=0;i<len;++i) G[i]=g[i]; for(reg int i=l;i<mid;++i) F[i-l]=f[i]; for(reg int i=mid;i<r;++i) F[i-l]=0;
NTT(F,len,1); NTT(G,len,1); for(reg int i=0;i<len;++i) F[i]=mul(F[i],G[i]); NTT(F,len,-1);
for(reg int i=mid;i<r;++i) f[i]=add(f[i],F[i-l]);
return solve(mid,r);
}
任意模数NTT
考虑 \(FFT\) 的精度如果在系数是 \(10^9\) 的情况会爆炸,那么考虑优化这个直接 \(FFT\) 再取模的做法
把当前多项式的每个 \(A_i\) 拆成 \(A_0[i]\times k+A_1[i]\) 这里 \(A_0[i]=A[i]/k,A_1[i]=A[i]\%k\),\(B\) 也拆一下系数
\(DFT\) 求点值的时候 考虑把两个并成一个,因为每个数一开始的虚部均为 \(0\)
构造 \(P_i=A_i+iB_i,Q_i=A_i-iB_i\),发现这两个点值显然共轭,所以可以直接反过来得到
那么解方程组可以的到 \(DFT\) 的结果(复数科技妙)
得到 \(DFT\) 的结果之后,考虑怎么 \(IDFT\) 四个多项式得到最终答案
接着构造 \(P_i=A_0[i]B_0[i]+iA_1[i]B_0[i],Q_i=A_1[i]B_1[i]+A_0[i]B_1[i]\)
因为这些个多项式都是实数,那么虚不会变成实,实不会变成虚,所以直接 \(IDFT\) 回去得到答案
很厉害的东西
inline void MTT(){
for(reg int i=0;i<=n;++i) A0[i].x=a[i]>>15,A0[i].y=a[i]&32767;
for(reg int i=0;i<=m;++i) B0[i].x=b[i]>>15,B0[i].y=b[i]&32767;
int lim=1; while(lim<=(n+m)) lim<<=1;
for(reg int i=0;i<lim;++i) r[i]=r[i>>1]>>1|((i&1)?(lim>>1):0);
FFT(A0,lim,1); FFT(B0,lim,1);
for(reg int i=0;i<lim;++i) A1[i]=A0[i?lim-i:i].conj(),B1[i]=B0[i?lim-i:i].conj();
node p,q;
for(reg int i=0;i<lim;++i) p=A0[i],q=A1[i],A0[i]=(p+q)*node(0.5,0),A1[i]=(q-p)*node(0,0.5);
for(reg int i=0;i<lim;++i) p=B0[i],q=B1[i],B0[i]=(p+q)*node(0.5,0),B1[i]=(q-p)*node(0,0.5);
for(reg int i=0;i<lim;++i) P[i]=A0[i]*B0[i]*node(1,0)+A1[i]*B0[i]*node(0,1),Q[i]=A1[i]*B1[i]*node(0,1)+A0[i]*B1[i]*node(1,0);
FFT(P,lim,-1); FFT(Q,lim,-1);
for(reg int i=0;i<=n+m;++i) printf("%lld ",((int)(P[i].x/lim+0.5)%mod*(1ll<<30)%mod+((int)(P[i].y/lim+0.5)%mod+(int)(Q[i].x/lim+0.5)%mod)%mod*(1ll<<15)%mod+(int)(Q[i].y/lim+0.5)%mod)%mod); puts("");
return ;
}
多项式求逆
已知 \(H(x)F(x)\equiv 1\mod x^{\frac n 2}\),求
\[G(x)F(x)=1 \mod x^n \]必有
\[G(x)F(x)=1 \mod x^{\frac n 2} \]所以做差,推一下式子有:
\[G(x)\equiv 2H(x)-F(x)H(x)^2\mod x^n \]分治即可,边界为 \(G(0)\equiv F(0) \mod 998244353\),复杂度 \(\Theta(n\log n)\)
任意模数写着 \(namespace\) 即可
inline void solve(int Len){
if(Len==1) return G[0]=ksm(F[0],mod-2),void(); solve((Len+1)>>1);
int lim=1; while(lim<=(Len<<1)) lim<<=1;
for(reg int i=0;i<Len;++i) H[i]=F[i]; for(reg int i=Len;i<lim;++i) H[i]=0;
P.MTT(H,G,Len,lim); P.MTT(H,G,Len,lim);
for(reg int i=0;i<Len;++i) G[i]=del(mul(G[i],2),H[i]);
for(reg int i=Len;i<lim;++i) G[i]=0;
return ;
}
城市规划
设 \(dp[i]\) 表示当前 \(i\) 个点构成的无项联通图的数目,\(dp[1]=1\)
转移考虑 \(i\) 个点每个点都可以向下一个点连边,所以是 \(dp[i+1]=dp[i] \times (2^i-1)\)
但是 \((5,7),(4,7)\) 的连边方式就可以保证 \((4,5)\)联通,并不需要 \(1\to 4\) 的哪个点和 \(5\) 连
设 \(n\) 个点的不要求联通的无向联通图数目为 \(g(n)=2^{\frac{n(n-1)}2-1}\)
设 \(f(n)\) 为联通的无向图,那么就有如下的式子:
\[g(n)=\sum_{i=1}^n \binom {n-1}{i-1} g(n-i) f(i) \]具体就是枚举第一个点所在的联通块里面有多少个点
设 \(F(x)=\frac{f(x)}{(x-1)!},G(x)=\frac{g(x)}{(x-1)!},H(x)=\frac{g(x)}{x!}\)
写成生成函数的式子发现就变成了卷积,多项式求逆即可
- 吃一堑涨一智:多项式求逆中的数组别用混
多项式 ln
现在一看,学了一年就是有点长进,至少会了个求导积分
\[G(x)=\ln(F(x)) \mod x^n \]两边同时对复合函数求导得到
\[G'(x)=\ln'(F(x))F'(x) \mod x^n \] \[G'(x)=\frac{F'(x)}{F(x)} \mod x^n \]多项式求逆完了积上去即可
inline void Dao(int *G,int *F,int len){
for(reg int i=0;i<len-1;++i) G[i]=mul(F[i+1],i+1);
return ;
}
inline void Ji(int *G,int *F,int len){
for(reg int i=0;i<len;++i) G[i+1]=mul(F[i],ksm(i+1,mod-2));
return ;
}
inline void Ln(){
P.Dao(G,F,n); P.Inv(F,H,n);
int lim=1; while(lim<(n<<1)) lim<<=1;
P.NTT(G,lim,1); P.NTT(H,lim,1);
for(reg int i=0;i<lim;++i) G[i]=mul(G[i],H[i]); P.NTT(G,lim,-1);
P.Ji(ans,G,n);
return ;
}
多项式 exp
先推式子:\(F(x)=e^{A(x)}\mod x^n\)
\[\ln F(x)-A(x)=0\mod x^n \]那么我们需要求出来的就是 \(G(F(x))=\ln F(x) -A(x)=0\)
这时候有个科技叫牛顿迭代,就是快速求函数零点的一个方法
问题是 \(G(F(x))=0\) 求 \(G(x)\),牛顿迭代的方法是 \(G(x)=G_0(x)-\frac{F(G_0(x))}{F'(G_0(x))}\)
这个式子比较好的地方是精度翻倍,也就是说对 \(\mod x^{\frac n2}\) 意义下做一次,精度就是 \(x^n\)
套回原来的式子,得到:
\[B(x)=B_0(x)(1-\ln B_0(x)+A(x))\mod x^n \] inline void exp(int *a,int *b,int n){
if(n==1) return b[0]=1,void(); exp(a,b,(n+1)>>1); int lim=1; while(lim<=(n<<1)) lim<<=1; Ln(b,Tln,n);
for(reg int i=0;i<n;++i) Tln[i]=del(a[i],Tln[i]); Tln[0]=add(Tln[0],1); for(reg int i=n;i<lim;++i) Tln[i]=0;
NTT(Tln,lim,1); NTT(b,lim,1); for(reg int i=0;i<lim;++i) b[i]=mul(b[i],Tln[i]); NTT(b,lim,-1);
for(reg int i=n;i<lim;++i) b[i]=0; return ;
}
这个板子易错点是 \(\ln\) 和 \(Inv\)
多项式快速幂
求 \(F(x)=G^k(x) \mod x^n\),两边同时 \(\ln\) 则
\[k \ln G(x)=\ln F(x) \mod x^n \]再来一把 \(exp\) 则答案就有了:
\[e^{k \ln G(x)}=F(x)\mod x^n \]又是大码量预警,我写一次 \(exp\) 要写十几分钟
不过对于 \(F[0]\neq 0\) 的情况,平移乱搞即可,注意特判各种情况
付公主的背包
涨教训了,生成函数有封闭形式
那么考虑构造出来封闭形式就是这个式子了
\[\frac{1}{e^{-\sum_{i=1}\ln (1+x^v)}} \]我一开始太年轻了,以为直接 \(\ln\),\(exp\),\(inv\)
但是每个都是 \(\Theta(n\log n)\) 所以不行,接着推式子得到:
\[\ln(1-x^{v_i})=\sum_{i\ge 1} \frac{x^{v_i\times i}}{i} \]直接做就行了
如何卡常?
\((1)\) 预处理单位根
\((2)\) 加法取模优化的函数直接挪到对应位置
\((3)\) 长度的处理