zoj 1670 Jewels from Heaven
题意:三个人,在给定正方形内,求第一个人拿到珠宝的概率。珠宝随机出现在正方形内。
思路:中垂线+半平面相交。
1 #include<cstdio> 2 #include<cstring> 3 #include<algorithm> 4 #include<iostream> 5 #include<cstdlib> 6 #include<string> 7 #include<cmath> 8 #include<vector> 9 using namespace std; 10 const int maxn=1e5+7; 11 const double eps=1e-8; 12 const double pi=acos(-1); 13 14 double dcmp(double x) 15 { 16 if(fabs(x) < eps) return 0; 17 else return x < 0 ? -1 : 1; 18 } 19 20 struct Point 21 { 22 double x, y; 23 Point(double x=0, double y=0):x(x),y(y) { } 24 }; 25 26 typedef Point Vector; 27 28 Vector operator + (const Point& A, const Point& B) 29 { 30 return Vector(A.x+B.x, A.y+B.y); 31 } 32 33 Vector operator - (const Point& A, const Point& B) 34 { 35 return Vector(A.x-B.x, A.y-B.y); 36 } 37 38 Vector operator * (const Point& A, double v) 39 { 40 return Vector(A.x*v, A.y*v); 41 } 42 43 Vector operator / (const Point& A, double v) 44 { 45 return Vector(A.x/v, A.y/v); 46 } 47 48 double Cross(const Vector& A, const Vector& B) 49 { 50 return A.x*B.y - A.y*B.x; 51 } 52 53 double Dot(const Vector& A, const Vector& B) 54 { 55 return A.x*B.x + A.y*B.y; 56 } 57 58 double Length(const Vector& A) 59 { 60 return sqrt(Dot(A,A)); 61 } 62 63 Vector Rotate(Vector A,double rad) 64 { 65 return Vector(A.x*cos(rad)-A.y*sin(rad),A.x*sin(rad)+A.y*cos(rad)); 66 } 67 68 bool operator < (const Point& p1, const Point& p2) 69 { 70 return p1.x < p2.x || (p1.x == p2.x && p1.y < p2.y); 71 } 72 73 bool operator == (const Point& p1, const Point& p2) 74 { 75 return p1.x == p2.x && p1.y == p2.y; 76 } 77 78 Vector Normal(Vector A) 79 { 80 double L=Length(A); 81 return Vector(-A.y/L,A.x/L); 82 } 83 struct Line 84 { 85 Point P; 86 Vector v; 87 double ang; 88 Line() {} 89 Line(Point P, Vector v):P(P),v(v) 90 { 91 ang = atan2(v.y, v.x); 92 } 93 bool operator < (const Line& L) const 94 { 95 return ang < L.ang; 96 } 97 }; 98 99 bool OnLeft(const Line& L, const Point& p) 100 { 101 return Cross(L.v, p-L.P) > 0; 102 } 103 104 105 Point GetLineIntersection(const Line& a, const Line& b) 106 { 107 Vector u = a.P-b.P; 108 double t = Cross(b.v, u) / Cross(a.v, b.v); 109 return a.P+a.v*t; 110 } 111 112 const double INF = 1e8; 113 114 Point ansPoly[maxn]; 115 int HalfplaneIntersection(vector<Line> L) //L为切割平面的直线集合,求半平面交,返回点的个数,点存在anspoly数组中 116 { 117 int n = L.size(); 118 sort(L.begin(), L.end()); // 按极角排序 119 int first, last; // 双端队列的第一个元素和最后一个元素的下标 120 vector<Point> p(n); // p[i]为q[i]和q[i+1]的交点 121 vector<Line> q(n); // 122 q[first=last=0] = L[0]; // 123 for(int i = 1; i < n; i++) 124 { 125 while(first < last && !OnLeft(L[i], p[last-1])) last--; 126 while(first < last && !OnLeft(L[i], p[first])) first++; 127 q[++last] = L[i]; 128 if(fabs(Cross(q[last].v, q[last-1].v)) < eps) // 129 { 130 last--; 131 if(OnLeft(q[last], L[i].P)) q[last] = L[i]; 132 } 133 if(first < last) p[last-1] = GetLineIntersection(q[last-1], q[last]); 134 } 135 while(first < last && !OnLeft(q[first], p[last-1])) last--; // 136 if(last - first <= 1) return 0; // 137 p[last] = GetLineIntersection(q[last], q[first]); // 138 // 从deque复制到输出中 139 int index=0; 140 for(int i = first; i <= last; i++) ansPoly[index++]=p[i]; 141 return index; 142 } 143 144 double PolygonArea(int n,Point *p) 145 { 146 double area=0; 147 for(int i=1; i<n-1; i++) 148 area+=Cross(p[i]-p[0],p[i+1]-p[0]); 149 return area/2; 150 } 151 152 Point p[5]; 153 int main() 154 { 155 // freopen("in.txt","r",stdin); 156 while(cin>>p[0].x>>p[0].y>>p[1].x>>p[1].y>>p[2].x>>p[2].y) 157 { 158 if(p[0].x==p[1].x&&p[0].y==p[1].y&&p[0].x==0)break; 159 // Point zd1,zd2; 160 vector<Line> vec; 161 vec.push_back(Line(Point(0,0),Point(1,0))); 162 vec.push_back(Line(Point(10000,0),Point(0,1))); 163 vec.push_back(Line(Point(10000,10000),Point(-1,0))); 164 vec.push_back(Line(Point(0,10000),Point(0,-1))); 165 Vector v=(p[1]-p[0]); 166 vec.push_back(Line((p[1]+p[0])*0.5,Normal(v))); 167 v=(p[2]-p[0]); 168 vec.push_back(Line((p[2]+p[0])*0.5,Normal(v))); 169 int m=HalfplaneIntersection(vec); 170 double ans=PolygonArea(m,ansPoly); 171 printf("%.3f\n",ans/(1.0*10000*10000)); 172 } 173 return 0; 174 }