多项式求逆 学习总结

感觉蒟蒻现在学多项式求逆貌似有些晚了

不过真的很有意思了(然而省选的时候自己还在玩泥巴什么也不会

多项式求逆是基于倍增的

假设我们知道

h(x)f(x)=1(mod x^n)

移项得

(h(x)*f(x)-1)=0(mod x^n)

两边同时求平方得

h(x)^2*f(x)^2 - 2*h(x)*f(x) +1 = 0 (mod x^2n)

设g(x)*f(x)=1(mod x^2n)

两边同时乘以g(x)可以得

h(x)^2*f(x) -2*h(x) + g(x) =0 （mod x^2n)

我们移项可以得到

g(x) = h(x) *(2- f(x)*h(x))

倍增即可，时间复杂度O(nlogn)

首先是picks的blog里提到的伯努利数

由于并没有找到相应的题目，于是就在cojs上自己出了一发

疯狂的求和问题 80分

我们知道伯努利数的生成函数是x/(e^x-1)

我们把下面做泰勒展开之后进行多项式求逆即可

#include<cstdio>
#include<cstring>
#include<iostream>
#include<algorithm>
#include<cstdlib>
#define G 3
using namespace std;

typedef long long LL;
const int mod=998244353;
const int maxn=500010;
int n,k,N,len,ans;
int jc[maxn],inv[maxn];
int rev[maxn];
int A[maxn],B[maxn],C[maxn];

void read(int &num){
	num=0;char ch=getchar();
	while(ch<‘!‘)ch=getchar();
	while(ch>=‘0‘&&ch<=‘9‘){
		num=(10LL*num+ch-‘0‘)%mod;
		ch=getchar();
	}return;
}
int pow_mod(int v,int p){
	int tmp=1;
	while(p){
		if(p&1)tmp=1LL*tmp*v%mod;
		v=1LL*v*v%mod;p>>=1;
	}return tmp;
}
void FFT(int *A,int n,int flag){
	for(len=0;(1<<len)<n;++len);
	for(int i=0;i<n;++i){
		rev[i]=rev[i>>1]>>1|((i&1)<<(len-1));
		C[i]=A[rev[i]];
	}
	for(int i=0;i<n;++i)A[i]=C[i];
	for(int k=2;k<=n;k<<=1){
		int mk=(k>>1);
		int wn=pow_mod(G,flag==1?(mod-1)/k:(mod-1)-(mod-1)/k);
		for(int i=0;i<n;i+=k){
			int w=1;
			for(int j=0;j<mk;++j){
				int a=i+j,b=i+j+mk;
				int x=A[a],y=1LL*A[b]*w%mod;
				A[a]=x+y;if(A[a]>=mod)A[a]-=mod;
				A[b]=x-y;if(A[b]<0)A[b]+=mod;
				w=1LL*w*wn%mod;
			}
		}
	}
	if(flag==-1){
		int c=pow_mod(n,mod-2);
		for(int i=0;i<n;++i)A[i]=1LL*A[i]*c%mod;
	}return;
}
void Get_inv(int n){
	if(n==1){
		B[0]=pow_mod(A[0],mod-2);
		return;
	}
	int now=(n>>1);Get_inv(now);
	static int tmp[maxn];
	for(int i=0;i<n;++i)tmp[i]=A[i];
	now=(n<<1);
	FFT(B,now,1);FFT(tmp,now,1);
	for(int i=0;i<now;++i)B[i]=1LL*B[i]*(2LL-1LL*tmp[i]*B[i]%mod+mod)%mod;
	FFT(B,now,-1);
	memset(B+n,0,sizeof(int)*n);
}
int Get_C(int n,int m){return 1LL*jc[n]*inv[m]%mod*inv[n-m]%mod;}

int main(){
	freopen("Crazy_Sum.in","r",stdin);
	freopen("Crazy_Sum.out","w",stdout);
	read(n);scanf("%d",&k);
	for(N=1;N<=k;N<<=1);
	jc[0]=1;
	for(int i=1;i<=N;++i)jc[i]=1LL*jc[i-1]*i%mod;
	inv[N]=pow_mod(jc[N],mod-2);
	for(int i=N-1;i>=0;--i)inv[i]=1LL*inv[i+1]*(i+1)%mod;
	for(int i=0;i<N;++i)A[i]=inv[i+1];
	Get_inv(N);
	for(int i=0;i<N;++i)B[i]=1LL*B[i]*jc[i]%mod;
	int now=n+1;
	for(int i=1;i<=k+1;++i){
		ans=ans+1LL*Get_C(k+1,i)*B[k+1-i]%mod*now%mod;
		if(ans>=mod)ans-=mod;
		now=1LL*now*(n+1)%mod;
	}ans=1LL*ans*pow_mod(k+1,mod-2)%mod;
	printf("%d\n",ans);
	return 0;
}

BZOJ 3456

很经典的题目啦，设fi表示i个点的联通图的个数，考虑1号点所在联通块的大小我们有

fn=2^C(n,2)-sigma(C(n-1,i-1)*fi*2^C(n-i,2))

很容易发现这是个CDQ+FFT的式子，直接上就可以了

#include<cstdio>
#include<cstring>
#include<iostream>
#include<algorithm>
#include<cstdlib>
#define G 3
using namespace std;

const int mod=1004535809;
const int maxn=300010;
int n,N,len;
int jc[maxn],inv[maxn],rev[maxn];
int h[maxn],f[maxn];
int A[maxn],B[maxn],C[maxn];

