題解P3711：【倉鼠的數學題】

這題黑的丫！怎麼會掉紫呢！

noteskey

伯努利數... 這裏 有介紹喲~ 寫的很是詳細呢~

反正這題就是推柿子...

另外就是黈力算法的運用 QWQ

咱們令 \(ANS(x)\) 爲答案多項式，那麼這個多項式能夠這麼求：

（下面咱們定義 \(S(n,k)\) 爲天然冪和函數（不是第二類斯特林數！），即 \(\sum_{i=0}^{ n} i^k\)）

\[\begin{aligned}ANS(x)=& \sum_{k=0}^{n} S_k (x) a_k \\ =& \sum_{k=0}^{n} (S(x,k)+x^k) a_k \\ =& \sum_{k=0}^n a_k \Big( x^k +{1\over k+1} \sum_{i=0}^k \begin{pmatrix} k+1\\i \end{pmatrix} B_i x^{k+1-i}\Big) \\ =& \sum_{k=0}^n a_k \Big( x^k + k! \sum_{i=0}^k {B_i\over i!} {x^{k+1-i} \over (k+1-i)! }\Big) \\ =& \sum_{k=0}^n a_k x^k + \sum_{k=0}^n a_k k! \sum_{i=0}^k {B_i\over i!} {x^{k+1-i} \over (k+1-i)! } \\ =& \sum_{i=0}^n a_i x^i + \sum_{i=1}^{n+1} x^i i! \sum_{k=i-1}^{n} a_k k!{B_{k+1-i}\over (k+1-i)!} \end{aligned}\]

而後感受作不下去了呢...後面雖然說像是卷積的形式然鵝根本就不是卷積呢QWQ

首先，咱們看着表達式太長了，因而乎把表達式換成一個函數：

\[\begin{aligned}f(i)=&~a_ii! \\g(i) =&~{B_i \over i!} \end{aligned}\]

那麼原來的式子就是：

\[ANS=\sum_{i=0}^n a_i x^i + \sum_{i=1}^{n+1} x^i i! \sum_{k=i-1}^{n} f(k)g(k+1-i) \]

這樣的話，咱們就更加清晰的發現後面的式子不是卷積的形式了 XD

那麼咱們令 \(gr(i)=g(n-i)\) （也就是翻轉多項式）

而後原來的式子就是：

\[ANS=\sum_{i=0}^n a_i x^i+ \sum_{i=1}^{n+1} x^i i! \sum_{k=i-1}^{n} f(k)gr(n+i-k-1) \]

發現這裏就是最高項改爲 \(n+i\) 項的卷積了，因而咱們發現 要計算後面的 \(\sum_{k=i-1}^n\) 什麼的算出來的是 \(f* gr\) 的第 \((n+i-1)\) 項

因而咱們把 \(f* gr\) 的多項式求出來，拿第 \((n-i+1)\) 項乘上 \(i!\) ，再加上 \(a_i\) 就是最後多項式的第 i 項答案了

好像沒有黑題難度？ 前提你得知道伯努利數這玩意兒丫！ 【霧

並且這個推導過程...好像比 shadowice 巨巨的短不少丫...【小聲

code

不知道爲何不開 Ofast 跑得比開了快... 600ms +

