MinimumCostPath (AOJ2445)

なんとなくそうっぽいなあ というのはすぐわかるんだけど, ちゃんとそれでいいかどうかというのがなかなか難しかった.

基本的に,あんまりでかいと,x+y=200と(N-x)+(N-y)=200とかで線を引いた時に(線a,bとする),a上の点からb上の点に行く時に遠回りをするやつは絶対に最短経路ではない というのがわかる(a上のある点pの右上とかに障害物が並んでたりするとそっからはそういう経路は取れないんだけど,そういう場合はそもそもそこを通る奴は最短経路にならない というのがわかるっぽい 示すべきことはpとの距離以下の距離で,そのpの右上の邪魔をスルーできるような点が存在することを示せばOKなんだけど,そんなに自明ではないよね・・・ イメージとしてはx+y=200の直前でいろいろ障害物を置かないとpのまわりは割と自由になるからOKで,おいた場合はそもそもそこら辺一帯を避けて行くような経路が存在する という感じ)

怪しかったのであとで解説見たんだけど,今の点が全然説明されてない気がする・・・ いざ証明書こうとすると難しいので飛ばしたのかなあ

設定も解法もすっきりしてて好き、な割には証明や実装がまあまあめんどかった.(というか証明はわかってない)

#include <bits/stdc++.h>
using namespace std;
#define rep(i,n) for(int i=0;i<(n);++i)
#define rep1(i,n) for(int i=1;i<=(n);++i)
#define all(c) (c).begin(),(c).end()
#define pb push_back
#define fs first
#define sc second
#define show(x) cout << #x << " = " << x << endl
#define chmin(x,y) x=min(x,y)
typedef long long ll;
typedef pair<int,ll> P;
ll mod=1e9+9;
void add(ll &x,ll y){
	x+=y;
	x%=mod;
}
int inf=1e9;
int dx[4]={0,1,0,-1},dy[4]={1,0,-1,0};
int N,M,x[50],y[50];
set<P> st;
P ds[100][100],dt[100][100];
P operator+ (const P& a,const P& b){
	if(a.fs<b.fs) return a;
	if(a.fs==b.fs) return P(a.fs,(a.sc+b.sc)%mod);
	return b;
}
bool is(int x,int y){
	return 0<=x&&x<N&&0<=y&&y<N&&st.find(P(x,y))==st.end();
}
void solve1(){
	rep(i,M) st.insert(P(x[i],y[i]));
	P d[300][300];
	rep(i,N) rep(j,N) d[i][j]=P(inf,0);
	queue<int> qx,qy;
	d[0][0]=P(0,1);
	qx.push(0),qy.push(0);
	while(!qx.empty()){
		int x=qx.front(),y=qy.front();qx.pop(),qy.pop();
		rep(di,4){
			int nx=x+dx[di],ny=y+dy[di];
			if(is(nx,ny)){
				if(d[nx][ny].fs==inf){qx.push(nx),qy.push(ny);}
				P q=P(d[x][y].fs+1,d[x][y].sc);
				d[nx][ny]=d[nx][ny]+q;
			}
		}
	}
	cout<<d[N-1][N-1].sc<<endl;
}
const int X=2000000;
ll f[X+1],inv[X+1],g[X+1];
ll Co(int x,int y){
	if(x<0||y<0) return 0;
	return f[x+y]*g[x]%mod*g[y]%mod;
}
void solve2(){
	rep(i,M) st.insert(P(x[i],y[i]));
	{
		rep(i,100) rep(j,100) ds[i][j]=P(inf,0);
		queue<int> qx,qy;
		ds[0][0]=P(0,1);
		qx.push(0),qy.push(0);
		while(!qx.empty()){
			int x=qx.front(),y=qy.front();qx.pop(),qy.pop();
			rep(di,4){
				int nx=x+dx[di],ny=y+dy[di];
				if(is(nx,ny)&&nx+ny<100){
					if(ds[nx][ny].fs==inf){qx.push(nx),qy.push(ny);}
					P q=P(ds[x][y].fs+1,ds[x][y].sc);
					ds[nx][ny]=ds[nx][ny]+q;
				}
			}
		}
	}
	{
		rep(i,100) rep(j,100) dt[i][j]=P(inf,0);
		queue<int> qx,qy;
		dt[0][0]=P(0,1);
		qx.push(0),qy.push(0);
		while(!qx.empty()){
			int x=qx.front(),y=qy.front();qx.pop(),qy.pop();
			rep(di,4){
				int nx=x+dx[di],ny=y+dy[di];
				if(is(N-1-nx,N-1-ny)&&nx+ny<100){
					if(dt[nx][ny].fs==inf){qx.push(nx),qy.push(ny);}
					P q=P(dt[x][y].fs+1,dt[x][y].sc);
					dt[nx][ny]=dt[nx][ny]+q;
				}
			}
		}
	}
	ll d[250]={};
	vector<P> ps;
	int mn=inf,mn2=inf;
	rep(i,100) chmin(mn,ds[i][99-i].fs);
	rep(i,100) chmin(mn2,dt[i][99-i].fs);
	if(mn==inf||mn2==inf){
		puts("0");
		return;
	}
	rep(i,100) if(dt[i][99-i].fs==mn2) ps.pb(P(N-1-i,N-1-(99-i)));
	rep(i,M){
		int t=x[i]+y[i];
		if(100<=t&&100<=(N-1)*2-t) ps.pb(P(x[i],y[i]));
	}
	sort(all(ps));
	ps.erase(unique(all(ps)),ps.end());
	int K=ps.size();
	ll ans=0;
	rep(i,K){
		int x=ps[i].fs,y=ps[i].sc;
		rep(j,100) if(ds[j][99-j].fs==mn){
			add(d[i],ds[j][99-j].sc*Co(x-j,y-(99-j)));
		}
		rep(j,i){
			add(d[i],-d[j]*Co(x-ps[j].fs,y-ps[j].sc));
		}
		if((N-1)*2-x-y<100){
			add(ans,d[i]*dt[N-1-x][N-1-y].sc);
		}
	}
	cout<<(ans+mod)%mod<<endl;
}
int main(){
	f[0]=1;
	rep1(i,X) f[i]=f[i-1]*i%mod;
	inv[1]=1;
	for(ll i=2;i<=X;i++) inv[i]=mod-mod/i*inv[mod%i]%mod;
	g[0]=1;
	rep1(i,X) g[i]=g[i-1]*inv[i]%mod;
	cin>>N>>M;
	rep(i,M){
		cin>>x[i]>>y[i];
		x[i]--,y[i]--;
	}
	if(N<=300) solve1();
	else solve2();
}