线性代数

掌握矩阵运算、线性系统求解和高级线性代数技术

概述

线性代数是现代算法的数学基础，在图形学、机器学习、网络分析等领域有广泛应用。本章节深入探讨矩阵运算、高斯消元法、特征值计算等核心技术，为解决复杂的线性系统问题提供强大工具。

核心内容

1. 基础矩阵运算

矩阵运算是线性代数的基础，包括加法、乘法、转置等基本操作。 

// 基础矩阵运算实现
#include <iostream>
#include <vector>
#include <iomanip>
using namespace std;

class Matrix {
private:
    vector<vector<double>> data;
    int rows, cols;

public:
    Matrix(int r, int c) : rows(r), cols(c) {
        data.resize(r, vector<double>(c, 0));
    }
    
    Matrix(vector<vector<double>>& mat) {
        data = mat;
        rows = mat.size();
        cols = mat[0].size();
    }
    
    // 矩阵加法
    Matrix operator+(const Matrix& other) const {
        if (rows != other.rows || cols != other.cols) {
            throw invalid_argument("矩阵维度不匹配");
        }
        
        Matrix result(rows, cols);
        for (int i = 0; i < rows; i++) {
            for (int j = 0; j < cols; j++) {
                result.data[i][j] = data[i][j] + other.data[i][j];
            }
        }
        return result;
    }
    
    // 矩阵乘法
    Matrix operator*(const Matrix& other) const {
        if (cols != other.rows) {
            throw invalid_argument("矩阵维度不匹配");
        }
        
        Matrix result(rows, other.cols);
        for (int i = 0; i < rows; i++) {
            for (int j = 0; j < other.cols; j++) {
                for (int k = 0; k < cols; k++) {
                    result.data[i][j] += data[i][k] * other.data[k][j];
                }
            }
        }
        return result;
    }
    
    // 矩阵转置
    Matrix transpose() const {
        Matrix result(cols, rows);
        for (int i = 0; i < rows; i++) {
            for (int j = 0; j < cols; j++) {
                result.data[j][i] = data[i][j];
            }
        }
        return result;
    }
    
    // 矩阵快速幂（用于递推关系）
    Matrix power(long long n) const {
        if (rows != cols) {
            throw invalid_argument("非方阵无法求幂");
        }
        
        Matrix result = identity(rows);
        Matrix base = *this;
        
        while (n > 0) {
            if (n % 2 == 1) {
                result = result * base;
            }
            base = base * base;
            n /= 2;
        }
        return result;
    }
    
    // 创建单位矩阵
    static Matrix identity(int n) {
        Matrix result(n, n);
        for (int i = 0; i < n; i++) {
            result.data[i][i] = 1;
        }
        return result;
    }
    
    // 打印矩阵
    void print() const {
        for (int i = 0; i < rows; i++) {
            for (int j = 0; j < cols; j++) {
                cout << setw(8) << fixed << setprecision(2) << data[i][j] << " ";
            }
            cout << endl;
        }
    }
    
    double& operator()(int i, int j) { return data[i][j]; }
    const double& operator()(int i, int j) const { return data[i][j]; }
    int getRows() const { return rows; }
    int getCols() const { return cols; }
};

// 应用：斐波那契数列的矩阵快速计算
long long fibonacciMatrix(int n) {
    if (n <= 1) return n;
    
    // 转移矩阵 [[1,1],[1,0]]
    vector<vector<double>> trans = {{1, 1}, {1, 0}};
    Matrix transition(trans);
    
    Matrix result = transition.power(n - 1);
    return (long long)(result(0, 0) + 0.5); // 四舍五入
}

int main() {
    // 测试基础运算
    vector<vector<double>> a = {{1, 2}, {3, 4}};
    vector<vector<double>> b = {{5, 6}, {7, 8}};
    
    Matrix A(a), B(b);
    
    cout << "矩阵A：" << endl;
    A.print();
    
    cout << "\n矩阵B：" << endl;
    B.print();
    
    cout << "\nA + B：" << endl;
    (A + B).print();
    
    cout << "\nA * B：" << endl;
    (A * B).print();
    
    cout << "\nA的转置：" << endl;
    A.transpose().print();
    
    // 测试矩阵快速幂
    cout << "\n斐波那契数列前10项：" << endl;
    for (int i = 0; i <= 10; i++) {
        cout << "F(" << i << ") = " << fibonacciMatrix(i) << endl;
    }
    
    return 0;
}

