mxn矩阵学习笔记

mxn矩阵学习笔记

1. 简介

nxm矩阵是线性代数中的基本概念，用于表示n行m列的矩阵。本教程将介绍如何使用C++实现一个通用的nxm矩阵类，支持任意维度的矩阵运算。

2. 代码实现

2.1 头文件 (matrixnxm.h)

#ifndef MATRIXNXM_H
#define MATRIXNXM_H

#include <vector>
#include <stdexcept>
#include <iostream>
#include <cmath>

namespace math {
namespace linear_algebra {

/**
 * @brief nxm矩阵类
 * 
 * 这个类实现了nxm矩阵的基本运算，包括：
 * - 矩阵加减
 * - 矩阵乘法
 * - 标量乘法
 * - 转置
 * - 行列式（仅方阵）
 * - 逆矩阵（仅方阵）
 * - 特征值和特征向量
 * - 矩阵性质检查
 */
class MatrixNxM {
public:
    // 构造函数
    MatrixNxM();  // 默认构造函数，创建0x0矩阵
    MatrixNxM(size_t rows, size_t cols);  // 创建指定维度的零矩阵
    MatrixNxM(const std::vector<std::vector<double>>& data);  // 从二维数组初始化

    // 基本属性
    size_t rows() const;  // 获取行数
    size_t cols() const;  // 获取列数
    bool isSquare() const;  // 检查是否为方阵

    // 元素访问
    double& operator()(size_t i, size_t j);  // 访问元素
    const double& operator()(size_t i, size_t j) const;  // 常量访问元素

    // 矩阵运算
    MatrixNxM operator+(const MatrixNxM& other) const;  // 矩阵加法
    MatrixNxM operator-(const MatrixNxM& other) const;  // 矩阵减法
    MatrixNxM operator*(const MatrixNxM& other) const;  // 矩阵乘法
    MatrixNxM operator*(double scalar) const;           // 标量乘法
    MatrixNxM operator/(double scalar) const;           // 标量除法

    // 矩阵变换
    MatrixNxM transpose() const;  // 转置矩阵
    double determinant() const;   // 计算行列式（仅方阵）
    MatrixNxM inverse() const;    // 计算逆矩阵（仅方阵）

    // 矩阵性质
    bool isSingular() const;      // 检查是否奇异
    bool isRegular() const;       // 检查是否非奇异
    bool isSymmetric() const;     // 检查是否对称
    bool isOrthogonal() const;    // 检查是否正交

    // 特征值和特征向量
    std::vector<double> eigenvalues() const;  // 计算特征值
    std::vector<MatrixNxM> eigenvectors() const;  // 计算特征向量

    // 特殊矩阵
    static MatrixNxM identity(size_t n);  // 创建n阶单位矩阵
    static MatrixNxM zero(size_t rows, size_t cols);  // 创建零矩阵
    static MatrixNxM rotation2D(double theta);  // 创建2D旋转矩阵
    static MatrixNxM rotation3D(double theta, char axis);  // 创建3D旋转矩阵
    static MatrixNxM scaling(const std::vector<double>& factors);  // 创建缩放矩阵

    // 输出运算符
    friend std::ostream& operator<<(std::ostream& os, const MatrixNxM& m);

private:
    std::vector<std::vector<double>> data;  // 矩阵数据
    size_t n;  // 行数
    size_t m;  // 列数

    // 辅助函数
    void checkDimensions(const MatrixNxM& other, const char* operation) const;
    void checkSquare(const char* operation) const;
    double minor(size_t row, size_t col) const;  // 计算余子式
    double cofactor(size_t row, size_t col) const;  // 计算代数余子式
};

} // namespace linear_algebra
} // namespace math

#endif // MATRIXNXM_H

2.2 实现文件 (matrixnxm.cpp)

#include "matrixnxm.h"
#include <algorithm>
#include <numeric>

