RcppArmadillo：复杂矩阵求逆编译错误

Question

RcppArmadillo：复杂矩阵求逆编译错误

3

我正在尝试使用RcppArmadillo反转一个复数方阵：

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

// [[Rcpp::export]]
cx_mat fn(cx_mat x) {
  return x.i();
}

然而，当我使用sourceCpp时会出现错误："undefined reference to zgetri_'"。如果我只是将cx_mat替换为mat，那么它就可以编译并正常工作，但这只适用于实矩阵。使用inv()会导致相同的编译错误。有趣的是，使用伪逆pinv()可以通过编译，但与R中的solve()结果略有不同。

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

// [[Rcpp::export]]
cx_mat fn(cx_mat x) {
  return pinv(x);
}

> a<-matrix(c(3+0.1i,7+2i,4+0i,2+0.5i),ncol=2)
> a
       [,1]   [,2]
[1,] 3+0.1i 4+0.0i
[2,] 7+2.0i 2+0.5i
> identical(solve(a), fn(a))
[1] FALSE
> solve(a) - fn(a)
                           [,1]                       [,2]
[1,] -6.938894e-17-7.80626e-18i  8.326673e-17-6.93889e-18i
[2,]  1.665335e-16+1.95156e-17i -1.665335e-16+4.16334e-17i

我知道在这种情况下，差异是机器精度级别的调整，但我仍然想知道是否有办法使inv()或.i()适用于复杂矩阵。谢谢。

- aenima

你确定不想对mat和t(conjugate(mat))的乘积进行求逆吗？我在想是否提供一个更详细、警示性的错误信息会更有帮助。 - IRTFM

这是编译/链接失败，而提问者没有展示。 - Dirk Eddelbuettel

1个回答

网页内容由stack overflow 提供, 点击上面的

可以查看英文原文，
原文链接

- Dirk Eddelbuettel · Accepted Answer

如果您在使用Rlapack.so安装的R时遇到问题，那么这是一个已知的问题，特别是在Windows或某些RedHat系统上使用RcppArmadillo时。

最好的解决方法是：

不使用未包含在Rlapack中的复杂值函数（其中一些可能没有选择）
重新配置您的R安装以使用完整的外部liblapack.so

我们在RcppArmadillo repo上有三个未解决的问题（实际上，我今天还添加了一个列表{{link3:缺少的十二个复杂函数，它们被Armadillo使用但在Rlapack.so中缺失）}，并且我刚刚恳求R Core添加更多的复杂值函数到Rlapack中。

只是为了强调第二点，您的示例在这里可以正常工作，因为我在Debian/Ubuntu构建中没有使用Rlapack:

R> library(Rcpp)
R> sourceCpp("/tmp/aenima.cpp")

R> a <- matrix(c(3+0.1i,7+2i,4+0i,2+0.5i),ncol=2)

R> fn(a)
                      [,1]                [,2]
[1,] -0.0898473+0.0029949i  0.167715-0.047919i
[2,]  0.3174603+0.0000000i -0.126984+0.031746i
R>

使用稍微修改过的文件，结尾处带有示例：

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

// [[Rcpp::export]]
cx_mat fn(cx_mat x) {
  return x.i();
}

/*** R
a <- matrix(c(3+0.1i,7+2i,4+0i,2+0.5i),ncol=2)
fn(a)
*/