hdu4273 Rescue（三维凸包）

题目连接

分析：
三维凸包重心到表面的最短距离

求出凸包重心后，O(n)求出重心到凸包各面的距离，取最小值即可

一开始不大理解：向凸包上一面做垂线，垂足可能不在这个面上，怎么破？
实际上，如果出现这样的状况，那么这个面上的垂线一定不是最优解

#include<bits/stdc++.h>

using namespace std;

const int N=110;
const double eps=1e-8;

struct node{
    double x,y,z;
    node (double xx=0,double yy=0,double zz=0)
    {
        x=xx;y=yy;z=zz;
    }
};

node operator +(const node &a,const node &b) {return node(a.x+b.x,a.y+b.y,a.z+b.z);}
node operator -(const node &a,const node &b) {return node(a.x-b.x,a.y-b.y,a.z-b.z);}
node operator *(const node &a,const double &b) {return node(a.x*b,a.y*b,a.z*b);}
node operator /(const node &a,const double &b) {return node(a.x/b,a.y/b,a.z/b);}

double Dot(node a,node b) {return a.x*b.x+a.y*b.y+a.z*b.z;}
node Cross(node a,node b) {return node(a.y*b.z-a.z*b.y,a.z*b.x-a.x*b.z,a.x*b.y-a.y*b.x);}

struct CH3D{
    struct face{
        int a,b,c;
        bool ff;
    };

    int n;
    node P[N];
    int num;
    face F[8*N];
    int vis[N][N];

    double lenth(node a) {
        return sqrt(Dot(a,a));
    }
    double S(node a,node b,node c) {
        return lenth(Cross(b-a,c-a)); 
    }
    double V(node a,node b,node c,node d) {
        return Dot(d-a,Cross(b-a,c-a));
    }
    int cansee(node p,face f) {
        double x=Dot(p-P[f.a],Cross(P[f.b]-P[f.a],P[f.c]-P[f.a]));
        if (x>eps) return 1;
        else return 0;
    }

    void deal(int p,int a,int b) {
        int f=vis[a][b];
        face add;
        if (F[f].ff) {
            if (cansee(P[p],F[f]))
                dfs(p,f);
            else {
                add.a=b;
                add.b=a;
                add.c=p;
                add.ff=1;
                vis[b][a]=vis[a][p]=vis[p][b]=num;
                F[num++]=add;
            }
        }
    }
    void dfs(int p,int now) {
        F[now].ff=0;
        deal(p,F[now].b,F[now].a);   //反着转一圈 
        deal(p,F[now].c,F[now].b);
        deal(p,F[now].a,F[now].c);
    }
    void create() {
        face add;

        num=0;
        if (n<4) return;

        bool flag=1;
        for (int i=1;i<n;i++)
            if (lenth(P[0]-P[i])>eps) {
                swap(P[1],P[i]);
                flag=0; break;
            }
        if (flag) return;
        flag=1;
        for (int i=2;i<n;i++)
            if (S(P[0],P[1],P[i])>eps) {
                swap(P[2],P[i]);
                flag=0; break;
            }
        if (flag) return;
        flag=1;
        for (int i=3;i<n;i++)
            if (fabs(V(P[0],P[1],P[2],P[i]))>eps) {
                swap(P[3],P[i]);
                flag=0; break;
            }
        if (flag) return;

        for (int i=0;i<4;i++)
        {
            add.a=(i+1)%4;
            add.b=(i+2)%4;
            add.c=(i+3)%4;
            add.ff=1;
            if (cansee(P[i],add)) swap(add.b,add.c);
            vis[add.a][add.b]=vis[add.b][add.c]=vis[add.c][add.a]=num;
            F[num++]=add;
        }
        for (int i=4;i<n;i++)
            for (int j=0;j<num;j++)
                if (F[j].ff&&cansee(P[i],F[j])) {
                    dfs(i,j);
                    break;
                }
        int tmp=num; num=0;
        for (int i=0;i<tmp;i++)
            if (F[i].ff) 
                F[num++]=F[i];
    }

    node centre()
    {
        node ans(0,0,0),o(0,0,0);
        double all=0;
        for (int i=0;i<num;i++) {
            double vol=V(o,P[F[i].a],P[F[i].b],P[F[i].c]);
            ans=ans+(o+P[F[i].a]+P[F[i].b]+P[F[i].c])/4.0*vol;
            all+=vol;
        }
        ans=ans/all;
        return ans;
    }
    double dis_face(node p,int i)
    {
        return fabs(V(P[F[i].a],P[F[i].b],P[F[i].c],p)/lenth(Cross(P[F[i].b]-P[F[i].a],P[F[i].c]-P[F[i].a])));
    }
};

CH3D hull;

int main()
{
    while (scanf("%d",&hull.n)!=EOF)
    {
        for (int i=0;i<hull.n;i++)
            scanf("%lf%lf%lf",&hull.P[i].x,&hull.P[i].y,&hull.P[i].z);
        hull.create();
        node p=hull.centre();
        double ans=1e30;
        for (int i=0;i<hull.num;i++)
            ans=min(ans,hull.dis_face(p,i));
        printf("%0.3lf\n",ans);
    }
}