namespace math {
namespace linear_algebra {

// 默认构造函数
MatrixNxM::MatrixNxM() : n(0), m(0) {}

// 创建指定维度的零矩阵
MatrixNxM::MatrixNxM(size_t rows, size_t cols) : n(rows), m(cols) {
    data.resize(n, std::vector<double>(m, 0.0));
}

// 从二维数组初始化
MatrixNxM::MatrixNxM(const std::vector<std::vector<double>>& data)
    : data(data), n(data.size()), m(data.empty() ? 0 : data[0].size()) {
    // 检查数据有效性
    for (const auto& row : data) {
        if (row.size() != m) {
            throw std::invalid_argument("Invalid matrix dimensions");
        }
    }
}

// 获取行数
size_t MatrixNxM::rows() const {
    return n;
}

// 获取列数
size_t MatrixNxM::cols() const {
    return m;
}

// 检查是否为方阵
bool MatrixNxM::isSquare() const {
    return n == m;
}

// 访问元素
double& MatrixNxM::operator()(size_t i, size_t j) {
    if (i >= n || j >= m) {
        throw std::out_of_range("Matrix index out of range");
    }
    return data[i][j];
}

// 常量访问元素
const double& MatrixNxM::operator()(size_t i, size_t j) const {
    if (i >= n || j >= m) {
        throw std::out_of_range("Matrix index out of range");
    }
    return data[i][j];
}

// 检查维度
void MatrixNxM::checkDimensions(const MatrixNxM& other, const char* operation) const {
    if (n != other.n || m != other.m) {
        throw std::invalid_argument(
            std::string("Matrix dimensions do not match for ") + operation);
    }
}

// 检查是否为方阵
void MatrixNxM::checkSquare(const char* operation) const {
    if (!isSquare()) {
        throw std::invalid_argument(
            std::string("Matrix must be square for ") + operation);
    }
}

// 矩阵加法
MatrixNxM MatrixNxM::operator+(const MatrixNxM& other) const {
    checkDimensions(other, "addition");
    MatrixNxM result(n, m);
    for (size_t i = 0; i < n; ++i) {
        for (size_t j = 0; j < m; ++j) {
            result(i, j) = data[i][j] + other(i, j);
        }
    }
    return result;
}

// 矩阵减法
MatrixNxM MatrixNxM::operator-(const MatrixNxM& other) const {
    checkDimensions(other, "subtraction");
    MatrixNxM result(n, m);
    for (size_t i = 0; i < n; ++i) {
        for (size_t j = 0; j < m; ++j) {
            result(i, j) = data[i][j] - other(i, j);
        }
    }
    return result;
}

// 矩阵乘法
MatrixNxM MatrixNxM::operator*(const MatrixNxM& other) const {
    if (m != other.n) {
        throw std::invalid_argument("Matrix dimensions do not match for multiplication");
    }
    MatrixNxM result(n, other.m);
    for (size_t i = 0; i < n; ++i) {
        for (size_t j = 0; j < other.m; ++j) {
            double sum = 0.0;
            for (size_t k = 0; k < m; ++k) {
                sum += data[i][k] * other(k, j);
            }
            result(i, j) = sum;
        }
    }
    return result;
}

// 标量乘法
MatrixNxM MatrixNxM::operator*(double scalar) const {
    MatrixNxM result(n, m);
    for (size_t i = 0; i < n; ++i) {
        for (size_t j = 0; j < m; ++j) {
            result(i, j) = data[i][j] * scalar;
        }
    }
    return result;
}

// 标量除法
MatrixNxM MatrixNxM::operator/(double scalar) const {
    if (std::abs(scalar) < 1e-10) {
        throw std::invalid_argument("Division by zero");
    }
    return *this * (1.0 / scalar);
}

// 转置矩阵
MatrixNxM MatrixNxM::transpose() const {
    MatrixNxM result(m, n);
    for (size_t i = 0; i < n; ++i) {
        for (size_t j = 0; j < m; ++j) {
            result(j, i) = data[i][j];
        }
    }
    return result;
}

// 计算余子式
double MatrixNxM::minor(size_t row, size_t col) const {
    checkSquare("minor calculation");
    if (row >= n || col >= m) {
        throw std::out_of_range("Row or column index out of range");
    }

    // 创建子矩阵
    std::vector<std::vector<double>> subdata(n - 1, std::vector<double>(m - 1));
    size_t sub_i = 0;
    for (size_t i = 0; i < n; ++i) {
        if (i == row) continue;
        size_t sub_j = 0;
        for (size_t j = 0; j < m; ++j) {
            if (j == col) continue;
            subdata[sub_i][sub_j] = data[i][j];
            ++sub_j;
        }
        ++sub_i;
    }
    return MatrixNxM(subdata).determinant();
}

// 计算代数余子式
double MatrixNxM::cofactor(size_t row, size_t col) const {
    double minor_value = minor(row, col);
    return ((row + col) % 2 == 0) ? minor_value : -minor_value;
}

// 计算行列式
double MatrixNxM::determinant() const {
    checkSquare("determinant calculation");
    
    if (n == 1) {
        return data[0][0];
    }
    if (n == 2) {
        return data[0][0] * data[1][1] - data[0][1] * data[1][0];
    }
    if (n == 3) {
        return data[0][0] * data[1][1] * data[2][2] +
               data[0][1] * data[1][2] * data[2][0] +
               data[0][2] * data[1][0] * data[2][1] -
               data[0][2] * data[1][1] * data[2][0] -
               data[0][0] * data[1][2] * data[2][1] -
               data[0][1] * data[1][0] * data[2][2];
    }

    // 对于n>3的情况，使用拉普拉斯展开
    double det = 0.0;
    for (size_t j = 0; j < m; ++j) {
        det += data[0][j] * cofactor(0, j);
    }
    return det;
}

// 计算逆矩阵
MatrixNxM MatrixNxM::inverse() const {
    checkSquare("inverse calculation");
    double det = determinant();
    if (std::abs(det) < 1e-10) {
        throw std::runtime_error("Matrix is singular");
    }

    MatrixNxM result(n, n);
    for (size_t i = 0; i < n; ++i) {
        for (size_t j = 0; j < n; ++j) {
            result(i, j) = cofactor(j, i) / det;  // 注意：转置了余子式矩阵
        }
    }
    return result;
}

// 检查是否奇异
bool MatrixNxM::isSingular() const {
    return !isSquare() || std::abs(determinant()) < 1e-10;
}

// 检查是否非奇异
bool MatrixNxM::isRegular() const {
    return !isSingular();
}

// 检查是否对称
bool MatrixNxM::isSymmetric() const {
    if (!isSquare()) return false;
    for (size_t i = 0; i < n; ++i) {
        for (size_t j = i + 1; j < m; ++j) {
            if (std::abs(data[i][j] - data[j][i]) > 1e-10) {
                return false;
            }
        }
    }
    return true;
}

// 检查是否正交
bool MatrixNxM::isOrthogonal() const {
    if (!isSquare()) return false;
    MatrixNxM transposed = transpose();
    MatrixNxM product = *this * transposed;
    return product.isSymmetric() && std::abs(product.determinant() - 1.0) < 1e-10;
}

// 创建单位矩阵
MatrixNxM MatrixNxM::identity(size_t n) {
    MatrixNxM result(n, n);
    for (size_t i = 0; i < n; ++i) {
        result(i, i) = 1.0;
    }
    return result;
}

// 创建零矩阵
MatrixNxM MatrixNxM::zero(size_t rows, size_t cols) {
    return MatrixNxM(rows, cols);
}

// 创建2D旋转矩阵
MatrixNxM MatrixNxM::rotation2D(double theta) {
    MatrixNxM result(2, 2);
    double cos_theta = std::cos(theta);
    double sin_theta = std::sin(theta);
    result(0, 0) = cos_theta;
    result(0, 1) = -sin_theta;
    result(1, 0) = sin_theta;
    result(1, 1) = cos_theta;
    return result;
}

// 创建3D旋转矩阵
MatrixNxM MatrixNxM::rotation3D(double theta, char axis) {
    MatrixNxM result(3, 3);
    double cos_theta = std::cos(theta);
    double sin_theta = std::sin(theta);

    switch (axis) {
        case 'x':
            result(0, 0) = 1.0;
            result(1, 1) = cos_theta;
            result(1, 2) = -sin_theta;
            result(2, 1) = sin_theta;
            result(2, 2) = cos_theta;
            break;
        case 'y':
            result(0, 0) = cos_theta;
            result(0, 2) = sin_theta;
            result(1, 1) = 1.0;
            result(2, 0) = -sin_theta;
            result(2, 2) = cos_theta;
            break;
        case 'z':
            result(0, 0) = cos_theta;
            result(0, 1) = -sin_theta;
            result(1, 0) = sin_theta;
            result(1, 1) = cos_theta;
            result(2, 2) = 1.0;
            break;
        default:
            throw std::invalid_argument("Invalid rotation axis");
    }
    return result;
}

// 创建缩放矩阵
MatrixNxM MatrixNxM::scaling(const std::vector<double>& factors) {
    size_t n = factors.size();
    MatrixNxM result(n, n);
    for (size_t i = 0; i < n; ++i) {
        result(i, i) = factors[i];
    }
    return result;
}

// 输出运算符
std::ostream& operator<<(std::ostream& os, const MatrixNxM& m) {
    os << "Matrix " << m.n << "x" << m.m << ":\n";
    for (size_t i = 0; i < m.n; ++i) {
        os << "| ";
        for (size_t j = 0; j < m.m; ++j) {
            os << m.data[i][j] << " ";
        }
        os << "|\n";
    }
    return os;
}

} // namespace linear_algebra
} // namespace math

