COO格式稀疏矩阵进行Permutation置换操作
2026/7/22 5:44:12 网站建设 项目流程

文章目录

      • ✅ 步骤 1:应用对称置换
      • ✅ 步骤 2:构造逆排列
      • ✅ 步骤 3:应用置换并排序 COO
      • 🔁 如何“逆置换”?
      • ✅ 补充说明

要对一个COO 格式(Coordinate Format)的稀疏方阵A ∈ R n × n A \in \mathbb{R}^{n \times n}ARn×n应用对称置换(symmetric permutation),即计算:

A ′ = P A P T A' = P A P^TA=PAPT

其中 ( P ) 是由给定的一维排列数组perm(长度为 ( n ))所定义的置换矩阵,其作用为:

  • P [ i , perm [ i ] ] = 1 P[i, \text{perm}[i]] = 1P[i,perm[i]]=1
  • 等价地,对向量x xx,有( P x ) i = x perm [ i ] (Px)_i = x_{\text{perm}[i]}(Px)i=xperm[i]

那么对矩阵的对称置换等价于:

  • 行重排:行i ii变为原来的行perm [ i ] \text{perm}[i]perm[i]
  • 列重排:列j jj变为原来的列perm [ j ] \text{perm}[j]perm[j]

对 COO 格式的三元组(row, col, val),操作非常直接:


✅ 步骤 1:应用对称置换

对每一个非零元素(r, c, v)

  • 新行号 =inv_perm[r]
  • 新列号 =inv_perm[c]

为什么是inv_perm

因为perm[i] = j表示“新位置 i 来自旧位置 j”,即:

  • 新矩阵的第 i 行 = 原矩阵的第perm[i]
    => 原矩阵的第 r 行 → 出现在新矩阵的inv_perm[r]

所以要将原(r, c)映射到新位置(inv_perm[r], inv_perm[c])

结论:对称置换P A P T P A P^TPAPT对应于将每个 COO 元素的行/列索引通过逆排列映射。


✅ 步骤 2:构造逆排列

给定perm(0-based),其逆排列inv_perm满足:

inv_perm[perm[i]]=i;

代码:

std::vector<int>inv_perm(n);for(inti=0;i<n;++i){inv_perm[perm[i]]=i;}

✅ 步骤 3:应用置换并排序 COO

COO 格式通常要求按(row, col)字典序排序(便于后续转换为 CSR 或其他格式)。

C++ 实现如下(假设使用 0-based 索引):

#include<vector>#include<algorithm>#include<tuple>structCOOMatrix{intn;// dimensionstd::vector<int>row,col;// row/col indicesstd::vector<double>val;// values (or use template)};COOMatrixapply_symmetric_permutation(constCOOMatrix&A,conststd::vector<int>&perm){intnnz=A.row.size();intn=A.n;// Step 1: build inverse permutationstd::vector<int>inv_perm(n);for(inti=0;i<n;++i){inv_perm[perm[i]]=i;}// Step 2: apply permutation to each entrystd::vector<int>new_row(nnz),new_col(nnz);std::vector<double>new_val=A.val;// copy valuesfor(intk=0;k<nnz;++k){new_row[k]=inv_perm[A.row[k]];new_col[k]=inv_perm[A.col[k]];}// Step 3: sort by (row, col)std::vector<std::tuple<int,int,double>>entries(nnz);for(intk=0;k<nnz;++k){entries[k]=std::make_tuple(new_row[k],new_col[k],new_val[k]);}std::sort(entries.begin(),entries.end());// Step 4: unpack backCOOMatrix A_perm;A_perm.n=n;A_perm.row.resize(nnz);A_perm.col.resize(nnz);A_perm.val.resize(nnz);for(intk=0;k<nnz;++k){std::tie(A_perm.row[k],A_perm.col[k],A_perm.val[k])=entries[k];}returnA_perm;}

🔁 如何“逆置换”?

如果你已经对A AA应用了对称置换得到A ′ = P A P T A' = P A P^TA=PAPT,那么恢复原矩阵只需应用逆排列的对称置换:

  • 原排列为perm,则逆操作的排列是inv_perm
  • 即:A = P^T A' P = apply_symmetric_permutation(A_prime, inv_perm)

或者更简单:A'再次用perm作为排列调用上述函数?❌ 不行!

正确做法:要撤销P A P T P A P^TPAPT,应使用P T A ′ P = ( P − 1 ) A ′ ( P − 1 ) T P^T A' P = (P^{-1}) A' (P^{-1})^TPTAP=(P1)A(P1)T,而P − 1 P^{-1}P1对应的排列就是inv_perm

所以:

COOMatrix A_restored=apply_symmetric_permutation(A_perm,inv_perm_of_perm);

但注意:inv_perm_of_perm就是原始的perm!因为inv_perm的逆是perm

因此:

  • 正向:A1 = apply(..., perm)
  • 逆向:A0 = apply(..., inv_perm),其中inv_permperm构造

✅ 补充说明

  • 该方法适用于任意稀疏结构(包括非对称矩阵),但对称置换常用于对称矩阵(如刚度矩阵)以改善数值性质。
  • COO 输出已按(row, col)排序,可直接用于Eigen::SparseMatrix::setFromTriplets或其他库。
  • 若你使用的是size_tint64_t索引,请相应调整类型。

需要专业的网站建设服务?

联系我们获取免费的网站建设咨询和方案报价,让我们帮助您实现业务目标

立即咨询