Eigen: 如果innerIndexPtr向量排序不好,selfadjointView就不工作。

问题描述 投票:0回答:1

在一个实时系统中,我使用低级别的API来操作Eigen稀疏矩阵,以提高特定问题集的性能。cs_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

selfadjointView 不起作用。

y = A.selfadjointView<Upper>() * x;
w = B.selfadjointView<Upper>() * x;
assert( y( 0 ) == w( 0 ) ); // <-- Fail!

Eigen文档中的例子(https:/eigen.tuxfamily.orgdoxgroup__TutorialSparse.html。)显示有序数据,但没有明确的指示。

不幸的是,我不能用Eigen得到Ap,因为时间对象的动态分配。有什么办法吗?

已经使用Eigen v3.3.7进行了测试。

eigen eigen3
1个回答
2
投票

要对一个矩阵的条目进行排序,你可以将它转置两次。

B = B.transpose(); B = B.transpose();

或者你可以将它转置一次,然后使用 selfadjointView<Lower>()或者你可以将它分配给一个行大矩阵(这也是隐式的转置)。

SparseMatrix<double, RowMajor> C = B;

w = C.selfadjointView<Upper>() * x;
© www.soinside.com 2019 - 2024. All rights reserved.