我正在通过 RcppArmadillo 使用 Armadillo 和 R。
我有一个稀疏矩阵和一个行号作为输入。我想搜索矩阵的相应行并返回非零值的位置。
到目前为止我的函数看起来像
// [[Rcpp::export]]
arma::uvec findAdjacentStates(arma::sp_mat adjacency, int row) {
arma::uvec out(adjacency.n_cols);
arma::SpSubview<double>::const_iterator start = adjacency.row(row).begin();
arma::SpSubview<double>::const_iterator end = adjacency.row(row).end();
for(arma::SpSubview<double>::const_iterator it = start; it != end; ++it)
{
Rcout << "location: " << it.row() << "," << it.col() << " ";
Rcout << "value: " << (*it) << std::endl;
}
return out;
}
这是基于之前的 SO answer .
函数使 R 崩溃。
require(Matrix)
x = rsparsematrix(10, 10, .2)
x = x > 1
x = as(x, "dgCMatrix")
findAdjacentStates(x, 1)
我不清楚的一件事是如何具体迭代行 vector ;之前的 SO 答案是针对稀疏矩阵进行迭代。
更新:根据 valgrind 的说法,问题出在 operator++ (SpSubview_iterators_meat.hpp:319) 中,因此这似乎不是迭代稀疏行 vector 的正确方法
最佳答案
迭代稀疏矩阵行的方法是使用 sp_mat::row_iterator .遗憾的是,无法提前知道输出 vector 的大小,而且 uvec
对象不像常规 vector 那样具有 push_back
。这是我的建议:
#include <RcppArmadillo.h>
// [[Rcpp::depends(RcppArmadillo)]]
using namespace Rcpp;
using namespace arma;
// [[Rcpp::export]]
IntegerVector findAdjacentStates(sp_mat adjacency, int row) {
IntegerVector out;
sp_mat::const_row_iterator start = adjacency.begin_row(row);
sp_mat::const_row_iterator end = adjacency.end_row(row);
for ( sp_mat::const_row_iterator i = start; i != end; ++i )
{
out.push_back(i.col());
}
return out;
}
我们可以很容易地对其进行测试:
# We need Rcpp and Matrix:
library(Rcpp)
library(Matrix)
# This is the C++ script I wrote:
sourceCpp('SOans.cpp')
# Make example data (setting seed for reproducibility):
set.seed(123)
x = rsparsematrix(10, 10, .2)
# And test the function:
findAdjacentStates(x, 1)
[1] 4
x
10 x 10 sparse Matrix of class "dgCMatrix"
[1,] . . 0.84 . 0.40 0.7 . . . -0.56
[2,] . . . . -0.47 . . . . .
[3,] . . . . . . . . -2.00 .
[4,] 0.15 . . . . . . . . -0.73
[5,] 1.80 . . . . . . . . .
[6,] . . . . . . . . 0.11 .
[7,] . . -1.10 . . . . -1.70 -1.10 .
[8,] . . . . . . . 1.30 . -0.22
[9,] -0.63 . 1.20 . . . . 0.36 . .
[10,] . . . . 0.50 -1.0 . . . .
所以,我们可以看到这是可行的;第 1 行(在 C++ 术语中;在 R 术语中是第 2 行)只有一个非零元素,它在第 4 列中(在 C++ 术语中;在 R 术语中是第 5 列)。如果您想将输出返回到 R,这应该可以工作。如果您想在另一个 C++ 函数中使用输出,根据您正在做的事情,您可能更喜欢使用 uvec
而不是 IntegerVector
,但您可以将 IntegerVector
转换为 uvec
(可能不是最有效的解决方案,但我想到的最好的解决方案现在):
// [[Rcpp::export]]
uvec findAdjacentStates(sp_mat adjacency, int row) {
IntegerVector tmp;
sp_mat::const_row_iterator start = adjacency.begin_row(row);
sp_mat::const_row_iterator end = adjacency.end_row(row);
for ( sp_mat::const_row_iterator i = start; i != end; ++i )
{
tmp.push_back(i.col());
}
uvec out = as<uvec>(tmp);
return out;
}
关于c++ - Armadillo :从稀疏矩阵中获取稀疏行 vector 的非零位置,我们在Stack Overflow上找到一个类似的问题: https://stackoverflow.com/questions/47147845/