我需要根据四个角点的坐标计算四面体的体积。
假设您有 4 个顶点 a、b、c、d(3-D 向量)。
现在,问题归结为编写解决向量的叉积和点积的代码。如果您来自 python,您可以使用 NumPy,否则您可以自己编写代码。
维基百科链接肯定会对您有所帮助。 链接
计算此体积的一种方法是:
1 [ax bx cx dx]
V = --- det [ay by cy dy]
6 [az bz cz dz]
[ 1 1 1 1]
这涉及到 4×4 行列式的评估。它可以很好地推广到更高维度的单纯形,其中 6 是 n!(维度的阶乘)的特例。生成的体积将是定向,即根据点的顺序可能为负。如果您不想这样,请取结果的绝对值。
如果您手头有数学库,则上述公式可能是最容易写下的公式之一,并且软件可以从那里获取它。如果没有,您可以首先通过从 a 到 c 减去 d 坐标来简化事情。这不会改变音量,但会将最右边的列变成
(0, 0, 0, 1)
。因此,您可以简单地将矩阵的值计算为左上角 3×3 子矩阵的行列式。并使用方程
det(a, b, c) = a · (b × c)
你最终得到了Surya的答案中的公式。
如果您没有点的坐标,而只有点之间的距离,请查看塔塔利亚公式,它本质上是上述公式的平方版本,尽管它并不像乍一看那么简单。
Ivan Seidel 的示例,Python 语言(答案是 1.3333...)
def determinant_3x3(m):
return (m[0][0] * (m[1][1] * m[2][2] - m[1][2] * m[2][1]) -
m[1][0] * (m[0][1] * m[2][2] - m[0][2] * m[2][1]) +
m[2][0] * (m[0][1] * m[1][2] - m[0][2] * m[1][1]))
def subtract(a, b):
return (a[0] - b[0],
a[1] - b[1],
a[2] - b[2])
def tetrahedron_calc_volume(a, b, c, d):
return (abs(determinant_3x3((subtract(a, b),
subtract(b, c),
subtract(c, d),
))) / 6.0)
a = [0.0, 0.0, 0.0]
d = [2.0, 0.0, 0.0]
c = [0.0, 2.0, 0.0]
b = [0.0, 0.0, 2.0]
print(tetrahedron_calc_volume(a, b, c, d))
以下是 PHP 代码,用于计算给定 4 个点的任何四面体的体积:
class Math{
public static function determinant(array $vals){
$s = sizeof($vals);
$sum = 0.0;
for( $i=0; $i < $s ; $i++ ){
$mult = 1.0;
for($j=0; $j < $s ; $j++ ){
$mult *= $vals[$j][ ($i+$j)%$s ];
}
$sum += $mult;
}
for( $i=0; $i < $s ; $i++ ){
$mult = 1;
for($j=0; $j < $s ; $j++ ){
$mult *= $vals[$j][ ($i-$j < 0? $s - ($j-$i) :($i-$j)) ];
}
$sum -= $mult;
}
return $sum;
}
public static function subtract(array $a, array $b){
for($i = 0; $i < sizeof($a); $i++)
$a[$i] -= $b[$i];
return $a;
}
}
// TEST CASE
$a = array(0,0,0);
$d = array(2,0,0);
$c = array(0,2,0);
$b = array(0,0,2);
echo abs(Math::determinant(array(
Math::subtract($a, $b),
Math::subtract($b, $c),
Math::subtract($c, $d),
)))/6;
对于
python
,numpy
,你可以这样编写代码。
import numpy as np
def tetrahedron_volume(points):
p1, p2, p3, p4 = list(map(np.array, points))
d0, d1, d2 = p4 - p1, p3 - p1, p2 - p1
normal = np.cross(d1, d2)
return 1/6 * np.abs(np.dot(d0, normal))