BZOJ 3782: 上学路线

3782: 上学路线

Time Limit: 10 Sec  Memory Limit: 128 MB
Submit: 211  Solved: 82
[Submit][Status][Discuss]

Description

小C所在的城市的道路构成了一个方形网格，它的西南角为(0,0)，东北角为(N,M)。小C家住在西南角，学校在东北角。现在有T个路口进行施工，小C不能通过这些路口。小C喜欢走最短的路径到达目的地，因此他每天上学时都只会向东或北行走；而小C又喜欢走不同的路径，因此他问你按照他走最短路径的规则，他可以选择的不同的上学路线有多少条。由于答案可能很大，所以小C只需要让你求出路径数mod P的值。

 

Input

第一行，四个整数N、M、T、P。
接下来的T行，每行两个整数，表示施工的路口的坐标。

 

Output

一行，一个整数，路径数mod P的值。
 
 

Sample Input

3 4 3 1019663265
3 0
1 1
2 2

Sample Output

8

HINT

 

1<=N,M<=10^10

0<=T<=200

p=1000003或p=1019663265

 

 

Source

 
[Submit][Status][Discuss]

 

 

Lucas定理 + 中国剩余定理CRT

 

  1 #include <cstdio>
  2 #include <algorithm>
  3 
  4 #define int long long
  5 
  6 const int mxt = 205;
  7 
  8 int n, m, t, p, f[mxt];
  9 
 10 struct point
 11 {
 12     int x, y;
 13 }pt[mxt];
 14 
 15 inline bool cmp(const point &a, const point &b)
 16 {
 17     if (a.x != b.x)
 18         return a.x < b.x;
 19     else
 20         return a.y < b.y;
 21 }
 22 
 23 namespace case1    // p = 1000003
 24 {
 25     const int siz = 1000005;
 26     
 27     int fac[siz], inv[siz];
 28     
 29     int lucas(int a, int b)
 30     {
 31         if (a < b)
 32             return 0;
 33         
 34         if (a < p && b < p)
 35             return fac[a] * inv[b] % p * inv[a - b] % p;
 36         
 37         return lucas(a / p, b / p) * lucas(a % p, b % p) % p;
 38     }
 39     
 40     inline void prework(void)
 41     {
 42         fac[0] = fac[1] = 1;
 43         inv[0] = inv[1] = 1;
 44         
 45         for (int i = 2; i <= p; ++i)
 46         {
 47             fac[i] = fac[i - 1] * i % p;
 48             inv[i] = inv[p % i] * (p - p/i) % p;
 49         }
 50         
 51         for (int i = 2; i <= p; ++i)
 52             inv[i] = inv[i] * inv[i - 1] % p;
 53     }
 54     
 55     inline int calc(int a, int b)
 56     {
 57         return lucas(a + b, a);
 58     }
 59     
 60     inline void main(void)
 61     {
 62         prework();
 63         
 64         int ans = calc(n, m);
 65         
 66         for (int i = 1; i <= t; ++i)
 67         {
 68             f[i] = calc(pt[i].x, pt[i].y);
 69             
 70             for (int j = 1; j < i; ++j)
 71                 if (pt[j].x <= pt[i].x && pt[j].y <= pt[i].y)
 72                     f[i] = (f[i] - f[j] * calc(pt[i].x - pt[j].x, pt[i].y - pt[j].y) % p + p) % p;
 73         }
 74         
 75         for (int i = 1; i <= t; ++i)
 76             if (pt[i].x <= n && pt[i].y <= m)
 77                 ans = (ans - f[i] * calc(n - pt[i].x, m - pt[i].y) % p + p) % p;
 78         
 79         printf("%lld\n", ans);
 80     }
 81 }
 82 
 83 namespace case2    // p = 3 * 5 * 6793 * 10007
 84 {
 85     const int siz = 10010;
 86     
 87     const int mod[4] = 
 88     {
 89         3, 5, 6793, 10007
 90     };
 91     
 92     int fac[4][siz], inv[4][siz];
 93     
 94     inline void prework(void)
 95     {
 96         for (int i = 0; i < 4; ++i)
 97         {
 98             fac[i][0] = fac[i][1] = 1;
 99             inv[i][0] = inv[i][1] = 1;
100             
101             for (int j = 2; j <= mod[i]; ++j)
102             {
103                 fac[i][j] = fac[i][j - 1] * j % mod[i];
104                 inv[i][j] = inv[i][mod[i] % j] * (mod[i] - mod[i]/j) % mod[i];
105             }
106             
107             for (int j = 2; j <= mod[i]; ++j)
108                 inv[i][j] = inv[i][j - 1] * inv[i][j] % mod[i];
109         }
110     }
111     
112     int lucas(int a, int b, int k)
113     {
114         if (a < b)
115             return 0;
116             
117         if (a < mod[k] && b < mod[k])
118             return fac[k][a] * inv[k][b] % mod[k] * inv[k][a - b] % mod[k];
119         
120         return lucas(a / mod[k], b / mod[k], k) * lucas(a % mod[k], b % mod[k], k) % mod[k];
121     }
122     
123     inline int calc(int a, int b)
124     {
125         return (
126             lucas(a + b, a, 0) * 339887755 % p
127         +    lucas(a + b, a, 1) * 407865306 % p
128         +    lucas(a + b, a, 2) * 673070820 % p
129         +    lucas(a + b, a, 3) * 618502650 % p
130         ) % p;
131     }
132     
133     inline void main(void)
134     {
135         prework();
136         
137         int ans = calc(n, m);
138         
139         for (int i = 1; i <= t; ++i)
140         {
141             f[i] = calc(pt[i].x, pt[i].y);
142             
143             for (int j = 1; j < i; ++j)
144                 if (pt[j].x <= pt[i].x && pt[j].y <= pt[i].y)
145                     f[i] = (f[i] - f[j] * calc(pt[i].x - pt[j].x, pt[i].y - pt[j].y) % p + p) % p;
146         }
147         
148         for (int i = 1; i <= t; ++i)
149             if (pt[i].x <= n && pt[i].y <= m)
150                 ans = (ans - f[i] * calc(n - pt[i].x, m - pt[i].y) % p + p) % p;
151         
152         printf("%lld\n", ans);
153     }
154 }
155 
156 signed main(void)
157 {
158     scanf("%lld%lld%lld%lld", &n, &m, &t, &p);
159     
160     for (int i = 1; i <= t; ++i)
161         scanf("%lld%lld", &pt[i].x, &pt[i].y);
162     
163     std::sort(pt + 1, pt + t + 1, cmp);
164     
165     if (p == 1000003)
166         case1::main();
167     else
168         case2::main();
169 }

 

@Author: YouSiki

 

posted @ 2017-02-24 15:41  YouSiki  阅读(242)  评论(0)  编辑  收藏  举报