P4724 【模板】三维凸包

题目连接

题面:
在这里插入图片描述

代码:

#include<iostream>
#include<cstdlib>
#include<cstdio>
#include<algorithm>
#include<cstring>
#include<string>
#include<cmath>
#include<queue>
#include<map>
#include<vector>
#define ll long long
#define llu unsigned ll
using namespace std;

const double eps=1e-8;
const double dnf=1e20;
const double pi=acos(-1.0);
const int maxn=2050;

//浮点型数值是否为0
int sgn(double x)
{
    if(abs(x)<eps) return 0;
    if(x<0) return -1;
    return 1;
}

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

    void input(void)
    {
        scanf("%lf%lf%lf",&x,&y,&z);
    }
    void output(void)
    {
        printf("%.2f %.2f %.2f\n",x,y,z);
    }

    double len(void)
    {
        return sqrt(x*x+y*y+z*z);
    }

    double len2(void)
    {
        return x*x+y*y+z*z;
    }

    double dis(const Point3 &b) const
    {
        return sqrt((x-b.x)*(x-b.x)+(y-b.y)*(y-b.y)+(z-b.z)*(z-b.z));
    }

    bool operator ==(const Point3 &b) const
    {
        return sgn(x-b.x)==0&&sgn(y-b.y)==0&&sgn(z-b.z)==0;
    }

    bool operator <(const Point3 &b) const
    {
        if(sgn(x-b.x)!=0) return x<b.x;
        if(sgn(y-b.y)!=0) return y<b.y;
        return sgn(z-b.z)<0;
    }

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

    Point3 operator ^ (const Point3 &b) const
    {
        return Point3(y*b.z-z*b.y,z*b.x-x*b.z,x*b.y-y*b.x);
    }
    double operator * (const Point3 &b) const
    {
        return x*b.x+y*b.y+z*b.z;
    }
};

struct CH3D
{
    struct face
    {
        int a,b,c;//表示凸包一个面上的三个点的编号
        bool ok;//表示该面是否属于最终的凸包上的面
    };

    int n;//初始顶点数
    Point3 P[maxn];
    int num;//凸包表面的三角形数
    face F[maxn<<3];//凸包表面的三角形
    int g[maxn][maxn];

    Point3 cross(const Point3 &a,const Point3 &b,const Point3 &c)
    {
        return (b-a)^(c-a);
    }

    //三角形面积*2
    double area_triangle(Point3 a,Point3 b,Point3 c)
    {
        return ((b-a)^(c-a)).len();
    }

    //四面体有向体积*6
    double volume_four(Point3 a,Point3 b,Point3 c,Point3 d)
    {
        return ((b-a)^(c-a))*(d-a);
    }

    //正:点在面同向
    double dblcmp(Point3 &p,face &f)
    {
        Point3 p1=P[f.b]-P[f.a];
        Point3 p2=P[f.c]-P[f.a];
        Point3 p3=p-P[f.a];
        return (p1^p2)*p3;
    }