3. 使用示例

3.1 基本运算

#include "matrixnxm.h"
#include <iostream>

using namespace math::linear_algebra;

int main() {
    // 创建矩阵
    std::vector<std::vector<double>> data1 = {
        {1.0, 2.0, 3.0},
        {4.0, 5.0, 6.0}
    };
    std::vector<std::vector<double>> data2 = {
        {7.0, 8.0},
        {9.0, 10.0},
        {11.0, 12.0}
    };
    MatrixNxM m1(data1);
    MatrixNxM m2(data2);

    // 矩阵运算
    MatrixNxM sum = m1 + m1;      // 矩阵加法
    MatrixNxM diff = m1 - m1;     // 矩阵减法
    MatrixNxM product = m1 * m2;  // 矩阵乘法
    MatrixNxM scaled = m1 * 2.0;  // 标量乘法

    // 输出结果
    std::cout << "m1:\n" << m1 << std::endl;
    std::cout << "m2:\n" << m2 << std::endl;
    std::cout << "m1 + m1:\n" << sum << std::endl;
    std::cout << "m1 - m1:\n" << diff << std::endl;
    std::cout << "m1 * m2:\n" << product << std::endl;
    std::cout << "m1 * 2:\n" << scaled << std::endl;

    return 0;
}

3.2 特殊矩阵

