1

我正在尝试移植一个函数，该函数使用 RcppEigen 将稀疏矩阵中的所有非零条目转换为与 RcppParallel 并行运行，但我无法让它工作。

即使用此来源：

#include <RcppEigen.h>
#include <RcppParallel.h>

using namespace Rcpp;

typedef Eigen::SparseMatrix<double> SpMat;
typedef Eigen::Map<SpMat> MSpMat;

// [[Rcpp::export]]
S4 serial_test( S4 m ) {
    SpMat src = as<MSpMat>( m );
    SpMat tgt = as<MSpMat>( clone( m ) );
    for( int i = 0; i < src.rows(); ++i ) {
        for( SpMat::InnerIterator srcIt( src, i ), tgtIt( tgt, i ); srcIt; ++srcIt, ++tgtIt ) {
            tgtIt.valueRef() = -1;
        }
    }
    S4 out = wrap( tgt );
    return out;
}

struct MyWorker : public RcppParallel::Worker {
    SpMat src;
    SpMat tgt;

    MyWorker( const SpMat& src, SpMat& tgt )
        : src( src ), tgt( tgt ) {}

    void operator()( std::size_t begin, std::size_t end ) {
        for( int i = begin; i < end; ++i ) {
            for( SpMat::InnerIterator srcIt( src, i ), tgtIt( tgt, i ); srcIt; ++srcIt, ++tgtIt ) {
                tgtIt.valueRef() = -1;
            }
        }
    }
};

// [[Rcpp::export]]
S4 para_test( S4 m ) {
    SpMat src = as<MSpMat>( m );
    SpMat tgt = as<MSpMat>( clone( m ) );
    MyWorker wrkr( src, tgt );
    RcppParallel::parallelFor( 0, src.outerSize(), wrkr );
    S4 out = wrap( tgt );
    return out;
}

这个输入：

> m <- Matrix::rsparsematrix( 1000, 1000, 0.1 )
> m[ 1:10, 1:10 ]
10 x 10 sparse Matrix of class "dgCMatrix"

 [1,] . . . . 0.52 . . .     . .
 [2,] . . . . 0.59 . . .     . .
 [3,] . . . . .    . . 0.47  . .
 [4,] . . . . .    . . 0.25  . .
 [5,] . . . . .    . . .     . .
 [6,] . . . . .    . . .     . .
 [7,] . . . . .    . . .     . .
 [8,] . . . . .    . . .     . .
 [9,] . . . . .    . . .    -1 .
[10,] . . . . .    . . .     . .

调用 serial_test 有效：

> m1 <- serial_test( m )
> m1[ 1:10, 1:10 ]
10 x 10 sparse Matrix of class "dgCMatrix"

 [1,] . . . . -1 . .  .  . .
 [2,] . . . . -1 . .  .  . .
 [3,] . . . .  . . . -1  . .
 [4,] . . . .  . . . -1  . .
 [5,] . . . .  . . .  .  . .
 [6,] . . . .  . . .  .  . .
 [7,] . . . .  . . .  .  . .
 [8,] . . . .  . . .  .  . .
 [9,] . . . .  . . .  . -1 .
[10,] . . . .  . . .  .  . .

但是 para_test 什么也没做，即使它运行良好：

> m2 <- para_test( m )
> m2[ 1:10, 1:10 ]
10 x 10 sparse Matrix of class "dgCMatrix"

 [1,] . . . . 0.52 . . .     . .
 [2,] . . . . 0.59 . . .     . .
 [3,] . . . . .    . . 0.47  . .
 [4,] . . . . .    . . 0.25  . .
 [5,] . . . . .    . . .     . .
 [6,] . . . . .    . . .     . .
 [7,] . . . . .    . . .     . .
 [8,] . . . . .    . . .     . .
 [9,] . . . . .    . . .    -1 .
[10,] . . . . .    . . .     . .

是否可以将 RcppParallel 与 Eigen 的 SparseMatrix 类一起使用？

4

0 回答 0