int pow_mod(int v,int p){
    int tmp=1;
    while(p){
        if(p&1)tmp=1LL*tmp*v%mod;
        v=1LL*v*v%mod;p>>=1;
    }return tmp;
}
void FFT(int *A,int n,int flag){
    for(int i=0;i<n;++i)if(i<rev[i])swap(A[i],A[rev[i]]);
    for(int k=2;k<=n;k<<=1){
        int mk=(k>>1);
        int wn=pow_mod(G,flag==1?(mod-1)/k:(mod-1)-(mod-1)/k);
        for(int i=0;i<n;i+=k){
            int w=1;
            for(int j=0;j<mk;++j){
                int a=i+j,b=i+j+mk;
                int x=A[a],y=1LL*A[b]*w%mod;
                A[a]=x+y;if(A[a]>=mod)A[a]-=mod;
                A[b]=x-y;if(A[b]<0)A[b]+=mod;
                w=1LL*w*wn%mod;
            }
        }
    }
    if(flag==-1){
        int c=pow_mod(n,mod-2);
        for(int i=0;i<n;++i)A[i]=1LL*A[i]*c%mod;
    }return;
}
void Solve(int L,int R){
    if(L==R)return;
    int mid=(L+R)>>1;
    Solve(L,mid);
    for(N=1,len=0;N<=(R-L+1);N<<=1,len++);
    for(int i=0;i<N;++i)A[i]=0,B[i]=h[i];
    for(int i=L;i<=mid;++i)A[i-L]=1LL*f[i]*inv[i-1]%mod;
    if(R-L<=500){
        int lim=R-L+1;
        for(int i=0;i<lim;++i){
            C[i]=0;
            for(int j=0;j<=i;++j){
                C[i]=C[i]+1LL*A[j]*B[i-j]%mod;
                if(C[i]>=mod)C[i]-=mod;
            }
        }
    }else{
        for(int i=0;i<N;++i)rev[i]=rev[i>>1]>>1|((i&1)<<(len-1));
        FFT(A,N,1);FFT(B,N,1);
        for(int i=0;i<N;++i)C[i]=1LL*A[i]*B[i]%mod;
        FFT(C,N,-1);
    }
    for(int i=mid+1;i<=R;++i){
        f[i]=f[i]-1LL*jc[i-1]*C[i-L]%mod;
        if(f[i]<0)f[i]+=mod;
    }
    Solve(mid+1,R);
}

int main(){
    scanf("%d",&n);
    jc[0]=1;
    for(int i=1;i<=n;++i)jc[i]=1LL*jc[i-1]*i%mod;
    inv[n]=pow_mod(jc[n],mod-2);
    for(int i=n-1;i>=0;--i)inv[i]=1LL*inv[i+1]*(i+1)%mod;
    h[1]=1;f[1]=1;
    for(int i=2;i<=n;++i){
        int tmp=(1LL*i*(i-1)/2)%(mod-1);
        h[i]=pow_mod(2,tmp);f[i]=h[i];
        h[i]=1LL*h[i]*inv[i]%mod;
    }
    Solve(1,n);
    printf("%d\n",f[n]);
    return 0;
}

但是这样我们是两个log的

我们注意到可以把上面的式子等价转化成

2^C(n,2)/(n-1)!=sigma( 2^C(n-i,2)/(n-i)! *fi/(i-1)! )

可以构造出h=g*f的形式

之后我们已知h和g，求f，求h*g^(-1)即可,也就是多项式求逆

#include<cstdio>
#include<cstring>
#include<cstdlib>
#include<iostream>
#include<algorithm>
#define G 3
using namespace std;

typedef long long LL;
const int maxn=500010;
const int mod=1004535809;
int n,N,len,t;
int jc[maxn],inv[maxn],rev[maxn];
int A[maxn],B[maxn];

int pow_mod(int v,int p){
    int tmp=1;
    while(p){
        if(p&1)tmp=1LL*tmp*v%mod;
        v=1LL*v*v%mod;p>>=1;
    }return tmp;
}
void FFT(int *A,int n,int flag){
    for(int i=0;i<n;++i){
        if(i<rev[i]){t=A[i];A[i]=A[rev[i]];A[rev[i]]=t;}
    }
    for(int k=2;k<=n;k<<=1){
        int mk=(k>>1);
        int wn=pow_mod(G,flag==1?(mod-1)/k:(mod-1)-(mod-1)/k);
        for(int i=0;i<n;i+=k){
            int w=1;
            for(int j=0;j<mk;++j){
                int a=i+j,b=i+j+mk;
                int x=A[a],y=1LL*A[b]*w%mod;
                A[a]=x+y;if(A[a]>=mod)A[a]-=mod;
                A[b]=x-y;if(A[b]<0)A[b]+=mod;
                w=1LL*w*wn%mod;
            }
        }
    }
    if(flag==-1){
        int c=pow_mod(n,mod-2);
        for(int i=0;i<n;++i)A[i]=1LL*A[i]*c%mod;
    }return;
}
void Get_inv(int n){
    if(n==1){
        B[0]=pow_mod(A[0],mod-2);
        return;
    }
    Get_inv(n>>1);
    int now=(n<<1);
    for(len=0;(1<<len)<now;++len);
    for(int i=0;i<now;++i)rev[i]=rev[i>>1]>>1|((i&1)<<(len-1));
    static int tmp[maxn];
    for(int i=0;i<n;++i)tmp[i]=A[i];
    FFT(tmp,now,1);FFT(B,now,1);
    for(int i=0;i<now;++i)B[i]=1LL*B[i]*(2LL-1LL*tmp[i]*B[i]%mod+mod)%mod;
    FFT(B,now,-1);
    memset(B+n,0,sizeof(int)*n);
}

