在实时系统上工作时,我使用低级API来操作特征稀疏矩阵,以提高特定问题集的性能。 特别是我将特征与 CSparsecs_symperm算法的变体混合,以获得 A 矩阵的排列:Ap = PAP'。cs_symperm对稀疏矩阵使用相同的低级结构,但对于固定列,行索引可能没有很好的排序。
这是一个简单的示例,由手工构建,说明了可能发生的情况:
SparseMatrix<double> A( 2, 2 );
A.insert( 0, 0 ) = 1.0;
A.insert( 0, 1 ) = 2.0;
A.insert( 1, 1 ) = 3.0;
A.makeCompressed();
SparseMatrix<double> B = A;
B.innerIndexPtr()[0] = 0;
B.innerIndexPtr()[1] = 1; // <-- Not ordered
B.innerIndexPtr()[2] = 0; // <-- Not ordered
B.valuePtr()[0] = 1.0;
B.valuePtr()[1] = 3.0; // <-- Not ordered
B.valuePtr()[2] = 2.0; // <-- Not ordered
这里 A 和 B 是相同的矩阵。唯一的区别是数据顺序。
矩阵向量积正常工作:
VectorXd x( 2 );
x << 1.0, 2.0;
VectorXd y = A * x;
VectorXd w = B * x;
assert( y( 0 ) == w( 0 ) ); // <-- OK
assert( y( 1 ) == w( 1 ) ); // <-- OK
自伴随视图不起作用:
y = A.selfadjointView<Upper>() * x;
w = B.selfadjointView<Upper>() * x;
assert( y( 0 ) == w( 0 ) ); // <-- Fail!
特征文档 (https://eigen.tuxfamily.org/dox/group__TutorialSparse.html( 中的示例显示了有序数据,但没有明确的指示。
不幸的是,由于临时对象的动态分配,我无法使用特征获取 Ap。知道吗?
已使用 Eigen v3.3.7 执行测试。
要对矩阵的条目进行排序,您可以将其转置两次:
B = B.transpose(); B = B.transpose();
或者你可以转置它一次并使用selfadjointView<Lower>()
,或者你可以把它分配给一个行主矩阵(这也隐式转置它(:
SparseMatrix<double, RowMajor> C = B;
w = C.selfadjointView<Upper>() * x;