d
p
(
i
,
j
)
dp(i,j)
dp(i,j) 表示以
a
a
a 为起点、
i
i
i 为终点,路径长度为
j
j
j 的路径数。
容易得到状态转移方程:
d
p
(
u
,
j
)
=
∑
(
v
,
u
)
d
p
(
v
,
j
−
1
)
dp(u,j)=sum_{(v,u)}dp(v,j-1)
dp(u,j)=∑(v,u)dp(v,j−1)。
但是我们发现一个问题:题目中要求“人物不能立刻沿着刚刚走来的路走回”,但显然这种 dp 方式是不符合要求的,所以想怎么解决。
这里用到一个常用的小trick:我们先把双向边拆成两条单向边,共
2
m
2m
2m 条,并编号。然后设
d
p
(
i
,
j
)
dp(i,j)
dp(i,j) 表示以
a
a
a 为起点、刚好走完编号为
i
i
i 的这条单向边,路径长度为
j
j
j 的路径数。你也可以理解成把边当成点。
状态转移方程同样也容易得到:

如图,设单向边
(
u
,
v
)
(u,v)
(u,v) 的编号为
i
i
i,所有终点为
u
u
u 的单向边的边集为
S
i
S_i
Si(也就是红圈圈出来的部分)。
显然有:
d
p
(
i
,
j
)
=
∑
k
∈
S
i
d
p
(
k
,
j
−
1
)
dp(i,j)=sum_{kin S_i}dp(k,j-1)
dp(i,j)=∑k∈Sidp(k,j−1)。
把 dp 的第二维滚动一下,即可
O
(
m
t
)
O(mt)
O(mt) 获得 30pts。
发现问题的瓶颈主要在
t
t
t 上,思考如何优化。
仔细观察 dp 式,发现和 【NOI Online #3 提高组】魔法值 的做法很相似。
于是联想到用矩阵快速幂来加速 dp 过程。
设计矩阵:(第二维被滚动掉了)
[
d
p
(
1
)
d
p
(
2
)
⋯
d
p
(
2
m
)
]
egin{bmatrix} dp(1)\ dp(2)\ cdots\ dp(2m) end{bmatrix}
⎣⎢⎢⎡dp(1)dp(2)⋯dp(2m)⎦⎥⎥⎤
然后根据上面推出的式子:
d
p
(
i
,
j
)
=
∑
k
∈
S
d
p
(
k
,
j
−
1
)
dp(i,j)=sum_{kin S}dp(k,j-1)
dp(i,j)=∑k∈Sdp(k,j−1),可以设计出转移矩阵:(这里的
c
i
,
j
c_{i,j}
ci,j 表示
j
j
j 是否属于集合
S
i
S_i
Si,属于则
c
i
,
j
=
1
c_{i,j}=1
ci,j=1,否则
c
i
,
j
=
0
c_{i,j}=0
ci,j=0)
B
=
[
c
1
,
1
,
c
1
,
2
,
…
,
c
1
,
2
m
c
2
,
1
,
c
2
,
2
,
…
,
c
2
,
2
m
⋯
⋯
⋯
⋯
c
2
m
,
1
,
c
2
m
,
2
,
…
,
c
2
m
,
2
m
]
B=egin{bmatrix} c_{1,1},&c_{1,2},&dots,&c_{1,2m}\ c_{2,1},&c_{2,2},&dots,&c_{2,2m}\ cdots&cdots&cdots&cdots\ c_{2m,1},&c_{2m,2},&dots,&c_{2m,2m}\ end{bmatrix}
B=⎣⎢⎢⎡c1,1,c2,1,⋯c2m,1,c1,2,c2,2,⋯c2m,2,…,…,⋯…,c1,2mc2,2m⋯c2m,2m⎦⎥⎥⎤
设初始矩阵是
A
A
A(表示所有
d
p
(
i
,
1
)
dp(i,1)
dp(i,1) 的值),那么答案矩阵即为:
A
n
s
=
A
×
B
t
−
1
Ans=A imes B^{t-1}
Ans=A×Bt−1,表示的是所有
d
p
(
i
,
t
)
dp(i,t)
dp(i,t) 的值。
那么答案就是终点为
b
b
b 的每一条单向边的
d
p
dp
dp 值的总和。
总时间复杂度为
O
(
m
3
log
t
)
O(m^3log t)
O(m3logt),可以跑过。
顺带提一下我一开始的错误想法:设
d
p
(
i
,
j
)
dp(i,j)
dp(i,j) 表示以
a
a
a 为起点、
i
i
i 为终点,路径长度为
j
j
j 的路径数,然后找到与
i
i
i 距离为
2
2
2 的点
v
v
v,用
d
p
(
v
,
j
−
2
)
dp(v,j-2)
dp(v,j−2) 来更新
d
p
(
i
,
j
)
dp(i,j)
dp(i,j),这样就能保证不会出现
i
i
i 走到一个点又走回头路的情况了。但是下面这种图就能卡掉这个算法:

询问
t
=
3
t=3
t=3,
a
=
2
a=2
a=2,
b
=
3
b=3
b=3。显然答案为
0
0
0,但是上面那种做法的答案是
1
1
1。
至于为什么是错的我不详细阐述了,如果也是像我一样相同想法的可以结合这个例子自己回去模拟一下并想一想。
正解的代码及注释如下:
#include<bits/stdc++.h>
#define N 25
#define M 65
#define mod 45989
using namespace std;
int n,m,t,a,b;
int cnt=1,head[N],to[M<<1],nxt[M<<1];
struct Matrix
{
int a[M<<1][M<<1];
Matrix(){memset(a,0,sizeof(a));}
void init(){for(int i=1;i<=cnt;i++) a[i][i]=1;}
}st,ch;
Matrix operator * (Matrix a,Matrix b)
{
Matrix c;
for(int i=1;i<=cnt;i++)
for(int j=1;j<=cnt;j++)
for(int k=1;k<=cnt;k++)
c.a[i][j]=(c.a[i][j]+a.a[i][k]*b.a[k][j])%mod;
return c;
}
Matrix poww(Matrix a,int b)
{
Matrix ans;
ans.init();
while(b)
{
if(b&1) ans=ans*a;
a=a*a;
b>>=1;
}
return ans;
}
void adde(int u,int v)
{
to[++cnt]=v;
nxt[cnt]=head[u];
head[u]=cnt;
}
void init()
{
for(int i=head[a];i;i=nxt[i])
st.a[1][i]=1;
for(int i=2;i<=cnt;i++)
{
int u=to[i];
for(int j=head[u];j;j=nxt[j])
if(j!=(i^1))
ch.a[i][j]++;
}
}
int main()
{
scanf("%d%d%d%d%d",&n,&m,&t,&a,&b);
a++,b++;
for(int i=1;i<=m;i++)
{
int u,v;
scanf("%d%d",&u,&v);
u++,v++;
adde(u,v),adde(v,u);
}
init();
Matrix tmp=st*poww(ch,t-1);
int ans=0;
for(int i=2;i<=cnt;i++)
if(to[i]==b)
ans=(ans+tmp.a[1][i])%mod;
printf("%d
",ans);
return 0;
}