P4478 [BJWC2018] 上学路线

思路

看到有 个障碍点并且障碍点非常少,所以想到暴力容斥,用总方案减去不合法的方案数。

只要经过任意一个障碍点,就是一个不合法方案。所以令 表示只经过障碍点 的方案数。

设当前障碍点为 ,如果不管途中是否经过其他障碍点,那么走到该点的总方案数为
但是其中还包含了经过其他障碍点的方案数。

设还经过了障碍点 ,满足 ,那么多统计的方案数为 减掉即可。
一个点能从其左下角的每个点转移过来,所以先把点按从左到右、从上到下排好序,才能保证转移是正确的。

我们发现 都非常大,所以要用 Lucas 定理,
然而 Lucas 定理只适用于 为质数,如果 不是质数怎么办呢。
我们可以把 分解质因数,对于每个质因子都列出一个方程,得到同余方程组,然后用中国剩余定理合并求解就可以了。

我们发现 ,所以不用繁琐的对组合数拆开取模,直接做就可以了。

代码

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
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
#include<iostream>
#include<cstdio>
#include<algorithm>
#include<cstring>
#include<climits>
#include<cstdlib>
using namespace std;
typedef long long ll;
typedef unsigned long long ull;
template <typename T>
inline void read(T&x){// 快读
int w = 0;x = 0;
char ch = getchar();
while(ch<'0' || ch>'9'){
if(ch=='-') w = 1;
ch = getchar();
}
while(ch>='0' && ch<='9'){
x = (x<<1)+(x<<3)+(ch^48);
ch = getchar();
}
if(w) x = ~x+1;
}
template <typename T,typename...Args>
inline void read(T&t,Args&...args){
read(t);read(args...);
}
template <typename T>
inline T Min(T x,T y){ return (x < y ? x : y); }
template <typename T>
inline T Max(T x,T y){ return (x > y ? x : y); }
template <typename T>
inline T Abs(T x){ return (x < 0 ? ~x+1 : x); }
const int N = 210,M = 1e6+10;
ll prime[4] = {3,5,6793,10007};
ll T,ans[N],P; int cnt;
struct node{
ll x,y;
inline int operator < (const node&G) const{
if(x^G.x) return x < G.x;
return y < G.y;
}
}g[N];
ll mul[4][M],inv[4][M];
inline ll quick_pow(ll x,ll y,int t){
ll res = 1;
while(y){
if(y&1) (res *= x)%=prime[t];
(x *= x)%=prime[t];
y >>= 1;
}
return res;
}
inline ll C(ll n,ll m,int t){// 组合数
if(m>n) return 0;
return mul[t][n]*inv[t][n-m]%prime[t]*inv[t][m]%prime[t];
}
inline ll Lucas(ll n,ll m,int t){// Lucas
if(!m) return 1;
return C(n%prime[t],m%prime[t],t)*Lucas(n/prime[t],m/prime[t],t)%prime[t];
}
inline ll Gcd(ll x,ll y){
if(!x || !y) return x|y;
ll xz = __builtin_ctzll(x);
ll yz = __builtin_ctzll(y);
ll z = Min(xz,yz);
ll diff;
y >>= yz;
while(x){
x >>= xz;
diff = x-y;
xz = __builtin_ctzll(diff);
y = Min(x,y);
x = Abs(diff);
}
return y << z;
}
inline void Exgcd(ll x,ll y,ll &Ex,ll &Ey){// 扩欧
if(!y){
Ex = 1;
Ey = 0;
return ;
}
Exgcd(y,x%y,Ey,Ex);
Ey = Ey-x/y*Ex;
}
inline ll ExCRT(ll n,ll m){// ExCRT 合并方程
ll r1 = Lucas(n,m,0);
ll m1 = prime[0];
ll r2,m2,Ex,Ey;
for(int t=1;t<=cnt;++t){
r2 = Lucas(n,m,t);
m2 = prime[t];
ll d = Gcd(m1,m2);
// if((r2-r1)%d) exit(0);
Exgcd(m1/d,m2/d,Ex,Ey);

Ex = ((r2-r1)/d*Ex%(m2/d)+(m2/d))%(m2/d);
r1 = ((r1+m1*Ex)%(m1/d*m2)+(m1/d*m2))%(m1/d*m2);
m1 = m1/d*m2;
}
return r1;
}
int main(){
// freopen("in.in","r",stdin);
// freopen("out.out","w",stdout);

ll n,m;read(n,m,T,P);
for(int i=1;i<=T;++i) read(g[i].x,g[i].y);
g[++T] = {n,m};
sort(g+1,g+1+T);// 排好序
if(P==1000003) prime[0] = 1000003,cnt = 0;
else cnt = 3;
for(int t=0;t<=cnt;++t){
mul[t][0] = 1;
for(ll i=1;i<prime[t];++i) mul[t][i] = mul[t][i-1]*i%prime[t];
inv[t][prime[t]-1] = quick_pow(mul[t][prime[t]-1],prime[t]-2,t);
for(ll i=prime[t]-1;i;--i) inv[t][i-1] = inv[t][i]*i%prime[t];
}
for(int i=1;i<=T;++i){
ans[i] = ExCRT(g[i].x+g[i].y,g[i].x);
for(int j=1;j<i;++j){
if(g[j].x<=g[i].x && g[j].y<=g[i].y){// 左下角的点
ans[i] = (ans[i]-ans[j]*ExCRT(g[i].x+g[i].y-g[j].x-g[j].y,g[i].x-g[j].x)%P+P)%P;
}
}
}
printf("%lld",ans[T]);

fclose(stdin);
fclose(stdout);
return 0;
}