3D数学矩阵的行列式、逆矩阵、正交矩阵、齐次矩阵详解与实例
1. 行列式
行列式是与方阵相关联的一个标量值，它在3D图形学中有着广泛的应用，包括判断变换的特性、计算逆矩阵、计算体积缩放等。
1.1 行列式的定义与性质
行列式是一个将方阵映射到标量的函数，通常表示为 det(A) 或 |A|。
cpp
// 计算2x2矩阵的行列式
float determinant2x2(const Matrix2x2& matrix) {
return matrix(0, 0) * matrix(1, 1) - matrix(0, 1) * matrix(1, 0);
}
// 计算3x3矩阵的行列式
float determinant3x3(const Matrix3x3& matrix) {
return matrix(0, 0) * (matrix(1, 1) * matrix(2, 2) - matrix(1, 2) * matrix(2, 1)) -
matrix(0, 1) * (matrix(1, 0) * matrix(2, 2) - matrix(1, 2) * matrix(2, 0)) +
matrix(0, 2) * (matrix(1, 0) * matrix(2, 1) - matrix(1, 1) * matrix(2, 0));
}
// 计算4x4矩阵的行列式（使用余子式展开）
float determinant4x4(const Matrix4x4& matrix) {
// 对第一行展开
float det = 0.0f;
for (int j = 0; j < 4; ++j) {
// 提取余子式矩阵
Matrix3x3 minor;
for (int row = 1; row < 4; ++row) {
for (int col = 0; col < 4; ++col) {
if (col < j) {
minor(row - 1, col) = matrix(row, col);
} else if (col > j) {
minor(row - 1, col - 1) = matrix(row, col);
}
}
}
// 交替加减
float sign = (j % 2 == 0) ? 1.0f : -1.0f;
det += sign * matrix(0, j) * determinant3x3(minor);
}
return det;
}
1.2 行列式的几何意义
在3D空间中，行列式有着重要的几何意义：
cpp
// 行列式作为体积缩放因子
float calculateVolumeScaleFactor(const Matrix3x3& transform) {
float det = determinant3x3(transform);
return fabs(det); // 取绝对值，因为我们关心的是体积变化
}
// 计算变换后的平行六面体体积
float calculateParallelepipedVolume(const Vector3& v1, const Vector3& v2, const Vector3& v3) {
// 构建一个3x3矩阵，其列为给定向量
Matrix3x3 matrix;
matrix(0, 0) = v1.x; matrix(1, 0) = v1.y; matrix(2, 0) = v1.z;
matrix(0, 1) = v2.x; matrix(1, 1) = v2.y; matrix(2, 1) = v2.z;
matrix(0, 2) = v3.x; matrix(1, 2) = v3.y; matrix(2, 2) = v3.z;
// 体积是行列式的绝对值
return fabs(determinant3x3(matrix));
}
// 判断三点是否共线（面积为零）
bool arePointsCollinear(const Vector3& p1, const Vector3& p2, const Vector3& p3) {
// 计算两个边向量
Vector3 v1 = p2 - p1;
Vector3 v2 = p3 - p1;
// 计算叉积的模长（即平行四边形的面积）
Vector3 cross = Vector3::cross(v1, v2);
float area = cross.length();
// 如果面积接近零，则三点共线
return area < 1e-6f;
}
// 判断四点是否共面（体积为零）
bool arePointsCoplanar(const Vector3& p1, const Vector3& p2, const Vector3& p3, const Vector3& p4) {
// 计算三个边向量
Vector3 v1 = p2 - p1;
Vector3 v2 = p3 - p1;
Vector3 v3 = p4 - p1;
// 创建包含这三个向量的矩阵
Matrix3x3 matrix;
matrix(0, 0) = v1.x; matrix(1, 0) = v1.y; matrix(2, 0) = v1.z;
matrix(0, 1) = v2.x; matrix(1, 1) = v2.y; matrix(2, 1) = v2.z;
matrix(0, 2) = v3.x; matrix(1, 2) = v3.y; matrix(2, 2) = v3.z;
// 计算行列式（即平行六面体的体积）
float volume = fabs(determinant3x3(matrix));
// 如果体积接近零，则四点共面
return volume < 1e-6f;
}
1.3 行列式的计算优化
对于特殊类型的矩阵，行列式计算可以简化：
cpp
// 计算对角矩阵的行列式
float determinantDiagonal(const Matrix4x4& matrix) {
// 对角矩阵的行列式是对角元素的乘积
return matrix(0, 0) * matrix(1, 1) * matrix(2, 2) * matrix(3, 3);
}
// 计算上三角/下三角矩阵的行列式
float determinantTriangular(const Matrix4x4& matrix, bool isUpper) {
// 三角矩阵的行列式也是对角元素的乘积
return matrix(0, 0) * matrix(1, 1) * matrix(2, 2) * matrix(3, 3);
}
// 使用LU分解计算行列式
float determinantLU(const Matrix4x4& matrix) {
// 使用LU分解可以高效计算大矩阵的行列式
// 这里提供一个简化版本
// 复制原矩阵避免修改它
Matrix4x4 LU = matrix;
int n = 4; // 矩阵大小
float det = 1.0f;
for (int i = 0; i < n; ++i) {
// 如果对角元素接近零，行列式接近零
if (fabs(LU(i, i)) < 1e-10f) {
return 0.0f;
}
det *= LU(i, i);
for (int j = i + 1; j < n; ++j) {
float factor = LU(j, i) / LU(i, i);
for (int k = i + 1; k < n; ++k) {
LU(j, k) -= factor * LU(i, k);
}
}
}
return det;
}
2. 逆矩阵
逆矩阵是矩阵代数中的一个基本概念，它在变换的反向操作、解线性方程组和许多图形学算法中都有应用。
2.1 逆矩阵的定义与性质
如果一个n×n方阵A的逆矩阵A^(-1)满足A * A^(-1) = A^(-1) * A = I（其中I是单位矩阵），则称A为可逆矩阵。
cpp
// 检查矩阵是否可逆
bool isInvertible(const Matrix4x4& matrix) {
// 一个矩阵可逆当且仅当其行列式不为零
return fabs(determinant4x4(matrix)) > 1e-10f;
}
// 检验逆矩阵的正确性
bool validateInverse(const Matrix4x4& matrix, const Matrix4x4& inverse) {
// 计算A * A^(-1)和A^(-1) * A
Matrix4x4 product1 = matrix * inverse;
Matrix4x4 product2 = inverse * matrix;
// 检查结果是否接近单位矩阵
Matrix4x4 identity = Matrix4x4::identity();
bool isValid = true;
for (int i = 0; i < 4; ++i) {
for (int j = 0; j < 4; ++j) {
isValid = isValid && (fabs(product1(i, j) - identity(i, j)) < 1e-5f);
isValid = isValid && (fabs(product2(i, j) - identity(i, j)) < 1e-5f);
}
}
return isValid;
}
2.2 计算逆矩阵的方法
有多种方法可以计算逆矩阵，每种方法有其适用场景：
2.2.1 使用伴随矩阵计算逆矩阵
cpp
// 计算2x2矩阵的逆
Matrix2x2 inverse2x2(const Matrix2x2& matrix) {
float det = determinant2x2(matrix);
if (fabs(det) < 1e-10f) {
throw std::runtime_error("Matrix is not invertible");
}
float invDet = 1.0f / det;
Matrix2x2 result;
result(0, 0) = matrix(1, 1) * invDet;
result(0, 1) = -matrix(0, 1) * invDet;
result(1, 0) = -matrix(1, 0) * invDet;
result(1, 1) = matrix(0, 0) * invDet;
return result;
}
// 计算3x3矩阵的伴随矩阵（代数余子式矩阵的转置）
Matrix3x3 adjugate3x3(const Matrix3x3& matrix) {
Matrix3x3 result;
// 计算余子式
result(0, 0) = matrix(1, 1) * matrix(2, 2) - matrix(1, 2) * matrix(2, 1);
result(0, 1) = -(matrix(0, 1) * matrix(2, 2) - matrix(0, 2) * matrix(2, 1));
result(0, 2) = matrix(0, 1) * matrix(1, 2) - matrix(0, 2) * matrix(1, 1);
result(1, 0) = -(matrix(1, 0) * matrix(2, 2) - matrix(1, 2) * matrix(2, 0));
result(1, 1) = matrix(0, 0) * matrix(2, 2) - matrix(0, 2) * matrix(2, 0);
result(1, 2) = -(matrix(0, 0) * matrix(1, 2) - matrix(0, 2) * matrix(1, 0));
result(2, 0) = matrix(1, 0) * matrix(2, 1) - matrix(1, 1) * matrix(2, 0);
result(2, 1) = -(matrix(0, 0) * matrix(2, 1) - matrix(0, 1) * matrix(2, 0));
result(2, 2) = matrix(0, 0) * matrix(1, 1) - matrix(0, 1) * matrix(1, 0);
return result;
}
// 计算3x3矩阵的逆
Matrix3x3 inverse3x3(const Matrix3x3& matrix) {
float det = determinant3x3(matrix);
if (fabs(det) < 1e-10f) {
throw std::runtime_error("Matrix is not invertible");
}
Matrix3x3 adj = adjugate3x3(matrix);
float invDet = 1.0f / det;
Matrix3x3 result;
for (int i = 0; i < 3; ++i) {
for (int j = 0; j < 3; ++j) {
result(i, j) = adj(i, j) * invDet;
}
}
return result;
}
2.2.2 高斯-约旦消元法（Gauss-Jordan Elimination）
cpp
// 使用高斯-约旦消元法计算4x4矩阵的逆
Matrix4x4 inverse4x4GaussJordan(const Matrix4x4& matrix) {
// 创建增广矩阵[A|I]
float augmented[4][8];
// 初始化增广矩阵
for (int i = 0; i < 4; ++i) {
for (int j = 0; j < 4; ++j) {
augmented[i][j] = matrix(i, j);
augmented[i][j + 4] = (i == j) ? 1.0f : 0.0f;
}
}
// 高斯-约旦消元
for (int i = 0; i < 4; ++i) {
// 寻找主元
int maxRow = i;
float maxVal = fabs(augmented[i][i]);
for (int k = i + 1; k < 4; ++k) {
if (fabs(augmented[k][i]) > maxVal) {
maxRow = k;
maxVal = fabs(augmented[k][i]);
}
}
// 如果主元接近零，矩阵不可逆
if (maxVal < 1e-10f) {
throw std::runtime_error("Matrix is not invertible");
}
// 交换行
if (maxRow != i) {
for (int j = 0; j < 8; ++j) {
std::swap(augmented[i][j], augmented[maxRow][j]);
}
}
// 缩放行使主元为1
float scale = 1.0f / augmented[i][i];
for (int j = 0; j < 8; ++j) {
augmented[i][j] *= scale;
}
// 消元，使其他行在当前列上为0
for (int k = 0; k < 4; ++k) {
if (k != i) {
float factor = augmented[k][i];
for (int j = 0; j < 8; ++j) {
augmented[k][j] -= factor * augmented[i][j];
}
}
}
}
// 提取结果（右半部分为逆矩阵）
Matrix4x4 result;
for (int i = 0; i < 4; ++i) {
for (int j = 0; j < 4; ++j) {
result(i, j) = augmented[i][j + 4];
}
}
return result;
}
2.2.3 特殊矩阵的逆
对于特殊类型的矩阵，可以更高效地计算逆矩阵：
cpp
// 计算对角矩阵的逆
Matrix4x4 inverseDiagonal(const Matrix4x4& matrix) {
Matrix4x4 result = Matrix4x4::zero();
for (int i = 0; i < 4; ++i) {
if (fabs(matrix(i, i)) < 1e-10f) {
throw std::runtime_error("Diagonal matrix is not invertible");
}
result(i, i) = 1.0f / matrix(i, i);
}
return result;
}
// 计算正交矩阵的逆（等于其转置）
Matrix4x4 inverseOrthogona
最低0.47元/天 解锁文章
扫一扫
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
