【问题标题】:RcppArmadillo: Complex matrix inverse compilation errorRcppArmadillo:复矩阵逆编译错误
【发布时间】:2016-01-03 05:19:23
【问题描述】:

我正在尝试使用 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 它时,这会引发错误: “未定义对 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() 处理复杂的矩阵。谢谢。

【问题讨论】:

  • 您确定不想用 t(conjugate(mat)) 反转 mat 的乘积吗?我确实想知道是否提供更多信息、警告性的错误消息可能会有所帮助。
  • 这是 OP 未显示的编译/链接失败。

标签: r matrix rcpp


【解决方案1】:

这是一个已知问题如果您在使用 Rlapack.so 的 R 安装中使用 RcppArmadillo,例如在 Windows 或某些 RedHat 系统上。

最好的答案是任何一个

  • 不使用 Rlapack has some 中未包含的复值函数,但遗憾的是并非全部(但您可能别无选择)
  • 重新配置您的 R 安装以使用完整的外部 liblapack.so

我们有三张公开票 at the RcppArmadillo repo(事实上,我今天什至添加了一张 the twelve missing complex functions 供 Armadillo 使用但在 Rlapack.so 中丢失),我只是恳求 R Core 添加更多复值函数拉拉包。

为了强调第二点,您的示例在这里运行良好,因为我不在 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)
*/

【讨论】:

  • 感谢您的解释。我在默认 R 安装的 Windows 8.1 上运行它。我想我现在必须和pinv 一起生活。我希望错误传播在我的应用程序中保持微不足道,因为这个矩阵逆只是我的应用程序的起点。
猜你喜欢
  • 2014-07-18
  • 2017-08-27
  • 1970-01-01
  • 1970-01-01
  • 1970-01-01
  • 1970-01-01
  • 2016-12-13
  • 1970-01-01
  • 1970-01-01
相关资源
最近更新 更多