luoguP4102 [HEOI2014]林中路径

luoguP4102 [HEOI2014]林中路径

原题传送门:>Here<

当路径长度$=L$,且只要求路径条数时,我们是很容易做的,只要用邻接矩阵乘一下就好了。

然而此题需要路径长度$\le L$的路径的长度平方和,感觉些微复杂。

我们先考虑$\le L$怎么做。

定义一个矩阵二元组$(A,B)$,其中$A$表示$=L$的路径条数,$B$表示$\le L$的路径条数。

定义二元组的乘法为:

$$(A,B)\times (C,D)=(A\times C,B+A\times D)$$

正确性比较显然:

$A\times C$显然正确。

设二元组$(A,B)$表示的长度为$L_0$,二元组$(C,D)$表示的长度为$L_1$。

则$\le L_0$的路径可以用$B$表示。

对于长度在$[L_0+1,L_1]$中的路径条数,很明显就是长度$=L_0$的矩阵乘上长度$\le L_1$的矩阵,即$B+A\times D$。

至于求平方和,我们不妨设矩阵的每一个元素为一个三元组$(A,B,C)$,其中$A$表示路径条数,$B$表示路径长度和,$C$表示路径长度的平方和。这也就分别是路径长度的$0,1,2$次方和。

我们考虑这样一个三元组的乘法:

对于两个三元组$F,G$,首先$(F\times G)_A=F_A\times G_A$显然。

对于$(F\times G)_B$,我们发现$F$中的所有路径会与$G$中的所有路径一一匹配,即$F_B$会被算$G_A$次,$G_B$同理。

所以$(F\times G)_B=F_A\times G_B+G_A\times F_B$。

对于$(F\times G)_C$,我们不妨推推式子:

设$F$中有 $num_F$ 条路径,长度分别为 $f_1,f_2,\cdots,f_{num_F}$ 。 $G$ 同理。

$$
(F\times G)C \
= \sum
{i=1}^{num_F}\sum_{j=1}^{num_G}(f_i+g_j)^2 \
= \sum_{i=1}^{num_F}\sum_{j=1}^{num_G}f_i^2+2f_ig_j+g_j^2 \
= F_C\times G_A+2F_B\times G_B+F_A\times G_C
$$

代码:

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
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
// luogu-judger-enable-o2
#include <bits/stdc++.h>

typedef long long ll;

namespace IO
{
const unsigned int Buffsize=1<<25,Output=1<<25;
static char Ch[Buffsize],*St=Ch,*T=Ch;
inline char getc()
{
return((St==T)&&(T=(St=Ch)+fread(Ch,1,Buffsize,stdin),St==T)?0:*St++);
}
static char Out[Output],*nowps=Out;

inline void flush(){fwrite(Out,1,nowps-Out,stdout);nowps=Out;}

inline int read()
{
int x=0;static char ch;int f=1;
for(ch=getc();!isdigit(ch);ch=getc())if(ch=='-')f=-1;
for(;isdigit(ch);ch=getc())x=x*10+(ch^48);
return x*f;
}

template<typename T>inline void write(T x,char ch='\n')
{
if(!x)*nowps++=48;
if(x<0)*nowps++='-',x=-x;
static unsigned int sta[111],tp;
for(tp=0;x;x/=10)sta[++tp]=x%10;
for(;tp;*nowps++=sta[tp--]^48);
*nowps++=ch;
}
}
using namespace IO;

int P=1000000007,n,m,k,p,s,t;
struct matrix{
ll A[100][100],B[100][100],C[100][100];
inline matrix(){memset(A,0,sizeof A);memset(B,0,sizeof B);memset(C,0,sizeof C);}
inline void set(const int &x,const int &y){++A[x][y],++B[x][y],++C[x][y];}
inline matrix operator*(matrix a){
static matrix ans;
ans=matrix();
for(register int k=0;k<n;++k)
for(register int i=0;i<n;++i)
for(register int j=0;j<n;++j){
ans.A[i][j]=(ans.A[i][j]+1ll*A[i][k]*a.A[k][j]);if(ans.A[i][j]>6e18)ans.A[i][j]%=P;//玄学
ans.B[i][j]=(ans.B[i][j]+1ll*A[i][k]*a.B[k][j]+1ll*B[i][k]*a.A[k][j]);if(ans.B[i][j]>4e18)ans.B[i][j]%=P;//玄学
ans.C[i][j]=(ans.C[i][j]+1ll*A[i][k]*a.C[k][j]+2ll*B[i][k]*a.B[k][j]+1ll*C[i][k]*a.A[k][j]);if(ans.C[i][j]>3.5e18)ans.C[i][j]%=P;//玄学
}
for(int i=0;i<n;i++)
for(int j=0;j<n;j++){
ans.A[i][j]%=P;
ans.B[i][j]%=P;
ans.C[i][j]%=P;
}
return ans;
}
inline matrix operator+(matrix a){
static matrix ans;
for(register int i=0;i<n;++i)
for(register int j=0;j<n;++j){
ans.A[i][j]=(A[i][j]+a.A[i][j])%P;
ans.B[i][j]=(B[i][j]+a.B[i][j])%P;
ans.C[i][j]=(C[i][j]+a.C[i][j])%P;
}
return ans;
}
}A,B,C;
int main(){
n=read(),m=read(),k=read(),p=read();
for(int i=1;i<=m;++i){
s=read()-1,t=read()-1;
++A.A[s][t];
++A.B[s][t];
++A.C[s][t];
}
B=A;
if(k)--k,C=A;
while(k){
if(k&1)C=C*A+B;
B=B*A+B,A=A*A;
k>>=1;
}
for(register int i=1;i<=p;++i){
s=read()-1,t=read()-1;
write(C.C[s][t]);
}
return flush(),0;
}

(常数太不优秀,只能靠玄学过去了)

Comments

Your browser is out-of-date!

Update your browser to view this website correctly. Update my browser now

×