    void deal(int p,int a,int b)
    {
        int f=g[a][b];
        face add;
        if(F[f].ok)
        {
            if(sgn(dblcmp(P[p],F[f]))>0)
                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)
    {
        Point3 &a=P[F[s].a];
        Point3 &b=P[F[s].b];
        Point3 &c=P[F[s].c];
        int d1=sgn(volume_four(a,b,c,P[F[t].a]));
        int d2=sgn(volume_four(a,b,c,P[F[t].b]));
        int d3=sgn(volume_four(a,b,c,P[F[t].c]));
        return (d1==0)&&(d2==0)&&(d3==0);
    }

    //构建三维凸包
    void create(void)
    {
        num=0;
        face add;
        //保证前四个点不共面
        //*********************************
        bool flag=true;
        for(int i=1;i<n;i++)
        {
            if(!(P[0]==P[i]))
            {
                swap(P[1],P[i]);
                flag=false;
                break;
            }
        }
        if(flag) return ;
        flag=true;
        for(int i=2;i<n;i++)
        {
            if(sgn(((P[1]-P[0])^(P[i]-P[0])).len())>0)
            {
                swap(P[2],P[i]);
                flag=false;
                break;
            }
        }
        if(flag) return ;
        flag=true;
        for(int i=3;i<n;i++)
        {
            if(sgn(((P[1]-P[0])^(P[2]-P[0]))*(P[i]-P[0]))!=0)
            {
                swap(P[3],P[i]);
                flag=false;
                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.ok=true;
            if(sgn(dblcmp(P[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(int i=4;i<n;i++)
        {
            for(int j=0;j<num;j++)
            {
                if(F[j].ok&&sgn(dblcmp(P[i],F[j]))>0)
                {
                    dfs(i,j);
                    break;
                }
            }
        }

        int tmp=num;
        num=0;
        for(int i=0;i<tmp;i++)
            if(F[i].ok)
                F[num++]=F[i];
    }

    //表面积
    double area(void)
    {
        double res=0;
        if(n==3)
        {
            Point3 p=cross(P[0],P[1],P[2]);
            return p.len()/2;
        }
        for(int i=0;i<num;i++)
            res+=area_triangle(P[F[i].a],P[F[i].b],P[F[i].c]);
        return res/2.0;
    }

    //体积
    double volume(void)
    {
        double res=0;
        Point3 tmp=Point3(0,0,0);
        for(int i=0;i<num;i++)
            res+=volume_four(tmp,P[F[i].a],P[F[i].b],P[F[i].c]);
        return abs(res/6.0);
    }

    //表面三角形个数
    int sum_of_triangle(void)
    {
        return num;
    }

    //表面多边形个数
    int sum_of_polygon(void)
    {
        int res=0;
        for(int i=0;i<num;i++)
        {
            bool flag=true;
            for(int j=0;j<i;j++)
            {
                if(same(i,j))
                {
                    flag=false;
                    break;
                }
            }
            res+=flag;
        }
        return res;
    }

    //重心
    Point3 bary_center(void)
    {
        Point3 ans=Point3(0,0,0);
        Point3 o=Point3(0,0,0);
        double all=0;
        for(int i=0;i<num;i++)
        {
            double vol=volume_four(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_point_to_face(Point3 p,int i)
    {
        double tmp1=abs(volume_four(P[F[i].a],P[F[i].b],P[F[i].c],p));
        double tmp2=((P[F[i].b]-P[F[i].a])^(P[F[i].c]-P[F[i].a])).len();
        return tmp1/tmp2;
    }
};

CH3D hull;

int main(void)
{
    while(scanf("%d",&hull.n)!=EOF)
    {
        for(int i=0;i<hull.n;i++) hull.P[i].input();
        hull.create();
        printf("%.3f\n",hull.area());
    }
    return 0;
}

  • 0
    点赞
  • 0
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
经导师精心指导并认可、获 98 分的毕业设计项目!【项目资源】:微信小程序。【项目说明】:聚焦计算机相关专业毕设及实战操练,可作课程设计与期末大作业,含全部源码,能直用于毕设,经严格调试,运行有保障!【项目服务】:有任何使用上的问题,欢迎随时与博主沟通,博主会及时解答。 经导师精心指导并认可、获 98 分的毕业设计项目!【项目资源】:微信小程序。【项目说明】:聚焦计算机相关专业毕设及实战操练,可作课程设计与期末大作业,含全部源码,能直用于毕设,经严格调试,运行有保障!【项目服务】:有任何使用上的问题,欢迎随时与博主沟通,博主会及时解答。 经导师精心指导并认可、获 98 分的毕业设计项目!【项目资源】:微信小程序。【项目说明】:聚焦计算机相关专业毕设及实战操练,可作课程设计与期末大作业,含全部源码,能直用于毕设,经严格调试,运行有保障!【项目服务】:有任何使用上的问题,欢迎随时与博主沟通,博主会及时解答。 经导师精心指导并认可、获 98 分的毕业设计项目!【项目资源】:微信小程序。【项目说明】:聚焦计算机相关专业毕设及实战操练,可作课程设计与期末大作业,含全部源码,能直用于毕设,经严格调试,运行有保障!【项目服务】:有任何使用上的问题,欢迎随时与博主沟通,博主会及时解答。
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

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

抵扣说明:

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

余额充值