c++ - Eigen 为函数的就地版本和非就地版本给出不同的结果

Question

我遇到了一个奇怪的问题，两个应该给出相同结果的函数却不一致。我已经包含了下面的代码。我知道 test1 的结果是正确的，而 test2 是错误的。

#include <RcppEigen.h>
using namespace Rcpp;
// [[Rcpp::depends(RcppEigen)]]

// [[Rcpp::export]]
Eigen::MatrixXd test1(Eigen::MatrixXd A){
  int p = A.rows();
  return A.triangularView<Eigen::Lower>().solve(Eigen::MatrixXd::Identity(p,p)).transpose();
}

// [[Rcpp::export]]
Eigen::MatrixXd test2(Eigen::MatrixXd A){
  int p = A.rows();
  Eigen::MatrixXd I = Eigen::MatrixXd::Identity(p,p);
  A.triangularView<Eigen::Lower>().solveInPlace(I);
  A.transposeInPlace();
  return A;
}


/*** R
A <- rWishart(1, 10, diag(4))[,,1]
A <- t(chol(A))

test1(A)
test2(A)
*/

这是输出

> test1(A)
          [,1]        [,2]        [,3]        [,4]
[1,] 0.2251857 -0.01455544 -0.20205410 -0.08993337
[2,] 0.0000000  0.32498583 -0.06486972 -0.14006616
[3,] 0.0000000  0.00000000  0.60379357  0.27294390
[4,] 0.0000000  0.00000000  0.00000000  0.37409978

> test2(A)
         [,1]      [,2]      [,3]        [,4]
[1,] 4.440779 0.1988932 1.5074352  0.04220045
[2,] 0.000000 3.0770572 0.3305895  0.91087781
[3,] 0.000000 0.0000000 1.6561952 -1.20836313
[4,] 0.000000 0.0000000 0.0000000  2.67308367

我的问题是如何编写不正确的 test1 的就地版本？还有为什么 test2 不正确？

Answer 1

行:

A.triangularView<Eigen::Lower>().solveInPlace(I);

修改 I 而不是 A。所以你需要结束 test2 :

I.transposeInPlace();
return I;