015-3D数学矩阵的行列式、逆矩阵、正交矩阵、齐次矩阵详解与实例

最新推荐文章于 2025-06-10 21:26:52 发布

原创

最新推荐文章于 2025-06-10 21:26:52 发布

· 480 阅读

5 ·

版权

文章标签：

#3d #矩阵 #线性代数

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