2. 高斯消元法

高斯消元法是求解线性方程组的经典算法，也是计算矩阵逆、行列式的基础。 

// 高斯消元法实现
#include <iostream>
#include <vector>
#include <cmath>
#include <iomanip>
using namespace std;

const double EPS = 1e-9;

class LinearSystem {
private:
    vector<vector<double>> augmented;
    int n, m;

public:
    LinearSystem(vector<vector<double>>& A, vector<double>& b) {
        n = A.size();
        m = A[0].size();
        augmented.resize(n, vector<double>(m + 1));
        
        for (int i = 0; i < n; i++) {
            for (int j = 0; j < m; j++) {
                augmented[i][j] = A[i][j];
            }
            augmented[i][m] = b[i];
        }
    }
    
    // 高斯消元求解线性方程组
    vector<double> solve() {
        int rank = 0;
        vector<int> pivot_col(n, -1);
        
        for (int col = 0; col < m && rank < n; col++) {
            // 寻找主元
            int pivot_row = -1;
            for (int row = rank; row < n; row++) {
                if (abs(augmented[row][col]) > EPS) {
                    if (pivot_row == -1 || abs(augmented[row][col]) > abs(augmented[pivot_row][col])) {
                        pivot_row = row;
                    }
                }
            }
            
            if (pivot_row == -1) continue; // 当前列全为0
            
            // 交换行
            if (pivot_row != rank) {
                swap(augmented[pivot_row], augmented[rank]);
            }
            
            pivot_col[rank] = col;
            
            // 消元
            for (int row = 0; row < n; row++) {
                if (row != rank && abs(augmented[row][col]) > EPS) {
                    double factor = augmented[row][col] / augmented[rank][col];
                    for (int j = col; j <= m; j++) {
                        augmented[row][j] -= factor * augmented[rank][j];
                    }
                }
            }
            rank++;
        }
        
        // 检查解的存在性
        for (int i = rank; i < n; i++) {
            if (abs(augmented[i][m]) > EPS) {
                throw runtime_error("方程组无解");
            }
        }
        
        vector<double> solution(m, 0);
        for (int i = 0; i < rank; i++) {
            if (pivot_col[i] != -1) {
                solution[pivot_col[i]] = augmented[i][m] / augmented[i][pivot_col[i]];
            }
        }
        
        return solution;
    }
    
    // 计算矩阵的秩
    int getRank() {
        int rank = 0;
        for (int col = 0; col < m && rank < n; col++) {
            int pivot_row = -1;
            for (int row = rank; row < n; row++) {
                if (abs(augmented[row][col]) > EPS) {
                    pivot_row = row;
                    break;
                }
            }
            
            if (pivot_row == -1) continue;
            
            if (pivot_row != rank) {
                swap(augmented[pivot_row], augmented[rank]);
            }
            
            for (int row = rank + 1; row < n; row++) {
                if (abs(augmented[row][col]) > EPS) {
                    double factor = augmented[row][col] / augmented[rank][col];
                    for (int j = col; j <= m; j++) {
                        augmented[row][j] -= factor * augmented[rank][j];
                    }
                }
            }
            rank++;
        }
        return rank;
    }
};

// 计算行列式
double determinant(vector<vector<double>> matrix) {
    int n = matrix.size();
    double det = 1;
    
    for (int i = 0; i < n; i++) {
        // 寻找主元
        int pivot = i;
        for (int j = i + 1; j < n; j++) {
            if (abs(matrix[j][i]) > abs(matrix[pivot][i])) {
                pivot = j;
            }
        }
        
        if (abs(matrix[pivot][i]) < EPS) return 0; // 矩阵奇异
        
        if (pivot != i) {
            swap(matrix[i], matrix[pivot]);
            det = -det; // 交换行改变符号
        }
        
        det *= matrix[i][i];
        
        // 消元
        for (int j = i + 1; j < n; j++) {
            if (abs(matrix[j][i]) > EPS) {
                double factor = matrix[j][i] / matrix[i][i];
                for (int k = i; k < n; k++) {
                    matrix[j][k] -= factor * matrix[i][k];
                }
            }
        }
    }
    
    return det;
}