int main(){
    scanf("%d",&n);
    for(N=1;N<=n;N<<=1);
    jc[0]=1;
    for(int i=1;i<N;++i)jc[i]=1LL*jc[i-1]*i%mod;
    inv[N-1]=pow_mod(jc[N-1],mod-2);
    for(int i=N-2;i>=0;--i)inv[i]=1LL*inv[i+1]*(i+1)%mod;
    for(int i=0;i<N;++i){
        int tmp=(1LL*i*(i-1)/2)%(mod-1);
        A[i]=1LL*pow_mod(2,tmp)*inv[i]%mod;
    }
    Get_inv(N);
    memset(B+n,0,sizeof(int)*n);
    for(len=0;(1<<len)<N;++len);
    for(int i=0;i<N;++i)rev[i]=rev[i>>1]>>1|((i&1)<<(len-1));
    for(int i=0;i<N;++i){
        int tmp=(1LL*i*(i-1)/2)%(mod-1);
        A[i]=1LL*pow_mod(2,tmp)*inv[i-1]%mod;
    }
    for(len=0;(1<<len)<N;++len);
    for(int i=0;i<N;++i)rev[i]=rev[i>>1]>>1|((i&1)<<(len-1));
    FFT(A,N,1);FFT(B,N,1);
    for(int i=0;i<N;++i)A[i]=1LL*A[i]*B[i]%mod;
    FFT(A,N,-1);
    printf("%d\n",1LL*A[n]*jc[n-1]%mod);
    return 0;
}

cojs 疯狂的欧拉图

搞完了城市规划之后自己把原来的考试题扩大了数据范围出到了cojs上(原来n^2可过)

式子还是一样的推导，之后也可以转换成多项式求逆

但是由于模数不兹瓷，所以我们要写模任意质数的FFT QAQ

CDQ+FFT

#include<cstdio>
#include<cstring>
#include<iostream>
#include<cstdlib>
#include<algorithm>
#include<cmath>
#define G 3
using namespace std;

const int maxn=400010;
typedef long long LL;
int n,N,len;
const int mod=1e9+7;
const int M=30000;
const long double pi=acos(-1.0);
int jc[maxn],inv[maxn];
int h[maxn],f[maxn];
int rev[maxn];
int a[maxn],b[maxn],c[maxn];
int a0[maxn],b0[maxn],a1[maxn],b1[maxn];
struct cpx{
    long double r,i;
    cpx(long double r=0,long double i=0):r(r),i(i){}
}A[maxn],B[maxn],C[maxn];
cpx operator +(const cpx &A,const cpx &B){return cpx(A.r+B.r,A.i+B.i);}
cpx operator -(const cpx &A,const cpx &B){return cpx(A.r-B.r,A.i-B.i);}
cpx operator *(const cpx &A,const cpx &B){return cpx(A.r*B.r-A.i*B.i,A.r*B.i+A.i*B.r);}
void FFT(cpx *A,int n,int type){
    for(int i=0;i<n;++i)if(i<rev[i])swap(A[i],A[rev[i]]);
    for(int k=0;(1<<k)<n;++k){
        int m=(1<<k),m2=(m<<1);
        long double o=2*pi/m2*type;
        cpx wn(cos(o),sin(o));
        for(int i=0;i<n;i+=m2){
            cpx w(1,0);
            for(int j=0;j<m;++j){
                cpx x=A[i+j],y=A[i+j+m]*w;
                A[i+j]=x+y;A[i+j+m]=x-y;
                w=w*wn;
            }
        }
    }
    if(type==-1)for(int i=0;i<n;++i)A[i].r/=n;
}
void mul(int *a,int *b,int *c){
    for(int i=0;i<N;++i)A[i]=cpx(a[i],0),B[i]=cpx(b[i],0);
    FFT(A,N,1);FFT(B,N,1);
    for(int i=0;i<N;++i)A[i]=A[i]*B[i];
    FFT(A,N,-1);
    for(int i=0;i<N;++i)c[i]=((LL)(A[i].r+0.5))%mod;
}
void mul_mod(int *a,int *b,int *c){
    for(int i=0;i<N;++i)a0[i]=a[i]/M,b0[i]=b[i]/M;
    mul(a0,b0,a0);
    for(int i=0;i<N;++i){
        c[i]=1LL*a0[i]*M*M%mod;
        a1[i]=a[i]%M;b1[i]=b[i]%M;
    }
    mul(a1,b1,a1);
    for(int i=0;i<N;++i){
        c[i]=c[i]+a1[i];
        if(c[i]>=mod)c[i]-=mod;
        a1[i]=a1[i]+a0[i];
        if(a1[i]>=mod)a1[i]-=mod;
        a0[i]=a[i]/M+a[i]%M;
        b0[i]=b[i]/M+b[i]%M;
    }
    mul(a0,b0,a0);
    for(int i=0;i<N;++i){
        c[i]=c[i]+1LL*M*(a0[i]-a1[i]+mod)%mod;
        if(c[i]>=mod)c[i]-=mod;
    }return;
}
LL pow_mod(LL v,int p){
    LL tmp=1;
    while(p){
        if(p&1)tmp=tmp*v%mod;
        v=v*v%mod;p>>=1;
    }return tmp;
}
void Solve(int L,int R){
    if(L==R)return;
    int mid=(L+R)>>1;
    Solve(L,mid);
    for(N=1,len=0;N<(R-L+1);N<<=1,len++);
    for(int i=0;i<N;++i)a[i]=b[i]=0;
    for(int i=L;i<=mid;++i)a[i-L]=1LL*f[i]*inv[i-1]%mod;
    for(int i=0;i<N;++i)b[i]=h[i];
    if(R-L<=1000){
    	int lim=R-L+1;
    	for(int i=0;i<lim;++i){
    		c[i]=0;
    		for(int j=0;j<=i;++j){
    			c[i]=c[i]+1LL*a[j]*b[i-j]%mod;
    			if(c[i]>=mod)c[i]-=mod;
    		}
    	}
    }else{
    	for(int i=0;i<N;++i)rev[i]=rev[i>>1]>>1|((i&1)<<(len-1));
    	mul_mod(a,b,c);
    }
    for(int i=mid+1;i<=R;++i){
        f[i]=f[i]-1LL*jc[i-1]*c[i-L]%mod;
        if(f[i]<0)f[i]+=mod;
    }
    Solve(mid+1,R);
}

