HDU 4273 Rescue(三维凸包重心到表面最短距离)

Rescue

Time Limit: 2000/1000 MS (Java/Others)    Memory Limit: 32768/32768 K (Java/Others)
Total Submission(s): 393    Accepted Submission(s): 291


Problem Description
I work at NASA outer space rescue team which needs much courage and patient. In daily life, I always receive a lot of mission, and I must complete it right now.
Today, team leader announced me that there is a huge spaceship dropping anchor in the out space, and we should reach there for rescue. As a working principle, at first, we should check whether there are persons living in the spaceship. So we carry a kind of machine called life sensor which can sense the life phenomenon when the distance between the machine and the living is not farther than the sense radius.
I have read the designing paper of the spaceship in advance. It has a form of a convex polyhedron, and we can assume it is isodense. For best control, control center of the whole ship is located at the center of the mass. It is sure that if someone is still alive, he will stay at the control center.
It's unfortunately that I find the door is stocked when I try to enter into the spaceship, so I can only sense the living out of the space ship. Now I have opened the machine and it's time to set the sense radius of it. I wonder the minimal radius of the machine which can allowe me to check whether there are persons living in the spaceship.
 

Input
There are multiple test cases.
The first line contains an integer n indicating the number of vertices of the polyhedron. (4 <= n <= 100)
Each of the next n lines contains three integers xi, yi, zi, the coordinates of the polyhedron vertices (-10,000 <= xi, yi, zi <= 10,000).
It guaranteed that the given points are vertices of the convex polyhedron, and the polyhedron is non-degenerate.
 

Output
For each test case, output a float number indicating the minimal radius of the machine. Your answer should accurate up to 0.001.
 

Sample Input
  
  
4 0 0 0 1 0 0 0 1 0 0 0 1 8 0 0 0 0 0 2 0 2 0 0 2 2 2 0 0 2 0 2 2 2 0 2 2 2
 

