大SpMat對象與RcppArmadillo

Question

我想學習和使用稀疏線性代數例程Rcpp和RcppArmadillo。下面

代碼是這裏的例子改編：http://gallery.rcpp.org/articles/armadillo-sparse-matrix/

code <- ' 
    S4 matx(x); 
    IntegerVector Xd = matx.slot("Dim"); 
    IntegerVector Xi = matx.slot("i"); 
    IntegerVector Xp = matx.slot("p"); 
    NumericVector Xx = matx.slot("x"); 

    arma::sp_mat Xsp(Xd[0], Xd[1]); 

    // create space for values, and copy 
    arma::access::rw(Xsp.values) = arma::memory::acquire_chunked<double>(Xx.size() + 1); 
    arma::arrayops::copy(arma::access::rwp(Xsp.values), 
      Xx.begin(), 
      Xx.size() + 1); 

    // create space for row_indices, and copy -- so far in a lame loop 
    arma::access::rw(Xsp.row_indices) = arma::memory::acquire_chunked<arma::uword>(Xx.size() + 1); 
    for (int j=0; j<Xi.size(); j++) 
    arma::access::rwp(Xsp.row_indices)[j] = Xi[j]; 

    // create space for col_ptrs, and copy -- so far in a lame loop 
    arma::access::rw(Xsp.col_ptrs) = arma::memory::acquire_chunked<arma::uword>(Xp.size() + 1); 
    for (int j=0; j<Xp.size(); j++) 
     arma::access::rwp(Xsp.col_ptrs)[j] = Xp[j]; 

    // important: set the sentinel as well 
    arma::access::rwp(Xsp.col_ptrs)[Xp.size()+1] = std::numeric_limits<arma::uword>::max(); 

    // set the number of non-zero elements 
    arma::access::rw(Xsp.n_nonzero) = Xx.size(); 

    Rcout << "SpMat Xsp:\\n" << arma::dot(Xsp,Xsp) << std::endl; 
' 

norm2 <- cxxfunction(signature(x="Matrix"), 
        code,plugin="RcppArmadillo") 

當我使用1E4的載體，東西很好地工作：

> p <- 10000 
> X <- Matrix(rnorm(p),sparse=TRUE) 
> norm2(X) 
SpMat Xsp: 
9997.14 
NULL 

然而，當我使用長度1E5的矢量，錯誤產生

> p <- 100000 
> X <- Matrix(rnorm(p),sparse=TRUE) 
> norm2(X) 

error: SpMat::init(): requested size is too large 

Error: 
> 

我似乎無法弄清楚我是什麼做錯了。 任何指針，將不勝感激。

==============更多信息==============

這個問題似乎是與具有尺寸> = 2^16 = 65536

以下工作：

> m <- 1000 
> n <- 65535 
> nnz <- 10000 
> iind <- sample.int(m,nnz,replace=TRUE) 
> jind <- sample.int(n,nnz,replace=TRUE) 
> xval <- rnorm(nnz) 
> X <- sparseMatrix(i=iind,j=jind,x=xval,dims=c(m,n)) 
> norm2(X) 
SpMat Xsp: 
10029.8 
NULL 

繼不起作用：

> m <- 1000 
> n <- 65536 
> nnz <- 10000 
> iind <- sample.int(m,nnz,replace=TRUE) 
> jind <- sample.int(n,nnz,replace=TRUE) 
> xval <- rnorm(nnz) 
> X <- sparseMatrix(i=iind,j=jind,x=xval,dims=c(m,n)) 
> norm2(X) 

error: SpMat::init(): requested size is too large 

Error: 
> 

爲什麼會是這種情況？

Answer 1

您的矩陣看起來很奇怪。通過說

Matrix(rnorm(p),sparse=TRUE) 

你最終得到p×1矩陣，雖然稀疏。如果我只分配10行或colums的東西 只是工作。

R> p <- 100000 
R> X <- Matrix(rnorm(p),nrow=10,sparse=TRUE) 
R> dim(X) 
[1] 10 10000 
R> norm2(X) 
SpMat Xsp: 
100832 
NULL 
R> 

因此，我認爲你只需要一個更好的稀疏矩陣一起工作 - 轉換代碼，犰狳的稀疏矩陣類型，都很好。

更新於2013-04-30：這實際上是一個犰狳錯誤，它剛剛固定在上游。一個新的RcppArmadillo版本0.3.810.2現在在SVN中，很快應該很快遷移到CRAN。您不再需要定義ARMA_64BIT_WORD。

Answer 2

我遇到了這個：http://arma.sourceforge.net/docs.html#config_hpp。

該溶液是通過將一行中的頭文件上述犰狳到64位整數選項設置：

#define ARMA_64BIT_WORD 

我把這個在[R根]/LIB/R /庫/ RcppArmadillo /包括/RcppArmadilloConfig.h

使用我曾試圖「包括」爲cxxfunction參數，但是我覺得這個定義語句是包括上述聲明，

#include RcppArmadillo.h 

在CPP代碼

。我不知道內聯軟件包的cxxfunction是否允許這樣做。

修改RcppArmadilloConfig.h後，我現在可以聲明一個大於2^16的矩陣。

> m <- 1000 
> n <- 65536+1000 
> nnz <- 10000 
> iind <- sample.int(m,nnz,replace=TRUE) 
> jind <- sample.int(n,nnz,replace=TRUE) 
> xval <- rnorm(nnz) 
> X <- sparseMatrix(i=iind,j=jind,x=xval,dims=c(m,n)) 
> norm2(X) 
SpMat Xsp: 
10218.8 
NULL 
>