正紧地多项式求exp 牛客练习赛24F题解
题目:https://www.nowcoder.com/acm/contest/157/F
题解:https://www.nowcoder.com/discuss/93531?type=101&order=0&pos=2&page=0
#include<cstdio> #include<algorithm> #include<cstring> #include<cmath> #define re register #define rep(i,s,t) for(re int i=s;i<=t;++i) #define _rep(i,s,t) for(re int i=s;i>=t;--i) #define Rep(i,s,t) for(re int i=s;i<t;++i) #define go(x) for(re int e=las[x];e;e=nxt[e]) #define re register #define fi first #define se second #define mp make_pair #define pb push_back #define pii pair<int,int> #define pi acos(-1) #define gi(x) read(x) #define gii(x,y) read(x),read(y) #define giii(x,y,z) read(x),read(y),read(z) #define ms(f,x) memset(f,x,sizeof f) #define open(x) freopen(#x".in","r",stdin),freopen(#x".out","w",stdout) namespace IO{ #define gc getchar() #define pc(x) putchar(x) template<typename T>inline void read(T &x){ x=0;int f=1;char ch=gc;while(ch>'9'||ch<'0'){if(ch=='-')f=-1;ch=gc;} while(ch>='0'&&ch<='9')x=(x<<3)+(x<<1)+ch-'0',ch=gc;x*=f;return; } template<typename T>inline void write(T x=0){ T wr[51];wr[0]=0;if(x<0)pc('-'),x=-x;if(!x)pc(48); while(x)wr[++wr[0]]=x%10,x/=10;while(wr[0])pc(48+wr[wr[0]--]);return; } } using IO::read; using IO::write; using namespace std; typedef long long ll; const int N=1e6+11,mod=19260817; int n,m; int a[N],b[N],p[N],c[N], d[N],e[N],f[N],inv[N]; int lg2; typedef double db; struct cp{ db r,i; cp(db a=0.,db b=0.){r=a,i=b;} inline cp operator +(cp A)const{return cp(r+A.r,i+A.i);} inline cp operator -(cp A)const{return cp(r-A.r,i-A.i);} inline cp operator *(cp A)const{return cp(r*A.r-i*A.i,r*A.i+i*A.r);} inline cp operator *(int A)const{return cp(A,0);} }a1[N],a2[N],b1[N],b2[N],c1[N],c2[N],c3[N]; inline int inc(int x,int y){ re int res=x+y; if(res>=mod)res-=mod; if(res<0)res+=mod; return res; } inline int fp(int a,int b){ if(b>=mod-1)b-=mod-1; if(b<0)b+=mod-1; re int res=1; for(;b;b>>=1,a=1ll*a*a%mod) if(b&1) res=1ll*res*a%mod; return res; } inline void fft(cp *a,int n,int f){ re int l=0,d=1; for(;d<n;d<<=1)++l; Rep(i,0,n) p[i]=(p[i>>1]>>1)^((i&1)<<(l-1)); Rep(i,0,n) if(i<p[i]) swap(a[i],a[p[i]]); for(re int i=1;i<n;i<<=1){ cp gn=cp(cos(1.*pi/i),sin(f*1.*pi/i)),w; for(re int j=0;w=cp(1,0),j<n;j+=(i<<1)) for(re int k=j;k<i+j;++k,w=w*gn){ cp x=a[k],y=w*a[i+k]; a[k]=x+y,a[i+k]=x-y; } } if(f==-1) Rep(i,0,n) a[i].r=1.*a[i].r/n; } inline void Mul(int *a,int *b,int *c,int len){ re int sqr=sqrt(mod); re cp temp; Rep(i,0,len){ a1[i]=a[i]/sqr,b1[i]=a[i]%sqr; a2[i]=b[i]/sqr,b2[i]=b[i]%sqr; } fft(a1,len,1),fft(b1,len,1),fft(a2,len,1),fft(b2,len,1); Rep(i,0,len){ c1[i]=a1[i]*a2[i]; c2[i]=a1[i]*b2[i]+a2[i]*b1[i]; c3[i]=b1[i]*b2[i]; } fft(c1,len,-1),fft(c2,len,-1),fft(c3,len,-1); Rep(i,0,len) c[i]=((ll)(round(c1[i].r))%mod*sqr%mod*sqr%mod+(ll)(round(c2[i].r))%mod*sqr%mod+(ll)(round(c3[i].r))%mod)%mod; } inline void Det(int *a,int *b,int len){ b[len-1]=0; Rep(i,1,len) b[i-1]=1ll*a[i]*i%mod; } inline void Area(int *a,int *b,int len){ b[0]=0; Rep(i,1,len) b[i]=1ll*a[i-1]*inv[i]%mod; } inline void Inv(int *a,int *b,int len){ if(len==1){ b[0]=fp(a[0],mod-2); return; } Inv(a,b,len>>1); Rep(i,0,len)d[i]=a[i],e[i]=b[i]; /* fft(d,tmp,1),fft(e,tmp,1); Rep(i,0,tmp) d[i]=1ll*d[i]*e[i]%mod*e[i]%mod; fft(d,tmp,-1); */ re int tmp=len<<1; Mul(d,e,d,tmp),Mul(d,e,d,tmp); Rep(i,0,len) b[i]=inc(b[i],b[i]), b[i]=inc(b[i],mod-d[i]); Rep(i,0,tmp) d[i]=e[i]=0; } inline void Ln(int *a,int *b,int len){ Inv(a,f,len),Det(a,d,len); re int tmp=len<<1; /* fft(f,tmp,1),fft(d,tmp,1); Rep(i,0,tmp) f[i]=1ll*f[i]*d[i]%mod; fft(f,tmp,-1), */ Mul(f,d,f,tmp); Area(f,b,len); Rep(i,0,tmp) f[i]=d[i]=0; } inline void Exp(int *a,int *b,int len){ if(len==1){ b[0]=1; return ; } Exp(a,b,len>>1),Ln(b,c,len); Rep(i,0,len) c[i]=inc(a[i],mod-c[i]),f[i]=b[i]; c[0]=inc(c[0],1); re int tmp=len<<1; /* fft(c,tmp,1),fft(f,tmp,1); Rep(i,0,tmp)c[i]=1ll*c[i]*f[i]%mod; fft(c,tmp,-1); */ Mul(c,f,c,tmp); Rep(i,0,len) b[i]=c[i]; Rep(i,0,tmp) c[i]=f[i]=0; } int main(){ gii(n,m); inv[1]=1; rep(i,2,1e5) inv[i]=mod-1ll*mod/i*inv[mod%i]%mod; re int d=1,v,ans=0; rep(i,1,n){ gi(v); for(re int j=v,k=1;j<=m;j+=v,++k) a[j]=inc(a[j],inv[k]); } while(d<=m)d<<=1; Exp(a,b,d); rep(i,1,m) ans=inc(ans,b[i]); printf("%d\n",ans); return 0; }