BZOJ #3453. tyvj 1858 XLkxc
阿新 • • 發佈:2020-08-09
題意:
設
\[f(i)=1^k+2^k+\cdots + i^k\\ g(i)=f(1)+f(2)+\cdots+ f(i) \]
給定 \(k,a,i,d\) ,求
\[h(i)=g(a)+g(a+d)+g(a+2d)+\cdots+ g(a+id) \]
引理:若一個多項式差分 \(k+1\) 次後為 \(0\) ,則這是一個 \(k\) 次多項式(若 \(g(x)=f(x+1)-f(x)\),則 \(g(x)\) 為差分一次後的多項式)。
根據之前的知識,我們知道 \(f(i)\) 是一個關於 \(i\) 的 \(k+1\) 次多項式,而 \(g(i)\) 是一個差分一次後是 \(f(i+1)\)
\(g(i)\) 可以用拉格朗日插值在 \(O(k)\) 的時間求出來,那麼我們可以暴力求出 \(h(1)+h(2)+\cdots+h(k+4)\) 的值,再用拉格朗日插值求出要求的 \(h(i)\)。這個故事告訴我們插值可以巢狀。
程式碼:
#include<cstdio> #include<cstring> #include<algorithm> #include<iostream> #define int long long #include<assert.h> using namespace std; typedef long long LL; const int K=200,M=1234567891; int a[K+9],p[K],cnt; int k,A,n,D; LL d[K+9],pre[K+9],suf[K+9],fac[K+9],inv_fac[K+9],h[K+9]; int ksm(int a,int b) { int res=1; while(b) { if(b&1) res=1ll*res*a%M; b>>=1,a=1ll*a*a%M; } return res; } void prework() { d[1]=1; for (int i=2;i<=K;i++) a[i]=0; cnt=0; for (int i=2;i<=K;i++) { if(!a[i]) a[i]=i,p[++cnt]=i,d[i]=ksm(i,k); for (int j=1;j<=cnt;j++) { if(p[j]>a[i]||p[j]>K/i) break; a[p[j]*i]=p[j]; d[p[j]*i]=d[p[j]]*d[i]%M; } } for (int i=1;i<=K;i++) d[i]=(d[i]+d[i-1])%M; for (int i=1;i<=K;i++) d[i]=(d[i]+d[i-1])%M; } void init() { scanf("%lld %lld %lld %lld",&k,&A,&n,&D); prework(); } LL Get_Sum(LL n,LL d[],int k) { pre[0]=1; for (int i=1;i<=k+3;i++) pre[i]=pre[i-1]*((n-i)%M)%M; suf[k+4]=1; for (int i=k+3;i>=1;i--) suf[i]=suf[i+1]*((n-i)%M)%M; fac[0]=1; for (int i=1;i<=k+3;i++) fac[i]=fac[i-1]*i%M; inv_fac[k+3]=ksm(fac[k+3],M-2); for (int i=k+2;i>=0;i--) inv_fac[i]=inv_fac[i+1]*(i+1)%M; for (int i=1;i<=k+3;i++) assert(fac[i]*inv_fac[i]%M==1); LL ans=0; for (int i=1;i<=k+3;i++) ans=(ans+d[i]*pre[i-1]%M*suf[i+1]%M*inv_fac[i-1]%M*inv_fac[k+3-i]%M*(k+3-i&1?-1:1))%M; return ans; } void work() { h[0]=Get_Sum(A,d,k); for (int i=1;i<=k+4;i++) h[i]=(h[i-1]+Get_Sum(A+1ll*i*D,d,k))%M; printf("%lld\n",(Get_Sum(n,h,k+1)+M)%M); } signed main() { int T; scanf("%lld",&T); while(T--) { init(); work(); } return 0; }