多项式求卡常/wq
查看原帖
多项式求卡常/wq
482049
Alex_wcq楼主2023/8/18 22:20

在这里能过,但跑得很慢,要 900+ ms;在其他OJ上跑得更慢,过不了。所以求帮忙卡常/wq

#include<bits/stdc++.h>
using namespace std;
#define mod 167772161
#define G 3
#define N (1<<21)
#define poly vector<int>
#define minus MMinus
int qpow(int x,int y){
	if(y==0) return 1;
	int t=qpow(x,y/2);
	if(y%2==0) return 1ll*t*t%mod;
	return 1ll*t*t%mod*x%mod;
}
int dim(poly x){
	return x.size()-1;
}
void write(poly x){
	int n=dim(x);
	for(int i=0;i<=n;++i) printf("%d ",x[i]);
	puts("");
}
poly read(int n){
	poly ret; ret.resize(n+1);
	for(int i=0;i<=n;++i) scanf("%d",&ret[i]);
	return ret;
}
poly read(){
	int n; scanf("%d",&n);
	return read(n);
}
int A[N+5],B[N+5],r[N+5];
void ntt(int *x, int lim, int opt) {
  	for(int i=0;i<lim;++i)
    	if(r[i]<i) swap(x[i],x[r[i]]);
  	for(int m=2;m<=lim;m*=2) {
    	int k=m/2;
    	int gn=qpow(G,(mod-1)/m);
    	for(int i=0;i<lim;i+=m){
      		int g=1;
      		for(int j=0;j<k;++j){
        		int tmp=1ll*x[i+j+k]*g%mod;
        		x[i+j+k]=(x[i+j]-tmp+mod)%mod;
       		 	x[i+j]=(x[i+j]+tmp)%mod;
       		 	g=1ll*g*gn%mod;
      		}
    	}
  	}
  	if(opt==-1){
    	reverse(x+1,x+lim);
    	int inv=qpow(lim,mod-2);
    	for(int i=0;i<lim;++i) 
    		x[i]=1ll*x[i]*inv%mod;
  	}
}
poly multi(poly a,poly b){
	int d1=dim(a),d2=dim(b),lim=1;
	r[0]=0;
  	while(lim<d1+d2+1) lim*=2;
  	for(int i=0;i<lim;++i) r[i]=(i&1)*(lim/2)+(r[i/2]/2);
	for(int i=0;i<lim;++i) A[i]=B[i]=0;
	for(int i=0;i<=d1;++i) A[i]=a[i];
	for(int i=0;i<=d2;++i) B[i]=b[i];
	ntt(A,lim,1); ntt(B,lim,1);
	for(int i=0;i<lim;++i) A[i]=1ll*A[i]*B[i]%mod;
	ntt(A,lim,-1); 
	poly ret; ret.resize(d1+d2+1);
	for(int i=0;i<=d1+d2;++i) ret[i]=A[i];
	return ret;
}
poly add(poly a,poly b){
	int d1=dim(a),d2=dim(b);
	int d=max(d1,d2);
	poly ret; ret.resize(d+1,0);
	for(int i=0;i<=d1;++i) ret[i]=a[i];
	for(int i=0;i<=d2;++i) ret[i]=(ret[i]+b[i])%mod;
	return ret;	
}
poly multi(poly a,int y){
	int d=dim(a);
	poly ret; ret.resize(d+1);
	for(int i=0;i<=d;++i) ret[i]=1ll*a[i]*y%mod;
	return ret;	
}
poly minus(poly x,poly y){
	return add(x,multi(y,mod-1));	
}
poly cut(poly x,int t){
	x.resize(t+1);
	return x;
}
poly inv(poly x,int l=-1){
	int d=dim(x),t=1;
	poly ret; ret.resize(1);
	ret[0]=qpow(x[0],mod-2);
	if(l==-1) l=d;
	while(t<=l){
		t*=2;
		ret=minus(multi(ret,2),multi(multi(ret,ret),cut(x,t)));
		ret.resize(t);
	}
	ret.resize(l+1);
	return ret;
}
poly sqr(poly x,int l=-1){
	int d=dim(x),t=1;
	poly ret; ret.resize(1);
	map<int,int> mp;
	int kc=sqrt(mod);
	int T=1,Q=1,ans=-1;
	mp[1]=0;
	for(int i=0;i<kc;++i){
		int v=1ll*x[0]*T%mod;
		mp[v]=i; T=1ll*T*G%mod;
	}
	for(int i=1;1ll*i*kc<=mod;++i){
		Q=1ll*Q*T%mod;
		if(mp[Q]){
			ans=kc*i-mp[Q];
			break;
		}
	}
	ret[0]=qpow(G,ans/2);
	if(ret[0]>mod/2) ret[0]=mod-ret[0];
	if(l==-1) l=d;
	while(t<=l){
		t*=2;
		ret=multi(add(multi(ret,ret),cut(x,t)),inv(multi(ret,2),t));
		ret.resize(t);
	}
	ret.resize(l+1);
	return ret;	
}
poly deri(poly x){
	int d=dim(x);
	poly ret; ret.resize(d);
	for(int i=0;i<d;++i)
		ret[i]=1ll*x[i+1]*(i+1)%mod;
	return ret;
}
poly inte(poly x,int C=0){
	int d=dim(x);
	poly ret; ret.resize(d+2);
	for(int i=1;i<=d+1;++i)
		ret[i]=1ll*x[i-1]*qpow(i,mod-2)%mod;
	ret[0]=C;
	return ret;
}
poly ln(poly x,int l=-1){
	int d=dim(x);
	if(l==-1) l=d;
	return cut(inte(multi(deri(x),inv(x,l))),l);
}
poly exp(poly x,int l=-1){
	int d=dim(x),t=1;
	poly ret,one; ret.resize(1);
	ret[0]=1; one=ret;
	if(l==-1) l=d;
	while(t<=l){
		t*=2;
		ret=multi(ret,add(minus(one,ln(ret,t)),cut(x,t)));
		ret.resize(t);
	}
	ret.resize(l+1);
	return ret;
}
poly qpow(poly x,int k,int l=-1){
	int d=dim(x);
	if(l==-1) l=d;
    poly t=ln(x,l);
    t=multi(t,k);
    return exp(t,l);
}
int fac[N+5],ifac[N+5];
int main(){
    int n,k; scanf("%d%d",&n,&k);
	fac[0]=1;
	for(int i=1;i<=n;++i) fac[i]=1ll*fac[i-1]*i%mod;
	ifac[n]=qpow(fac[n],mod-2);
	for(int i=n-1;i>=0;--i) ifac[i]=1ll*ifac[i+1]*(i+1)%mod;
    poly f; f.resize(n);
    for(int i=0;i<n;++i) f[i]=1ll*fac[i]*ifac[i+1]%mod;
    f=qpow(f,k,n);
    for(int i=0;i<min(k,n+1);++i) printf("0 ");
    for(int i=k;i<=n;++i)
        printf("%lld ",1ll*f[i-k]*fac[i]%mod*ifac[k]%mod);
	return 0;
}
2023/8/18 22:20
加载中...