4

在我的实现中,我经常使用子矩阵和矩阵块。我想知道犰狳中是否有一种方法可以让我提取一个更大的矩阵块,并为这个子矩阵使用与原始矩阵中的块相同的内存。我的问题是我不知道该怎么做,因为原始矩阵中的位置不连续。

这是一个简单的示例,说明了当我的原始矩阵为 时我想要做什么A = [A1 A2]

#include <RcppArmadillo.h>
// [[Rcpp::depends(RcppArmadillo)]]
// [[Rcpp::export]]
arma::mat foo(arma::mat A, arma::uword block_start, arma::uword block_length) {
  arma::uword rows = A.n_rows;
  arma::mat B = arma::mat(A.begin_col(block_start), rows, block_length, false);
// fill B here for illustration; 
// in the real piece of code I do various manipulations, multiplications, etc using B
  B.fill(1.0);
  return A;
}

/*** R
A <- matrix(0, 4, 4)
foo(A, 0, 2)
> A <- matrix(0, 4, 4)
> foo(A, 0, 2)
     [,1] [,2] [,3] [,4]
[1,]    1    1    0    0
[2,]    1    1    0    0
[3,]    1    1    0    0
[4,]    1    1    0    0
*/

在这种情况下,子矩阵的位置是连续的,我可以使用高级构造函数来链接内存。

现在假设我希望子矩阵是A[1:2, 1:2]. 我可以B在犰狳中获得与原始元素使用相同内存的副本A吗?(理想情况下,这个问题的解决方案也可以推广到列也不连续的情况,例如A[c(1, 3), c(1, 3)]。)

编辑:为了澄清,我真的需要B上面函数中的矩阵独立存在。我不在fill我的真实代码中使用它,而是在各种矩阵乘法等中使用它(和多个其他子矩阵)。所以我追求的是一种创建B非连续子矩阵的方法,A同时确保它们使用相同的记忆。

4

1 回答 1

1

您可以使用子矩阵视图写入连续或非连续内存:

#include <RcppArmadillo.h>
// [[Rcpp::depends(RcppArmadillo)]]
// [[Rcpp::export]]
arma::mat foo(arma::mat A) {
    A.submat(0, 0, 1, 1).fill(1.0);
    return A;
}

// [[Rcpp::export]]
arma::mat bar(arma::mat A) {
    arma::uvec rows;
    rows << 0 << 2;
    arma::uvec cols;
    cols << 0 << 2;
    A.submat(rows, cols).fill(2.0);
    return A;
}

/*** R
A <- matrix(0, 4, 4)
foo(A)
bar(A)
*/

输出:

> A <- matrix(0, 4, 4)

> foo(A)
     [,1] [,2] [,3] [,4]
[1,]    1    1    0    0
[2,]    1    1    0    0
[3,]    0    0    0    0
[4,]    0    0    0    0

> bar(A)
     [,1] [,2] [,3] [,4]
[1,]    2    0    2    0
[2,]    0    0    0    0
[3,]    2    0    2    0
[4,]    0    0    0    0
于 2019-05-31T10:44:02.067 回答