BZOJ #4767. 两双手

题面传送门
首先我们设走到终点走了$T1$步第一种走法，$T2$步第二种走法，则可得方程组：

$$\begin{cases}T1A_x+T2B_x=E_x\\T1A_y+T2B_y=E_y\end{cases}$$

解得$\begin{cases}T1=\frac{E_yB_x-E_xB_y}{A_yB_x-A_xB_y}\\ T2=\frac{E_yA_x-E_xA_y}{B_yA_x-B_xA_y}\end{cases}$
所以我们可以得出起点到任意点的步数。如果不能整除到终点先特判掉。
然后就变成一个坐标系走路问题：你要走$T1$横向步，$T2$纵向步，且不能经过$O(n)$个形如$(x_i,y_i)$的障碍点，问你方案数。
看上去可以直接设$f_{i,j}$表示走到第$i$步，走了$j$步横向的方案数，但是这样复杂度达到$O(n^4)$不可取。考虑容斥，先将障碍点按照横纵坐标升序排序，然后设$dp_{i,j}$表示第$i$个障碍点，已经选了$j$个障碍点的方案数，两个障碍点之间的转移可以用组合意义$C_{\Delta x+\Delta y}^{\Delta x}$求解，最后统计答案即可，时间复杂度$O(n^3)$，进一步的，因为容斥系数为$-1$，可以看作选了一个就将贡献乘上$-1$，因此可以将后一维省去，时间复杂度降至$O(n^2)$。
所以题目只给500的原因是要预处理$4n^2$范围内的阶乘逆元
code:

#include<bits/stdc++.h>
#define I inline
#define max(a,b) ((a)>(b)?(a):(b))
#define min(a,b) ((a)<(b)?(a):(b))
#define abs(x) ((x)>0?(x):-(x))
#define ll long long
#define db double
#define lb long db
#define N (500+5)
#define M ((N<<5)+5)
#define K (350)
#define mod 1000000007
#define Mod (mod-1)
#define eps (1e-9)
#define ull unsigned ll
#define it iterator
#define Gc() getchar() 
#define Me(x,y) memset(x,y,sizeof(x))
#define Mc(x,y) memcpy(x,y,sizeof(x))
#define d(x,y) ((k+1)*(x)+(y))
#define R(n) (1ll*rand()*rand()%(n)+1)
#define Pc(x) putchar(x)
#define LB lower_bound
#define UB upper_bound
#define PB push_back
using namespace std;
int n,m,k,x,y,z,Ex,Ey,E,Px,Py,H,Ax,Ay,Bx,By,X[N],Y[N];ll dp[N],Ans,frc[N*N*4],Inv[N*N*4];
struct Node{int x,y;}S[N];I bool cmp(Node x,Node y){return x.x^y.x?x.x<y.x:x.y<y.y;}
I ll mpow(ll x,int y=mod-2){ll Ans=1;while(y) y&1&&(Ans=Ans*x%mod),y>>=1,x=x*x%mod;return Ans;}
I ll C(int x,int y){return frc[x+y]*Inv[y]%mod*Inv[x]%mod;}
int main(){
	freopen("1.in","r",stdin);
	int i,j,h;scanf("%d%d%d",&Ex,&Ey,&n);scanf("%d%d%d%d",&Ax,&Ay,&Bx,&By);for(i=1;i<=n;i++) scanf("%d%d",&X[i],&Y[i]);for(frc[0]=Inv[0]=i=1;i<=1000000;i++) frc[i]=frc[i-1]*i%mod;Inv[1000000]=mpow(frc[1000000]);for(i=999999;~i;i--) Inv[i]=Inv[i+1]*(i+1)%mod;
	E=Ay*Bx-Ax*By;if((Ey*Bx-Ex*By)%E||(Ey*Ax-Ex*Ay)%E) {puts("0");return 0;}for(i=1;i<=n;i++){if((Y[i]*Bx-X[i]*By)%E||(Y[i]*Ax-X[i]*Ay)%E) continue;S[++H]=(Node){(Y[i]*Bx-X[i]*By)/E,-(Y[i]*Ax-X[i]*Ay)/E};}
	Px=(Ey*Bx-Ex*By)/E;Py=-(Ey*Ax-Ex*Ay)/E;sort(S+1,S+H+1,cmp);dp[0]=1;for(i=0;i<=H;i++){
		for(j=0;j<i;j++){if(S[j].x>S[i].x||S[j].y>S[i].y) continue;dp[i]-=dp[j]*C(S[i].x-S[j].x,S[i].y-S[j].y)%mod;}dp[i]=(dp[i]%mod+mod)%mod;S[i].x<=Px&&S[i].y<=Py&&(Ans+=dp[i]*C(Px-S[i].x,Py-S[i].y)%mod);
	}printf("%lld\n",Ans%mod);
}