POJ 2187 Beauty Contest（凸包+旋转卡壳）

Description

Bessie, Farmer John's prize cow, has just won first place in a bovine beauty contest, earning the title 'Miss Cow World'. As a result, Bessie will make a tour of N (2 <= N <= 50,000) farms around the world in order to spread goodwill between farmers and their cows. For simplicity, the world will be represented as a two-dimensional plane, where each farm is located at a pair of integer coordinates (x,y), each having a value in the range -10,000 ... 10,000. No two farms share the same pair of coordinates.

Even though Bessie travels directly in a straight line between pairs of farms, the distance between some farms can be quite large, so she wants to bring a suitcase full of hay with her so she has enough food to eat on each leg of her journey. Since Bessie refills her suitcase at every farm she visits, she wants to determine the maximum possible distance she might need to travel so she knows the size of suitcase she must bring.Help Bessie by computing the maximum distance among all pairs of farms.

Input

* Line 1: A single integer, N

* Lines 2..N+1: Two space-separated integers x and y specifying coordinate of each farm

Output

* Line 1: A single integer that is the squared distance between the pair of farms that are farthest apart from each other.

题目大意：给一堆点，求最远点距离的平方。

思路：先求个凸包，再用旋转卡壳求最远点对（这题数据出得不好，凸包之后暴力求最远点对）

