一道三维凸包的题目，题目要求求出三维凸包的表面积。看懂了网上的三维凸包的代码以后，自己写的代码，跟网上的模板有所不同。调了一个晚上，结果发现错的只是数组开太小。_(:з」∠)_
代码如下：
1 #include <cstdio> 2 #include <iostream> 3 #include <algorithm> 4 #include <cstring> 5 #include <cmath> 6 7 using namespace std; 8 9 const double EPS = 1e-10; 10 const int N = 555; 11 inline int sgn(double x) { return (x > EPS) - (x < -EPS);} 12 struct Point { 13 double x, y, z; 14 Point() {} 15 Point(double x, double y, double z) : x(x), y(y), z(z) {} 16 bool operator < (Point a) const { 17 if (sgn(x - a.x)) return x < a.x; 18 if (sgn(y - a.y)) return y < a.y; 19 return z < a.z; 20 } 21 bool operator == (Point a) const { return !sgn(x - a.x) && !sgn(y - a.y) && !sgn(z - a.z);} 22 Point operator + (Point a) { return Point(x + a.x, y + a.y, z + a.z);} 23 Point operator - (Point a) { return Point(x - a.x, y - a.y, z - a.z);} 24 Point operator * (double p) { return Point(x * p, y * p, z * p);} 25 Point operator / (double p) { return Point(x / p, y / p, z / p);} 26 } ; 27 inline double dot(Point a, Point b) { return a.x * b.x + a.y * b.y + a.z * b.z;} 28 inline double dot(Point o, Point a, Point b) { return dot(a - o, b - o);} 29 inline Point cross(Point a, Point b) { return Point(a.y * b.z - a.z * b.y, a.z * b.x - a.x * b.z, a.x * b.y - a.y * b.x);} 30 inline Point cross(Point o, Point a, Point b) { return cross(a - o, b - o);} 31 inline double veclen(Point x) { return sqrt(dot(x, x));} 32 inline double area(Point o, Point a, Point b) { return veclen(cross(a - o, b - o)) / 2.0;} 33 // directed volume 34 inline double volume(Point o, Point a, Point b, Point c) { return dot(cross(a - o, b - o), c - o);} 35 36 struct _3DCH { 37 Point pt[N]; 38 int pcnt; 39 struct Face { // follow right-hand 40 int a, b, c; 41 bool ok; 42 Face() {} 43 Face(int a, int b, int c) : a(a), b(b), c(c) { ok = true;} 44 } ; 45 Face f[N << 3]; 46 int fcnt, fid[N][N]; 47 bool outside(int p, int a, int b, int c) { 48 int dir = sgn(volume(pt[a], pt[b], pt[c], pt[p])); 49 return dir < 0; 50 } 51 bool outside(int p, int x) { return outside(p, f[x].a, f[x].b, f[x].c);} 52 void addface(int p, int a, int b, int c) { 53 int dir = sgn(volume(pt[a], pt[b], pt[c], pt[p])); 54 if (dir == 0) return ; 55 if (dir > 0) { 56 f[fcnt] = Face(a, b, c); 57 fid[a][b] = fid[b][c] = fid[c][a] = fcnt++; 58 } else { 59 f[fcnt] = Face(c, b, a); 60 fid[c][b] = fid[b][a] = fid[a][c] = fcnt++; 61 } 62 } 63 bool dfs(int p, int x) { 64 if (!f[x].ok) return true; 65 if (outside(p, x)) { 66 f[x].ok = false; 67 if (!dfs(p, fid[f[x].b][f[x].a])) addface(f[x].c, f[x].b, f[x].a, p); 68 if (!dfs(p, fid[f[x].c][f[x].b])) addface(f[x].a, f[x].c, f[x].b, p); 69 if (!dfs(p, fid[f[x].a][f[x].c])) addface(f[x].b, f[x].a, f[x].c, p); 70 return true; 71 } 72 return false; 73 } 74 bool build() { 75 bool ok; 76 fcnt = 0; 77 sort(pt, pt + pcnt); 78 pcnt = unique(pt, pt + pcnt) - pt; 79 ok = false; 80 for (int i = 2; i < pcnt; i++) { 81 if (sgn(veclen(cross(pt[i], pt[0], pt[1])))) { 82 swap(pt[2], pt[i]); 83 ok = true; 84 break; 85 } 86 } 87 if (!ok) return false; 88 ok = false; 89 for (int i = 3; i < pcnt; i++) { 90 if (sgn(volume(pt[i], pt[0], pt[1], pt[2]))) { 91 swap(pt[3], pt[i]); 92 ok = true; 93 break; 94 } 95 } 96 if (!ok) return false; 97 addface(0, 1, 2, 3); 98 addface(1, 2, 3, 0); 99 addface(2, 3, 0, 1); 100 addface(3, 0, 1, 2); 101 for (int i = 4; i < pcnt; i++) { 102 for (int j = fcnt - 1; j >= 0; j--) { 103 if (f[j].ok && outside(i, j)) { 104 dfs(i, j); 105 break; 106 } 107 } 108 } 109 int sz = fcnt; 110 fcnt = 0; 111 for (int i = 0; i < sz; i++) if (f[i].ok) f[fcnt++] = f[i]; 112 return true; 113 } 114 double vol() { 115 Point O = Point(0.0, 0.0, 0.0); 116 double ret = 0.0; 117 for (int i = 0; i < fcnt; i++) ret += volume(O, pt[f[i].a], pt[f[i].b], pt[f[i].c]); 118 return fabs(ret) / 6.0; 119 } 120 double facearea() { 121 double ret = 0.0; 122 for (int i = 0; i < fcnt; i++) ret += area(pt[f[i].a], pt[f[i].b], pt[f[i].c]); 123 return ret; 124 } 125 } ch; 126 127 int main() { 128 // freopen("in", "r", stdin); 129 while (cin >> ch.pcnt) { 130 for (int i = 0; i < ch.pcnt; i++) scanf("%lf%lf%lf", &ch.pt[i].x, &ch.pt[i].y, &ch.pt[i].z); 131 ch.build(); 132 printf("%.3f ", ch.facearea()); 133 = } 134 return 0; 135 }
——written by Lyon