int main(){
	freopen("Crazy_Graph.in","r",stdin);
	freopen("Crazy_Graph.out","w",stdout);
    scanf("%d",&n);
    jc[0]=1;
    for(int i=1;i<=n;++i)jc[i]=1LL*jc[i-1]*i%mod;
    inv[n]=pow_mod(jc[n],mod-2);
    for(int i=n-1;i>=0;--i)inv[i]=1LL*inv[i+1]*(i+1)%mod;
    h[1]=1;f[1]=1;
    for(int i=2;i<=n;++i){
        h[i]=pow_mod(2LL,(1LL*(i-1)*(i-2)/2)%(mod-1));
        f[i]=h[i];
        h[i]=1LL*h[i]*inv[i]%mod;
    }
    Solve(1,n);
    f[n]=f[n]+1LL*n*(n-1)*pow_mod(2LL,mod-2)%mod*f[n]%mod;
    if(f[n]>=mod)f[n]-=mod;
    printf("%d\n",f[n]);
    return 0;
}

多项式求逆

#include<cstdio>
#include<cstring>
#include<iostream>
#include<algorithm>
#include<cstdlib>
#include<cmath>
using namespace std;

typedef long long LL;
const int mod=1e9+7;
const int M=32000;
const int maxn=400010;
const long double pi=acos(-1.0);
int n,len,N;
int jc[maxn],inv[maxn],rev[maxn];
int a[maxn],b[maxn],c[maxn];
int a0[maxn],b0[maxn],a1[maxn],b1[maxn];
struct cpx{
	long double r,i;
	cpx(long double r=0,long double i=0):r(r),i(i){}
}A[maxn],B[maxn];
cpx operator +(const cpx &A,const cpx &B){return cpx(A.r+B.r,A.i+B.i);}
cpx operator -(const cpx &A,const cpx &B){return cpx(A.r-B.r,A.i-B.i);}
cpx operator *(const cpx &A,const cpx &B){return cpx(A.r*B.r-A.i*B.i,A.r*B.i+A.i*B.r);}

int pow_mod(int v,int p){
	int tmp=1;
	while(p){
		if(p&1)tmp=1LL*tmp*v%mod;
		v=1LL*v*v%mod;p>>=1;
	}return tmp;
}
void FFT(cpx *A,int n,int flag){
	for(int i=0;i<n;++i)if(i<rev[i])swap(A[i],A[rev[i]]);
	for(int k=2;k<=n;k<<=1){
		int mk=(k>>1);
		long double o=2*pi/k*flag;
		cpx wn(cos(o),sin(o));
		for(int i=0;i<n;i+=k){
			cpx w(1,0);
			for(int j=0;j<mk;++j){
				int a=i+j,b=i+j+mk;
				cpx x=A[a],y=A[b]*w;
				A[a]=x+y;A[b]=x-y;
				w=w*wn;
			}
		}
	}
	if(flag==-1){for(int i=0;i<n;++i)A[i].r/=n;}
}
void mul(int *a,int *b,int *c,int N){
	for(int i=0;i<N;++i)A[i]=cpx(a[i],0),B[i]=cpx(b[i],0);
	FFT(A,N,1);FFT(B,N,1);
	for(int i=0;i<N;++i)A[i]=A[i]*B[i];
	FFT(A,N,-1);
	for(int i=0;i<N;++i)c[i]=(LL)(A[i].r+0.5)%mod;
}
void mul_mod(int *A,int *B,int *C,int N){
	for(int i=0;i<N;++i)a0[i]=A[i]%M,b0[i]=B[i]%M;
	mul(a0,b0,a0,N);
	for(int i=0;i<N;++i)C[i]=a0[i],a1[i]=A[i]/M,b1[i]=B[i]/M;
	mul(a1,b1,a1,N);
	for(int i=0;i<N;++i){
		C[i]=C[i]+1LL*M*M*a1[i]%mod;
		if(C[i]>=mod)C[i]-=mod;
		a1[i]=a1[i]+a0[i];
		if(a1[i]>=mod)a1[i]-=mod;
		a0[i]=A[i]/M+A[i]%M;
		b0[i]=B[i]/M+B[i]%M;
	}
	mul(a0,b0,a0,N);
	for(int i=0;i<N;++i){
		C[i]=C[i]+1LL*M*(a0[i]-a1[i]+mod)%mod;
		if(C[i]>=mod)C[i]-=mod;
	}return;
}
void Get_inv(int n){
	if(n==1){
		b[0]=pow_mod(a[0],mod-2);
		return;
	}
	Get_inv(n>>1);
	int now=(n<<1);
	for(len=0;(1<<len)<now;++len);
	for(int i=0;i<now;++i)rev[i]=rev[i>>1]>>1|((i&1)<<(len-1));
	static int tmp[maxn];
	for(int i=0;i<n;++i)tmp[i]=a[i];
	mul_mod(b,tmp,c,now);
	c[0]=2-c[0];
	if(c[0]<0)c[0]+=mod;
	for(int i=1;i<now;++i)c[i]=mod-c[i];
	mul_mod(b,c,tmp,now);
	for(int i=0;i<n;++i)b[i]=tmp[i];
	memset(b+n,0,sizeof(int)*n);
}

