r - 非连续矩阵的高级构造函数

Question

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

这是一个简单的示例，说明了当我的原始矩阵为 时我想要做什么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同时确保它们使用相同的记忆。

标签： rrcpparmadillorcpparmadillo

Answer 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