正经地多项式求exp 牛客练习赛24F题解

#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],cnt[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(){
    re int d=1,v,ans=0;
    gii(n,m);
    rep(i,1,n)
        gi(v),++cnt[v];
    inv[1]=1;
    rep(i,2,1e5)
        inv[i]=mod-1ll*mod/i*inv[mod%i]%mod;
    rep(i,1,m)
        if(cnt[i]){
            for(re int j=i,k=1;j<=m;j+=i,++k)
                a[j]=(a[j]+1ll*cnt[i]*inv[k])%mod;
        }
    while(d<=m)d<<=1;
    Exp(a,b,d);
    rep(i,1,m)
        ans=inc(ans,b[i]);
    printf("%d\n",ans);
    return 0;
}