poj 3862 Asteroids三维凸包➕重心
来源：互联网 发布：网络安全设备管理制度 编辑：程序博客网 时间：2024/05/21 15:00
POJ 3862 Asteroids （三维凸包，求两个凸包重心到表面的最短距离）
下边是队友写的，另附上kuangbin大神魔板
#include<stdio.h>#include<algorithm>#include<string.h>#include<math.h>#include<stdlib.h>#define eps 1e-8#define N 110using 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,k;int main(){ freopen("asteroids.in","r",stdin); freopen("asteroids.out","w",stdout); 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)); } scanf("%d",&k.n); for(int i=0;i<k.n;i++) { scanf("%lf%lf%lf",&k.tn[i].x,&k.tn[i].y,&k.tn[i].z); } k.create(); point pp = k.barycenter(); double answer2=1e20; for(int i=0;i<k.num;i++) { answer2=min(answer2,k.ptoface(pp,i)); } printf("%.5lf\n",answer+answer2); } return 0;}/*HDU 4273 Rescue给一个三维凸包，求重心到表面的最短距离模板题：三维凸包+多边形重心+点面距离*/#include<stdio.h>#include<algorithm>#include<string.h>#include<math.h>#include<stdlib.h>using namespace std;const int MAXN=550;const double eps=1e-8;struct Point{ double x,y,z; Point(){} Point(double xx,double yy,double zz):x(xx),y(yy),z(zz){} //两向量之差 Point operator -(const Point p1) { return Point(x-p1.x,y-p1.y,z-p1.z); } //两向量之和 Point operator +(const Point p1) { return Point(x+p1.x,y+p1.y,z+p1.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 d) { return Point(x*d,y*d,z*d); } Point operator / (double d) { return Point(x/d,y/d,z/d); } //点乘 double operator ^(Point p) { return (x*p.x+y*p.y+z*p.z); }};struct CH3D{ struct face { //表示凸包一个面上的三个点的编号 int a,b,c; //表示该面是否属于最终凸包上的面 bool ok; }; //初始顶点数 int n; //初始顶点 Point P[MAXN]; //凸包表面的三角形数 int num; //凸包表面的三角形 face F[8*MAXN]; //凸包表面的三角形 int g[MAXN][MAXN]; //向量长度 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) ); } //三角形面积*2 double area(Point a,Point b,Point c) { return vlen((b-a)*(c-a)); } //四面体有向体积*6 double volume(Point a,Point b,Point c,Point d) { return (b-a)*(c-a)^(d-a); } //正：点在面同向 double dblcmp(Point &p,face &f) { Point m=P[f.b]-P[f.a]; Point n=P[f.c]-P[f.a]; Point t=p-P[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(dblcmp(P[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=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); } bool same(int s,int t) { Point &a=P[F[s].a]; Point &b=P[F[s].b]; Point &c=P[F[s].c]; return fabs(volume(a,b,c,P[F[t].a]))<eps && fabs(volume(a,b,c,P[F[t].b]))<eps && fabs(volume(a,b,c,P[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(P[0]-P[i])>eps) { swap(P[1],P[i]); flag=false; break; } } if(flag)return; flag=true; //使前三个点不共线 for(i=2;i<n;i++) { if(vlen((P[0]-P[1])*(P[1]-P[i]))>eps) { swap(P[2],P[i]); flag=false; break; } } if(flag)return; flag=true; //使前四个点不共面 for(int i=3;i<n;i++) { if(fabs((P[0]-P[1])*(P[1]-P[2])^(P[0]-P[i]))>eps) { swap(P[3],P[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(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(i=4;i<n;i++) { for(j=0;j<num;j++) { if(F[j].ok&&dblcmp(P[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(P[0],P[1],P[2]); res=vlen(p)/2.0; return res; } for(int i=0;i<num;i++) res+=area(P[F[i].a],P[F[i].b],P[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,P[F[i].a],P[F[i].b],P[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,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 ptoface(Point p,int i) { return fabs(volume(P[F[i].a],P[F[i].b],P[F[i].c],p)/vlen((P[F[i].b]-P[F[i].a])*(P[F[i].c]-P[F[i].a]))); }};CH3D hull;int main(){ // freopen("in.txt","r",stdin); // freopen("out.txt","w",stdout); while(scanf("%d",&hull.n)==1) { 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(); Point p=hull.barycenter(); double ans1=1e20; for(int i=0;i<hull.num;i++) { ans1=min(ans1,hull.ptoface(p,i)); } scanf("%d",&hull.n); 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(); p=hull.barycenter(); double ans2=1e20; for(int i=0;i<hull.num;i++) { ans2=min(ans2,hull.ptoface(p,i)); } printf("%.5f\n",ans1+ans2); } return 0;} 0 0
- poj 3862 Asteroids三维凸包➕重心
- Asteroids (三维凸包+重心)
- POJ 3862 Asteroids（两个三维凸包的重心到表面最短距离和）
- Gym 101308A Asteroids 三维凸包 求重心
- UVA 1438 - Asteroids（三维凸包）
- LA 4589 - Asteroids 三维凸包
- hdu4273Rescue三维凸包的重心
- uva 1438 - Asteroids（几何重心+凸包）
- hdu 4273 Rescue（三维凸包重心）
- HDU 4273 Rescue（三维凸包+重心模板）
- HDU 4273 Rescue（三维凸包 + 重心）
- HDU 4273 Rescue（三维凸包 + 重心）模版
- 求三维凸包重心到表面最短距离
- 三维集合之三维凸包，和求密度均匀的凸多面体的重心
- POJ 3528 & POJ 2974 三维凸包
- hud4273--poj3862+三维凸包模板题+凸包重心到表面的最小距离
- 三维计算几何模板--表面三角形个数 表面多边形个数 三维凸包 表面积 凸包重心 点到面的距离
- HDU4273(求三维凸包重心到表面的最短距离)
- NYOJ—水池数目（DFS）
- 查找一条数据的临近数据（上一条、下一条）
- 小工具URL编码
- HCE基础知识普及
- C# WPF 小工具HTML转JS
- poj 3862 Asteroids三维凸包➕重心
- jeesite <meta name="decorator" content="default"/>
- 理解javascript中的MVVM开发模式
- Jquery File Upload 入门
- 【DRP】——JSTL之格式化标签与函数
- 蓝桥杯——带分数
- Youme语音格式及流量
- Java定时任务
- 补码疑云
