Rcpp NumericMatrix 在 left/right 上乘以标量时的奇怪行为
Bizarre behavior of Rcpp NumericMatrix when multiplied by scalar on left/right
任何人都可以解释以下行为吗?
当声明一个新的NumericMatrix
,y
作为原始矩阵,x
,乘以一个标量,c
,[=41]的顺序=] 乘法很重要。如果我将左边的标量和右边的矩阵相乘(例如 NumericMatrix y = c * x;
),那么我会得到奇怪的行为。原来的矩阵,x
,变了!
然而,如果我将原始矩阵放在左边并在右边乘以标量(例如NumericMatrix y = x * c;
),x
保持不变。
这似乎不会影响其他数据类型。我已经用 int
和 NumericVector
.
进行了测试
示例:使用NumericMatrix
时出现问题
#include <Rcpp.h>
using namespace Rcpp;
// [[Rcpp::export]]
NumericMatrix testfun(NumericMatrix x) {
NumericMatrix y(x.rows(), x.cols());
y = x * 2;
std::cout << x; // x is unmodified
y = 2 * x;
std::cout << x; // x is now modified
return x;
}
/*** R
x <- matrix(2, nrow = 3, ncol = 3)
print(x)
y <- testfun(x = x)
print(y)
print(x)
*/
输出结果如下
> x <- matrix(2, nrow = 3, ncol = 3)
> print(x)
[,1] [,2] [,3]
[1,] 2 2 2
[2,] 2 2 2
[3,] 2 2 2
> y <- testfun(x = x)
2.00000 2.00000 2.00000
2.00000 2.00000 2.00000
2.00000 2.00000 2.00000
4.00000 4.00000 4.00000
4.00000 4.00000 4.00000
4.00000 4.00000 4.00000
> print(y)
[,1] [,2] [,3]
[1,] 4 4 4
[2,] 4 4 4
[3,] 4 4 4
> print(x)
[,1] [,2] [,3]
[1,] 4 4 4
[2,] 4 4 4
[3,] 4 4 4
这是我的会话信息
> sessionInfo()
R version 3.6.1 (2019-07-05)
Platform: x86_64-apple-darwin15.6.0 (64-bit)
Running under: macOS Catalina 10.15.1
Matrix products: default
BLAS: /System/Library/Frameworks/Accelerate.framework/Versions/A/Frameworks/vecLib.framework/Versions/A/libBLAS.dylib
LAPACK: /Library/Frameworks/R.framework/Versions/3.6/Resources/lib/libRlapack.dylib
locale:
[1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
attached base packages:
[1] stats graphics grDevices utils datasets methods base
loaded via a namespace (and not attached):
[1] compiler_3.6.1 tools_3.6.1 RcppArmadillo_0.9.800.1.0
[4] Rcpp_1.0.2 RcppProgress_0.4.1 packrat_0.5.0
[7] RcppParallel_4.4.4
.
发布的问题太长了,隐藏了它的重点。我们不需要第二个和第三个例子。我们只需要这段代码:
代码
#include <Rcpp.h>
using namespace Rcpp;
// [[Rcpp::export]]
NumericMatrix testfun(NumericMatrix x) {
NumericMatrix y(x.rows(), x.cols());
y = x * 2;
std::cout << x; // x is unmodified
y = 2 * x; // contrast with x * 2
std::cout << x; // x is now modified
return x;
}
/*** R
print(x <- matrix(2, nrow = 2, ncol = 2))
print(y <- testfun(x = x))
print(x)
*/
输出
R> Rcpp::sourceCpp("~/git/Whosebug/59515517/question.cpp")
R> print(x <- matrix(2, nrow = 2, ncol = 2))
[,1] [,2]
[1,] 2 2
[2,] 2 2
R> print(y <- testfun(x = x))
2.00000 2.00000
2.00000 2.00000
4.00000 4.00000
4.00000 4.00000
[,1] [,2]
[1,] 4 4
[2,] 4 4
R> print(x)
[,1] [,2]
[1,] 4 4
[2,] 4 4
R>
问题
x * 2
的行为不同于 2 * x
。后者有副作用。这可能是一个错误,
更大的问题
用 Rcpp 做矩阵代数真的没有意义。实施是初步的和不完整的。如果你想做 "math",请使用 RcppArmadillo 或 RcppEigen。
Arma 实施
#include <RcppArmadillo.h>
// [[Rcpp::depends(RcppArmadillo)]]
// [[Rcpp::export]]
arma::mat testfun(arma::mat x) {
arma::mat y(x.n_rows, x.n_cols);
y = x * 2;
std::cout << x; // x is unmodified
y = 2 * x; // contrast with x * 2
std::cout << x; // x is now modified
return x;
}
/*** R
print(x <- matrix(2, nrow = 2, ncol = 2))
print(y <- testfun(x = x))
print(x)
*/
武装输出
R> Rcpp::sourceCpp("~/git/Whosebug/59515517/answer.cpp")
R> print(x <- matrix(2, nrow = 2, ncol = 2))
[,1] [,2]
[1,] 2 2
[2,] 2 2
R> print(y <- testfun(x = x))
2.0000 2.0000
2.0000 2.0000
2.0000 2.0000
2.0000 2.0000
[,1] [,2]
[1,] 2 2
[2,] 2 2
R> print(x)
[,1] [,2]
[1,] 2 2
[2,] 2 2
R>
我会看看能否解决你在这里发现的 Rcpp 错误。
编辑: 这是一个错误,我还不太明白,但我在代码来源的问题上添加了一些内容:https://github.com/RcppCore/Rcpp/issues/365
编辑 2: 修复现在在 master 中。感谢 KK 的 PR,感谢在评论中暗示这背后可能是什么的其他人,特别是感谢 Ralf 让我们所有人都尝试修复。
任何人都可以解释以下行为吗?
当声明一个新的NumericMatrix
,y
作为原始矩阵,x
,乘以一个标量,c
,[=41]的顺序=] 乘法很重要。如果我将左边的标量和右边的矩阵相乘(例如 NumericMatrix y = c * x;
),那么我会得到奇怪的行为。原来的矩阵,x
,变了!
然而,如果我将原始矩阵放在左边并在右边乘以标量(例如NumericMatrix y = x * c;
),x
保持不变。
这似乎不会影响其他数据类型。我已经用 int
和 NumericVector
.
示例:使用NumericMatrix
时出现问题
#include <Rcpp.h>
using namespace Rcpp;
// [[Rcpp::export]]
NumericMatrix testfun(NumericMatrix x) {
NumericMatrix y(x.rows(), x.cols());
y = x * 2;
std::cout << x; // x is unmodified
y = 2 * x;
std::cout << x; // x is now modified
return x;
}
/*** R
x <- matrix(2, nrow = 3, ncol = 3)
print(x)
y <- testfun(x = x)
print(y)
print(x)
*/
输出结果如下
> x <- matrix(2, nrow = 3, ncol = 3)
> print(x)
[,1] [,2] [,3]
[1,] 2 2 2
[2,] 2 2 2
[3,] 2 2 2
> y <- testfun(x = x)
2.00000 2.00000 2.00000
2.00000 2.00000 2.00000
2.00000 2.00000 2.00000
4.00000 4.00000 4.00000
4.00000 4.00000 4.00000
4.00000 4.00000 4.00000
> print(y)
[,1] [,2] [,3]
[1,] 4 4 4
[2,] 4 4 4
[3,] 4 4 4
> print(x)
[,1] [,2] [,3]
[1,] 4 4 4
[2,] 4 4 4
[3,] 4 4 4
这是我的会话信息
> sessionInfo()
R version 3.6.1 (2019-07-05)
Platform: x86_64-apple-darwin15.6.0 (64-bit)
Running under: macOS Catalina 10.15.1
Matrix products: default
BLAS: /System/Library/Frameworks/Accelerate.framework/Versions/A/Frameworks/vecLib.framework/Versions/A/libBLAS.dylib
LAPACK: /Library/Frameworks/R.framework/Versions/3.6/Resources/lib/libRlapack.dylib
locale:
[1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
attached base packages:
[1] stats graphics grDevices utils datasets methods base
loaded via a namespace (and not attached):
[1] compiler_3.6.1 tools_3.6.1 RcppArmadillo_0.9.800.1.0
[4] Rcpp_1.0.2 RcppProgress_0.4.1 packrat_0.5.0
[7] RcppParallel_4.4.4
.
发布的问题太长了,隐藏了它的重点。我们不需要第二个和第三个例子。我们只需要这段代码:
代码
#include <Rcpp.h>
using namespace Rcpp;
// [[Rcpp::export]]
NumericMatrix testfun(NumericMatrix x) {
NumericMatrix y(x.rows(), x.cols());
y = x * 2;
std::cout << x; // x is unmodified
y = 2 * x; // contrast with x * 2
std::cout << x; // x is now modified
return x;
}
/*** R
print(x <- matrix(2, nrow = 2, ncol = 2))
print(y <- testfun(x = x))
print(x)
*/
输出
R> Rcpp::sourceCpp("~/git/Whosebug/59515517/question.cpp")
R> print(x <- matrix(2, nrow = 2, ncol = 2))
[,1] [,2]
[1,] 2 2
[2,] 2 2
R> print(y <- testfun(x = x))
2.00000 2.00000
2.00000 2.00000
4.00000 4.00000
4.00000 4.00000
[,1] [,2]
[1,] 4 4
[2,] 4 4
R> print(x)
[,1] [,2]
[1,] 4 4
[2,] 4 4
R>
问题
x * 2
的行为不同于 2 * x
。后者有副作用。这可能是一个错误,
更大的问题
用 Rcpp 做矩阵代数真的没有意义。实施是初步的和不完整的。如果你想做 "math",请使用 RcppArmadillo 或 RcppEigen。
Arma 实施
#include <RcppArmadillo.h>
// [[Rcpp::depends(RcppArmadillo)]]
// [[Rcpp::export]]
arma::mat testfun(arma::mat x) {
arma::mat y(x.n_rows, x.n_cols);
y = x * 2;
std::cout << x; // x is unmodified
y = 2 * x; // contrast with x * 2
std::cout << x; // x is now modified
return x;
}
/*** R
print(x <- matrix(2, nrow = 2, ncol = 2))
print(y <- testfun(x = x))
print(x)
*/
武装输出
R> Rcpp::sourceCpp("~/git/Whosebug/59515517/answer.cpp")
R> print(x <- matrix(2, nrow = 2, ncol = 2))
[,1] [,2]
[1,] 2 2
[2,] 2 2
R> print(y <- testfun(x = x))
2.0000 2.0000
2.0000 2.0000
2.0000 2.0000
2.0000 2.0000
[,1] [,2]
[1,] 2 2
[2,] 2 2
R> print(x)
[,1] [,2]
[1,] 2 2
[2,] 2 2
R>
我会看看能否解决你在这里发现的 Rcpp 错误。
编辑: 这是一个错误,我还不太明白,但我在代码来源的问题上添加了一些内容:https://github.com/RcppCore/Rcpp/issues/365
编辑 2: 修复现在在 master 中。感谢 KK 的 PR,感谢在评论中暗示这背后可能是什么的其他人,特别是感谢 Ralf 让我们所有人都尝试修复。