三维凸包

三维凸包（已经忘光了

#include <bits/stdc++.h>
using namespace std;
typedef double db;
const int N = 2100;
const db eps = 1e-9;
int n;
int vis[N][N];
db Ran() { return rand() / (double)RAND_MAX; }
db reps() { return (Ran() - 0.5) * eps; }
struct vec
{
db x, y, z;
void in()
  {
scanf("%lf%lf%lf", &x, &y, &z);
return;
  }
vec(db a = 0, db b = 0, db c = 0) : x(a), y(b), z(c) {}
vec operator+(const vec r)
  {
return vec(x + r.x, y + r.y, z + r.z);
  }
vec operator-(const vec r)
  {
return vec(x - r.x, y - r.y, z - r.z);
  }
vec operator*(const vec r)
  {
return vec(y * r.z - z * r.y, z * r.x - x * r.z, x * r.y - y * r.x);
  }
db operator&(const vec r)
  {
return x * r.x + y * r.y + z * r.z;
  }
} A[N];
typedef vec pon;
db dis(vec t)
{
return sqrt(t.x * t.x + t.y * t.y + t.z * t.z);
}
void shake(pon &t)
{
t.x += reps(), t.y += reps(), t.z += reps();
return;
}
struct Face
{
int v[3];
Face() { v[0] = v[1] = v[2] = 0; }
Face(int a, int b, int c) { v[0] = a, v[1] = b, v[2] = c; }
vec norm()
  {
return (A[v[1]] - A[v[0]]) * (A[v[2]] - A[v[0]]);
  }
} f[N], c[N];
int see(Face a, pon b)
{
return ((b - A[a.v[0]]) & a.norm()) > 0;
}
db dist(Face a, pon b) //有向
{
vec nor = a.norm();
return (nor & (b - A[a.v[0]])) / dis(nor);
}

db area(Face a)
{
return dis(a.norm()) / 2;
}
db volu(Face a)
{
db d = dist(a, pon(0, 0, 0));
return area(a) * d / 3;
}
int cnt = 0;
void Convex3()
{
cnt = 0;
f[++cnt] = Face(1, 2, 3);
f[++cnt] = Face(3, 2, 1);
for (int i = 4, cc = 0; i <= n; i++)
  {
for (int j = 1, v; j <= cnt; j++)
    {
if (!(v = see(f[j], A[i])))
c[++cc] = f[j];
for (int k = 0; k < 3; k++)
vis[f[j].v[k]][f[j].v[(k + 1) % 3]] = v;
    }
for (int j = 1; j <= cnt; j++)
    {
for (int k = 0; k < 3; k++)
      {
int x = f[j].v[k], y = f[j].v[(k + 1) % 3];
if (vis[x][y] && !vis[y][x])
c[++cc] = Face(x, y, i);
      }
    }
for (int j = 1; j <= cc; j++)
f[j] = c[j];
cnt = cc, cc = 0;
  }
return;
}
int main()
{
srand(time(NULL));
n = 1;
while (scanf("%lf%lf%lf", &A[n].x, &A[n].y, &A[n].z) == 3)
shake(A[n]), n++;
n--;
Convex3();
return 0;
}