Sample Output
  
  
0.144 1.000
先构建一个三维凸包 求出其重心 ,然后求出重心到每个面的距离,找出最小值
#include <iostream>
#include <cstdio>
#include <cstring>
#include <stdlib.h>
#include <algorithm>
#include <queue>
#include <map>
#include <cmath>
#define N 101
#define ll long long int
#define eps 1e-8
#define PI acos(-1.0)
#define inf 0x3f3f3f3f
using namespace std;
struct point
{
    double x,y,z;
    point() {}
    point(double xx,double yy,double zz): x(xx),y(yy),z(zz) {}
    point operator -(const point p)
    {
        return point(x-p.x, y-p.y, z-p.z);
    }
    point operator +(const point p)
    {
        return point(x+p.x, y+p.y, z+p.z);
    }
    point operator *(const point p)
    {
        return point(y*p.z-z*p.y, z*p.x-x*p.z, x*p.y-y*p.x);
    }
    point operator *(double p)
    {
        return point(x*p, y*p,z*p);
    }
    point operator /(double p)
    {
        return point(x/p, y/p,z/p);
    }
    double operator ^( point p)
    {
        return x*p.x+y*p.y+z*p.z;
    }

};
struct node
{
    struct face
    {
        int a,b,c;
        bool ok;
    };
    int n;
    point tn[N];
    int num;
    face F[8*N];
    int g[N][N];
    double vlen(point a)
    {
        return sqrt(a.x*a.x+a.y*a.y+a.z*a.z);
    }
    point cross(const point &a,const point &b,const point &c)
    {
        return point ( (b.y-a.y)*(c.z-a.z) - (b.z-a.z)*(c.y-a.y),
                       (b.z-a.z)*(c.x-a.x) - (b.x-a.x)*(c.z-a.z),
                       (b.x-a.x)*(c.y-a.y) - (b.y-a.y)*(c.x-a.x)

                     );
    }
    double area(point a,point b,point c)
    {
        return vlen((b*a)*(c-a));
    }
    double volume(point a,point b,point c,point d)
    {
        return (b-a)*(c-a)^(d-a);
    }
    double dbcmp(point &p,face &f)
    {
        point m = tn[f.b] - tn[f.a];
        point n = tn[f.c] - tn[f.a];
        point t = p - tn[f.a];
        return (m*n)^t;
    }
    void deal(int p,int a,int b)
    {
        int f = g[a][b];
        face add;
        if(F[f].ok)
        {
            if(dbcmp(tn[p],F[f]) > eps)
                dfs(p,f);
            else
            {
                add.a=b;
                add.b=a;
                add.c=p;
                add.ok=true;
                g[p][b] = g[a][p] = g[b][a] =num;
                F[num++] = add;
            }
        }
    }
    void dfs(int p,int now)
    {
        F[now].ok = false;
        deal(p,F[now].b,F[now].a);
        deal(p,F[now].c,F[now].b);
        deal(p,F[now].a,F[now].c);

    }
    bool same(int s,int t)
    {
        point &a = tn[F[s].a];
        point &b = tn[F[s].b];
        point &c = tn[F[s].c];
        return fabs( volume(a,b,c,tn[F[t].a])) < eps &&
               fabs( volume(a,b,c,tn[F[t].b])) < eps &&
               fabs( volume(a,b,c,tn[F[t].c])) < eps ;
    }
    void create()
    {
        int i,j,tmp;
        face add;
        num = 0;
        if(n<4) return ;
        bool flag =true;
        for(i=1; i<n; i++)
        {
            if(vlen(tn[0]-tn[i]) > eps)
            {
                swap(tn[1],tn[i]);
                flag = false;
                break;
            }
        }
        if(flag) return ;
        flag =true;
        for(i=2; i<n; i++)
        {
            if(vlen((tn[0]-tn[1])*(tn[1]-tn[i])) > eps )
            {
                swap(tn[2],tn[i]);
                flag =false;
                break;
            }
        }
        if(flag) return ;
        flag =true;
        for(i =3; i<n; i++)
        {
            if(fabs((tn[0]-tn[1])*(tn[1]-tn[2])^(tn[0]-tn[i])) >eps    )
            {
                swap(tn[3],tn[i]);
                flag =false;
                break;
            }
        }
        if(flag) return ;
        for(i=0; i<4; i++)
        {
            add.a = (i+1)%4;
            add.b = (i+2)%4;
            add.c = (i+3)%4;
            add.ok = true;
            if(dbcmp(tn[i],add) > 0)
            {
                swap(add.b,add.c);
            }
            g[add.a][add.b] =   g[add.b][add.c] =  g[add.c][add.a] =num;
            F[num++]=add;
        }
        for(i =4; i<n; i++)
        {
            for(j =0; j<num; j++)
            {
                if(F[j].ok && dbcmp(tn[i],F[j])>eps)
                {
                    dfs(i,j);
                    break;
                }
            }
        }
        tmp =num;
        for(i=num=0; i<tmp; i++)
        {
            if(F[i].ok)
                F[num++] = F[i];
        }
    }
    double area()
    {
        double res=0;
        if(n==3)
        {
            point p = cross(tn[0],tn[1],tn[2]);
            res = vlen(p)/2.0;
            return res;
        }
        for(int i=0; i<num; i++)
        {
            res += area(tn[F[i].a],tn[F[i].b],tn[F[i].c]);
        }
        return res/2.0;
    }
    double volume()
    {
        double res=0;
        point tmp(0,0,0);
        for(int i = 0; i<num; i++)
            res += volume(tmp,tn[F[i].a],tn[F[i].b],tn[F[i].c]);
        return fabs(res/6.0);
    }
    int triangle()
    {
        return num;
    }
    int polygon()
    {
        int i,j,res,flag;
        for(i=res=0; i<num; i++)
        {
            flag=1;
            for(j=0; j<i; j++)
            {
                if(same(i,j))
                {
                    flag=0;
                    break;
                }
            }
            res+=flag;
        }
        return res;
    }
    point barycenter()
    {
        point ans(0,0,0),o(0,0,0);
        double all=0;
        for(int i=0; i<num; i++)
        {
            double vol = volume(o,tn[F[i].a],tn[F[i].b],tn[F[i].c]);
            ans = ans+(o+tn[F[i].a]+tn[F[i].b]+tn[F[i].c])/4.0*vol;
            all+=vol;
        }
        ans = ans/all;
        return ans;
    }
    double ptoface(point p,int i)
    {
        return fabs(volume(tn[F[i].a],tn[F[i].b],tn[F[i].c],p) /
                    vlen( (tn[F[i].b]-tn[F[i].a])*(tn[F[i].c]-tn[F[i].a]) ));
    }
};
node ans;
int main()
{
    while(~scanf("%d",&ans.n))
    {
        for(int i=0; i<ans.n; i++)
        {
            scanf("%lf%lf%lf",&ans.tn[i].x,&ans.tn[i].y,&ans.tn[i].z);
        }
        ans.create();
        point p=ans.barycenter();
        double answer=1e20;
        for(int i=0; i<ans.num; i++)
        {
            answer=min(answer,ans.ptoface(p,i));
        }
        printf("%.3lf\n",answer);
    }
    return 0;
}


  • 0
    点赞
  • 0
    收藏
    觉得还不错? 一键收藏
  • 0
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值