zoukankan      html  css  js  c++  java
  • 【洛谷 P5409】 【模板】—第一类斯特林数·列(多项式Ln+多项式快速幂)

    传送门

    斯特林数学习笔记


    考虑这样一个展开

    (x+1)n=i=0xi(ni)=i=0xii!ni (x+1)^n=sum_{i=0}^{infty}x^i{nchoose i}\=sum_{i=0}^{infty}frac{x^i}{i!}n^{underline i}


    xn=i=0n(1)nis(n,i)xix^{underline n}=sum_{i=0}^{n}(-1)^{n-i}s(n,i)x^i

    所以(x+1)n=i=0xii!j=0i(1)ijs(i,j)nj=j=0nji=j(1)ijs(i,j)xii!=i=0nij=i(1)jis(j,i)xjj!(x+1)^n=sum_{i=0}^{infty}frac{x^i}{i!}sum_{j=0}^{i}(-1)^{i-j}s(i,j)n^j \=sum_{j=0}^{infty}n^jsum_{i=j}^{infty}(-1)^{i-j}s(i,j)frac{x^i}{i!}\ =sum_{i=0}^{infty}n^isum_{j=i}^{infty}(-1)^{j-i}s(j,i)frac {x^j}{j!}

    又由于(x+1)n=enLn(x+1)=i=0niLn(x+1)ii!(x+1)^n=e^{n*Ln(x+1)}=sum_{i=0}^{infty}n^ifrac{Ln(x+1)^i}{i!}

    所以i=0niLn(x+1)ii!=i=0nij=i(1)jis(j,i)xjj!sum_{i=0}^{infty}n^ifrac{Ln(x+1)^i}{i!}=sum_{i=0}^{infty}n^isum_{j=i}^{infty}(-1)^{j-i}s(j,i)frac {x^j}{j!}

    Ln(x+1)kk!=j=k(1)jks(j,k)zjj!frac{Ln(x+1)^k}{k!}=sum_{j=k}^{infty}(-1)^{j-k}s(j,k)frac{z^j}{j!}

    对左边做一个LnLn再做个快速幂就可以了
    可以用能处理a0≠1a_0= ot1的快速幂解决

    复杂度O(nlogn)O(nlogn)

    #include<bits/stdc++.h>
    using namespace std;
    #define gc getchar
    inline int read(){
    	char ch=gc();
    	int res=0,f=1;
    	while(!isdigit(ch))f^=ch=='-',ch=gc();
    	while(isdigit(ch))res=(res+(res<<2)<<1)+(ch^48),ch=gc();
    	return f?res:-res;
    }
    #define re register
    #define pb push_back
    #define cs const
    #define pii pair<int,int>
    #define fi first
    #define se second
    #define ll long long
    #define poly vector<int>
    #define bg begin
    cs int mod=167772161,G=3;
    inline int add(int a,int b){return (a+=b)>=mod?a-mod:a;}
    inline void Add(int &a,int b){(a+=b)>=mod?(a-=mod):0;}
    inline int dec(int a,int b){return (a-=b)<0?a+mod:a;}
    inline void Dec(int &a,int b){(a-=b)<0?(a+=mod):0;}
    inline int mul(int a,int b){return 1ll*a*b>=mod?1ll*a*b%mod:a*b;}
    inline void Mul(int &a,int b){a=mul(a,b);}
    inline int ksm(int a,int b,int res=1){
    	for(;b;b>>=1,a=mul(a,a))(b&1)&&(res=mul(res,a));return res;
    }
    inline void chemx(int &a,int b){a<b?a=b:0;}
    inline void chemn(int &a,int b){a>b?a=b:0;}
    cs int N=(1<<20)|5,C=20;
    poly w[C+1];
    int rev[N],fac[N],ifac[N],inv[N];
    inline void init(cs int len=N-5){
    	fac[0]=ifac[0]=inv[0]=inv[1]=1;
    	for(int i=1;i<=len;i++)fac[i]=mul(fac[i-1],i);
    	ifac[len]=ksm(fac[len],mod-2);
    	for(int i=len-1;i;i--)ifac[i]=mul(ifac[i+1],i+1);
    	for(int i=2;i<=len;i++)inv[i]=mul(mod-mod/i,inv[mod%i]);
    }
    inline void init_w(){
    	for(int i=1;i<=C;i++)w[i].resize(1<<(i-1));
    	int wn=ksm(G,(mod-1)/(1<<C));
    	w[C][0]=1;
    	for(int i=1;i<(1<<(C-1));i++)w[C][i]=mul(w[C][i-1],wn);
    	for(int i=C-1;i;i--)
    	for(int j=0;j<(1<<(i-1));j++)
    	w[i][j]=w[i+1][j<<1];
    }
    inline void init_rev(int lim){
    	for(int i=0;i<lim;i++)rev[i]=(rev[i>>1]>>1)|((i&1)*(lim>>1));
    }
    inline void ntt(poly &f,int lim,int kd){
    	for(int i=0;i<lim;i++)if(i>rev[i])swap(f[i],f[rev[i]]);
    	for(int a0,a1,l=1,mid=1;mid<lim;mid<<=1,l++)
    	for(int i=0;i<lim;i+=(mid<<1))
    	for(int j=0;j<mid;j++)
    	a0=f[i+j],a1=mul(w[l][j],f[i+j+mid]),f[i+j]=add(a0,a1),f[i+j+mid]=dec(a0,a1);
    	if(kd==-1){
    		reverse(f.bg()+1,f.bg()+lim);
    		for(int i=0;i<lim;i++)Mul(f[i],inv[lim]);
    	}
    }
    inline poly operator +(poly a,poly b){
    	if(a.size()<b.size())a.resize(b.size());
    	for(int i=0;i<b.size();i++)Add(a[i],b[i]);
    	return a;
    }
    inline poly operator -(poly a,poly b){
    	if(a.size()<b.size())a.resize(b.size());
    	for(int i=0;i<b.size();i++)Dec(a[i],b[i]);
    	return a;
    }
    inline poly operator *(poly a,int b){
    	for(int i=0;i<a.size();i++)Mul(a[i],b);
    	return a;
    }
    inline poly operator *(poly a,poly b){
    	int deg=a.size()+b.size()-1,lim=1;
    	if(deg<=64){
    		poly c(deg,0);
    		for(int i=0;i<a.size();i++)
    		for(int j=0;j<b.size();j++)
    		Add(c[i+j],mul(a[i],b[j]));
    		return c;
    	}
    	while(lim<deg)lim<<=1;
    	init_rev(lim);
    	a.resize(lim),ntt(a,lim,1);
    	b.resize(lim),ntt(b,lim,1);
    	for(int i=0;i<lim;i++)Mul(a[i],b[i]);
    	ntt(a,lim,-1),a.resize(deg);
    	return a;
    }
    inline poly Inv(poly a,int deg){
    	poly b(1,ksm(a[0],mod-2)),c;
    	for(int lim=4;lim<(deg<<2);lim<<=1){
    		c=a,c.resize(lim>>1);
    		init_rev(lim);
    		c.resize(lim),ntt(c,lim,1);
    		b.resize(lim),ntt(b,lim,1);
    		for(int i=0;i<lim;i++)Mul(b[i],dec(2,mul(b[i],c[i])));
    		ntt(b,lim,-1),b.resize(lim>>1);
    	}b.resize(deg);return b;
    }
    inline poly deriv(poly a){
    	for(int i=0;i<a.size()-1;i++)a[i]=mul(a[i+1],i+1);
    	a.pop_back();return a;
    }
    inline poly integ(poly a){
    	a.pb(0);
    	for(int i=a.size()-1;i;i--)a[i]=mul(a[i-1],inv[i]);
    	a[0]=0;return a;
    }
    inline poly Ln(poly a,int deg){
    	a=integ(Inv(a,deg)*deriv(a)),a.resize(deg);return a;
    }
    inline poly ksm(poly a,int b,int deg){
    	poly res(1,1);
    	for(;b;b>>=1){
    		if(b&1){
    			res=res*a;if(res.size()>deg)res.resize(deg);
    		}
    		a=a*a;
    		if(a.size()>deg)a.resize(deg);
    	}
    	res.resize(deg);
    	return res;
    }
    poly f;
    int n,k;
    int main(){
    	n=read(),k=read();
    	if(n<k){for(int i=0;i<=n;i++)cout<<0<<" ";return 0;}
    	init_w(),init();
    	f.pb(1),f.pb(1);
    	f=Ln(f,n+1),f=ksm(f,k,n+1);
    	f=f*ifac[k];
    	for(int i=0;i<k;i++)cout<<0<<" ";
    	for(int i=k;i<=n;i++){
    		int res=f[i];
    		if((i-k)&1)Mul(res,mod-fac[i]);
    		else Mul(res,fac[i]);
    		cout<<res<<" ";
    	}
    }
    
  • 相关阅读:
    editplus 快捷键
    python 堆栈
    python 矩阵转置
    python 单向链表
    python 读空的json文件
    c++ 结构体
    手把手教你如何利用Meterpreter渗透Windows系统
    vuejs npm chromedriver 报错
    强大的开源网络侦查工具:IVRE
    在vue 中使用Stylus
  • 原文地址:https://www.cnblogs.com/stargazer-cyk/p/12328649.html
Copyright © 2011-2022 走看看