BZOJ4555 NTT或多项式求逆或CDQ+NTT-CSDN博客

本文链接：https://blog.csdn.net/animalcoder/article/details/81435884

求这个东西i,j,k<=n,NTT板子题

1884ms


//如果 r* 2^k +1 是个素数，那么在modr*2^k+1意义下，可以处理 2^k 以内规模的数据，

//2281701377=17*2^27+1 是一个挺好的数，平方刚好不会爆 long long

//1004535809=479*2^21+1 加起来刚好不会爆 int 也不错

//还有就是 998244353=119*2^23+1

//2^21=2e6 


#include<bits/stdc++.h> 
#define ll long long
using namespace std;

inline int gi(){//快读 
    int x=0,f=1;char ch=getchar();
    while(ch<'0'||ch>'9'){if(ch=='-')f=-1;ch=getchar();}
    while(ch>='0'&&ch<='9'){x=(x<<3)+(x<<1)+ch-'0';ch=getchar();}
    return f?x:-x;
}
const int _ = 400005;
const int mod = 998244353;//模数 
int n,m,s,E,lim,ans;int jc[_];int p2[_];
int inv[_],a[_],b[_],rev[_],l,og[_],len;//NTT本身用到的数组 
int qmod(int a,int b){//快速幂 
    int res=1;
    while (b) {if (b&1) res=1ll*res*a%mod;a=1ll*a*a%mod;b>>=1;}
    return res;
}
void ntt(int *P,int opt){//基本上不用改 
    for (int i=0;i<len;++i) if (i<rev[i]) swap(P[i],P[rev[i]]);
    for (int i=1;i<len;i<<=1){
        int W=qmod(3,(mod-1)/(i<<1));//3 根据模数的原根  可能要改 
        if (opt==-1) W=qmod(W,mod-2);
        og[0]=1;
        for (int j=1;j<i;++j) og[j]=1ll*og[j-1]*W%mod;
        for (int p=i<<1,j=0;j<len;j+=p)
            for (int k=0;k<i;++k){
                int x=P[j+k],y=1ll*og[k]*P[j+k+i]%mod;
                P[j+k]=(x+y)%mod,P[j+k+i]=(x-y+mod)%mod;
            }
    }
    if (opt==-1) for (int i=0,Inv=qmod(len,mod-2);i<len;++i) P[i]=1ll*P[i]*Inv%mod;
}
int main(){
    n=gi();//快读
	 
	E=n;        //求和表达式上限 (3.会用到) 
	lim=n;     //预处理上限(非NTT)
	
	//预处理阶乘 跟逆元阶乘 (非NTT)
    jc[0]=p2[0]=1; for (int i=1;i<=lim;++i) jc[i]=1ll*jc[i-1]*i%mod,p2[i]=1ll*p2[i-1]*2%mod;
    inv[lim]=qmod(jc[lim],mod-2);
    for (int i=lim;i>=1;--i) inv[i-1]=1ll*inv[i]*i%mod;
    
    //1.NTT   固定操作，处理数据的两倍 对2^k向上取整=len   len卷积上限 
    for (len=1;len<=(E<<1);len<<=1) ++l;--l;
    
    //2.NTT	  固定操作，预处理rev 蝶形变换 
    for (int i=0;i<len;++i) rev[i]=(rev[i>>1]>>1)|((i&1)<<l);
    
    //3.NTT   灵活操作  ai bi的表达式 
    for (int i=0;i<=E;++i) a[i]=i&1?mod-inv[i]:inv[i];
    for (int i=2;i<=E;++i) b[i]=1ll*(qmod(i,n+1)-1)%mod*qmod(i-1,mod-2)%mod*inv[i]%mod;
    b[0]=1;b[1]=n+1;//这题的bi是等比数列求和 
    
    //4.NTT   固定操作  ai卷积bi  也就是(a0+a1x+a2x^2...)*(b0+b1x+b2x^2....)=sum_{i=0...n}sum_{j=0...i}ai*b(j-i)
    ntt(a,1);ntt(b,1); 
    for (int i=0;i<len;++i) a[i]=1ll*a[i]*b[i]%mod;
    ntt(a,-1);
    
    //5.NTT   非固定操作  计算求和表达式的答案ans  
    for (int i=0;i<=E;++i) (ans+=1ll*p2[i]*jc[i]%mod*a[i]%mod)%=mod;
    //ans=1ll*jc[n]*jc[m]%mod*ans%mod;
    printf("%d\n",ans);
    return 0;
}

$\small \frac{dp[n]}{n!}=\sum_{i=1}^{n}\frac{dp[i]}{i!}*\frac{2}{(n-i)!},ans=\sum_{i=0}^{n}dp[i]$

F=F*G+1,,,,,1是常数项,,,,也就是dp[0]=1

F=1/(1-G) 对G求逆即可


//如果 r* 2^k +1 是个素数，那么在modr*2^k+1意义下，可以处理 2^k 以内规模的数据，

//2281701377=17*2^27+1 是一个挺好的数，平方刚好不会爆 long long