代码（188MS）：

  1 #include <cstdio>
  2 #include <cstring>
  3 #include <iostream>
  4 #include <algorithm>
  5 #include <cmath>
  6 using namespace std;
  7 
  8 const int MAXN = 50010;
  9 const double EPS = 1e-10;
 10 const double PI = acos(-1.0);//3.14159265358979323846
 11 const double INF = 1;
 12 
 13 inline int sgn(double x) {
 14     return (x > EPS) - (x < -EPS);
 15 }
 16 
 17 struct Point {
 18     double x, y, ag;
 19     Point() {}
 20     Point(double x, double y): x(x), y(y) {}
 21     void read() {
 22         scanf("%lf%lf", &x, &y);
 23     }
 24     bool operator == (const Point &rhs) const {
 25         return sgn(x - rhs.x) == 0 && sgn(y - rhs.y) == 0;
 26     }
 27     bool operator < (const Point &rhs) const {
 28         if(y != rhs.y) return y < rhs.y;
 29         return x < rhs.x;
 30     }
 31     Point operator + (const Point &rhs) const {
 32         return Point(x + rhs.x, y + rhs.y);
 33     }
 34     Point operator - (const Point &rhs) const {
 35         return Point(x - rhs.x, y - rhs.y);
 36     }
 37     Point operator * (const double &b) const {
 38         return Point(x * b, y * b);
 39     }
 40     Point operator / (const double &b) const {
 41         return Point(x / b, y / b);
 42     }
 43     double length() {
 44         return sqrt(x * x + y * y);
 45     }
 46     Point unit() {
 47         return *this / length();
 48     }
 49     void makeAg() {
 50         ag = atan2(y, x);
 51     }
 52     void print() {
 53         printf("%.10f %.10f\n", x, y);
 54     }
 55 };
 56 typedef Point Vector;
 57 
 58 double dist(const Point &a, const Point &b) {
 59     return (a - b).length();
 60 }
 61 
 62 double cross(const Point &a, const Point &b) {
 63     return a.x * b.y - a.y * b.x;
 64 }
 65 //ret >= 0 means turn left
 66 double cross(const Point &sp, const Point &ed, const Point &op) {
 67     return sgn(cross(sp - op, ed - op));
 68 }
 69 
 70 double area(const Point& a, const Point &b, const Point &c) {
 71     return fabs(cross(a - c, b - c)) / 2;
 72 }
 73 //counter-clockwise
 74 Point rotate(const Point &p, double angle, const Point &o = Point(0, 0)) {
 75     Point t = p - o;
 76     double x = t.x * cos(angle) - t.y * sin(angle);
 77     double y = t.y * cos(angle) + t.x * sin(angle);
 78     return Point(x, y) + o;
 79 }
 80 
 81 struct Seg {
 82     Point st, ed;
 83     double ag;
 84     Seg() {}
 85     Seg(Point st, Point ed): st(st), ed(ed) {}
 86     void read() {
 87         st.read(); ed.read();
 88     }
 89     void makeAg() {
 90         ag = atan2(ed.y - st.y, ed.x - st.x);
 91     }
 92 };
 93 typedef Seg Line;
 94 
 95 //ax + by + c > 0
 96 Line buildLine(double a, double b, double c) {
 97     if(sgn(a) == 0 && sgn(b) == 0) return Line(Point(sgn(c) > 0 ? -1 : 1, INF), Point(0, INF));
 98     if(sgn(a) == 0) return Line(Point(sgn(b), -c/b), Point(0, -c/b));
 99     if(sgn(b) == 0) return Line(Point(-c/a, 0), Point(-c/a, sgn(a)));
100     if(b < 0) return Line(Point(0, -c/b), Point(1, -(a + c) / b));
101     else return Line(Point(1, -(a + c) / b), Point(0, -c/b));
102 }
103 
104 void moveRight(Line &v, double r) {
105     double dx = v.ed.x - v.st.x, dy = v.ed.y - v.st.y;
106     dx = dx / dist(v.st, v.ed) * r;
107     dy = dy / dist(v.st, v.ed) * r;
108     v.st.x += dy; v.ed.x += dy;
109     v.st.y -= dx; v.ed.y -= dx;
110 }
111 
112 bool isOnSeg(const Seg &s, const Point &p) {
113     return (p == s.st || p == s.ed) ||
114         (((p.x - s.st.x) * (p.x - s.ed.x) < 0 ||
115           (p.y - s.st.y) * (p.y - s.ed.y) < 0) &&
116          sgn(cross(s.ed, p, s.st) == 0));
117 }
118 
119 bool isIntersected(const Point &s1, const Point &e1, const Point &s2, const Point &e2) {
120     return (max(s1.x, e1.x) >= min(s2.x, e2.x)) &&
121         (max(s2.x, e2.x) >= min(s1.x, e1.x)) &&
122         (max(s1.y, e1.y) >= min(s2.y, e2.y)) &&
123         (max(s2.y, e2.y) >= min(s1.y, e1.y)) &&
124         (cross(s2, e1, s1) * cross(e1, e2, s1) >= 0) &&
125         (cross(s1, e2, s2) * cross(e2, e1, s2) >= 0);
126 }
127 
128 bool isIntersected(const Seg &a, const Seg &b) {
129     return isIntersected(a.st, a.ed, b.st, b.ed);
130 }
131 
132 bool isParallel(const Seg &a, const Seg &b) {
133     return sgn(cross(a.ed - a.st, b.ed - b.st)) == 0;
134 }
135 
136 //return Ax + By + C =0 's A, B, C
137 void Coefficient(const Line &L, double &A, double &B, double &C) {
138     A = L.ed.y - L.st.y;
139     B = L.st.x - L.ed.x;
140     C = L.ed.x * L.st.y - L.st.x * L.ed.y;
141 }
142 //point of intersection
143 Point operator * (const Line &a, const Line &b) {
144     double A1, B1, C1;
145     double A2, B2, C2;
146     Coefficient(a, A1, B1, C1);
147     Coefficient(b, A2, B2, C2);
148     Point I;
149     I.x = - (B2 * C1 - B1 * C2) / (A1 * B2 - A2 * B1);
150     I.y =   (A2 * C1 - A1 * C2) / (A1 * B2 - A2 * B1);
151     return I;
152 }
153 
154 bool isEqual(const Line &a, const Line &b) {
155     double A1, B1, C1;
156     double A2, B2, C2;
157     Coefficient(a, A1, B1, C1);
158     Coefficient(b, A2, B2, C2);
159     return sgn(A1 * B2 - A2 * B1) == 0 && sgn(A1 * C2 - A2 * C1) == 0 && sgn(B1 * C2 - B2 * C1) == 0;
160 }
161 
162 struct Poly {
163     int n;
164     Point p[MAXN];//p[n] = p[0]
165     void init(Point *pp, int nn) {
166         n = nn;
167         for(int i = 0; i < n; ++i) p[i] = pp[i];
168         p[n] = p[0];
169     }
170     double area() {
171         if(n < 3) return 0;
172         double s = p[0].y * (p[n - 1].x - p[1].x);
173         for(int i = 1; i < n; ++i)
174             s += p[i].y * (p[i - 1].x - p[i + 1].x);
175         return s / 2;
176     }
177 };
178 
179 void Graham_scan(Point *p, int n, int *stk, int &top) {//stk[0] = stk[top]
180     sort(p, p + n);
181     top = 1;
182     stk[0] = 0; stk[1] = 1;
183     for(int i = 2; i < n; ++i) {
184         while(top && cross(p[i], p[stk[top]], p[stk[top - 1]]) >= 0) --top;
185         stk[++top] = i;
186     }
187     int len = top;
188     stk[++top] = n - 2;
189     for(int i = n - 3; i >= 0; --i) {
190         while(top != len && cross(p[i], p[stk[top]], p[stk[top - 1]]) >= 0) --top;
191         stk[++top] = i;
192     }
193 }
194 //use for half_planes_cross
195 bool cmpAg(const Line &a, const Line &b) {
196     if(sgn(a.ag - b.ag) == 0)
197         return sgn(cross(b.ed, a.st, b.st)) < 0;
198     return a.ag < b.ag;
199 }
200 //clockwise, plane is on the right
201 bool half_planes_cross(Line *v, int vn, Poly &res, Line *deq) {
202     int i, n;
203     sort(v, v + vn, cmpAg);
204     for(i = n = 1; i < vn; ++i) {
205         if(sgn(v[i].ag - v[i-1].ag) == 0) continue;
206         v[n++] = v[i];
207     }
208     int head = 0, tail = 1;
209     deq[0] = v[0], deq[1] = v[1];
210     for(i = 2; i < n; ++i) {
211         if(isParallel(deq[tail - 1], deq[tail]) || isParallel(deq[head], deq[head + 1]))
212             return false;
213         while(head < tail && sgn(cross(v[i].ed, deq[tail - 1] * deq[tail], v[i].st)) > 0)
214             --tail;
215         while(head < tail && sgn(cross(v[i].ed, deq[head] * deq[head + 1], v[i].st)) > 0)
216             ++head;
217         deq[++tail] = v[i];
218     }
219     while(head < tail && sgn(cross(deq[head].ed, deq[tail - 1] * deq[tail], deq[head].st)) > 0)
220         --tail;
221     while(head < tail && sgn(cross(deq[tail].ed, deq[head] * deq[head + 1], deq[tail].st)) > 0)
222         ++head;
223     if(tail <= head + 1) return false;
224     res.n = 0;
225     for(i = head; i < tail; ++i)
226         res.p[res.n++] = deq[i] * deq[i + 1];
227     res.p[res.n++] = deq[head] * deq[tail];
228     res.n = unique(res.p, res.p + res.n) - res.p;
229     res.p[res.n] = res.p[0];
230     return true;
231 }
232 
233 //ix and jx is the points whose distance is return, res.p[n - 1] = res.p[0]
234 double dia_roataing_calipers(Poly &res, int &ix, int &jx) {
235     double dia = 0;
236     int q = 1;
237     for(int i = 0; i < res.n - 1; ++i) {
238         while(sgn(area(res.p[i], res.p[q + 1], res.p[i + 1]) - area(res.p[i], res.p[q], res.p[i + 1])) > 0)
239             q = (q + 1) % (res.n - 1);
240         if(sgn(dist(res.p[i], res.p[q]) - dia) > 0) {
241             dia = dist(res.p[i], res.p[q]);
242             ix = i; jx = q;
243         }
244         if(sgn(dist(res.p[i + 1], res.p[q]) - dia) > 0) {
245             dia = dist(res.p[i + 1], res.p[q]);
246             ix = i + 1; jx = q;
247         }
248     }
249     return dia;
250 }
251 
252 /*******************************************************************************************/
253 
254 Point p[MAXN];
255 int stk[MAXN], top, n;
256 Poly res;
257 
258 int dist2(const Point &a, const Point &b) {
259     int x = a.x - b.x, y = a.y - b.y;
260     return x * x + y * y;
261 }
262 
263 int main() {
264     scanf("%d", &n);
265     for(int i = 0; i < n; ++i) p[i].read();
266     Graham_scan(p, n, stk, top);
267     for(int i = 0; i <= top; ++i) res.p[res.n++] = p[stk[i]];
268     int ans1, ans2;
269     dia_roataing_calipers(res, ans1, ans2);
270     printf("%d\n", dist2(res.p[ans1], res.p[ans2]));
271 }

