Sichuan State Programming Contest 2011 —— D.Division

http://acm.uestc.edu.cn/problem.php?pid=1560

會升讓我看的題。給你一個凸包，和凸包外一點，讓找到一條射線把這個凸包分成面積相等的兩塊，求這個射線的機關向量。

話說，以前見過類似的，求切割成兩個相等面積的。會升說枚舉弧度，我的第一反應也是這個，寫了好長好長，我用二分找的。結果一直WA。後來會升說可以枚舉線段長，就是從射線的端點的凸包張角最大的那兩個點，然後枚舉線段上的點，二分即可，因為枚舉弧度的話，會用到三角函數，而這個是很很損失精度的行為 = =。然後晚上吃飯回來就寫差不多了，但是一直WA。想到一個BUG,因為我找兩個端點是用向量找的，是以，如果兩端的兩個點正好一個在第二象限，一個在第三象限，那麼用向量排（atan值）是不對的！改了下，用極角排序，就是用叉積排，然後就過了！！！！激動死了都！！這題比賽貌似就6個人過也。。。

至于找線段上的點使得射線切割成兩個相等面積，直接找到點後利用半平面交，二分查找這個點。如果這條射線左邊的凸包面積比凸包面積一半大，那麼肯定要右移，就這麼二分查找。

附自己做的資料兩組

3
7
0 -2
1 -1
2 0
2 2
0 5
-1 5
-1 0
3 0
7
0 -2
1 -1
2 0
2 2
0 5
-1 5
-1 0
-3 0
ans：
-0.8575 0.5145
0.0513 0.9987
           

#include <cmath>
#include <cstdio>
#include <cstdlib>
#include <iostream>
#include <cstring>
#include <string>
#include <algorithm>
#include <queue>
#include <stack>
#include <climits>
#include <iomanip>

using namespace std;