int main(){
	freopen("Crazy_Graph.in","r",stdin);
	freopen("Crazy_Graph.out","w",stdout);
	scanf("%d",&n);
	for(N=1;N<=n;N<<=1);
	jc[0]=1;
	for(int i=1;i<N;++i)jc[i]=1LL*jc[i-1]*i%mod;
	inv[N-1]=pow_mod(jc[N-1],mod-2);
	for(int i=N-2;i>=0;--i)inv[i]=1LL*inv[i+1]*(i+1)%mod;
	a[0]=1;
	for(int i=1;i<N;++i){
		int tmp=(1LL*(i-1)*(i-2)/2)%(mod-1);
		a[i]=1LL*pow_mod(2,tmp)*inv[i]%mod;
	}
	Get_inv(N);
	memset(b+n,0,sizeof(int)*n);
	for(len=0;(1<<len)<N;++len);
	for(int i=0;i<N;++i)rev[i]=rev[i>>1]>>1|((i&1)<<(len-1));
	a[0]=1;
	for(int i=1;i<N;++i){
		int tmp=(1LL*(i-1)*(i-2)/2)%(mod-1);
		a[i]=1LL*pow_mod(2,tmp)*inv[i-1]%mod;
	}
	mul_mod(a,b,c,N);
	int ans=1LL*jc[n-1]*c[n]%mod;
	ans=ans+1LL*n*(n-1)*pow_mod(2LL,mod-2)%mod*ans%mod;
	if(ans>=mod)ans-=mod;
	printf("%d\n",ans);
	return 0;
}

cojs 异化多肽

貌似是多项式求逆的裸题，Nescafe的题目

构造生成函数f，不难发现我们要求的是f^1+f^2+……

之后求和得 1/(1-f)

多项式求逆即可

#include<cstdio>
#include<cstring>
#include<iostream>
#include<cstdlib>
#include<algorithm>
#define G 5
using namespace std;

typedef long long LL;
const int maxn=300010;
const int mod=1005060097;
int n,m,k,N,L;
int A[maxn],B[maxn];
int rev[maxn];

void read(int &num){
	num=0;char ch=getchar();
	while(ch<‘!‘)ch=getchar();
	while(ch>=‘0‘&&ch<=‘9‘)num=num*10+ch-‘0‘,ch=getchar();
}
int pow_mod(int v,int p){
	int tmp=1;
	while(p){
		if(p&1)tmp=1LL*tmp*v%mod;
		v=1LL*v*v%mod;p>>=1;
	}return tmp;
}
void FFT(int *A,int n,int flag){
	for(L=0;(1<<L)<n;++L);
	for(int i=0;i<n;++i)rev[i]=rev[i>>1]>>1|((i&1)<<(L-1));
	for(int i=0;i<n;++i)if(i<rev[i])swap(A[i],A[rev[i]]);
	for(int k=2;k<=n;k<<=1){
		int mk=(k>>1);
		int wn=pow_mod(G,flag==1?(mod-1)/k:(mod-1)-(mod-1)/k);
		for(int i=0;i<n;i+=k){
			int w=1;
			for(int j=0;j<mk;++j){
				int x=A[i+j],y=1LL*A[i+j+mk]*w%mod;
				A[i+j]=x+y;if(A[i+j]>=mod)A[i+j]-=mod;
				A[i+j+mk]=x-y;if(A[i+j+mk]<0)A[i+j+mk]+=mod;
				w=1LL*w*wn%mod;
			}
		}
	}
	if(flag==-1){
		int inv=pow_mod(n,mod-2);
		for(int i=0;i<n;++i)A[i]=1LL*A[i]*inv%mod;
	}return;
}
void Get_inv(int n){
	if(n==1){
		B[0]=pow_mod(A[0],mod-2);
		return;
	}
	int now=(n>>1);
	Get_inv(now);
	static int tmp[maxn];
	for(int i=0;i<now;++i)tmp[i]=A[i];
	FFT(tmp,n,1);FFT(B,n,1);
	for(int i=0;i<n;++i)B[i]=1LL*B[i]*(2-1LL*tmp[i]*B[i]%mod+mod)%mod;
	FFT(B,n,-1);
	memset(B+now,0,sizeof(int)*now);
}

int main(){
	freopen("polypeptide.in","r",stdin);
	freopen("polypeptide.out","w",stdout);
	read(n);read(m);
	for(int i=1;i<=m;++i){
		read(k);
		if(k<=n)A[k]++;
	}A[0]--;if(A[0]<0)A[0]+=mod;
	for(N=1;N<=n;N<<=1);N<<=1;
	Get_inv(N);
	int ans=-B[n];
	ans=(ans%mod+mod)%mod;
	printf("%d\n",ans);
	return 0;
}

HEOI 2016 求和

QAQ 当前考场上没推出式子来

其实当时并不是很会FFT，现在也不是很会

我们不难发现如果没有2^j的话求的是n个数划分成若干个有序集合的方案数

