1. 程式人生 > 其它 >增量法(三維凸包)

增量法(三維凸包)

題意:告訴你n個三維的點,讓你求出這些點構成的最小表面積

#include<iostream>
#include<cmath>
#include<algorithm>
using namespace std;
const int maxn = 210;
const double eps = 1e-12;    //精度要求高一點

int n,m;    //n是總點數,m是凸包中平面的數量
bool g[maxn][maxn];    //g用來判斷一條邊被照到幾次

double rand_eps()    //用rand函式來生成一個非常小的隨機數
{
    return
((double)rand() / RAND_MAX - 0.5) * eps; //用rand生成一個-0.5到0.5之間的數,再乘eps,就得到了一個非常小的隨機數 } struct Point{ //定義點的結構體 double x,y,z; //xyz三個座標 void shake() //微小擾動,給每個座標都加一個極小的隨機數 { x += rand_eps(), y += rand_eps(), z += rand_eps(); } Point operator-(Point t) //過載一下減號運算子 {
return {x-t.x, y-t.y, z-t.z}; } Point operator*(Point t) //叉乘 { return {y*t.z - t.y*z, t.x*z - x*t.z, x*t.y - y*t.x}; } double operator&(Point t) //點積 { return x*t.x + y*t.y + z*t.z; } double len() //求向量的模長,三個座標也能存向量 { return sqrt(x*x + y*y + z*z); } }p[maxn];
//p來存所有點 struct Plane{ //定義平面的結構體 int v[3]; //三個頂點 Point norm() //求法向量 { return (p[v[1]] - p[v[0]]) * (p[v[2]] - p[v[0]]); //平面中兩向量的叉積 } bool above(Point t ) //判斷一個點是否在平面上方 { return ((t-p[v[0]]) & norm()) >= 0; //用向量和法向量的點積判斷 } double area() //求一平面的面積 { return norm().len() / 2; //法向量的模長除以2 } }plane[maxn],tp[maxn]; //plane存凸包上的平面,tp用來更新凸包,每次凸包上要留的平面和新加的平面都存進tp,要刪的平面不存,最後將tp在複製給plane,就實現了凸包的更新 void convex() { plane[m++] = {0,1,2}; //初始化凸包,隨便三個點存入,確定最開始的一個平面,這裡取得是前三個點 plane[m++] = {2,1,0}; //因為不知道第一個平面怎麼樣是逆時針,所以都存一遍,順時針存的一會會被刪掉 for(int i = 3;i < n;i++) //從第四個點開始迴圈每個點 { int cnt = 0; for(int j = 0;j < m;j++) //迴圈每個平面 { bool fg = plane[j].above(p[i]); //判斷這個點是否在該平面上方 if(!fg) //如果是下方的話,說明照不到 tp[cnt++] = plane[j]; //存進tp陣列 for(int k = 0;k < 3;k++) //迴圈該平面的三條邊 g[plane[j].v[k]][plane[j].v[(k+1)%3]] = fg; //ab邊照不照得到情況賦值給g[a][b] } for(int j = 0;j < m;j++) //然後就迴圈每個平面的每條邊 { for(int k = 0;k < 3;k++) { int a = plane[j].v[k],b = plane[j].v[(k+1)%3]; if(g[a][b] && !g[b][a]) //判斷該邊是否被照到了一次,即是否是交界線的邊 tp[cnt++] = {a,b,i}; //若是,加新平面abi,ab一定是逆時針的,i在後面 } } m = cnt; //將tp再賦值給plane for(int j = 0;j < m;j++) plane[j] = tp[j]; } } int main() { cin >> n; for(int i = 0;i < n;i++) { cin >> p[i].x >> p[i].y >> p[i].z; //輸入n個點 p[i].shake(); //微小擾動 } convex(); //求三維凸包 double ans = 0; //求面積和 for(int i = 0;i < m;i++) //迴圈最終凸包上的m個平面 ans += plane[i].area(); //將平面的面積加和 printf("%.6lf\n",ans); //輸出答案 return 0; }