poj 3743 LL’s cake (PSLG,Accepted)
搞了好久都过不了,看了下题解是用PSLG来做的。POJ 2164 && LA 3218 Find the Border (Geometry, PSLG 平面直线图) - LyonLys - 博客园 这篇里面写过一下,就是把点都提取出来,然后模拟沿着边界移动,找到多边形并计算面积。
而我的做法是直接模拟多边形切割,各种超时爆内存。先留着,看以后能不能用这个来过。
没过的代码:
1 #include <cstdio> 2 #include <cstring> 3 #include <iostream> 4 #include <algorithm> 5 #include <cmath> 6 #include <vector> 7 #include <queue> 8 9 using namespace std; 10 11 const double EPS = 1e-4; 12 inline int sgn(double x) { return (x > EPS) - (x < -EPS);} 13 struct Point { 14 double x, y; 15 Point() {} 16 Point(double x, double y) : x(x), y(y) {} 17 Point operator + (Point a) { return Point(x + a.x, y + a.y);} 18 Point operator - (Point a) { return Point(x - a.x, y - a.y);} 19 Point operator * (double p) { return Point(x * p, y * p);} 20 Point operator / (double p) { return Point(x / p, y / p);} 21 bool operator < (Point a) const { return sgn(x - a.x) < 0 || sgn(x - a.x) == 0 && y < a.y;} 22 bool operator == (Point a) const { return sgn(x - a.x) == 0 && sgn(y - a.y) == 0;} 23 } ; 24 25 inline double cross(Point a, Point b) { return a.x * b.y - a.y * b.x;} 26 inline double dot(Point a, Point b) { return a.x * b.x + a.y * b.y;} 27 inline double veclen(Point x) { return sqrt(dot(x, x));} 28 inline Point normal(Point x) { return Point(-x.y, x.x) / veclen(x);} 29 inline Point vecunit(Point x) { return x / veclen(x);} 30 31 struct Line { 32 Point s, t; 33 Line() {} 34 Line(Point s, Point t) : s(s), t(t) {} 35 Point vec() { return t - s;} 36 Point point(double x) { return s + vec() * x;} 37 } ; 38 inline Point llint(Line a, Line b) { return a.point(cross(b.vec(), a.s - b.s) / cross(a.vec(), b.vec()));} 39 inline bool onseg(Point x, Point s, Point t) { return sgn(cross(s - x, t - x)) == 0 && sgn(dot(s - x, t - x)) < 0;} 40 inline bool onseg(Point x, Line a) { return onseg(x, a.s, a.t);} 41 42 struct Circle { 43 Point c; 44 double r; 45 Circle() {} 46 Circle(Point c, double r) : c(c), r(r) {} 47 bool in(Point x) { return sgn(veclen(x - c) - r) <= 0;} 48 Point point(double x) { return Point(c.x + cos(x) * r, c.y + sin(x) * r);} 49 } ; 50 const double R = 10.0; 51 Circle cake = Circle(Point(0.0, 0.0), R); 52 const double PI = acos(-1.0); 53 template<class T> T sqr(T x) { return x * x;} 54 inline double angle(Point x) { return atan2(x.y, x.x);} 55 56 int clint(Line s, Point *sol) { 57 Point nor = normal(s.vec()), ip = llint(s, Line(cake.c, cake.c + nor)); 58 double dis = veclen(cake.c - ip); 59 if (sgn(dis - cake.r) >= 0) return 0; 60 Point dxy = vecunit(s.vec()) * sqrt(sqr(cake.r) - sqr(dis)); 61 int ret = 0; 62 sol[ret] = ip + dxy; 63 if (onseg(sol[ret], s)) ret++; 64 sol[ret] = ip - dxy; 65 if (onseg(sol[ret], s)) ret++; 66 return ret; 67 } 68 69 double getsec(Point a, Point b) { 70 double a1 = angle(a - cake.c); 71 double a2 = angle(b - cake.c); 72 double da = fabs(a1 - a2); 73 if (da > PI) da = PI * 2.0 - da; 74 return sqr(cake.r) * da * sgn(cross(a - cake.c, b - cake.c)) / 2.0; 75 } 76 77 inline double gettri(Point a, Point b) { return cross(a - cake.c, b - cake.c) / 2.0;} 78 //typedef vector<Point> VP; 79 const int N = 111; 80 struct VP { 81 Point vex[N]; 82 int n; 83 void clear() { n = 0;} 84 void push_back(Point x) { vex[n++] = x;} 85 void pop_back() { n--;} 86 int size() { return n;} 87 } ; 88 89 double cpint(VP pt) { 90 double ret = 0.0; 91 int n = pt.size(); 92 Point tmp[2]; 93 pt.vex[n] = pt.vex[0]; 94 for (int i = 0; i < n; i++) { 95 int ic = clint(Line(pt.vex[i], pt.vex[i + 1]), tmp); 96 if (ic == 0) { 97 if (!cake.in(pt.vex[i]) || !cake.in(pt.vex[i + 1])) ret += getsec(pt.vex[i], pt.vex[i + 1]); 98 else ret += gettri(pt.vex[i], pt.vex[i + 1]); 99 } else if (ic == 1) { 100 if (cake.in(pt.vex[i])) ret += gettri(pt.vex[i], tmp[0]), ret += getsec(tmp[0], pt.vex[i + 1]); 101 else ret += getsec(pt.vex[i], tmp[0]), ret += gettri(tmp[0], pt.vex[i + 1]); 102 } else { 103 if (pt.vex[i] < pt.vex[i + 1] ^ tmp[0] < tmp[1]) swap(tmp[0], tmp[1]); 104 ret += getsec(pt.vex[i], tmp[0]); 105 ret += gettri(tmp[0], tmp[1]); 106 ret += getsec(tmp[1], pt.vex[i + 1]); 107 } 108 // cout << "~~ic " << ic << ' ' << ret << endl; 109 } 110 return fabs(ret); 111 } 112 113 bool fixpoly(VP &poly) { 114 double sum = 0.0; 115 int n = poly.size(); 116 poly.vex[n] = poly.vex[0]; 117 for (int i = 0; i < n; i++) sum += cross(poly.vex[i], poly.vex[i + 1]); 118 if (sgn(sum) == 0) return false; 119 if (sgn(sum) < 0) reverse(poly.vex, poly.vex + n); 120 return true; 121 } 122 123 void cutpoly(VP &poly, Line l, VP &ret) { 124 ret.clear(); 125 int n = poly.size(); 126 // cout << n << endl; 127 poly.vex[n] = poly.vex[0]; 128 for (int i = 0; i < n; i++) { 129 if (sgn(cross(l.vec(), poly.vex[i] - l.s)) >= 0) ret.push_back(poly.vex[i]); 130 if (sgn(cross(l.vec(), poly.vex[i] - poly.vex[i + 1]))) { 131 Point ip = llint(l, Line(poly.vex[i], poly.vex[i + 1])); 132 // cout << "ip " << ip.x << ' ' << ip.y << endl; 133 if (onseg(ip, poly.vex[i], poly.vex[i + 1]) || poly.vex[i] == ip) ret.push_back(ip); 134 } 135 } 136 // cout << "cp sz " << ret.size() << endl; 137 } 138 139 const int M = 22222; 140 int q[1111111], qh, qt, nu; 141 VP rec[M]; 142 queue<int> recycle; 143 144 int getID() { 145 int ret; 146 if (nu >= M) { 147 if (recycle.empty()) { puts("shit!"); while (1) ;} 148 ret = recycle.front(); 149 recycle.pop(); 150 } else ret = nu++; 151 return ret; 152 } 153 154 void retID(int x) { recycle.push(x);} 155 156 int main() { 157 // freopen("in", "r", stdin); 158 // freopen("out", "w", stdout); 159 int T, n, tmp; 160 double x, y; 161 cin >> T; 162 while (T-- && cin >> n) { 163 while (!recycle.empty()) recycle.pop(); 164 qh = qt = nu = 0; 165 tmp = getID(); 166 rec[tmp].clear(); 167 rec[tmp].push_back(Point(-R * 2.0, -R * 2.0)); 168 rec[tmp].push_back(Point(R * 2.0, -R * 2.0)); 169 rec[tmp].push_back(Point(R * 2.0, R * 2.0)); 170 rec[tmp].push_back(Point(-R * 2.0, R * 2.0)); 171 fixpoly(rec[tmp]); 172 q[qt++] = tmp; 173 for (int i = 0; i < n; i++) { 174 cin >> x >> y; 175 int sz = qt - qh; 176 Line t = Line(cake.point(x), cake.point(y)); 177 // cout << cake.point(x).x << '=' << cake.point(x).y << endl; 178 // cout << cake.point(y).x << '~' << cake.point(y).y << endl; 179 for (int j = 0; j < sz; j++) { 180 tmp = getID(); 181 // cout << "qh ?? " << qh << ' ' << q[qh] << ' ' << rec[q[qh]].size() << endl; 182 cutpoly(rec[q[qh]], t, rec[tmp]); 183 if (fixpoly(rec[tmp])) { 184 // cout << j << "~~1 " << rec[tmp].size() << endl; 185 // for (int k = 0; k < rec[tmp].size(); k++) cout << rec[tmp].vex[k].x << ' ' << rec[tmp].vex[k].y << endl; 186 q[qt++] = tmp; 187 } 188 swap(t.s, t.t); 189 tmp = getID(); 190 cutpoly(rec[q[qh]], t, rec[tmp]); 191 if (fixpoly(rec[tmp])) { 192 // cout << j << "~~2 " << rec[tmp].size() << endl; 193 // for (int k = 0; k < rec[tmp].size(); k++) cout << rec[tmp].vex[k].x << ' ' << rec[tmp].vex[k].y << endl; 194 q[qt++] = tmp; 195 } 196 retID(q[qh++]); 197 } 198 // cout << "sz~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ " << qt - qh << endl; 199 } 200 double mx = 0.0; 201 while (qh < qt) { 202 mx = max(mx, cpint(rec[q[qh++]])); 203 // cout << ".. " << mx << endl; 204 } 205 printf("%.2f\n", mx); 206 } 207 return 0; 208 } 209 210 /* 211 6 212 2 213 -3.140000 0.000000 214 -1.000000 1.000000 215 2 216 -3.141592 0.000000 217 -1.570796 1.570796 218 3 219 -3.000000 3.000000 220 -2.000000 2.000000 221 -1.000000 1.000000 222 4 223 -3.140000 0.000000 224 -1.000000 1.000000 225 -3.140000 -1.000000 226 1.000000 0.000000 227 6 228 -3.140000 0.000000 229 -1.000000 1.000000 230 -3.140000 -1.000000 231 1.000000 0.000000 232 -3.140000 -1.000000 233 1.000000 0.000000 234 6 235 -3.141592 0.000000 236 -1.570796 1.570796 237 -3.141592 -1.570796 238 0.000000 1.570796 239 -3.141592 1.570796 240 0.000000 -1.570796 241 */
PSLG的方法将尽快更新上来!
UPD:
模拟遍历边界,985ms压线过,因为有圆弧,所以有几个特判。比较好奇别人那些稳稳的不超时除了少用了STL还做了些什么?
代码如下:
1 #include <cstdio> 2 #include <iostream> 3 #include <algorithm> 4 #include <cstring> 5 #include <cmath> 6 #include <vector> 7 #include <map> 8 #include <set> 9 10 using namespace std; 11 12 const double EPS = 1e-8; 13 inline int sgn(double x) { return (x > EPS) - (x < -EPS);} 14 15 struct Point { 16 double x, y; 17 int id; 18 Point() {} 19 Point(double x, double y) : x(x), y(y) {} 20 bool operator < (Point a) const { return sgn(x - a.x) < 0 || sgn(x - a.x) == 0 && y < a.y;} 21 bool operator == (Point a) const { return sgn(x - a.x) == 0 && sgn(y - a.y) == 0;} 22 Point operator + (Point a) { return Point(x + a.x, y + a.y);} 23 Point operator - (Point a) { return Point(x - a.x, y - a.y);} 24 Point operator * (double p) { return Point(x * p, y * p);} 25 Point operator / (double p) { return Point(x / p, y / p);} 26 } ; 27 28 inline double cross(Point a, Point b) { return a.x * b.y - a.y * b.x;} 29 inline double dot(Point a, Point b) { return a.x * b.x + a.y * b.y;} 30 inline double veclen(Point x) { return sqrt(dot(x, x));} 31 inline Point vecunit(Point x) { return x / veclen(x);} 32 inline Point normal(Point x) { return Point(-x.y, x.x) / veclen(x);} 33 34 const int N = 111; 35 Point pts[N * N]; 36 int ptcnt; 37 38 struct Line { 39 Point s, t; 40 Line() {} 41 Line(Point s, Point t) : s(s), t(t) {} 42 Point vec() { return t - s;} 43 Point point(double p) { return s + vec() * p;} 44 } ; 45 46 inline Point llint(Line a, Line b) { return a.point(cross(b.vec(), a.s - b.s) / cross(a.vec(), b.vec()));} 47 inline bool onseg(Point x, Point a, Point b) { return sgn(cross(a - x, b - x)) == 0 && sgn(dot(a - x, b - x)) <= 0;} 48 inline bool onseg(Point x, Line l) { return onseg(x, l.s, l.t);} 49 50 const double R = 10.0; 51 inline bool oncircle(Point x) { return sgn(veclen(x) - R) == 0;} 52 inline Point getpt(double p) { return Point(cos(p) * R, sin(p) * R);} 53 inline double angle(Point x) { return atan2(x.y, x.x);} 54 55 struct Node { 56 double ang; 57 int id; 58 bool arc; 59 Node() {} 60 Node(double ang, int id) : ang(ang), id(id) { arc = false;} 61 bool operator < (Node x) const { return sgn(ang - x.ang) < 0 || sgn(ang - x.ang) == 0 && arc > x.arc;} 62 } ; 63 Line cut[N]; 64 65 const double PI = acos(-1.0); 66 template<class T> T sqr(T x) { return x * x;} 67 Point ori, tmp[N << 1]; 68 vector<Node> nb[N * N], oc; 69 typedef pair<int, int> PII; 70 typedef pair<int, bool> PIB; 71 set<PII> used; 72 map<int, PIB> nx[N * N], anx[N * N]; 73 74 75 inline double caltri(Point a, Point b) { return cross(a, b) / 2.0;} 76 double calsec(Point a, Point b) { 77 double da = atan2(b.y, b.x) - atan2(a.y, a.x); 78 da += da < 0 ? PI * 2.0 : 0.0; 79 return sqr(R) * da / 2.0; 80 } 81 82 int main() { 83 // freopen("in", "r", stdin); 84 // freopen("out", "w", stdout); 85 int T, n; 86 double s, t; 87 Point ip; 88 scanf("%d", &T); 89 while (T-- && ~scanf("%d", &n)) { 90 ptcnt = 0; 91 used.clear(); 92 for (int i = 0; i < n; i++) { 93 scanf("%lf%lf", &s, &t); 94 cut[i] = Line(getpt(s), getpt(t)); 95 pts[ptcnt++] = getpt(s); 96 pts[ptcnt++] = getpt(t); 97 for (int j = 0; j < i; j++) { 98 if (sgn(cross(cut[i].vec(), cut[j].vec()))) { 99 ip = llint(cut[i], cut[j]); 100 // cout << "ip " << ip.x << ' ' << ip.y << endl; 101 if (onseg(ip, cut[i])) pts[ptcnt++] = ip; 102 } 103 } 104 } 105 // cout << "pt " << ptcnt << endl; 106 sort(pts, pts + ptcnt); 107 ptcnt = unique(pts, pts + ptcnt) - pts; 108 // cout << "npt " << ptcnt << endl; 109 for (int i = 1; i <= ptcnt; i++) nb[pts[i - 1].id = i].clear(), nx[i].clear(), anx[i].clear(); 110 int ptn; 111 for (int i = 0; i < n; i++) { 112 ptn = 0; 113 for (int j = 1; j <= ptcnt; j++) { 114 if (onseg(pts[j - 1], cut[i])) tmp[ptn++] = pts[j - 1]; 115 } 116 sort(tmp, tmp + ptn); 117 for (int j = 1; j < ptn; j++) { 118 nb[tmp[j].id].push_back(Node(angle(tmp[j - 1] - tmp[j]), tmp[j - 1].id)); 119 nb[tmp[j - 1].id].push_back(Node(angle(tmp[j] - tmp[j - 1]), tmp[j].id)); 120 } 121 } 122 oc.clear(); 123 for (int i = 1; i <= ptcnt; i++) if (oncircle(pts[i - 1])) oc.push_back(Node(angle(pts[i - 1]), i)); 124 sort(oc.begin(), oc.end()); 125 // for (int i = 0; i < oc.size(); i++) cout << oc[i].id << ' '; cout << endl; 126 oc.push_back(oc[0]); 127 for (int i = 1, sz = oc.size(); i < sz; i++) { 128 nb[oc[i].id].push_back(Node(angle(pts[oc[i - 1].id - 1] - pts[oc[i].id - 1]), oc[i - 1].id)); 129 nb[oc[i].id][nb[oc[i].id].size() - 1].arc = true; 130 nb[oc[i - 1].id].push_back(Node(angle(pts[oc[i].id - 1] - pts[oc[i - 1].id - 1]), -oc[i].id)); 131 nb[oc[i - 1].id][nb[oc[i - 1].id].size() - 1].arc = true; 132 } 133 for (int i = 1; i <= ptcnt; i++) { 134 sort(nb[i].begin(), nb[i].end()); 135 // cout << i << " : " << pts[i - 1].x << ' ' << pts[i - 1].y << endl; 136 nb[i].push_back(nb[i][0]); 137 // for (int j = 0; j < nb[i].size(); j++) cout << nb[i][j].id << '-' << nb[i][j].arc << ' '; cout << endl; 138 for (int j = 1, sz = nb[i].size(); j < sz; j++) { 139 if (nb[i][j].id < 0) continue; 140 if (nb[i][j].arc) { 141 if (nb[i][j - 1].arc) { if (j < nb[i].size() - 1) nx[nb[i][j + 1].id][i] = PIB(abs(nb[i][j - 1].id), true), j++; else nx[nb[i][0].id][i] = PIB(abs(nb[i][j - 1].id), true);} 142 else anx[nb[i][j].id][i] = PIB(abs(nb[i][j - 1].id), false); 143 } else { 144 if (!nb[i][j - 1].arc || nb[i][j - 1].id < 0) nx[nb[i][j].id][i] = PIB(abs(nb[i][j - 1].id), nb[i][j - 1].arc); 145 } 146 } 147 nb[i].pop_back(); 148 } 149 // for (int i = 1; i <= ptcnt; i++) { 150 // if (anx[i].size()) cout << anx[i].size() << '~' << (*anx[i].begin()).first << '~' << (*anx[i].begin()).second.first << ' ' << nx[i].size() << endl; 151 // else puts("~~~"); 152 // } 153 double mx = 0.0, area; 154 int ls, cur; 155 PIB tt; 156 bool arc; 157 for (int i = 0; i < ptcnt; i++) { 158 for (int j = 0, sz = nb[i].size(); j < sz; j++) { 159 if (nb[i][j].arc) continue; 160 ls = i, cur = nb[i][j].id; 161 if (used.find(PII(ls, cur)) != used.end()) continue; 162 arc = false; 163 area = caltri(pts[ls - 1], pts[cur - 1]); 164 used.insert(PII(ls, cur)); 165 // cout << "start " << ls << ' '; 166 int cnt = 10; 167 while (cur != i && cnt--) { 168 // cout << cur << ' '; 169 if (arc) tt = anx[ls][cur]; 170 else tt = nx[ls][cur]; 171 ls = cur, cur = tt.first, arc = tt.second; 172 if (arc) area += calsec(pts[ls - 1], pts[cur - 1]); 173 else area += caltri(pts[ls - 1], pts[cur - 1]), used.insert(PII(ls, cur)); 174 } 175 // cout << area << endl; 176 mx = max(mx, fabs(area)); 177 } 178 } 179 printf("%.2f\n", mx); 180 } 181 return 0; 182 }
——written by Lyon