多项式求逆 学习总结

最新推荐文章于 2023-06-07 15:17:29 发布
最新推荐文章于 2023-06-07 15:17:29 发布
阅读量357 收藏 1
点赞数
原文链接：http://www.cnblogs.com/joyouth/p/5587648.html
版权

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

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

 

多项式求逆是基于倍增的

假设我们知道

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都能很轻松的转成多项式求逆啊

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

 

转载于:https://www.cnblogs.com/joyouth/p/5587648.html

weixin_30687051
关注
  • 0
    点赞
  • 踩
  • 1
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
shuzhijisuan.rar_多项式求根
09-19
包括 ：定积分计算、多项式求根(牛顿法)周期性方程(追赶法)等三种算法
多项式求逆
a_forever_dream的博客
10-10 7718
介绍 准确来说，是求多项式的逆元。 对于多项式 A(x)A(x)A(x)，如果存在 A(x)B(x)≡1(modxn)A(x)B(x)\equiv 1 \pmod {x^n}A(x)B(x)≡1(modxn)，那么称 B(x)B(x)B(x) 为 A(x)A(x)A(x) 在模 xnx^nxn 意义下的逆元。 注意，这里的模 xnx^nxn，是指舍弃含有 xnx^nxn 及更高次的项。 比如说有一...
【学习笔记】超简单的多项式求逆（含全套证明）
繁凡さん的博客
12-30 3131
整理的算法模板合集： ACM模板 目录多项式求逆一、分治FFT二、倍增法及其证明三、多项式求逆例题[P4238 【模板】多项式乘法逆](https://www.luogu.com.cn/problem/P4238) 点我看多项式全家桶(●^◡_◡◡​^●) 多项式求逆 一、分治FFT 二、倍增法及其证明 三、多项式求逆例题 P4238 【模板】多项式乘法逆 直接按照上述的思路实现一下模板即可。注意要中间 h(x)h(x)h(x) 要取的是模 x⌈n⌉2x^{\frac{\lceil n \rceil}
多项式算法5：多项式求逆
ha_ing的博客
06-03 2363
多项式算法5：多项式求逆
逆元模板
弱智就要多努力
10-10 593
扩展欧几里得求逆元：int extgcd(int a, int b, int &x, int &y) { int d = a; if(b != 0) { d = extgcd(b, a%b, y, x); y -= (a / b) * x; } else x = 1, y = 0; return d; } int mod
两类矩阵多项式的逆矩阵的求法 (2013年)
06-19
矩阵的求逆是矩阵论中研究的重要问题,尤其是一些矩阵多项式的求逆问题。在求矩阵多项式的逆矩阵过程中,研究发现一些特殊矩阵多项式与其逆...对丰富矩阵多项式的求逆理论具有重要意义,对学习求逆知识也具有借鉴作用。
jisuanqi.rar_求多项式系数
09-23
求两个一元多项式A（X）=a0+a1x+a2x^2+....+anx^n和B（X）=b0+b1x+b2x^2+....+bnx^n 的和C（X），差D（X），积（X）。 算法的输入：两个多项式（分别输入各项的系数和指数） 算法的输出：多项式的和,差，积。
Matlab多项式运算与方程求根
08-29
Matlab多项式运算与方程求根对各种多项式运算给出实列和方法
机器学习多项式拟合
01-11
哈工大机器学习多项式拟合，使用MATLAB实现，有图像对比
扩展欧几里得算法的matlab程序（求多项式的乘法逆元）
06-23
matlab的M函数文件，附带了函数的使用说明了
多项式求逆 与 多项式开方
XianHaoMing的博客
12-09 1169
介绍与基本概念多项式求逆元是多项式除法、多项式开根和多项式取模的必要过程，用快速傅里叶变换及倍增算法可以做到用OO(nn loglog nn)的时间复杂度求出一个多项式的逆元。对于一个多项式 A(x)A{(x)}，称其最高项的次数为这个多项式的度（degreedegree），记作 degAdeg_A 对于多项式 A(x)A{(x)},B(x)B{(x)}，存在唯一的Q(x)Q{(x)},R(x)R
多项式算法2：多项式求逆
alpha202的博客
06-02 1499
已知多项式F(x)F(x)F(x)，求满足条件的G(x)G(x)G(x)使得 F(x)G(x)≡1(modxn)F(x)G(x)\equiv1\pmod{x^n}F(x)G(x)≡1(modxn) 考虑迭代求解 令 F(x)G0(x)≡1 (modx⌈n2⌉)F(x)G_0(x)\equiv1\ \pmod{x^{\left\lceil\frac{n}{2}\right\rceil}}F...
GPU加速NTT算法
weixin_45825274的博客
06-07 934
快速数论变换-NTT算法：在离散数值域上进行快速多项式乘法运算的高效算法。NTT的理论推导可参考下面列出的两个视频。关于NTT的作用，简单来讲，有一个n次多项式A和一个m次多项式B，A和B相乘后得到一个m+n次多项式C，给定A的n+1个系数和B的m+1个系数求C的m+n+1个系数。这是一个完整的过程，需要对A和B分别进行NTT变换（得到两组长度为m+n+1的点值），然后把两组点值对应相乘得到C的m+n+1个点值，再对C的m+n+1个点值进行NTT的逆变换才能得到C的每一个系数。
【学习笔记】多项式全家桶（包含全套证明）
热门推荐
繁凡さん的博客
12-29 1万+
目录0x00.多项式基本概述0x01.多项式的基础概念0x02. 多项式的度0x03. 多项式的乘法0x04. 多项式的逆元0x10.拉格朗日插值0x20.多项式乘法0x21.快速傅里叶变换（FFT）0x22.快速数论变换（NTT）0x23.快速沃尔什变换（FWT）0x24.快速莫比乌斯变换（FMT）0x30.多项式求逆0x40.多项式开方0x50.多项式除法 / 取模0x60.多项式牛顿迭代0x70.多项式多点求值 / 快速插值0x80.多项式三角函数0x90.多项式反...
多项式求逆与多项式开根
某蒟蒻的博客
02-13 1098
闲着没事干研究些黑科技（雾） 多项式求逆 求 A(x)*B(x)==1(mod x^n) 其中n为A(x)，B(x)的度的较大值 已知A(x)求B(x)，B(x)=A(x)^(-1)(mod x^n) 假设n=1，则B(x)=A(x)常数项在mod p 意义下的的逆元 假设n&gt;1 已知 A(x)*B’(x)==1(mod x^(n/2)) A(x)*B(x)==1(mod...
python 多项式求根
最新发布
01-17
多项式求根是一个常见的数值计算问题。在Python中，有多种方法可以用来求解多项式的根。其中一种常用的方法是使用YRoots包。YRoots是一个用于数值求根的Python包，可以用来求解多元多项式的根。 以下是一个使用...

“相关推荐”对你有帮助么？

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论
添加红包

请填写红包祝福语或标题

个

红包个数最小为10个

元

红包金额最低5元

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

抵扣说明：

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

余额充值