Bzoj3782 上学路线
Submit: 239 Solved: 95
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
3 0
1 1
2 2
Sample Output
8
HINT
1<=N,M<=10^10
0<=T<=200
p=1000003或p=1019663265
Source
数学问题 容斥 动规 组合数 lucas定理 中国剩余定理 脑洞题
Bzoj4767 两双手 ←我还以为这道上学路线是Bzoj4767的简化版,欢脱地开了坑。
神特么是超超超超强化版?!
刚开始的时候情况看上去挺好,博主欢快地敲着lucas定理的板子
直到博主注意到1019663265是个合数
exm?合数怎么玩?
想到了中国剩余定理强行合并,但是不会真的这么麻烦吧……怂一波看了题解,真这么麻烦啊……
把数字丢到wolframalpha,得到了这个:
$1019663265=3*5*6793*10007$
只有四个质数,也就是说只要调用四次函数再合并即可,不用多写一个完全体的CRT板子。
于是就搞出了长长的代码:
1、如果模数是1000003,可以直接用lucas定理算容斥,解法同Bzoj4767
2、如果模数是1019663265,需要在那四个质因数的模意义下分别算一遍,用中国剩余定理合并答案。
合并的过程可以暴力写四行,学自popoQQQ的写法
3、基本原理就是这么简单,剩下的就是码农的工作了
刚开始写出来,编译器提示代码里有C++11限定的写法,丢到B站一测果然CE,然后换成了现在这样。交了一发WA了,心中无比绝望,感觉又要调一下午,好在只是有个地方没开LL
LL好啊!(意味深)
PS:写了那样的结构体名以表心情
1 /*by SilverN*/ 2 #include<algorithm> 3 #include<iostream> 4 #include<cstring> 5 #include<cstdio> 6 #include<cmath> 7 #include<vector> 8 #define LL long long 9 using namespace std; 10 const int mxn=110010; 11 LL read(){ 12 LL x=0,f=1;char ch=getchar(); 13 while(ch<'0' || ch>'9'){if(ch=='-')f=-1;ch=getchar();} 14 while(ch>='0' && ch<='9'){x=x*10+ch-'0';ch=getchar();} 15 return x*f; 16 } 17 struct point{ 18 LL x,y;bool ban; 19 bool operator < (const point b)const{ 20 return (x<b.x)||(x==b.x && y<=b.y); 21 } 22 }s[mxn]; 23 int n; 24 LL Ex,Ey,p; 25 // 26 struct SKY{ 27 LL fac[mxn*10]; 28 LL inv[mxn*10]; 29 int P; 30 void init(){ 31 P=1000003; 32 fac[1]=1;inv[1]=1; 33 fac[0]=1;inv[0]=1; 34 for(int i=2;i<p;i++){ 35 fac[i]=(LL)fac[i-1]*i%P; 36 inv[i]=((-P/i)*inv[P%i]%P+P)%P; 37 } 38 for(int i=2;i<p;i++) inv[i]=inv[i]*inv[i-1]%P; 39 return; 40 } 41 LL calc(LL n,LL m){return fac[n]*inv[m]%P*inv[n-m]%P;} 42 LL lucas(LL n,LL m){ 43 if(!m)return 1; 44 return calc(n%P,m%P)*lucas(n/P,m/P)%P; 45 } 46 }SK; 47 struct FXXKER{ 48 LL fac[mxn][4]; 49 LL inv[mxn][4]; 50 LL P[5]; 51 void init(){ 52 P[0]=3; 53 P[1]=5; 54 P[2]=6793; 55 P[3]=10007; 56 P[4]=1019663265; 57 for(int k=0;k<4;k++){ 58 fac[1][k]=1;inv[1][k]=1; 59 fac[0][k]=1;inv[0][k]=1; 60 for(int i=2;i<mxn;i++){ 61 fac[i][k]=(LL)fac[i-1][k]*i%P[k]; 62 inv[i][k]=((-P[k]/i)*inv[P[k]%i][k]%P[k]+P[k])%P[k]; 63 } 64 for(int i=2;i<mxn;i++) inv[i][k]=inv[i][k]*inv[i-1][k]%P[k]; 65 } 66 return; 67 } 68 LL C(LL n,LL m,int k){ 69 return fac[n][k]*inv[m][k]%P[k]*inv[n-m][k]%P[k]; 70 } 71 LL lucas(LL n,LL m,int k){ 72 if(!m)return 1; 73 return C(n%P[k],m%P[k],k)*lucas(n/P[k],m/P[k],k)%P[k]; 74 } 75 LL Crt(LL n,LL m){ 76 LL r0=lucas(n,m,0); LL r1=lucas(n,m,1); 77 LL r2=lucas(n,m,2); LL r3=lucas(n,m,3); 78 return (r0*339887755+r1*407865306+r2*673070820+r3*618502650)%P[4]; 79 } 80 }FK; 81 82 LL calc(LL n,LL m){ 83 if(p==1000003)return SK.lucas(n,m); 84 else return FK.Crt(n,m); 85 } 86 LL f[300]; 87 void solve(){ 88 int i,j; 89 for(i=1;i<=n;i++){ 90 f[i]=calc(s[i].x+s[i].y,s[i].x); 91 // printf("i:%d X+Y:%lld Y:%lld f:%lld\n ",i,s[i].x+s[i].y,s[i].x,f[i]); 92 for(j=1;j<i;j++){ 93 if(s[j].x<=s[i].x && s[j].y<=s[i].y) 94 f[i]=((f[i]-f[j]*calc(s[i].x-s[j].x+s[i].y-s[j].y,s[i].x-s[j].x)%p)+p)%p; 95 } 96 } 97 printf("%lld\n",f[n]); 98 return; 99 } 100 int main(){ 101 int i,j; 102 Ex=read();Ey=read();n=read();p=read(); 103 if(p==1000003) SK.init(); 104 else FK.init(); 105 for(i=1;i<=n;i++){ 106 s[i].x=read();s[i].y=read(); 107 } 108 s[++n].x=Ex;s[n].y=Ey; 109 sort(s+1,s+n+1); 110 solve(); 111 return 0; 112 }
本文为博主原创文章,转载请注明出处。