设fn为这个方案数，考虑枚举第一个集合的大小

我们有fn=sigma( C(n,i) * f(n-i) )

考虑2^j实际上是每个集合都要*2

那么我们的递推式只要改成 fn=sigma( C(n,i) * f(n-i) *2)就可以了

然后上CDQ+FFT

#include<cstdio>
#include<cstring>
#include<cstdlib>
#include<iostream>
#include<algorithm>
#define G 3
using namespace std;

typedef long long LL;
const int maxn=500010;
const int mod=998244353;
int n,N,len,ans;
int jc[maxn],inv[maxn],rev[maxn];
int f[maxn],A[maxn],B[maxn],C[maxn];

int pow_mod(int v,int p){
	int tmp=1;
	while(p){
		if(p&1)tmp=1LL*tmp*v%mod;
		v=1LL*v*v%mod;p>>=1;
	}return tmp;
}
void FFT(int *A,int n,int flag){
	for(int i=0;i<n;++i)if(i<rev[i])swap(A[i],A[rev[i]]);
	for(int k=2;k<=n;k<<=1){
		int mk=(k>>1);
		int wn=pow_mod(G,flag==1?(mod-1)/k:(mod-1)-(mod-1)/k);
		for(int i=0;i<n;i+=k){
			int w=1;
			for(int j=0;j<mk;++j){
				int a=i+j,b=i+j+mk;
				int x=A[a],y=1LL*A[b]*w%mod;
				A[a]=x+y;if(A[a]>=mod)A[a]-=mod;
				A[b]=x-y;if(A[b]<0)A[b]+=mod;
				w=1LL*w*wn%mod;
			}
		}
	}
	if(flag==-1){
		int c=pow_mod(n,mod-2);
		for(int i=0;i<n;++i)A[i]=1LL*A[i]*c%mod;
	}return;
}
void Solve(int L,int R){
	if(L==R)return;
	int mid=(L+R)>>1;
	Solve(L,mid);
	for(N=1,len=0;N<=(R-L+1);N<<=1,len++);
	for(int i=0;i<N;++i)A[i]=0,B[i]=inv[i];
	for(int i=L;i<=mid;++i)A[i-L]=f[i];
	if(R-L<=500){
		int lim=R-L+1;
		for(int i=0;i<lim;++i){
			C[i]=0;
			for(int j=0;j<=i;++j){
				C[i]=C[i]+1LL*A[j]*B[i-j]%mod;
				if(C[i]>=mod)C[i]-=mod;
			}
		}
	}else{
		for(int i=0;i<N;++i)rev[i]=rev[i>>1]>>1|((i&1)<<(len-1));
		FFT(A,N,1);FFT(B,N,1);
		for(int i=0;i<N;++i)C[i]=1LL*A[i]*B[i]%mod;
		FFT(C,N,-1);
	}
	for(int i=mid+1;i<=R;++i){
		f[i]=f[i]+2LL*C[i-L]%mod;
		if(f[i]>=mod)f[i]-=mod;
	}
	Solve(mid+1,R);
}

int main(){
	freopen("heoi2016_sum.in","r",stdin);
	freopen("heoi2016_sum.out","w",stdout);
	scanf("%d",&n);
	jc[0]=1;
	for(int i=1;i<=n;++i)jc[i]=1LL*jc[i-1]*i%mod;
	inv[n]=pow_mod(jc[n],mod-2);f[0]=1;
	for(int i=n-1;i>=0;--i)inv[i]=1LL*inv[i+1]*(i+1)%mod;
	Solve(0,n);
	for(int i=0;i<=n;++i){
		ans=ans+1LL*jc[i]*f[i]%mod;
		if(ans>=mod)ans-=mod;
	}printf("%d\n",ans);
	return 0;
}

之后我们考虑化简一下式子

我们有fn-sigma(C(n,i)*fi*2)=0

之后构造多项式可以得到f - g*f =1 (因为f0=1)

则f = 1/(1-g)

多项式求逆即可

#include<cstdio>
#include<cstring>
#include<iostream>
#include<cstdlib>
#include<algorithm>
#define G 3
using namespace std;

typedef long long LL;
const int mod=998244353;
const int maxn=300010;
int n,N,len,ans=0;
int jc[maxn],inv[maxn],rev[maxn];
int A[maxn],B[maxn];

int pow_mod(int v,int p){
	int tmp=1;
	while(p){
		if(p&1)tmp=1LL*tmp*v%mod;
		v=1LL*v*v%mod;p>>=1;
	}return tmp;
}
void FFT(int *A,int n,int flag){
	for(int i=0;i<n;++i)if(i<rev[i])swap(A[i],A[rev[i]]);
	for(int k=2;k<=n;k<<=1){
		int mk=(k>>1);
		int wn=pow_mod(G,flag==1?(mod-1)/k:(mod-1)-(mod-1)/k);
		for(int i=0;i<n;i+=k){
			int w=1;
			for(int j=0;j<mk;++j){
				int a=i+j,b=i+j+mk;
				int x=A[a],y=1LL*A[b]*w%mod;
				A[a]=x+y;if(A[a]>=mod)A[a]-=mod;
				A[b]=x-y;if(A[b]<0)A[b]+=mod;
				w=1LL*w*wn%mod;
			}
		}
	}
	if(flag==-1){
		int c=pow_mod(n,mod-2);
		for(int i=0;i<n;++i)A[i]=1LL*A[i]*c%mod;
	}return;
}
void Get_inv(int n){
	if(n==1){
		B[0]=pow_mod(A[0],mod-2);
		return;
	}
	Get_inv(n>>1);
	int now=(n<<1);
	for(len=0;(1<<len)<now;++len);
	for(int i=0;i<now;++i)rev[i]=rev[i>>1]>>1|((i&1)<<(len-1));
	static int tmp[maxn];
	for(int i=0;i<n;++i)tmp[i]=A[i];
	FFT(tmp,now,1);FFT(B,now,1);
	for(int i=0;i<now;++i)B[i]=1LL*B[i]*(2LL-1LL*tmp[i]*B[i]%mod+mod)%mod;
	FFT(B,now,-1);
	memset(B+n,0,sizeof(int)*n);
}