//1004535809=479*2^21+1 加起来刚好不会爆 int 也不错

//还有就是 998244353=119*2^23+1

//2^21=2e6 
//BZOJ4555 2196ms
#include<bits/stdc++.h> 
#define ll long long
using namespace std;

inline int gi(){//快读 
    int x=0,f=1;char ch=getchar();
    while(ch<'0'||ch>'9'){if(ch=='-')f=-1;ch=getchar();}
    while(ch>='0'&&ch<='9'){x=(x<<3)+(x<<1)+ch-'0';ch=getchar();}
    return f?x:-x;
}
const int _ = 800005;
const int mod = 998244353;//模数 
int n,m,s,E,lim,ans;int jc[_];
int inv[_],a[_],b[_],rev[_],l,og[_],len;//NTT本身用到的数组 
int qmod(int a,int b){//快速幂 
    int res=1;
    while (b) {if (b&1) res=1ll*res*a%mod;a=1ll*a*a%mod;b>>=1;}
    return res;
}
void ntt(int *P,int opt,int len){//基本上不用改  
    for (int i=0;i<len;++i) if (i<rev[i]) swap(P[i],P[rev[i]]);
    for (int i=1;i<len;i<<=1){
        int W=qmod(3,(mod-1)/(i<<1));//3 根据模数的原根  可能要改 
        if (opt==-1) W=qmod(W,mod-2);
        og[0]=1;
        for (int j=1;j<i;++j) og[j]=1ll*og[j-1]*W%mod;
        for (int p=i<<1,j=0;j<len;j+=p)
            for (int k=0;k<i;++k){
                int x=P[j+k],y=1ll*og[k]*P[j+k+i]%mod;
                P[j+k]=(x+y)%mod,P[j+k+i]=(x-y+mod)%mod;
            }
    }
    if (opt==-1) for (int i=0,Inv=qmod(len,mod-2);i<len;++i) P[i]=1ll*P[i]*Inv%mod;
}
void getInv(int *a,int *b,int n)//基本不改，套路 
{
	static int tmp[_];
	if(n==1)return void(b[0]=qmod(a[0],mod-2));
	getInv(a,b,n>>1);
	for (int i=0;i<n;i++) tmp[i]=a[i],tmp[n+i]=0;
	int L=0; while (!(n>>L&1)) L++;
	for (int i=1;i<(n<<1);i++) rev[i]=(rev[i>>1]>>1)|((i&1)<<L);
	ntt(tmp,1,n<<1); ntt(b,1,n<<1);
	for (int i=0;i<(n<<1);i++)
    	tmp[i]=(ll)b[i]*(2+mod-(ll)tmp[i]*b[i]%mod)%mod;//固定操作 
    ntt(tmp,-1,n<<1);
    for (int i=0;i<n;i++) b[i]=tmp[i],b[n+i]=0;
}


int main(){
    n=gi();//快读
	 
	E=n;        //求和表达式上限 (3.会用到) 
	lim=n; //预处理上限(非NTT)
	
	//预处理阶乘 跟逆元阶乘 (非NTT)
    jc[0]=1; for (int i=1;i<=lim;++i) jc[i]=1ll*jc[i-1]*i%mod;
    inv[lim]=qmod(jc[lim],mod-2);
    for (int i=lim;i;--i) inv[i-1]=1ll*inv[i]*i%mod;
    
    //1.NTT   固定操作，处理数据的两倍 对2^k向上取整=len   len卷积上限 
    for (len=1;len<=(E<<1);len<<=1) ++l;--l;//len<=(E<<1)
    
    //2.NTT	  固定操作，预处理rev 蝶形变换 
    for (int i=0;i<len;++i) rev[i]=(rev[i>>1]>>1)|((i&1)<<l);
    
    //3.NTT   灵活操作  ai bi的表达式 
    for (int i=1;i<=E;++i) a[i]=1ll*2*(mod-inv[i])%mod;a[0]=1;
    //for (int i=0;i<=E;++i) b[i]=i&1?mod-inv[i]:inv[i];
    
    //3.1 NTT求逆操作  对a求逆，a保持不变，B为求逆结果 
    getInv(a,b,len>>1);
    
    //4.NTT   固定操作  ai卷积bi  也就是(a0+a1x+a2x^2...)*(b0+b1x+b2x^2....)=sum_{i=0...n}sum_{j=0...i}ai*b(j-i)
    //ntt(a,1);ntt(b,1); 
    //for (int i=0;i<len;++i) a[i]=1ll*a[i]*b[i]%mod;
    //ntt(a,-1);
    
    
    // 5.NTT   非固定操作  计算求和表达式的答案ans  
	//这里是ans=bn*n!+b[n-1]*(n-1)!+..... 
    int ans=b[n];
    for (int i=n;i;i--) ans=((ll)ans*i+b[i-1])%mod;
 	printf("%d\n",ans);

    
    //

    return 0;
}

上面的柿子CDQ+NTT是通解：https://blog.csdn.net/animalcoder/article/details/81268201