// 矩阵求逆
vector<vector<double>> inverse(vector<vector<double>> A) {
    int n = A.size();
    vector<vector<double>> inv(n, vector<double>(n, 0));
    
    // 创建增广矩阵 [A|I]
    for (int i = 0; i < n; i++) {
        inv[i][i] = 1;
    }
    
    // 高斯-约旦消元
    for (int i = 0; i < n; i++) {
        // 寻找主元
        int pivot = i;
        for (int j = i + 1; j < n; j++) {
            if (abs(A[j][i]) > abs(A[pivot][i])) {
                pivot = j;
            }
        }
        
        if (abs(A[pivot][i]) < EPS) {
            throw runtime_error("矩阵不可逆");
        }
        
        if (pivot != i) {
            swap(A[i], A[pivot]);
            swap(inv[i], inv[pivot]);
        }
        
        // 归一化主元行
        double pivot_val = A[i][i];
        for (int j = 0; j < n; j++) {
            A[i][j] /= pivot_val;
            inv[i][j] /= pivot_val;
        }
        
        // 消元
        for (int j = 0; j < n; j++) {
            if (j != i && abs(A[j][i]) > EPS) {
                double factor = A[j][i];
                for (int k = 0; k < n; k++) {
                    A[j][k] -= factor * A[i][k];
                    inv[j][k] -= factor * inv[i][k];
                }
            }
        }
    }
    
    return inv;
}

int main() {
    // 测试线性方程组求解
    vector<vector<double>> A = {
        {2, 1, -1},
        {-3, -1, 2},
        {-2, 1, 2}
    };
    vector<double> b = {8, -11, -3};
    
    LinearSystem system(A, b);
    
    try {
        vector<double> solution = system.solve();
        cout << "方程组解：" << endl;
        for (int i = 0; i < solution.size(); i++) {
            cout << "x" << i + 1 << " = " << fixed << setprecision(4) << solution[i] << endl;
        }
    } catch (const exception& e) {
        cout << "错误： " << e.what() << endl;
    }
    
    // 测试行列式计算
    vector<vector<double>> B = {{1, 2, 3}, {4, 5, 6}, {7, 8, 10}};
    cout << "\n行列式： " << determinant(B) << endl;
    
    // 测试矩阵求逆
    try {
        vector<vector<double>> C = {{2, 1}, {1, 1}};
        vector<vector<double>> inv_C = inverse(C);
        cout << "\n逆矩阵：" << endl;
        for (auto& row : inv_C) {
            for (double val : row) {
                cout << setw(8) << fixed << setprecision(4) << val << " ";
            }
            cout << endl;
        }
    } catch (const exception& e) {
        cout << "错误： " << e.what() << endl;
    }
    
    return 0;
}

3. 高级应用

特征值与特征向量

特征值分解在图论、主成分分析、马尔可夫链等问题中有重要应用。 

// 幂迭代法求最大特征值
#include <iostream>
#include <vector>
#include <cmath>
#include <random>
using namespace std;

// 向量运算
double dotProduct(const vector<double>& a, const vector<double>& b) {
    double sum = 0;
    for (int i = 0; i < a.size(); i++) {
        sum += a[i] * b[i];
    }
    return sum;
}

double vectorNorm(const vector<double>& v) {
    return sqrt(dotProduct(v, v));
}

vector<double> normalize(const vector<double>& v) {
    double norm = vectorNorm(v);
    vector<double> result(v.size());
    for (int i = 0; i < v.size(); i++) {
        result[i] = v[i] / norm;
    }
    return result;
}

vector<double> matrixVectorMultiply(const vector<vector<double>>& A, const vector<double>& v) {
    int n = A.size();
    vector<double> result(n, 0);
    for (int i = 0; i < n; i++) {
        for (int j = 0; j < n; j++) {
            result[i] += A[i][j] * v[j];
        }
    }
    return result;
}