int main(){
	freopen("heoi2016_sum.in","r",stdin);
	freopen("heoi2016_sum.out","w",stdout);
	scanf("%d",&n);
	for(N=1;N<=n;N<<=1);
	jc[0]=1;
	for(int i=1;i<N;++i)jc[i]=1LL*jc[i-1]*i%mod;
	inv[N-1]=pow_mod(jc[N-1],mod-2);
	for(int i=N-2;i>=0;--i)inv[i]=1LL*inv[i+1]*(i+1)%mod;
	for(int i=1;i<N;++i){
		A[i]=(inv[i]<<1);
		if(A[i]>=mod)A[i]-=mod;
	}A[0]--;if(A[0]<0)A[0]+=mod;
	Get_inv(N);
	for(int i=0;i<=n;++i){
		ans=ans+1LL*(mod-B[i])*jc[i]%mod;
		if(ans>=mod)ans-=mod;
	}printf("%d\n",ans);
	return 0;
}

另外，求贝尔数的CDQ+FFT的程序

#include<cstdio>
#include<cstring>
#include<iostream>
#include<algorithm>
#include<cstdlib>
#define G 11
using namespace std;

typedef long long LL;
const int maxn=300010;
const int mod=786433;
int n,N,len;
int jc[maxn],inv[maxn],rev[maxn];
int f[maxn],A[maxn],B[maxn],C[maxn];

int pow_mod(int v,int p){
	int tmp=1;
	while(p){
		if(p&1)tmp=1LL*tmp*v%mod;
		v=1LL*v*v%mod;p>>=1;
	}return tmp;
}
void FFT(int *A,int n,int flag){
	for(int i=0;i<n;++i)if(i<rev[i])swap(A[i],A[rev[i]]);
	for(int k=2;k<=n;k<<=1){
		int mk=(k>>1);
		int wn=pow_mod(G,flag==1?(mod-1)/k:(mod-1)-(mod-1)/k);
		for(int i=0;i<n;i+=k){
			int w=1;
			for(int j=0;j<mk;++j){
				int a=i+j,b=i+j+mk;
				int x=A[a],y=1LL*A[b]*w%mod;
				A[a]=x+y;if(A[a]>=mod)A[a]-=mod;
				A[b]=x-y;if(A[b]<0)A[b]+=mod;
				w=1LL*w*wn%mod;
			}
		}
	}
	if(flag==-1){
		int c=pow_mod(n,mod-2);
		for(int i=0;i<n;++i)A[i]=1LL*A[i]*c%mod;
	}return;
}
void Solve(int L,int R){
	if(L==R)return;
	int mid=(L+R)>>1;
	Solve(L,mid);
	for(N=1,len=0;N<=(R-L+1);N<<=1,len++);
	A[0]=B[0]=0;
	for(int i=1;i<N;++i)A[i]=0,B[i]=inv[i-1];
	for(int i=L;i<=mid;++i)A[i-L]=1LL*f[i]*inv[i]%mod;
	if(R-L<=500){
		int lim=R-L+1;
		for(int i=0;i<lim;++i){
			C[i]=0;
			for(int j=0;j<=i;++j){
				C[i]=C[i]+1LL*A[j]*B[i-j]%mod;
				if(C[i]>=mod)C[i]-=mod;
			}
		}
	}else{
		for(int i=0;i<N;++i)rev[i]=rev[i>>1]>>1|((i&1)<<(len-1));
		FFT(A,N,1);FFT(B,N,1);
		for(int i=0;i<N;++i)C[i]=1LL*A[i]*B[i]%mod;
		FFT(C,N,-1);
	}
	for(int i=mid+1;i<=R;++i){
		f[i]=f[i]+1LL*jc[i-1]*C[i-L]%mod;
		if(f[i]>=mod)f[i]-=mod;
	}Solve(mid+1,R);
}

int main(){
	scanf("%d",&n);
	jc[0]=1;
	for(int i=1;i<=n;++i)jc[i]=1LL*jc[i-1]*i%mod;
	inv[n]=pow_mod(jc[n],mod-2);
	for(int i=n-1;i>=0;--i)inv[i]=1LL*inv[i+1]*(i+1)%mod;
	f[0]=1;Solve(0,n);
	printf("%d\n",f[n]);
	return 0;
}

但是我并不知道怎么转成多项式求逆，如果有老司机知道还请告诉我QAQ

留下的一些坑：

1、城市规划的多项式求ln的写法

2、多项式开根

3、贝尔数的多项式求exp的写法

感觉很多CDQ+FFT都能很轻松的转成多项式求逆啊

虽然转了之后并没有快多少

时间： 2024-08-14 06:54:21

多项式求逆 学习总结的相关文章

NTT+多项式求逆+多项式开方(BZOJ3625)

