[Luogu5320][BJOI2019]堪破神机(DP+斯特林数)

https://www.cnblogs.com/cjyyb/p/10747543.html

特征方程+斯特林反演化简式子，要注意在模998244353意义下5没有二次剩余，所以每个数都要用$a+b\sqrt{5}$的形式表示，运算类似复数。

斯特林反演的几个用法：

1.下降幂转幂：连续求和时可以通过等比数列求和公式加速。

2.幂转下降幂：类似自然数幂和地用有限微积分加速，或帮助设计DP状态。

 1 #include<cstdio>
 2 #include<algorithm>
 3 #define rep(i,l,r) for (int i=(l); i<=(r); i++)
 4 typedef long long ll;
 5 using namespace std;
 6 
 7 const int N=510,mod=998244353,i2=499122177,i5=598946612,i6=166374059;
 8 
 9 ll V,L,R;
10 int T,m,K,S[N][N],C[N][N];
11 
12 struct num{ int a,b; }z;
13 num operator +(const num &a,const num &b){ return (num){(a.a+b.a)%mod,(a.b+b.b)%mod}; }
14 num operator -(const num &a,const num &b){ return (num){(a.a-b.a+mod)%mod,(a.b-b.b+mod)%mod}; }
15 num operator *(const num &a,const num &b){ return (num){int((1ll*a.a*b.a+V*a.b*b.b)%mod),int((1ll*a.b*b.a+1ll*a.a*b.b)%mod)}; }
16 num operator *(const num &a,int b){ return (num){int(1ll*a.a*b%mod),int(1ll*a.b*b%mod)}; }
17 
18 int ksm(int a,ll b){
19     int res=1;
20     for (; b; a=1ll*a*a%mod,b>>=1)
21         if (b & 1) res=1ll*res*a%mod;
22     return res;
23 }
24 
25 num ksm(num a,ll b){
26     num res=(num){1,0};
27     for (; b; a=a*a,b>>=1)
28         if (b & 1) res=res*a;
29     return res;
30 }
31 
32 num Inv(num a){ return (num){a.a,(mod-a.b)%mod}*ksm((1ll*a.a*a.a-V*a.b*a.b%mod+mod)%mod,mod-2); }
33 
34 int main(){
35     freopen("calc.in","r",stdin);
36     freopen("calc.out","w",stdout);
37     scanf("%d%d",&T,&m);
38     if (m==2) V=5; else V=3;
39     while (T--){
40         scanf("%lld%lld%d",&L,&R,&K); S[0][0]=1;
41         rep(i,1,K) rep(j,1,i) S[i][j]=(S[i-1][j-1]-1ll*S[i-1][j]*(i-1)%mod+mod)%mod;
42         rep(i,0,K) C[i][0]=1;
43         rep(i,1,K) rep(j,1,i) C[i][j]=(C[i-1][j]+C[i-1][j-1])%mod;
44         num a,b,A,B; ll len=R-L+1;
45         if (m==2) L++,R++,a=(num){i2,i2},b=(num){i2,mod-i2},A=(num){0,i5},B=(num){0,mod-i5};
46         else R>>=1,L=(L+1)>>1,a=(num){2,1},b=(num){2,mod-1},A=(num){i2,i6},B=(num){i2,mod-i6};
47         int ans=0; ll t=R-L;
48         rep(i,0,K){
49             num res=z;
50             rep(j,0,i){
51                 num s=ksm(A,j)*ksm(B,i-j)*C[i][j];
52                 num p=ksm(ksm(a,L),j)*ksm(ksm(b,L),(i-j));
53                 num q=ksm(a,j)*ksm(b,i-j);
54                 num w=p*((ksm(q,t+1)-(num){1,0})*Inv(q-(num){1,0}));
55                 if (q.a==1 && q.b==0) w=p*((R-L+1)%mod);
56                 res=res+s*w;
57             }
58             ans=(ans+1ll*res.a*S[K][i])%mod;
59         }
60         int Ans=1ll*ans*ksm(len%mod,mod-2)%mod;
61         rep(i,1,K) Ans=1ll*Ans*ksm(i,mod-2)%mod;
62         printf("%d\n",Ans);
63     }
64     return 0;
65 }

 

转载于:https://www.cnblogs.com/HocRiser/p/10799332.html