const int MAX = 105;
struct point {double x,y;int ind;};
struct line{point a,b; double ang;};
point p[MAX],s[MAX],c;
line ln[MAX],deq[MAX],ltmp[MAX];;
const double eps = 1e-6;
bool dy(double x,double y)	{	return x > y + eps;}	// x > y 
bool xy(double x,double y)	{	return x < y - eps;}	// x < y 
bool dyd(double x,double y)	{ 	return x > y - eps;}	// x >= y 
bool xyd(double x,double y)	{	return x < y + eps;} 	// x <= y 
bool dd(double x,double y) 	{	return fabs( x - y ) < eps;}  // x == y
double crossProduct(point a,point b,point c)//向量 ac 在 ab 的方向 
{
	return (c.x - a.x)*(b.y - a.y) - (b.x - a.x)*(c.y - a.y);
}
void makeline_hp(point a,point b,line &l) // 線段(向量ab)左側區域有效 
{
	l.a = a; l.b = b;
	l.ang = atan2(a.y - b.y,a.x - b.x);	// 如果是右側區域，改成b.y - a.y,b.x - a.x 
}
double area_polygon(point p[],int n)
{
	if( n < 3 ) return 0.0; 
	double s = 0.0;
	for(int i=0; i<n; i++)
		s += p[(i+1)%n].y * p[i].x - p[(i+1)%n].x * p[i].y;
	return fabs(s)/2.0;
}
bool parallel(line u,line v)
{
	return dd( (u.a.x - u.b.x)*(v.a.y - v.b.y) - (v.a.x - v.b.x)*(u.a.y - u.b.y) , 0.0 );
}
point l2l_inst_p(line l1,line l2)
{
	point ans = l1.a;
	double t = ((l1.a.x - l2.a.x)*(l2.a.y - l2.b.y) - (l1.a.y - l2.a.y)*(l2.a.x - l2.b.x))/
			   ((l1.a.x - l1.b.x)*(l2.a.y - l2.b.y) - (l1.a.y - l1.b.y)*(l2.a.x - l2.b.x));
	ans.x += (l1.b.x - l1.a.x)*t;
	ans.y += (l1.b.y - l1.a.y)*t;
	return ans;
}
bool equal_ang(line a,line b)	// 第一次unique的比較函數 
{
	return dd(a.ang,b.ang);
}
bool cmphp(line a,line b)	// 排序的比較函數 
{
	if( dd(a.ang,b.ang) ) return xy(crossProduct(b.a,b.b,a.a),0.0);
	return xy(a.ang,b.ang);
}
bool equal_p(point a,point b)//第二次unique的比較函數 
{
	return dd(a.x,b.x) && dd(a.y,b.y);
}
void makeline_hp(double x1,double y1,double x2,double y2,line &l)
{
	l.a.x = x1; l.a.y = y1; l.b.x = x2; l.b.y = y2;
	l.ang = atan2(y2 - y1,x2 - x1);
}
void inst_hp_nlogn(line *ln,int n,point *s,int &len)
{
	len = 0;
	sort(ln,ln+n,cmphp);
	n = unique(ln,ln+n,equal_ang) - ln;
	int bot = 0,top = 1;
	deq[0] = ln[0]; deq[1] = ln[1];
	for(int i=2; i<n; i++)
	{
		if( parallel(deq[top],deq[top-1]) || parallel(deq[bot],deq[bot+1]) )
			return ;
		while( bot < top && dy(crossProduct(ln[i].a,ln[i].b,
			l2l_inst_p(deq[top],deq[top-1])),0.0) )
			top--;
		while( bot < top && dy(crossProduct(ln[i].a,ln[i].b,
			l2l_inst_p(deq[bot],deq[bot+1])),0.0) )
			bot++;
		deq[++top] = ln[i];
	}
	while( bot < top && dy(crossProduct(deq[bot].a,deq[bot].b,
		l2l_inst_p(deq[top],deq[top-1])),0.0) )	top--;
	while( bot < top && dy(crossProduct(deq[top].a,deq[top].b,
		l2l_inst_p(deq[bot],deq[bot+1])),0.0) )	bot++;
	if( top <= bot + 1 ) return ;
	
	for(int i=bot; i<top; i++)
		s[len++] = l2l_inst_p(deq[i],deq[i+1]);
	if( bot < top + 1 ) s[len++] = l2l_inst_p(deq[bot],deq[top]);
	len = unique(s,s+len,equal_p) - s;
} 
double disp2p(point a,point b) //  a b 兩點之間的距離 
{
	return sqrt( ( a.x - b.x ) * ( a.x - b.x ) + ( a.y - b.y ) * ( a.y - b.y ) );
}
bool cmp1(point a,point b)  // 排序   
{  
    double len = crossProduct(c,a,b);  
    if( dd(len,0.0) )  
        return xy(disp2p(c,a),disp2p(c,b));  
    return xy(len,0.0);  
}  
int main()
{
	int ncases,n,len;
	point d,b,e,mid;;
	scanf("%d",&ncases);
	int ind = 1;
	while( ncases-- )
	{
		scanf("%d",&n);
		for(int i=0; i<n; i++)
		{
			scanf("%lf%lf",&p[i].x,&p[i].y);
			p[i].ind = i;
		}
		p[n] = p[0];
		for(int i=0; i<n; i++)
			makeline_hp(p[i],p[i+1],ln[i]);
		scanf("%lf%lf",&c.x,&c.y);
		memcpy(s,p,sizeof(p));
		sort(s,s+n,cmp1);
		double asum = area_polygon(p,n)/2;
		memcpy(ltmp,ln,sizeof(ln));
		b = p[s[0].ind]; e = p[s[n-1].ind];
		while( !(dd(b.x,e.x) && dd(b.y,e.y)) )
		{
			mid.x = (b.x + e.x)/2;
			mid.y = (b.y + e.y)/2;
			memcpy(ln,ltmp,sizeof(ltmp));
			makeline_hp(c,mid,ln[n]);
			inst_hp_nlogn(ln,n+1,s,len);
			double al = area_polygon(s,len);
			if( dd(al,asum) )
				break;
			if( dy(al,asum) )
				b = mid;
			else
				e = mid;	
		}
		double x = atan2(mid.y - c.y,mid.x - c.x);
		point ans; 
		ans.x = cos(x);	ans.y = sin(x);
		printf("Case #%d: ",ind++);
		printf("%.4lf %.4lf\n",ans.x,ans.y);
	}
return 0;
}