计算几何-给定四面体四个点坐标求其体积

前阵子笔试题目遇到一个很有意思的题目，给定四面体四个点坐标求其体积，一开始我还以为是要使用蛮力法，将一个点投影到其余三个点组成的面上，然后求高与底面积，吭吭哧哧写了半天发现比较难调试。

过后我找了一下网上其他人的做法，实际上有更好的解决办法，那就是把这个四面体给补全成一个平行六面体，然后这个四面体的体积就是这个六面体体积的六分之一。

原理

假设四面体为$S_{ABCD}$，四个顶点的三维坐标$ A(x_1,y_1,z_1)$、$B(x_2,y_2,z_2)$、$C(x_3,y_3,z_3)$ 和 $D(x_4,y_4,z_4)$

然后我们可以把$AB$，$AC$和$AD$三条线段组成的一个角看成是平行六面体的一个角，即

如果我们把底面积给补上去，我们很快就会发现这个四面体的体积就是平行四面体体积的六分之一（锥体体积是柱体体积的三分之一、三角形面积是平行四边形面积的二分之一）：

而这个平行六面体可以用下面的公式计算： $$ V_6=\left| \det \begin{bmatrix} x_2-x_1 & y_2-y_1 & z_2-z_1\ x_3-x_1 & y_3-y_1 & z_3-z_1\ x_4-x_1 & y_4-y_1 & z_4-z_1 \end{bmatrix} \right| $$ 矩阵的行向量分别是从顶点 $A$ 到其他三个顶点 $B$、$C$、$D$ 的向量。值得注意的是，如果四点共面，则矩阵的行列式为零，对应体积为零。

编程实现

#include <iostream>
#include<cmath>
using namespace std;
struct Point3D {
    double x, y, z;
};
// 函数用于计算行列式的值
double calculateDeterminant(const Point3D& p1, const Point3D& p2,
    const Point3D& p3, const Point3D& p4) {
    // 计算矩阵中的元素
    double a = p2.x - p1.x;
    double b = p2.y - p1.y;
    double c = p2.z - p1.z;
    double d = p3.x - p1.x;
    double e = p3.y - p1.y;
    double f = p3.z - p1.z;
    double g = p4.x - p1.x;
    double h = p4.y - p1.y;
    double i = p4.z - p1.z;

    // 使用行列式公式计算结果
    double determinant = a * (e * i - f * h) -
        b * (d * i - f * g) +
        c * (d * h - e * g);

    return determinant;
}
double tetrahedronVolume(const Point3D& p1, const Point3D& p2,
    const Point3D& p3, const Point3D& p4) {
    // 计算行列式的值
    double determinant = calculateDeterminant(p1, p2, p3, p4);

    return abs(determinant) / 6.0;
}

int main() {
    // 输入四个点的坐标
    Point3D p1, p2, p3, p4;
    cout << "Please input A(x1, y1, z1), B(x2, y2, z2), C(x3, y3, z3), D(x4, y4, z4):" << endl;
    cin >> p1.x >> p1.y >> p1.z
        >> p2.x >> p2.y >> p2.z
        >> p3.x >> p3.y >> p3.z
        >> p4.x >> p4.y >> p4.z;

    cout << "volume: " << tetrahedronVolume(p1, p2, p3, p4);
    
    return 0;
}