// 幂迭代法
pair<double, vector<double>> powerIteration(const vector<vector<double>>& A, int maxIter = 1000, double tol = 1e-9) {
    int n = A.size();
    
    // 初始化随机向量
    random_device rd;
    mt19937 gen(rd());
    uniform_real_distribution<> dis(-1.0, 1.0);
    
    vector<double> v(n);
    for (int i = 0; i < n; i++) {
        v[i] = dis(gen);
    }
    v = normalize(v);
    
    double eigenvalue = 0;
    for (int iter = 0; iter < maxIter; iter++) {
        vector<double> Av = matrixVectorMultiply(A, v);
        double newEigenvalue = dotProduct(v, Av);
        
        if (abs(newEigenvalue - eigenvalue) < tol) {
            return {newEigenvalue, normalize(Av)};
        }
        
        eigenvalue = newEigenvalue;
        v = normalize(Av);
    }
    
    return {eigenvalue, v};
}

// 图的邻接矩阵分析
class GraphAnalysis {
private:
    vector<vector<double>> adj;
    int n;

public:
    GraphAnalysis(vector<vector<int>>& edges, int nodes) : n(nodes) {
        adj.resize(n, vector<double>(n, 0));
        
        // 构建邻接矩阵
        for (auto& edge : edges) {
            adj[edge[0]][edge[1]] = 1;
            adj[edge[1]][edge[0]] = 1; // 无向图
        }
    }
    
    // 计算图的谱半径（最大特征值）
    double spectralRadius() {
        auto result = powerIteration(adj);
        return result.first;
    }
    
    // PageRank算法（修改的幂迭代）
    vector<double> pageRank(double damping = 0.85, int maxIter = 1000, double tol = 1e-9) {
        // 构建转移矩阵
        vector<vector<double>> P(n, vector<double>(n, 0));
        vector<int> outDegree(n, 0);
        
        for (int i = 0; i < n; i++) {
            for (int j = 0; j < n; j++) {
                if (adj[i][j] > 0) outDegree[i]++;
            }
        }
        
        for (int i = 0; i < n; i++) {
            for (int j = 0; j < n; j++) {
                if (outDegree[i] > 0) {
                    P[j][i] = adj[i][j] / outDegree[i];
                } else {
                    P[j][i] = 1.0 / n; // 处理悬挂节点
                }
            }
        }
        
        // PageRank迭代
        vector<double> rank(n, 1.0 / n);
        
        for (int iter = 0; iter < maxIter; iter++) {
            vector<double> newRank(n, (1 - damping) / n);
            
            for (int i = 0; i < n; i++) {
                for (int j = 0; j < n; j++) {
                    newRank[i] += damping * P[i][j] * rank[j];
                }
            }
            
            // 检查收敛
            double diff = 0;
            for (int i = 0; i < n; i++) {
                diff += abs(newRank[i] - rank[i]);
            }
            
            if (diff < tol) break;
            rank = newRank;
        }
        
        return rank;
    }
};

int main() {
    // 测试特征值计算
    vector<vector<double>> matrix = {
        {4, 1, 0},
        {1, 3, 1},
        {0, 1, 2}
    };
    
    auto result = powerIteration(matrix);
    cout << "最大特征值： " << result.first << endl;
    cout << "对应特征向量：" << endl;
    for (double val : result.second) {
        cout << val << " ";
    }
    cout << endl;
    
    // 测试图分析
    vector<vector<int>> edges = {{0, 1}, {1, 2}, {2, 0}, {1, 3}, {3, 4}, {4, 1}};
    GraphAnalysis graph(edges, 5);
    
    cout << "\n图的谱半径： " << graph.spectralRadius() << endl;
    
    vector<double> pagerank = graph.pageRank();
    cout << "\nPageRank值：" << endl;
    for (int i = 0; i < pagerank.size(); i++) {
        cout << "节点 " << i << ": " << pagerank[i] << endl;
    }
    
    return 0;
}

解题技巧

⚡ 数值稳定性

使用部分主元选择避免数值不稳定，设置合适的精度阈值，注意浮点数比较的陷阱。

🎯 算法选择

根据问题特点选择算法：稠密矩阵用高斯消元，稀疏矩阵考虑迭代方法，大规模问题使用近似算法。

💡 优化技巧

利用矩阵的特殊性质（对称、三对角等）优化算法，使用块矩阵分解处理大规模问题。

🔍 应用识别

识别线性代数在不同领域的应用：图论中的邻接矩阵、机器学习中的主成分分析、网络分析中的PageRank。

📚 章节目录