//by Judge
#include<cstdio>
#include<cstring>
#include<iostream>
#include<algorithm>
#define Rg register
#define fp(i,a,b) for(Rg int i=(a),I=(b)+1;i<I;++i)
#define ll long long
using namespace std;
const int mod=998244353;
const int iG=332748118;
const int M=1<<19|3;
typedef int arr[M];
#ifndef Judge
#define getchar() (p1==p2&&(p2=(p1=buf)+fread(buf,1,1<<21,stdin),p1==p2)?EOF:*p1++)
#endif
char buf[1<<21],*p1=buf,*p2=buf;
inline int mul(int x,int y){return 1ll*x*y%mod;}
inline int dec(int x,int y){return x<y?x-y+mod:x-y;}
inline int inc(int x,int y){return x+y>=mod?x+y-mod:x+y;}
inline int read(){ int x=0,f=1; char c=getchar();
    for(;!isdigit(c);c=getchar()) if(c=='-') f=-1;
    for(;isdigit(c);c=getchar()) x=x*10+c-'0'; return x*f;
} char sr[1<<21],z[20];int CCF=-1,Z;
inline void Ot(){fwrite(sr,1,CCF+1,stdout),CCF=-1;}
inline void print(int x,char chr=' '){
    if(CCF>1<<20)Ot();if(x<0)sr[++CCF]=45,x=-x;
    while(z[++Z]=x%10+48,x/=10);
    while(sr[++CCF]=z[Z],--Z);sr[++CCF]=chr;
} int n,limit; arr fac,ifac,a,B,C,r;
inline int qpow(int x,int p=mod-2,int s=1){
    for(;p;p>>=1,x=mul(x,x)) if(p&1) s=mul(s,x); return s;
}
inline void init(int n){ int len=-1;
    for(limit=1;limit<=n;limit<<=1) ++len;
    fp(i,0,limit-1) r[i]=(r[i>>1]>>1)|((i&1)<<len); 
}
inline void NTT(int* a,int tp){ fp(i,0,limit-1) if(i<r[i]) swap(a[i],a[r[i]]);
    for(Rg int mid=1;mid<limit;mid<<=1){ Rg int I=mid<<1,Gn=qpow(tp?3:iG,(mod-1)/I);
        for(Rg int j=0,x;j<limit;j+=I) for(Rg int k=0,g=1;k<mid;++k,g=mul(g,Gn))
            x=mul(a[j+k+mid],g),a[j+k+mid]=dec(a[j+k],x),a[j+k]=inc(a[j+k],x);
    } if(tp) return ; int inv=qpow(limit); fp(i,0,limit-1) a[i]=mul(a[i],inv);
}
void get_inv(int* a,int* b,int n){
    if(n==1) return b[0]=qpow(a[0]),void(); static int c[M],d[M];
    get_inv(a,b,n>>1),init(n); for(int i=0;i<n;++i) c[i]=a[i],d[i]=b[i];
    for(int i=n;i<limit;++i) c[i]=d[i]=0; NTT(c,1),NTT(d,1);
    fp(i,0,limit-1) c[i]=mul(c[i],mul(d[i],d[i])); NTT(c,0);
    fp(i,0,n-1) b[i]=dec(inc(b[i],b[i]),c[i]);
}
inline void prep(int len){ B[0]=ifac[0]=ifac[1]=fac[0]=fac[1]=1;
    fp(i,2,len) ifac[i]=mul(mod-mod/i,ifac[mod%i]);
    fp(i,2,len) ifac[i]=mul(ifac[i-1],ifac[i]),fac[i]=mul(fac[i-1],i);
    get_inv(ifac+1,B,len);
}
int main(){ n=read();
    for(limit=1;limit<=n;limit<<=1); prep(limit);
    fp(i,0,n) a[i]=read(),C[i]=mul(a[i],fac[i]);
    reverse(B,B+1+n),init(n<<1|1);
    fp(i,n+1,limit-1) B[i]=C[i]=0; NTT(B,1),NTT(C,1);
    fp(i,0,limit-1) B[i]=mul(B[i],C[i]); NTT(B,0);
    fp(i,0,n+1) B[n+i-1]=mul(B[n+i-1],ifac[i]);
    fp(i,0,n) B[n+i-1]=inc(B[n+i-1],a[i]);
    print(a[0]); fp(i,1,n+1) print(B[n+i-1]);
    return sr[CCF]='\n',Ot(),0;
}

more...

其實咱們只須要用另一種伯努利數，就可讓推導更加簡潔（沒簡潔多少，代碼複雜度也基本沒變的說）

這個伯努利數就是以前提到的那篇博客裏面講的 \(B^+\)，這個伯努利數列知足的性質更加適合作這道題...

因而咱們令 \(B=B^+\)， 而後重推一遍：

（其實不必，給個代碼就 OJBK 了？ 【滑稽）

\[\begin{aligned}ANS(x)=& \sum_{k=0}^{n} S_k (x) a_k \\ =& \sum_{k=0}^{n} S^+(x,k) a_k \\ =& \sum_{k=0}^n a_k \Big( {1\over k+1} \sum_{i=0}^k \begin{pmatrix} k+1\\i \end{pmatrix} B_i x^{k+1-i} \Big) \\ =& \sum_{k=0}^n a_k k! \sum_{i=0}^k {B_i\over i!} {x^{k+1-i} \over (k+1-i)! } \\ =& \sum_{i=1}^{n+1} x^i i! \sum_{k=i-1}^{n} a_k k!{B_{k+1-i}\over (k+1-i)!} \end{aligned}\]

依然令：

\[\begin{aligned}f(i)=&~a_ii! \\g(i) =&~{B_i \over i!} \\g_r(i)=&~g(n-i)\end{aligned}\]

（注意這裏的 \(B_i\) 就是 \(B_i^+\) ）

\[ANS=\sum_{i=1}^{n+1} x^i i! \sum_{k=i-1}^{n} f(k)gr(n+i-k-1) \]

