bzoj 1209: [HNOI2004]最佳包裹 三维凸包

1209: [HNOI2004]最佳包裹

Time Limit: 10 Sec  Memory Limit: 162 MB
Submit: 160  Solved: 58
[Submit][Status][Discuss]

Description

H 公司生产了一种金属制品，是由一些笔直的金属条支撑起来的，金属条和别的金属条在交点上被焊接在了一起。现在由于美观需要，在这个产品用一层特殊的材料包 裹起来。公司为了节约成本，希望消耗的材料最少（不计裁剪时的边角料的损失）。你的程序需要根据给定的输入，给出符合题意的输出： 输入包括该产品的顶点的个数，以及所有顶点的坐标； 你需要根据输入的计算出包裹这个产品所需要的材料的最小面积。 结果要求精确到小数点后第六位。（四舍五入）

Input

第1行是一个整数n(4 <= n <= 100)，表示顶点的个数；第2行到第n+1行，每行是3个实数xi,yi,zi，表示第i个顶点的坐标。每个顶点的位置各不相同。

Output

输出只有一个实数，表示包裹一个该产品所需的材料面积的最小值。

Sample Input

4
0 0 0
1 0 0
0 1 0
0 0 1 说明：该输入示例中共有4个点，可参见后面的图示。

Sample Output

2.366025

HINT

　　这道题算是三维凸包的模板题吧，话说网上的三维凸包的教程中判三点共线，四点共面的方法确实有些复杂，经过idy博客的启发，其实在做凸包前对每个点加一遍噪音就行了。如果一个点可以看见当前一个面，那么这个面就会被删除，其表现就是有向体积为负数，删除一个面的同时，删除这个面的三条边（也是有向的）最有扫一遍观察哪些边的反向边不存在，然后将他们的反向边与当前加的点一同构成新的凸包。

#include<iostream>
#include<cstdio>
#include<algorithm>
#include<cstring>
#include<cmath>
#include<assert.h>
#include<vector>
using namespace std;
#define MAXN 120
#define eps 1e-11
typedef double real;
inline int sgn(real x)
{
        if (abs(x)<eps)return 0;
        return x<0?-1:1;
}
struct point
{
        real x,y,z;
        point(real x,real y,real z):x(x),y(y),z(z){}
        point(){}
        int read()
        {
                return scanf("%lf %lf %lf\n",&x,&y,&z);
        }
        void noise()
        {
                x+=(real)(rand()%10000-5000)/10000*eps;
                y+=(real)(rand()%10000-5000)/10000*eps;
                z+=(real)(rand()%10000-5000)/10000*eps;
        }
        real len()
        {
                return sqrt(x*x+y*y+z*z);
        }
};
typedef point vect;
struct line
{
        point ps;
        real x,y,z;
        line(){}
        line(point p1,point p2)
        {
                ps=p1;
                x=p2.x-p1.x;
                y=p2.y-p1.y;
                z=p2.z-p1.z;
        }
        line(point ps,real x,real y,real z):ps(ps),x(x),y(y),z(z){}
};
vect xmul(line l1,line l2)
{
        return vect(l1.y*l2.z-l1.z*l2.y,l1.z*l2.x-l1.x*l2.z,l1.x*l2.y-l1.y*l2.x);
}
real volume(line l1,line l2,line l3)
{
        return + l1.x*l2.y*l3.z - l1.x*l2.z*l3.y 
               - l1.y*l2.x*l3.z + l1.y*l2.z*l3.x
               + l1.z*l2.x*l3.y - l1.z*l2.y*l3.x;
}
//+1 2 3
//-1 3 2
//-2 1 3
//+2 3 1
//+3 1 2
//-3 2 1
struct surface
{
        point ps;
        real x1,y1,z1;
        real x2,y2,z2;
        surface(){}
        surface(point p1,point p2,point p3)
        {
                ps=p1;
                x1=p2.x-p1.x,y1=p2.y-p1.y,z1=p2.z-p1.z;
                x2=p3.x-p1.x,y2=p3.y-p1.y,z2=p3.z-p1.z;
        }
        real volume(point pt)
        {
                return ::volume(line(ps,pt),line(ps,x1,y1,z1),line(ps,x2,y2,z2));
        }
        vect nvect()
        {
                return xmul(line(ps,x1,y1,z1),line(ps,x2,y2,z2));
        }
        void reverse()
        {
                swap(x1,x2);
                swap(y1,y2);
                swap(z1,z2);
        }
};
point pl[MAXN];
struct face
{
        int pt[3];
        face(int x,int y,int z)
        {
                pt[0]=x;pt[1]=y;pt[2]=z;
        }
        surface ToSurface()
        {
                return surface(pl[pt[0]],pl[pt[1]],pl[pt[2]]);
        }
        void print()
        {
                printf("Face:%d %d %d\n",pt[0],pt[1],pt[2]);
        }
};
vector<face> cc;
vector<pair<int,int> > chs;
bool status[MAXN][MAXN];
int main()
{
        freopen("input.txt","r",stdin);
        int n;
        scanf("%d",&n);
        for (int i=0;i<n;i++)
                pl[i].read();
        for (int i=0;i<n;i++)
                pl[i].noise();
    /*    for (int i=0;i<n;i++)
                swap(pl[rand()%n],pl[rand()%n]);*/
        cc.push_back(face(0,1,2));
        cc.push_back(face(2,1,0));
        for (int i=0;i<3;i++)
                status[i][(i+1)%3]=true;
        for (int i=1;i<4;i++)
                status[i%3][i-1]=true;
        for (int i=3;i<n;i++)
        {
                //for (int j=0;j<cc.size();j++)cc[j].print();    printf("\n");
                chs.clear();
                for (int j=0;j<cc.size();j++)
                {
                        if (cc[j].ToSurface().volume(pl[i])<0)
                        {
                                for (int k=0;k<3;k++)
                                {
                                        status[cc[j].pt[k]][cc[j].pt[(k+1)%3]]=false;
                                        chs.push_back(make_pair(cc[j].pt[k],cc[j].pt[(k+1)%3]));
                                }
                                swap(cc[j],cc[cc.size()-1]);
                                cc.pop_back();
                                j--;
                        }
                }
                for (int j=0;j<chs.size();j++)
                {
                        if (!status[chs[j].first][chs[j].second] && status[chs[j].second][chs[j].first])continue;
                        chs[j]=chs[chs.size()-1];
                        j--;
                        chs.pop_back();
                }
                for (int j=0;j<chs.size();j++)
                {
                        cc.push_back(face(i,chs[j].first,chs[j].second));
                        status[i][chs[j].first]=status[chs[j].first][chs[j].second]=status[chs[j].second][i]=true;
                }
                for (int j=0;j<n;j++)
                        for (int k=0;k<n;k++)
                                assert(!(status[j][k]^status[k][j]));
        }
        //for (int j=0;j<cc.size();j++)cc[j].print();    printf("\n");
        real ans=0;
        for (int i=0;i<cc.size();i++)
                ans+=cc[i].ToSurface().nvect().len()/2;
        printf("%.6lf\n",abs(ans));
}