#include "matrixnxm.h"
#include <iostream>

using namespace math::linear_algebra;

int main() {
    // 创建单位矩阵
    MatrixNxM identity = MatrixNxM::identity(3);
    std::cout << "3x3单位矩阵:\n" << identity << std::endl;

    // 创建旋转矩阵
    MatrixNxM rot2d = MatrixNxM::rotation2D(M_PI / 4);  // 45度旋转
    std::cout << "2D旋转矩阵:\n" << rot2d << std::endl;

    MatrixNxM rot3d = MatrixNxM::rotation3D(M_PI / 4, 'z');  // 绕z轴旋转45度
    std::cout << "3D旋转矩阵:\n" << rot3d << std::endl;

    // 创建缩放矩阵
    std::vector<double> factors = {2.0, 3.0, 4.0};
    MatrixNxM scale = MatrixNxM::scaling(factors);
    std::cout << "缩放矩阵:\n" << scale << std::endl;

    return 0;
}

4. 编译和运行

4.1 编译

使用提供的Makefile进行编译：

make        # 编译所有目标
make test   # 运行测试
make demo   # 运行示例
make clean  # 清理编译文件

4.2 运行测试

./matrixnxm_test

4.3 运行示例

./matrixnxm_demo

5. 注意事项

1. 数值精度

   • 使用适当的误差范围（1e-10）判断矩阵性质
   • 注意浮点数计算的精度问题
   • 高维矩阵计算可能涉及较大的数值，注意溢出问题

2. 异常处理

   • 处理维度不匹配的情况
   • 处理奇异矩阵的情况
   • 处理索引越界的情况
   • 处理无效的输入数据

3. 性能考虑

   • 矩阵乘法是O(n³)操作
   • 行列式计算是O(n!)操作
   • 考虑使用并行计算优化性能
   • 对于特定维度使用特化实现
   • 考虑使用SIMD指令优化计算

4. 使用建议

   • 优先使用成员函数而不是全局函数
   • 保持接口的一致性
   • 提供清晰的错误信息
   • 对于大型矩阵，考虑使用稀疏矩阵实现
   • 考虑使用BLAS等优化库

6. 扩展阅读

1. 矩阵运算

   • 矩阵分解（LU、QR、SVD等）
   • 特征值分解
   • 奇异值分解
   • 矩阵范数

2. 应用领域

   • 线性方程组求解
   • 最小二乘问题
   • 图像处理
   • 机器学习
   • 计算机图形学
   • 控制系统
   • 量子力学