一樣的，把 \(f* g_r\) 的第 \((n-i+1)\) 項乘上 \(i!\) 就是最後的答案了 ，不一樣的就是這裏不須要再加上 \(a_i\) XD

code*2

代碼...可謂沒有什麼區別

//by Judge
#include<cstdio>
#include<cstring>
#include<iostream>
#include<algorithm>
#define Rg register
#define fp(i,a,b) for(Rg int i=(a),I=(b)+1;i<I;++i)
#define ll long long
using namespace std;
const int mod=998244353;
const int iG=332748118;
const int M=1<<19|3;
typedef int arr[M];
#ifndef Judge
#define getchar() (p1==p2&&(p2=(p1=buf)+fread(buf,1,1<<21,stdin),p1==p2)?EOF:*p1++)
#endif
char buf[1<<21],*p1=buf,*p2=buf;
inline int mul(int x,int y){return 1ll*x*y%mod;}
inline int dec(int x,int y){return x<y?x-y+mod:x-y;}
inline int inc(int x,int y){return x+y>=mod?x+y-mod:x+y;}
inline int read(){ int x=0,f=1; char c=getchar();
    for(;!isdigit(c);c=getchar()) if(c=='-') f=-1;
    for(;isdigit(c);c=getchar()) x=x*10+c-'0'; return x*f;
} char sr[1<<21],z[20];int CCF=-1,Z;
inline void Ot(){fwrite(sr,1,CCF+1,stdout),CCF=-1;}
inline void print(int x,char chr=' '){
    if(CCF>1<<20)Ot();if(x<0)sr[++CCF]=45,x=-x;
    while(z[++Z]=x%10+48,x/=10);
    while(sr[++CCF]=z[Z],--Z);sr[++CCF]=chr;
} int n,limit; arr fac,ifac,a,B,C,r;
inline int qpow(int x,int p=mod-2,int s=1){
    for(;p;p>>=1,x=mul(x,x)) if(p&1) s=mul(s,x); return s;
}
inline void init(int n){ int len=-1;
    for(limit=1;limit<=n;limit<<=1) ++len;
    fp(i,0,limit-1) r[i]=(r[i>>1]>>1)|((i&1)<<len); 
}
inline void NTT(int* a,int tp){ fp(i,0,limit-1) if(i<r[i]) swap(a[i],a[r[i]]);
    for(Rg int mid=1;mid<limit;mid<<=1){ Rg int I=mid<<1,Gn=qpow(tp?3:iG,(mod-1)/I);
        for(Rg int j=0,x;j<limit;j+=I) for(Rg int k=0,g=1;k<mid;++k,g=mul(g,Gn))
            x=mul(a[j+k+mid],g),a[j+k+mid]=dec(a[j+k],x),a[j+k]=inc(a[j+k],x);
    } if(tp) return ; int inv=qpow(limit); fp(i,0,limit-1) a[i]=mul(a[i],inv);
}
void get_inv(int* a,int* b,int n){
    if(n==1) return b[0]=qpow(a[0]),void(); static int c[M],d[M];
    get_inv(a,b,n>>1),init(n); for(int i=0;i<n;++i) c[i]=a[i],d[i]=b[i];
    for(int i=n;i<limit;++i) c[i]=d[i]=0; NTT(c,1),NTT(d,1);
    fp(i,0,limit-1) c[i]=mul(c[i],mul(d[i],d[i])); NTT(c,0);
    fp(i,0,n-1) b[i]=dec(inc(b[i],b[i]),c[i]);
}
inline void prep(int len){
    B[0]=ifac[0]=ifac[1]=fac[0]=fac[1]=1;
    fp(i,2,len) ifac[i]=mul(mod-mod/i,ifac[mod%i]);
    fp(i,2,len) ifac[i]=mul(ifac[i-1],ifac[i]);
    fp(i,2,len) fac[i]=mul(fac[i-1],i);
    get_inv(ifac+1,B,len); B[1]=mod-B[1];
}
int main(){ n=read();
    for(limit=1;limit<=n;limit<<=1); prep(limit);
    fp(i,0,n) a[i]=read(),C[i]=mul(a[i],fac[i]);
    reverse(B,B+1+n),init(n<<1|1);
    fp(i,n+1,limit-1) B[i]=C[i]=0; NTT(B,1),NTT(C,1);
    fp(i,0,limit-1) B[i]=mul(B[i],C[i]); NTT(B,0);
    fp(i,0,n+1) B[n+i-1]=mul(B[n+i-1],ifac[i]);
    print(a[0]); fp(i,1,n+1) print(B[n+i-1]);
    return sr[CCF]='\n',Ot(),0;
}