如何在c++中实现一个高效的稀疏矩阵? (压缩存储格式)

CSR格式是C++中平衡内存、随机访问和乘法性能的首选稀疏矩阵表示,适用于行遍历多、列访问少的场景;其核心为row_ptr、col_indices和values三个数组,构建需两趟扫描或先COO后转换,避免频繁单点插入。

用 CSR 格式实现稀疏矩阵最实用

CSR(Compressed Sparse Row)是 C++ 中平衡内存、随机访问和乘法性能的首选压缩格式,尤其适合行遍历多、列访问少的场景(如大多数线性代数库和图算法)。它不支持高效列插入,但对已知非零结构的矩阵初始化和矩阵-向量乘非常快。

  • row_ptr 是长度为 n_rows + 1 的数组,row_ptr[i] 表示第 i 行第一个非零元在 valuescol_indices 中的起始下标
  • col_indices 存储每个非零元的列号(按行优先顺序)
  • values 存储对应非零值,与 col_indices 一一对应
  • 构造时需先统计每行非零个数,再做两趟扫描——这是避免动态重分配的关键

别直接手写 insert(),先用 COO 构建再转 CSR

频繁单点插入会破坏 CSR 的紧凑性,导致 O(n) 查找或 O(nnz) 重排。正确做法是:收集三元组 (row, col, value)std::vector,排序后去重(可选),最后一次性构建 CSR。

  • 排序用 std::sortrow 主序、col 次序:
    std::sort(triplets.begin(), triplets.end(),
        [](const auto& a, const auto& b) {
            return a.row != b.row ? a.row < b.row : a.col < b.col;
        });
  • 去重合并(如相同位置多次赋值)需手动遍历合并,不能依赖 std::unique
  • 转换时用 row_ptr[0] = 0,然后对每行计数填 row_ptr[i+1] = row_ptr[i] + count_of_row_i

矩阵-向量乘 y = A·x 必须手写循环,别依赖 operator[]

CSR 的核心优势在 y[i] += values[k] * x[col_indices[k]] 这种无分支密集访存,而通用 operator()(i,j) 查找会退化成 O(nnz) 每次调用。

  • 典型乘法循环结构:
    for (int i = 0; i < n_rows; ++i) {
        y[i] = 0.0;
        for (int k = row_ptr[i]; k < row_ptr[i+1]; ++k) {
            y[i] += values[k] * x[col_indices[k]];
        }
    }
  • x 很大且稀疏,可考虑转置 CSR(CSC)或引入索引过滤,但通常不值得
  • 注意 col_indices[k] 必须在 [0, n_cols) 范围内,否则是数据损坏,不是越界检查问题

用 std::vector 管理内存,但避免 vector

std::vectorstd::vector 是 CSR 的标准选择;vector 是位压缩特化,不满足 RandomAccessIterator 要求,会导致 data() 失效、迭代器失效、无法传给 BLAS 接口。

  • 三个核心缓冲区声明应为:
    std::vector values;
    std::vector col_indices;
    std::vector row_ptr;
  • 如果需要只读共享(如多个算法共用同一矩阵),可包装为 struct 并提供 const 成员函数,但不要用 shared_ptr 包裹原始指针——这增加间接层且无必要
  • 对超大规模矩阵(>10M 非零元),考虑用 std::vectorres

    erve()
    预分配,避免多次 realloc

CSR 的真正复杂点不在结构本身,而在构建阶段的排序/去重逻辑和乘法循环中对 row_ptr 边界的精确把握——这两处出错不会崩溃,但结果全错,且难以调试。