立命館合宿2012 day3 問題G (AOJ 2368)

問題

日本語なので本文参照（http://judge.u-aizu.ac.jp/onlinejudge/description.jsp?id=2368）

制約条件

n≦50
wi≦100
m≦10
0≦xi,yi≦100

方針

まずは、線分の端点および、線分同士の交点を全列挙する。
この点を頂点とするグラフを作る。

それぞれの線分について、
自分の上に乗っている点を見つけてきて、
隣接する点同士に、対応する容量を加算する。

全部の加算が終わったら、容量をもとにグラフを作って、
スタートからゴールへ最大流を流す。

EPSを変えまくったら通った。
spaghetti sourceのintersectSPは不正確（？）なので使うと通らないっぽい。
ccwを使って判定するとよい。
のだが、ccwにEPSがないのでそのままでは通らず、ccwにEPSを入れたら通った。

ソースコード

const double EPS=1e-4, INF=1e12;
//最大流のコードは省略
void add_edge(Graph &g, int s,int t,int w){
	g[s].pb(Edge(s,t,w));
	g[t].pb(Edge(t,s,0));
}


typedef complex<double> P;
namespace std {
  bool operator < (const P& a, const P& b) {
    return real(a) != real(b) ? real(a) < real(b) : imag(a) < imag(b);
  }
}
double cross(const P& a, const P& b) {
  return imag(conj(a)*b);
}
double dot(const P& a, const P& b) {
  return real(conj(a)*b);
}
struct L : public vector<P> {
  L(const P &a, const P &b) {
    push_back(a); push_back(b);
  }
};
typedef vector<P> G;
struct C {
  P p; double r;
  C(const P &p, double r) : p(p), r(r) { }
};
int ccw(P a, P b, P c) {
  b -= a; c -= a;
  if (cross(b, c) > EPS)   return +1;       // counter clockwise
  if (cross(b, c) < -EPS)   return -1;       // clockwise
  if (dot(b, c) < EPS)     return +2;       // c--a--b on line
  if (norm(b) < norm(c)) return -2;       // a--b--c on line
  return 0;
}
bool intersectLS(const L &l, const L &s) {
  return cross(l[1]-l[0], s[0]-l[0])*       // s[0] is left of l
         cross(l[1]-l[0], s[1]-l[0]) < EPS; // s[1] is right of l
}
bool intersectSS(const L &s, const L &t) {
  return ccw(s[0],s[1],t[0])*ccw(s[0],s[1],t[1]) <= 0 &&
         ccw(t[0],t[1],s[0])*ccw(t[0],t[1],s[1]) <= 0;
}
bool intersectSP(const L &s, const P &p) {
	rep(i,2)if(abs(s[i]-p)<EPS)return 1;
	return ccw(p,s[1],s[0])==2;
}
P crosspoint(const L &l, const L &m) {
  double A = cross(l[1] - l[0], m[1] - m[0]);
  double B = cross(l[1] - l[0], l[1] - m[0]);
  if (abs(A) < EPS && abs(B) < EPS) return m[0]; // same line
  if (abs(A) < EPS) assert(false); // !!!PRECONDITION NOT SATISFIED!!!
  return m[0] + B / A * (m[1] - m[0]);
}

map<P,int> id;
void push(vector<P> &ps,const P &p){
	if(id.count(p))return;
	rep(i,ps.size())if(abs(ps[i]-p)<EPS){
		id[p]=i;
		return;
	}
	id[p]=ps.size();
	ps.pb(p);
}
struct cmp{
	P o;
	cmp(P p):o(p){}
	bool operator()(const P&a,const P&b)const{
		return abs(a-o)<abs(b-o);
	}
};


int n,m,w[500];

int main(){
	vector<L> ls;
	vector<P> ps;
	cin>>n;
	rep(i,n){
		double ax,ay,bx,by;
		cin>>ax>>ay>>bx>>by>>w[i];
		ls.pb(L(P(ax,ay),P(bx,by)));
		push(ps,P(ax,ay)); push(ps,P(bx,by));
	}
	cin>>m;
	vector<P> ss;
	double x,y;
	rep(i,m){
		cin>>x>>y;
		ss.pb(P(x,y));
		push(ps,P(x,y));
	}
	cin>>x>>y;
	P g(x,y);
	push(ps,g);
	
	rep(i,n)rep(j,i)if(intersectSS(ls[i],ls[j])){
		P cp=crosspoint(ls[i],ls[j]);
		push(ps,cp);
	}
	
	int V=ps.size();
	
	Graph graph(V+1);
	int s=V, t=id[g];
	rep(i,m){
		add_edge(graph,s,id[ss[i]],inf);
	}
	
	vector<vector<double> > width(V,vector<double>(V));
	rep(i,n){
		vector<P> online;
		rep(j,V)if(intersectSP(ls[i],ps[j]))online.pb(ps[j]);
		sort(all(online),cmp(ls[i][0]));
		const double dist=abs(ls[i][0]-ls[i][1]);
		rep(j,online.size()-1){
			double tw=w[i]*abs(online[j]-online[j+1])/dist;
			width[id[online[j]]][id[online[j+1]]]+=tw;
			width[id[online[j+1]]][id[online[j]]]+=tw;
		}
	}
	rep(i,V)rep(j,V)if(width[i][j]>EPS)add_edge(graph,i,j,(int)(width[i][j]+EPS));
	
	int ans=maximumFlow(graph,s,t);
	cout<<ans<<endl;
	
	return 0;
	
}