我在MATLAB中有以下代码,在尝试计算lamb
之前,它会进行LU分解,我将其包含在内以提供一些上下文。
P=[1,2,3;4,5,6;7,8,9];
U=[0;1;2];
[F,J]=lu(P);
Jlamda=F\U;
lamb=J\Jlamda;
F为:
0.142857142857143 1 0
0.571428571428571 0.500000000000000 1
1 0 0
U为:
7 8 9
0 0.857142857142857 1.71428571428571
0 0 1.11022302462516e-16
当我尝试在Eigen中用以下代码复制这个时:
MatrixXd P(3, 3);
P << 1, 2, 3, 4, 5, 6, 7, 8, 9;
MatrixXd U(3, 1);
U << 0, 1, 2;
PartialPivLU<MatrixXd> lu = PartialPivLU<MatrixXd>(P);
MatrixXd J = lu.matrixLU().triangularView<UpLoType::Upper>();
MatrixXd F = lu.matrixLU().triangularView<UpLoType::UnitLower>();
MatrixXd Jlamda = F.lu().solve(U);
MatrixXd l = J.lu().solve(Jlamda);
cout << F << endl;
cout << endl;
cout << J << endl;
其中打印:
7 8 9
0 0.857143 1.71429
0 0 1.11022e-16
虽然我显然可以手工制作一个矩阵,将C++中的F行转换为MATLAB中的F行,但我不确定如何动态地完成这一操作。PartialPivLU
是实现这一点的最佳方法吗?还是我遗漏了一些更琐碎的东西?
2条答案
按热度按时间6pp0gazn1#
通过调用
[F,J]=lu(P)
,得到的矩阵F
是一个 * 置换 * 下三角矩阵。您可以将函数调用为[F,J,perm]=lu(P)
,以将F
接收为真正的下三角矩阵,并将P
接收为单独的置换矩阵,以便F*J = perm*P
。MatrixXd F = lu.matrixLU().triangularView<UpLoType::UnitLower>();
如果你想要像Matlab返回的那样的置换后的下三角矩阵,那么你可以通过调用
permutationP
然后将这个矩阵乘以F
来将置换矩阵存储在特征值中。5rgfhyps2#
DCSmith给出的答案非常接近,但是似乎还必须对置换矩阵调用transpose()才能获得正确的结果:
图纸:
这与MatLab示例相同