使用 PartialPivLU 进行 LU 分解

Nic*_*ton 2 c++ matlab matrix linear-algebra eigen3

我在 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;
Run Code Online (Sandbox Code Playgroud)

F 是:

P=[1,2,3;4,5,6;7,8,9];
U=[0;1;2];
[F,J]=lu(P);
Jlamda=F\U;
lamb=J\Jlamda;
Run Code Online (Sandbox Code Playgroud)

U 是:

0.142857142857143   1                   0
0.571428571428571   0.500000000000000   1
1                   0                   0
Run Code Online (Sandbox Code Playgroud)

当我尝试使用以下代码在 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;
Run Code Online (Sandbox Code Playgroud)

哪个打印:

7                   8                   9
0                   0.857142857142857   1.71428571428571
0                   0                   1.11022302462516e-16
Run Code Online (Sandbox Code Playgroud)
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;
Run Code Online (Sandbox Code Playgroud)

虽然我显然可以手工制作一个矩阵,将 C++ 中的 F 行转换为 MATLAB 中的 F 行,但我不知道如何动态地执行此操作。

最好的方法是PartialPivLU解决这个问题,还是我错过了一些更琐碎的事情?

小智 5

通过调用[F,J]=lu(P),所得矩阵F置换下三角矩阵。您可以调用函数 as[F,J,perm]=lu(P)来接收F真正的下三角矩阵和P单独的置换矩阵,以便F*J = perm*P。默认情况下,该行

MatrixXd F = lu.matrixLU().triangularView<UpLoType::UnitLower>();

使用 Eigen 返回真正的下三角矩阵。如果您想要像 Matlab 返回的置换下三角矩阵,那么您可以通过调用permutationP然后将该矩阵乘以 来将置换矩阵存储在 Eigen 中F