View Code

代码（204MS）（修改了一下，不过是在VJ上交的）：

  1 #include <cstdio>
  2 #include <cstring>
  3 #include <iostream>
  4 #include <algorithm>
  5 #include <cmath>
  6 using namespace std;
  7 
  8 const int MAXN = 50010;
  9 const double EPS = 1e-10;
 10 const double PI = acos(-1.0);//3.14159265358979323846
 11 const double INF = 1;
 12 
 13 inline int sgn(double x) {
 14     return (x > EPS) - (x < -EPS);
 15 }
 16 
 17 struct Point {
 18     double x, y, ag;
 19     Point() {}
 20     Point(double x, double y): x(x), y(y) {}
 21     void read() {
 22         scanf("%lf%lf", &x, &y);
 23     }
 24     bool operator == (const Point &rhs) const {
 25         return sgn(x - rhs.x) == 0 && sgn(y - rhs.y) == 0;
 26     }
 27     bool operator < (const Point &rhs) const {
 28         if(y != rhs.y) return y < rhs.y;
 29         return x < rhs.x;
 30     }
 31     Point operator + (const Point &rhs) const {
 32         return Point(x + rhs.x, y + rhs.y);
 33     }
 34     Point operator - (const Point &rhs) const {
 35         return Point(x - rhs.x, y - rhs.y);
 36     }
 37     Point operator * (const double &b) const {
 38         return Point(x * b, y * b);
 39     }
 40     Point operator / (const double &b) const {
 41         return Point(x / b, y / b);
 42     }
 43     double length() {
 44         return sqrt(x * x + y * y);
 45     }
 46     Point unit() {
 47         return *this / length();
 48     }
 49     void makeAg() {
 50         ag = atan2(y, x);
 51     }
 52     void print() {
 53         printf("%.10f %.10f\n", x, y);
 54     }
 55 };
 56 typedef Point Vector;
 57 
 58 double dist(const Point &a, const Point &b) {
 59     return (a - b).length();
 60 }
 61 
 62 double cross(const Point &a, const Point &b) {
 63     return a.x * b.y - a.y * b.x;
 64 }
 65 //ret >= 0 means turn right
 66 double cross(const Point &sp, const Point &ed, const Point &op) {
 67     return cross(sp - op, ed - op);
 68 }
 69 
 70 double area(const Point& a, const Point &b, const Point &c) {
 71     return fabs(cross(a - c, b - c)) / 2;
 72 }
 73 //counter-clockwise
 74 Point rotate(const Point &p, double angle, const Point &o = Point(0, 0)) {
 75     Point t = p - o;
 76     double x = t.x * cos(angle) - t.y * sin(angle);
 77     double y = t.y * cos(angle) + t.x * sin(angle);
 78     return Point(x, y) + o;
 79 }
 80 
 81 struct Seg {
 82     Point st, ed;
 83     double ag;
 84     Seg() {}
 85     Seg(Point st, Point ed): st(st), ed(ed) {}
 86     void read() {
 87         st.read(); ed.read();
 88     }
 89     void makeAg() {
 90         ag = atan2(ed.y - st.y, ed.x - st.x);
 91     }
 92 };
 93 typedef Seg Line;
 94 
 95 //ax + by + c > 0
 96 Line buildLine(double a, double b, double c) {
 97     if(sgn(a) == 0 && sgn(b) == 0) return Line(Point(sgn(c) > 0 ? -1 : 1, INF), Point(0, INF));
 98     if(sgn(a) == 0) return Line(Point(sgn(b), -c/b), Point(0, -c/b));
 99     if(sgn(b) == 0) return Line(Point(-c/a, 0), Point(-c/a, sgn(a)));
100     if(b < 0) return Line(Point(0, -c/b), Point(1, -(a + c) / b));
101     else return Line(Point(1, -(a + c) / b), Point(0, -c/b));
102 }
103 
104 void moveRight(Line &v, double r) {
105     double dx = v.ed.x - v.st.x, dy = v.ed.y - v.st.y;
106     dx = dx / dist(v.st, v.ed) * r;
107     dy = dy / dist(v.st, v.ed) * r;
108     v.st.x += dy; v.ed.x += dy;
109     v.st.y -= dx; v.ed.y -= dx;
110 }
111 
112 bool isOnSeg(const Seg &s, const Point &p) {
113     return (p == s.st || p == s.ed) ||
114         (((p.x - s.st.x) * (p.x - s.ed.x) < 0 ||
115           (p.y - s.st.y) * (p.y - s.ed.y) < 0) &&
116          sgn(cross(s.ed, p, s.st) == 0));
117 }
118 
119 bool isIntersected(const Point &s1, const Point &e1, const Point &s2, const Point &e2) {
120     return (max(s1.x, e1.x) >= min(s2.x, e2.x)) &&
121         (max(s2.x, e2.x) >= min(s1.x, e1.x)) &&
122         (max(s1.y, e1.y) >= min(s2.y, e2.y)) &&
123         (max(s2.y, e2.y) >= min(s1.y, e1.y)) &&
124         (cross(s2, e1, s1) * cross(e1, e2, s1) >= 0) &&
125         (cross(s1, e2, s2) * cross(e2, e1, s2) >= 0);
126 }
127 
128 bool isIntersected(const Seg &a, const Seg &b) {
129     return isIntersected(a.st, a.ed, b.st, b.ed);
130 }
131 
132 bool isParallel(const Seg &a, const Seg &b) {
133     return sgn(cross(a.ed - a.st, b.ed - b.st)) == 0;
134 }
135 
136 //return Ax + By + C =0 's A, B, C
137 void Coefficient(const Line &L, double &A, double &B, double &C) {
138     A = L.ed.y - L.st.y;
139     B = L.st.x - L.ed.x;
140     C = L.ed.x * L.st.y - L.st.x * L.ed.y;
141 }
142 //point of intersection
143 Point operator * (const Line &a, const Line &b) {
144     double A1, B1, C1;
145     double A2, B2, C2;
146     Coefficient(a, A1, B1, C1);
147     Coefficient(b, A2, B2, C2);
148     Point I;
149     I.x = - (B2 * C1 - B1 * C2) / (A1 * B2 - A2 * B1);
150     I.y =   (A2 * C1 - A1 * C2) / (A1 * B2 - A2 * B1);
151     return I;
152 }
153 
154 bool isEqual(const Line &a, const Line &b) {
155     double A1, B1, C1;
156     double A2, B2, C2;
157     Coefficient(a, A1, B1, C1);
158     Coefficient(b, A2, B2, C2);
159     return sgn(A1 * B2 - A2 * B1) == 0 && sgn(A1 * C2 - A2 * C1) == 0 && sgn(B1 * C2 - B2 * C1) == 0;
160 }
161 
162 struct Poly {
163     int n;
164     Point p[MAXN];//p[n] = p[0]
165     void init(Point *pp, int nn) {
166         n = nn;
167         for(int i = 0; i < n; ++i) p[i] = pp[i];
168         p[n] = p[0];
169     }
170     double area() {
171         if(n < 3) return 0;
172         double s = p[0].y * (p[n - 1].x - p[1].x);
173         for(int i = 1; i < n; ++i)
174             s += p[i].y * (p[i - 1].x - p[i + 1].x);
175         return s / 2;
176     }
177 };
178 //the convex hull is clockwise
179 void Graham_scan(Point *p, int n, int *stk, int &top) {//stk[0] = stk[top]
180     sort(p, p + n);
181     top = 1;
182     stk[0] = 0; stk[1] = 1;
183     for(int i = 2; i < n; ++i) {
184         while(top && cross(p[i], p[stk[top]], p[stk[top - 1]]) <= 0) --top;
185         stk[++top] = i;
186     }
187     int len = top;
188     stk[++top] = n - 2;
189     for(int i = n - 3; i >= 0; --i) {
190         while(top != len && cross(p[i], p[stk[top]], p[stk[top - 1]]) <= 0) --top;
191         stk[++top] = i;
192     }
193 }
194 //use for half_planes_cross
195 bool cmpAg(const Line &a, const Line &b) {
196     if(sgn(a.ag - b.ag) == 0)
197         return sgn(cross(b.ed, a.st, b.st)) < 0;
198     return a.ag < b.ag;
199 }
200 //clockwise, plane is on the right
201 bool half_planes_cross(Line *v, int vn, Poly &res, Line *deq) {
202     int i, n;
203     sort(v, v + vn, cmpAg);
204     for(i = n = 1; i < vn; ++i) {
205         if(sgn(v[i].ag - v[i-1].ag) == 0) continue;
206         v[n++] = v[i];
207     }
208     int head = 0, tail = 1;
209     deq[0] = v[0], deq[1] = v[1];
210     for(i = 2; i < n; ++i) {
211         if(isParallel(deq[tail - 1], deq[tail]) || isParallel(deq[head], deq[head + 1]))
212             return false;
213         while(head < tail && sgn(cross(v[i].ed, deq[tail - 1] * deq[tail], v[i].st)) > 0)
214             --tail;
215         while(head < tail && sgn(cross(v[i].ed, deq[head] * deq[head + 1], v[i].st)) > 0)
216             ++head;
217         deq[++tail] = v[i];
218     }
219     while(head < tail && sgn(cross(deq[head].ed, deq[tail - 1] * deq[tail], deq[head].st)) > 0)
220         --tail;
221     while(head < tail && sgn(cross(deq[tail].ed, deq[head] * deq[head + 1], deq[tail].st)) > 0)
222         ++head;
223     if(tail <= head + 1) return false;
224     res.n = 0;
225     for(i = head; i < tail; ++i)
226         res.p[res.n++] = deq[i] * deq[i + 1];
227     res.p[res.n++] = deq[head] * deq[tail];
228     res.n = unique(res.p, res.p + res.n) - res.p;
229     res.p[res.n] = res.p[0];
230     return true;
231 }
232 
233 //ix and jx is the points whose distance is return, res.p[n - 1] = res.p[0], res must be clockwise
234 double dia_roataing_calipers(Poly &res, int &ix, int &jx) {
235     double dia = 0;
236     int q = 1;
237     for(int i = 0; i < res.n - 1; ++i) {
238         while(sgn(cross(res.p[i], res.p[q + 1], res.p[i + 1]) - cross(res.p[i], res.p[q], res.p[i + 1])) > 0)
239             q = (q + 1) % (res.n - 1);
240         if(sgn(dist(res.p[i], res.p[q]) - dia) > 0) {
241             dia = dist(res.p[i], res.p[q]);
242             ix = i; jx = q;
243         }
244         if(sgn(dist(res.p[i + 1], res.p[q]) - dia) > 0) {
245             dia = dist(res.p[i + 1], res.p[q]);
246             ix = i + 1; jx = q;
247         }
248     }
249     return dia;
250 }
251 
252 /*******************************************************************************************/
253 
254 Point p[MAXN];
255 int stk[MAXN], top, n;
256 Poly res;
257 
258 int dist2(const Point &a, const Point &b) {
259     int x = a.x - b.x, y = a.y - b.y;
260     return x * x + y * y;
261 }
262 
263 int main() {
264     scanf("%d", &n);
265     for(int i = 0; i < n; ++i) p[i].read();
266     Graham_scan(p, n, stk, top);
267     for(int i = 0; i <= top; ++i) res.p[res.n++] = p[stk[i]];
268     int ans1, ans2;
269     dia_roataing_calipers(res, ans1, ans2);
270     printf("%d\n", dist2(res.p[ans1], res.p[ans2]));
271 }

View Code

posted @ 2013-11-11 19:39 Oyking 阅读(395) 评论(0) 编辑收藏举报

刷新页面返回顶部

Oyking

POJ 2187 Beauty Contest（凸包+旋转卡壳）

公告