链接：https://www.luogu.com.cn/problem/P4007
题解：可以先想到一个\(dp\)，令\(dp_{i,j,k}\)表示\(i\)个一滴血的随从，\(j\)个两滴血的随从,\(k\)个三滴血的随从的概率，这样可以直接转移。期望不好算，可以再开一维设期望，由于每一次加\(1\),所以将概率相加即期望。
这个东西可以矩阵快速幂进行优化，但优化后\(O(Tn^3logn)\)过不了。
考虑一个倍增，令\(dp_{i}\)表示\(2^i\)个矩阵的相乘结果，那么
\(ANS_{n}=dp_{0}\times ...\times p\)(\(n\)二进制表示)
由于\(p\)只有第\(1\)行有值，那么我们可以只算出第\(1\)行，然后就是一个行向量与若干矩阵相乘了，此时就是\(O(Tn^2logn)\)了。
#include<iostream>
#include<cstdio>
#pragma GCC optimize("Ofast")
#pragma GCC optimize("unroll-loops")
#pragma GCC target("sse,sse2,sse3,ssse3,sse4,popcnt,abm,mmx,avx,tune=native")
#define mod 998244353
using namespace std;
struct node
{
int map[201][201];
};
node e,now[61],c,res;
long long T,m,K,n,M[9][9][9],inv[25],length;
inline long long read()
{
long long sum=0;
char c=0;
while (c<'0'||c>'9')
c=getchar();
while ('0'<=c&&c<='9')
{
sum=sum*10+c-'0';
c=getchar();
}
return sum;
}
inline node mul(node a,node b)
{
for (int i=0;i<=length;++i)
for (int j=0;j<=length;++j)
c.map[i][j]=0;
__int128 t;
for (int i=0;i<=length;++i)
for (int j=0;j<=length;++j)
{
t=0;
for (int k=0;k<=length;++k)
t=(t+1ll*a.map[i][k]*b.map[k][j]);
c.map[i][j]=t%mod;
}
return c;
}
inline node mult(node a,node b)
{
for (int i=0;i<=length;++i)
c.map[0][i]=0;
__int128 t;
for (int i=0;i<=length;++i)
{
t=0;
for (int j=0;j<=length;++j)
t=(t+1ll*a.map[0][j]*b.map[j][i]);
c.map[0][i]=t%mod;
}
return c;
}
inline long long fast_pow(long long a,int b)
{
if (b==0)
return 1;
if (b&1)
return fast_pow(a*a%mod,b/2)*a%mod;
else
return fast_pow(a*a%mod,b/2);
}
signed main()
{
T=read(),m=read(),K=read();
for (int i=0;i<=8;++i)
for (int j=0;j<=8;++j)
for (int k=0;k<=8;++k)
if (i+j+k<=K)
M[i][j][k]=++length;
for (int i=0;i<=length;++i)
e.map[i][i]=1;
if (m==1)
{
length=K+1;
for (int i=0;i<=length;++i)
e.map[i][i]=1;
now[0].map[K+1][K+1]=1;
for (int i=0;i<=K;++i)
{
now[0].map[K+1][i]=fast_pow(i+1,mod-2);
now[0].map[i][i]=fast_pow(i+1,mod-2);
if (i>=1)
now[0].map[i-1][i]=i*fast_pow(i+1,mod-2)%mod;
}
for (int i=1;i<=60;++i)
now[i]=mul(now[i-1],now[i-1]);
while (T--)
{
n=read();
res=e;
for (long long i=0,m=1;i<=60;++i,m*=2)
if (m&n)
res=mul(res,now[i]);
printf("%lld\n",(long long)(res.map[K+1][1]));
}
return 0;
}
if (m==2)
{
length=2*K+1;
for (int i=0;i<=length;++i)
e.map[i][i]=1;
now[0].map[0][0]=1;
for (int i=1;i<=K;++i)
{
now[0].map[0][i]=fast_pow(i+1,mod-2);
if (i+1<=K)
now[0].map[i+1][i]=fast_pow(i+1,mod-2);
else
now[0].map[K+K+1][i]=fast_pow(i+1,mod-2);
now[0].map[i][i]=fast_pow(i+1,mod-2)%mod;
if (i>=2)
now[0].map[i-1][i]=(i-1)*fast_pow(i+1,mod-2)%mod;
}
for (int i=0;i<=K;++i)
{
now[0].map[0][i+K+1]=fast_pow(i+1,mod-2);
now[0].map[i+K+1][i+K+1]=fast_pow(i+1,mod-2);
if (i>=1)
now[0].map[i+K][i+K+1]=i*fast_pow(i+1,mod-2)%mod;
}
for (int i=1;i<=60;++i)
now[i]=mul(now[i-1],now[i-1]);
while (T--)
{
n=read();
res=e;
for (long long i=0,m=1;i<=60;++i,m*=2)
if (m&n)
res=mul(res,now[i]);
printf("%lld\n",(long long)(res.map[0][1]));
}
return 0;
}
now[0].map[0][0]=1;
for (int i=0;i<=24;++i)
inv[i]=fast_pow(i,mod-2);
for (int i=0;i<=8;++i)
for (int j=0;j<=8;++j)
for (int k=0;k<=8;++k)
if (i+j+k<=K)
{
now[0].map[M[i][j][k]][M[i][j][k]]=inv[i+j+k+1];
now[0].map[0][M[i][j][k]]=inv[i+j+k+1];
if (i>=1)
now[0].map[M[i-1][j][k]][M[i][j][k]]=i*inv[i+j+k+1]%mod;
if (j>=1&&i+j+k+1<=K)
now[0].map[M[i+1][j-1][k+1]][M[i][j][k]]=j*inv[i+j+k+1]%mod;
if (j>=1&&i+j+k+1>K)
now[0].map[M[i+1][j-1][k]][M[i][j][k]]=j*inv[i+j+k+1]%mod;
if (i+j+k+1<=K)
now[0].map[M[i][j+1][k]][M[i][j][k]]=k*inv[i+j+k+1]%mod;
if (i+j+k+1>K&&k>=1)
now[0].map[M[i][j+1][k-1]][M[i][j][k]]=k*inv[i+j+k+1]%mod;
}
for (int i=1;i<=60;++i)
now[i]=mul(now[i-1],now[i-1]);
while (T--)
{
n=read();
for (int i=0;i<=length;++i)
for (int j=0;j<=length;++j)
res.map[i][j]=0;
res.map[0][0]=1;
for (long long i=0,m=1;i<=60;++i,m*=2)
if (m&n)
res=mult(res,now[i]);
printf("%lld\n",(long long)(res.map[0][M[0][0][1]]));
}
return 0;
}