BZOJ1209 最佳包裹 （三维凸包 增量法）

题意

求三维凸包的表面积。
$N\le100$

题解

暴力往当前的凸包里加点。$O(n^2)$。题解详见大佬博客

扰动函数$shake$是为了避免四点共面。

CODE

实测$eps$开到$1e-10$才过

#include <bits/stdc++.h>
using namespace std;
#define il inline
const double eps = 1e-10;
const int MAXN = 105;
il double Rand() { return rand()/(double)RAND_MAX; }
il double reps() { return (Rand()-0.5)*eps; };
int n;
struct point {
	double x, y, z;
	il void shake() { x+=reps(); y+=reps(); z+=reps(); }
	il double len() { return sqrt(x*x + y*y + z*z); }
	il point operator -(point o) { return (point){ x-o.x , y-o.y, z-o.z }; }
	il point operator *(point o) { return (point){ y*o.z-z*o.y, z*o.x-x*o.z, x*o.y-y*o.x }; }
	il double operator &(point o) { return x*o.x + y*o.y + z*o.z; }
}A[MAXN];
struct face {
	int v[3];
	il point normal() { return (A[v[1]]-A[v[0]]) * (A[v[2]]-A[v[0]]); }
	il double area() { return normal().len() / 2; }
}f[MAXN<<1], tmp[MAXN<<1];
bool see(face a, point b) { return ( (b-A[a.v[0]])&a.normal() ) > 0; }
int cnt;
bool vis[MAXN][MAXN];
void Convex3D() {
	f[++cnt] = (face) { {1, 2, 3} };
	f[++cnt] = (face) { {3, 2, 1} };
	for(int i = 4, cur; i <= n; ++i) {
		cur = 0;
		for(int j = 1, can; j <= cnt; ++j) {
			if(!(can = see(f[j], A[i]))) tmp[++cur] = f[j];
			for(int k = 0; k < 3; ++k) vis[f[j].v[k]][f[j].v[(k+1)%3]] = can;
		}
		for(int j = 1; j <= cnt; ++j)
			for(int k = 0; k < 3; ++k) {
				int u = f[j].v[k], v = f[j].v[(k+1)%3];
				if(vis[u][v] && !vis[v][u]) tmp[++cur] = (face){ {u, v, i} };
			}
		for(int j = 1; j <= cur; ++j) f[j] = tmp[j]; cnt = cur;
	}
}
int main() {
	srand(19260817);
	cin>>n;
	for(int i = 1; i <= n; ++i) cin>>A[i].x>>A[i].y>>A[i].z, A[i].shake();
	Convex3D();
	double S = 0;
	for(int i = 1; i <= cnt; ++i) S += f[i].area();
	printf("%.6f\n", S);
}