Можно ли использовать RcppParallel с Eigen SparseMatrix?

Я пытаюсь портировать функцию, которая преобразует все ненулевые записи в sparsematrix, используя 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?

0 ответов

Другие вопросы по тегам