UVa 1433 - Aerodynamics

这是我第一次使用多个namespace 写程序 ， 对于3D和2D要同时出现的题目完全可以把3D中的每一个变量后加一个3 ， 但是完全没有下面的写法方便。

训练指南的三维凸包写法非常实用 ， 进行些许扰动 ， 可以减少很多麻烦。 这个题指明了： no two points coincide , 所以不需要判重

提示：
1. 求凸包
2. 然后关于z的一个平面一个平面的枚举就好了

#include <iostream>
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <cmath>
#include <string>
#include <vector>
#include <deque>
#include <stack>
#include <algorithm>

using namespace std;

const double eps = 1e-10;
const int maxn = 110;

int dcmp(double a) { return fabs(a)<eps?0:(a<0?-1:1); }

namespace TwoD
{
    struct points
    {
        double x , y;
        void read(){ scanf("%lf%lf",&x,&y); }
        points(double x = 0 , double y = 0):x(x),y(y){}
        bool operator <(const points& b)const { return dcmp(x-b.x)==-1 || (dcmp(x-b.x)==0 && dcmp(y-b.y)==-1); }
        bool operator ==(const points& b)const { return dcmp(x-b.x)==0 && dcmp(y-b.y)==0; }
    };

    typedef points Vector;
    typedef vector<points> polygon;

    Vector operator +(Vector a , Vector b) { return Vector(a.x+b.x , a.y+b.y); }
    Vector operator -(Vector a , Vector b) { return Vector(a.x-b.x , a.y-b.y); }
    Vector operator *(Vector a , double b) { return Vector(a.x*b   , a.y*b  ); }
    Vector operator /(Vector a , double b) { return Vector(a.x/b   , a.y/b  ); }

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

    points ch[maxn];
    int convexHull(polygon p)
    {
        sort(p.begin(), p.end());
        p.erase(unique(p.begin(), p.end()), p.end());

        int n = p.size();
        int m = 0;
        for(int i=0;i<n;i++)
        {
            while(m>1 && dcmp(Cross(ch[m-1]-ch[m-2], p[i]-ch[m-2]))<=0) m--;
            ch[m++] = p[i];
        }

        int k = m;

        for(int i=n-2;i>=0;i--)
        {
            while(m>k && dcmp(Cross(ch[m-1]-ch[m-2], p[i]-ch[m-2]))<=0) m--;
            ch[m++] = p[i];
        }

        if(n>1) m--;
        return m;
    }

    double area(int m)
    {
        double res = 0;
        for(int i=1;i<m-1;i++) res+= Cross(ch[i]-ch[0], ch[i+1]-ch[0]);
        return res/2.0;
    }
};

namespace ThreeD
{
    struct points
    {
        double x , y , z;
        void read(){ scanf("%lf%lf%lf",&x,&y,&z); }
        points(double x =0  , double y =0  , double z =0 ):x(x),y(y),z(z){}
    };
    typedef points Vector;

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

    double rand01() { return rand() / (double)RAND_MAX; }
    double randeps() { return (rand01()-0.5)*eps; }
    points addnoise(points p) { return points(p.x+randeps() , p.y+randeps() , p.z+randeps()); }

    struct face
    {
        int v[3];

        face(int v1=0, int v2=0, int v3=0){ v[0] = v1; v[1] = v2; v[2] = v3; }
        Vector normal(points *p)
        {
            return Cross(p[v[1]]-p[v[0]], p[v[2]]-p[v[0]]);
        }

        int cansee(points* p , int i) 
        {
            return Dot(p[i]-p[v[0]] , normal(p))>0?1:0;
        }
    };

    points p[maxn];
    points t[maxn];
    int vis[maxn][maxn];

    vector<face> convexHull(points* p , int n)
    {
        vector<face> res;
        res.push_back(face(0,1,2));
        res.push_back(face(2,1,0));

        for(int i=3;i<n;i++)
        {
            vector<face> Next;
            for(int j=0;j<res.size();j++)
            {
                int hi = res[j].cansee(p, i);
                if(!hi) Next.push_back(res[j]);
                for(int k=0;k<3;k++) vis[res[j].v[k]][res[j].v[(k+1)%3]] = hi;
            }

            for(int j=0;j<res.size();j++) for(int k=0;k<3;k++)
            {
                int a = res[j].v[k] , b = res[j].v[(k+1)%3];
                if(vis[a][b]!=vis[b][a] && vis[a][b]) Next.push_back(face(a,b,i));
            }
            res = Next; 
        }
        return res; 
    }
};

using namespace ThreeD;


int main(int argc, char *argv[]) {

    int n , zmin , zmax , Case=0;
    while(cin>>n>>zmin>>zmax)
    {
        if(Case) cout<<endl; Case++;

        for(int i=0;i<n;i++) p[i].read();

        for(int i=0;i<n;i++) t[i] = addnoise(p[i]);
        vector<face> c = convexHull(t, n);

        for(int z=zmin;z<=zmax;z++)
        {
            TwoD::polygon poly;
            for(int i=0;i<c.size();i++) for(int j=0;j<3;j++)
            {
                int a = c[i].v[j] , b = c[i].v[(j+1)%3];
                double z1 = p[a].z , z2 = p[b].z , x1 = p[a].x , x2 = p[b].x , y1 = p[a].y , y2 = p[b].y;
                if(dcmp(z1-z2)==0) continue;
                if(dcmp(p[a].z-z)*dcmp(p[b].z-z)<=0) 
                    poly.push_back(TwoD::points(x1+(z-z1)/(z2-z1)*(x2-x1) , y1+(z-z1)/(z2-z1)*(y2-y1)));
            }
            int m = TwoD::convexHull(poly);
//          cout<<poly.size()<<" "<<m<<endl; 
            printf("%.5lf\n",TwoD::area(m)); 
        }
    }

    return 0;
}

程序的时间表现相当优秀:-)

详解：
1. 在枚举每一个平面的时候 ， 找到这个平面在凸包上的所有交点 ， 然后把这些交点排序 ， 这里我用了二维凸包 ， 这更加鲁棒

注意：
两个z值相等的线段不做考量 ， 那会不会漏掉些点呢？ 不会 ， 因为每个点会在多个平面中出现 ， 所以不用担心