P4102 [HEOI2014] 林中路径 题解

设 $A[i][j]$ 为从 $i$ 点到 $j$ 点的边的数量。

设:

$$ \begin{aligned} E_0[L] & = \sum_{i = 1}^LA^i\\ E_1[L] & = \sum_{i = 1}^Li\times A^i\\ E_2[L] & = \sum_{i = 1}^Li^2\times A^i\\ \end{aligned} $$

答案就是 $E_2[k]$ 。

考虑用分治求解这三个矩阵。

当 $L$ 是偶数时:

$$ \begin{aligned} E_0[L] & = \sum_{i = 1}^LA^i\\& = \sum_{i = 1}^{L/2}A^i +A^{L/2}\times \sum_{i = 1}^{L/2}A^i\\&=E_0[L/2]+A^{L/2}\times E_0[L/2]\\ E_1[L]& = \sum_{i = 1}^Li\times A^i\\ & = \sum_{i = 1}^{L/2}i\times A^i+A^{L/2}\times \sum_{i=1}^{L/2}[(i+L/2)\times A^i]\\ &=E_1[L/2]+A^{L/2}\times(\sum_{i=1}^{L/2}i\times A^i+L/2\times\sum_{i=1}^{L/2}A^i)\\ &=E_1[L/2]+A^{L/2}\times(E_1[L/2]+L/2\times E_0[L/2])\\ E_2[L]& = \sum_{i = 1}^Li^2\times A^i\\ & = \sum_{i = 1}^{L/2}i^2\times A^i+A^{L/2}\times \sum_{i=1}^{L/2}[(i+L/2)^2\times A^i]\\ & = \sum_{i = 1}^{L/2}i^2\times A^i+A^{L/2}\times \sum_{i=1}^{L/2}[(i^2+2i(L/2)+(L/2)^2\times A^i]\\ &=E_2[L/2]+A^{L/2}\times(\sum_{i=1}^{L/2}i^2\times A^i+L\times\sum_{i=1}^{L/2}i\times A^i+L^2/4\times\sum_{i=1}^{L/2}A^i)\\ &=E_2[L/2]+A^{L/2}\times(E_2[L/2]+L\times E_1[L/2]+L^2/4\times E_0[L/2])\\ \end{aligned} $$
当 $L$ 是奇数时:
$$ \begin{aligned} E_0[L] & = \sum_{i = 1}^LA^i\\ &=A+A\times\sum_{i = 1}^{L-1}A^i\\&=A+A\times E_0[L-1]\\ E_1[L] & = \sum_{i = 1}^Li\times A^i\\& =A+A\times\sum_{i = 1}^{L-1}((i+1)\times A^i)\\ &=A+A\times(\sum_{i = 1}^{L-1}i\times A^i+\sum_{i = 1}^{L-1} A^i)\\ &=A+A\times (E_1[L-1]+E_0[L-1])\\ E_2[L] & = \sum_{i = 1}^Li^2\times A^i\\ &=A+A\times\sum_{i = 1}^{L-1}((i+1)^2\times A^i)\\ &=A+A\times\sum_{i = 1}^{L-1}((i^2+2i+1)\times A^i)\\ &=A+A\times(\sum_{i = 1}^{L-1}i^2\times A^i+2\sum_{i = 1}^{L-1} i\times A^i+\sum_{i = 1}^{L-1} A^i)\\ &=A+A\times (E_2[L-1]+2E_1[L-1]+E_0[L-1])\\ \end{aligned} $$
时间复杂度 $O(n^3\log_k)$

代码:

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
#include<bits/stdc++.h>
#define per(x,y,z) for(int x=(y);(x)<=(z);++(x))
#define ll long long
using namespace std;
const int N=110,p=1e9+7;
int n,m,k,q;
struct mx{ll s[N][N];mx(){memset(s,0,sizeof(s));}}A,E0,E1,E2,s;
mx operator*(mx a,mx b){
    mx c;per(k,1,n)per(i,1,n)per(j,1,n) 
        c.s[i][j]=(c.s[i][j]+(a.s[i][k]*b.s[k][j])%p)%p;
    return c;
}
mx operator+(mx a,mx b){
    per(i,1,n) per(j,1,n) 
        a.s[i][j]=(a.s[i][j]+b.s[i][j])%p;
    return a;
}
mx operator*(mx a,ll b){
    per(i,1,n) per(j,1,n) 
        a.s[i][j]=(a.s[i][j]*b)%p;
    return a;
}
void solve(int L){
    if(L==1){s=E0=E1=E2=A;return;}
    if(L&1){
        solve(L-1);
        E2=A+(A*(E2+E1*2+E0));
        E1=A+(A*(E1+E0));
        E0=A+(A*E0);
        s=s*A;
    }else{
        solve(L/2);
        ll L2=L/2;
        E2=E2+s*(E2+E1*L+E0*(L2*L2%p));
        E1=E1+s*(E1+E0*(L2%p));
        E0=E0+(s*E0);
        s=s*s;
    }
}
signed main(){
    ios::sync_with_stdio(false);
    cin>>n>>m>>k>>q;
    per(i,1,m){
        int x,y;cin>>x>>y;
        A.s[x][y]++;
    }
    solve(k);
    while(q--){
        int s,t;cin>>s>>t;
        cout<<E2.s[s][t]<<'\n';
    }
}