定义多项式h(x)的每一项系数hi,为i在c[1]~c[n]中的出现次数. 定义多项式f(x)的每一项系数fi,为权值为i的方案数. 通过简单的分析我们可以发现:f(x)=2/(sqrt(1-4h(x))+1) 于是我们需要多项式开方和多项式求逆. 多项式求逆: 求B(x),使得A(x)*B(x)=1 (mod x^m) 考虑倍增. 假设我们已知A(x)*B(x)=1 (mod x^m),要求C(x),使得A(x)*C(x)=1 (mod x^(2m)) 简单分析可得C(x)=B(x)*(2-A

【BZOJ 3456】 3456: 城市规划 （NTT+多项式求逆）

3456: 城市规划 Time Limit: 40 Sec  Memory Limit: 256 MBSubmit: 658  Solved: 364 Description 刚刚解决完电力网络的问题, 阿狸又被领导的任务给难住了. 刚才说过, 阿狸的国家有n个城市, 现在国家需要在某些城市对之间建立一些贸易路线, 使得整个国家的任意两个城市都直接或间接的连通. 为了省钱, 每两个城市之间最多只能有一条直接的贸易路径. 对于两个建立路线的方案, 如果存在一个城市对, 在两个方案中是否建立路线不一

Nescafe41 ProblemA 异化多肽 多项式求逆

题目大意:. 思路: 搞出C的生成函数F(x),那么: 长度为1的答案为F(x) 长度为2的答案为F2(x) - 故最终的答案为 F(x)+F2(x)+F3(x)+... =1?F+∞(x)1?F(x) =11?F(x) 然后就是多项式求逆了= = 跪picks大毒瘤 #include <cstdio> #include <cstring> #include <iostream> #include <algorithm> #define M (263000&

【BZOJ】4555: [Tjoi2016&amp;Heoi2016]求和 排列组合+多项式求逆 或 斯特林数+NTT

[题意]给定n,求Σi=0~nΣj=1~i s(i,j)*2^j*j!,n<=10^5. [算法]生成函数+排列组合+多项式求逆 [题解]参考: [BZOJ4555][Tjoi2016&Heoi2016]求和-NTT-多项式求逆 $ans=\sum_{i=0}^{n}\sum_{j=0}^{i}s(i,j)*2^j*j!$ 令$g(n)=\sum_{j=0}^{n}s(n,j)*2^j*j!$ 则ans是Σg(i),只要计算出g(i)的生成函数就可以统计答案. g(n)可以理解为将n个数划分

【bzoj3456】城市规划 容斥原理+NTT+多项式求逆

题目描述 求出n个点的简单(无重边无自环)无向连通图数目mod 1004535809(479 * 2 ^ 21 + 1). 输入 仅一行一个整数n(<=130000) 输出 仅一行一个整数, 为方案数 mod 1004535809. 样例输入 3 样例输出 4 题解 容斥原理+NTT+多项式求逆 设 $f_i$ 表示 $i$ 个点的简单无向连通图的数目,$g_i$ 表示 $i$ 个点的简单无向图的数目. 根据定义得 $g_i=2^{\frac{n(n-1}2}$ . 对于 $f_i$ ,考虑容斥

多项式求逆

多项式求逆 求 \(A(x)\) 在 \(\%x^{n}\) 意义下的逆元 \(B(x)\) 首先求出 \(A(x)\) 在 \(\%x^{\lceil \frac{n}{2} \rceil}\) 意义下的逆元 \(C(x)\),即 $A(x)C(x)=1 $ \((\%x^{\lceil \frac{n}{2} \rceil})\) 移项得 \(A(x)C(x)-1 = 0\) \((\%x^{\lceil \frac{n}{2} \rceil})\) 两边平方 \(A^{2}(x)C^{2}

[Codeforces438E][bzoj3625] 小朋友和二叉树 [多项式求逆+多项式开根]

题面 传送门 思路 首先,我们把这个输入的点的生成函数搞出来: \(C=\sum_{i=0}^{lim}s_ix^i\) 其中\(lim\)为集合里面出现过的最大的数,\(s_i\)表示大小为\(i\)的数是否出现过 我们再设另外一个函数\(F\),定义\(F_k\)表示总权值为\(k\)的二叉树个数 那么,一个二叉树显然可以通过两个子树(可以权值为0,也就是空子树)和一个节点构成 那么有如下求\(F\)的式子 \(F_0=1\) \(F_k=\sum_{i=0}^k s_i \sum_{j=0

【bzoj3456】城市规划（多项式求逆+dp）

Description 求\(~n~\)个点组成的有标号无向连通图的个数.\(~1 \leq n \leq 13 \times 10 ^ 4~\). Solution 这道题的弱化版是poj1737, 其中\(n \leq 50\), 先来解决这个弱化版的题.考虑\(~dp~\),直接统计答案难以入手,于是考虑容斥.显然有,符合条件的方案数\(=\)所有方案数\(-\)不符合条件的方案数,而这个不符合条件的方案数就是图没有完全联通的情况.设\(~dp_i~\)表示\(~i~\)个点组成的合法方案

CF438E The Child and Binary Tree（生成函数+多项式开根+多项式求逆）

传送门 可以……这很多项式开根模板……而且也完全不知道大佬们怎么把这题的式子推出来的…… 首先,这题需要多项式开根和多项式求逆.多项式求逆看这里->这里,这里讲一讲多项式开根 多项式开方:已知多项式$A$,求多项式$B$满足$A^2\equiv B\pmod{x^n}$(和多项式求逆一样这里需要取模,否则$A$可能会有无数项) 假设我们已经求出$A'^2\equiv B\pmod{x^n}$,考虑如何计算出$A^2\equiv B\pmod{x^{2n}}$ 首先肯定存在$A^2\equiv B