题解:P2287 [HNOI2004] 最佳包裹

· · 题解

题意简述:求 n 个点的三维凸包的表面积。(n \le 100

先对所有点的位置进行微小扰动,以防止出现四点共面的情况。

我们可以枚举这些点能够组成的所有三角形,如果剩余所有点全部在该三角形所在平面的同一侧,那这个三角形就在这个三维凸包上,直接让答案加上它的面积即可。

怎样求若干个点是否在三角形所在平面的同一侧呢?设三角形为 \triangle ABC,求出 \overrightarrow{AB} \times \overrightarrow{AC}(叉乘),也就是该平面的法向量,记为向量 \boldsymbol{v} 。枚举剩余的点,设当前点为 D,求出 \overrightarrow{AD} \cdot \boldsymbol{v}(点乘)。如果存在两个点求出来的值一正一负,就说明这两个点在平面的两侧。三角形的面积即为 \frac{|\boldsymbol{v}|}{2}

算法的时间复杂度为 O(n^4)

代码示例:

#include <iostream>
#include <random>
#include <ctime>
#include <limits> 
#include <cmath>

using namespace std;

constexpr int maxn = 100;
constexpr double eps = 1e-10;
mt19937 rd(time(nullptr));

struct point
{
    double x, y, z;
} a[maxn + 5];

struct vec
{
    double x, y, z;
    vec(double x, double y, double z) : x{x}, y{y}, z{z} {}
};

inline double get_delta()    // 微小扰动的值
{ return (1. * rd() / numeric_limits<unsigned>::max() - 0.5) * eps; }

inline vec operator-(const point& a, const point& b)      // 求出向量
{ return {b.x - a.x, b.y - a.y, b.z - a.z}; }

inline double operator*(const vec& a, const vec& b)      // 点乘
{ return a.x * b.x + a.y * b.y + a.z * b.z; }

inline vec operator^(const vec& a, const vec& b)     // 叉乘
{ return {a.y * b.z - a.z * b.y, a.z * b.x - a.x * b.z, a.x * b.y - a.y * b.x}; }

inline double get_mod(const vec& a)     // 求向量的模长
{ return sqrt(a.x * a.x + a.y * a.y + a.z * a.z); }

int main()
{
    ios_base::sync_with_stdio(false);
    cin.tie(nullptr);
    (cout << fixed).precision(6);

    int n; cin >> n;
    for (int i = 1; i <= n; i++)
        cin >> a[i].x >> a[i].y >> a[i].z;
    for (int i = 1; i <= n; i++)    // 微小扰动
    {
        a[i].x += get_delta();
        a[i].y += get_delta();
        a[i].z += get_delta();
    }

    double ans = 0;
    for (int i = 1; i <= n; i++)
        for (int j = i + 1; j <= n; j++)
            for (int k = j + 1; k <= n; k++)
            {
                vec cur = (a[j] - a[i]) ^ (a[k] - a[i]);
                bool flag1 = false, flag2 = false;
                for (int l = 1; l <= n; l++)
                    if (i != l && j != l && k != l)
                    {
                        if (cur * (a[l] - a[i]) > 0)
                            flag1 = true;
                        else
                            flag2 = true;
                    }
                if (!flag1 || !flag2)
                    ans += get_mod(cur) / 2;   // 加上面积
            }
    cout << ans << endl;

    return 0;
}