Come chiamare la funzione definita dall'utente in RcppParallel?

Question

Ispirato allo snodo http://gallery.rcpp.org/articles/parallel-distance-matrix/, provo a utilizzare RcppParallel per eseguire la ricerca di forza bruta nello spazio parametrico ad alta dimensionalità per il backtesting utilizzando i multithread. Sono bloccato in come chiamare una funzione auto-definita nella parte struct. L'idea è così:

Innanzitutto, creare una matrice parametrica NumericMatrix params_mat in R prima, e utilizzare i dati backtesting con List, NumericVector, CharacterVector tipo di dati, come List Data_1, NumericVector Data_2, CharacterVector Data_3, ..., che sono statici per ogni scenario parametrico params_vec (Si noti che è la fila di params_mat).

Successivamente, definire la funzione di backtesting che genera un vettore che consiste di 3 variabili chiave per valutare le prestazioni della strategia.

Ecco un esempio del mio params_mat e Backtesting_Fun che può essere eseguito in R e Rcpp, rispettivamente.

//[[Rcpp::export]] 
NumericMatrix data_frame_rcpp(const Rcpp::List& list_params) 
{ 
    NumericMatrix res = list_params[0]; 
    return res; 
} 

# R codes to generate params_mat 
params <- expand.grid (x_1=seq(1,100,1), x_2=seq(3,100,2), ..., x_n=seq(4,200,1));       
list_params = list(ts(params)); 
tmp_params_data = data_frame_rcpp(list_params);            
params_mat = matrix(tmp_params_data, ncol = ncol(tmp_params_data), dimnames = NULL); 
params_vec = params_mat[ii,]; 

# User-defined Rcpp codes for backtesting 
NumericVector Backtesting_Fun (List Data_1, NumericVector Data_2, CharacterVector Data_3, ..., NumericVector params_vec) 
{ 
    // Main function parts to run backtesting for each params_vec scenario. 
    ... etc 

    // save 3 key result variables together with each params_vec (just a simple illustration). 
    NumericVector res = NumericVector::create(params_vec[0],...,params_vec[size-1], 
              key_1, key_2, key_3); 
    return res; 
} 

Certamente abbiamo bisogno di riscrivere/modificare l'originale Rcpp Backtesting_Fun con i tipi RVector/RMatrix, e quindi utilizzare i seguenti RcppParallel codici di chiamare Backtesting_Fun in struct Backtest_parallel:

// [[Rcpp::depends(RcppParallel)]] 
#include <RcppParallel.h> 
using namespace RcppParallel; 

RVector<double> Backtesting_Fun (const RVector<double> Data_1, const RVector<double> Data_2, 
           const RVector<string> Data_3,..., const RVector<double> params_vec) 
{ 
    // Main function parts to run backtesting for each params_vec scenario. 
    ... etc; 

    // save 3 key result variables together with each params_vec 
    ... etc; 

    return res; 
} 

struct Backtest_parallel : public Worker 
{  
    // input matrix to read from 
    const RVector<List> Data_1; 
    const RVector<double> Data_2; 
    const RVector<string> Data_3; 
    ... 
    const RMatrix<double> params_mat; 

    // output matrix to write to 
    RMatrix<double> rmat; 

    // initialize from Rcpp input and output matrixes (the RMatrix class 
    // can be automatically converted to from the Rcpp matrix type) 
    Backtest_parallel(const List Data_1, const NumericVector Data_2, 
    const CharacterVector Data_3, ..., const NumericMatrix params_mat) 
     : Data_1(Data_1), Data_2(Data_2), Data_3(Data_3), ..., params_mat(params_mat) {} 

    // function call operator that work for the specified range (begin/end) 
    void operator()(std::size_t begin, std::size_t end) 
    { 
     for (std::size_t ii = begin; ii < end; i++) 
     { 
     // params rows that we will operate on 
     RMatrix<double>::Row params_row = params_mat.row(ii); 

     // Run the backtesting function defined above 
     RVector<double> res = Backtesting_Fun(Data_1, Data_2, ..., params_row) 
     for (std::size_t jj = 0; jj < res.length(); jj++) 
     { 
      // write to output matrix 
      rmat(ii,jj) = res[jj]; 
     } 
     } 
    } 
}; 

// [[Rcpp::export]] 
NumericMatrix rcpp_parallel_backtest(List Data_1, NumericVector Data_2, CharacterVector Data_3, 
            ..., NumericMatrix params_mat) 
{  
    // allocate the matrix we will return 
    NumericMatrix rmat(params_mat.nrow(), params_mat.nrow()+3); 

    // create the worker 
    Backtest_parallel backtest_parallel(Data_1, Date_2, ..., params_mat); 

    // call it with parallelFor 
    parallelFor(0, rmat.nrow(), backtest_parallel); 

    return rmat; 
} 

Ecco le mie domande:

1. Can RVector contiene List tipo di dati o c'è un contenitore specifico in RcppParallel per contenere List;

2. Nel Backtesting_Fun, l'ingresso dovrebbe essere RVector/RMatrix tipi, significa che abbiamo davvero bisogno di convertire i codici orginal principali Rcpp con NumericVector in RVector?

Oppure c'è un modo migliore per eseguire il calcolo parallelo per il mio caso in RcppParallel? Grazie in anticipo.

EDIT:

1. Guardo gli altri esempi per quanto riguarda RcppPararrel in http://gallery.rcpp.org/articles/parallel-matrix-transform/, http://gallery.rcpp.org/articles/parallel-inner-product/, l'idea comune in struct operator() è quello di utilizzare i puntatori per manipolare l'input dei dati per operator(), quindi c'è qualche modo costruire una funzione definita dall'utente nel mio caso con gli input del puntatore?

2. Se il modo precedente non funziona, è possibile l'uso wrap convertire RVector/RMatrix nuovamente dentro Rcpp tipo di dati, cioè, in NumericVector..operator() modo che i tipi di ingresso della funzione definita dall'utente Backtesting_Fun può rimanere invariata.

Answer 1

Credo di poter trovare un modo alternativo per risolvere questo problema: I tasti sono di utilizzare un filo di accesso sicuri per contenere variabili all'interno del struct, e rimanere RVector/RMatrix nella funzione principale esterno in modo che il parallelFor può funzionare bene, che è la parte più importante in questo algo parallelo.Qui ci sono le mie vie:

1. Sbarazzarsi della List tipo di dati: invece, siamo in grado di convertire la variabile List utilizzando NumericVector/NumericMatrix contenitore e registrare il suo indice corrispondente in modo che il subvector/sottomatrice punterà gli stessi elementi come se fosse un elemento di una lista.

2. Converti RVector/RMatrix in arma::vec/arma::mat: Come accennato in RcppParallel Github, C++ Armadillo sono thread-safe in operatore di struct. Qui, modifico l'esempio dato in Parallel Distance Matrix Calculation with RcppParallel usando questa idea, che rimane quasi la stessa velocità di test.

   struct JsDistance : public Worker 
   { 
       const RMatrix<double> tmp_MAT; // input matrix to read from 
       RMatrix<double> tmp_rmat;  // output matrix to write to 
       std::size_t row_size, col_size; 
   
       // Convert global input/output into RMatrix/RVector type 
       JsDistance(const NumericMatrix& matrix_input, NumericMatrix& matrix_output, 
         std::size_t row_size, std::size_t col_size) 
       : tmp_MAT(matrix_input), tmp_rmat(matrix_output), row_size(row_size), col_size(col_size) {} 
   
       // convert RVector/RMatrix into arma type for Rcpp function 
       // and the follwing arma data will be shared in parallel computing 
       arma::mat convert() 
       { 
       RMatrix<double> tmp_mat = tmp_MAT; 
       arma::mat MAT(tmp_mat.begin(), row_size, col_size, false); 
       return MAT; 
       } 
   
   
       void operator()(std::size_t begin, std::size_t end) 
       { 
       for (std::size_t i = begin; i < end; i++) 
       { 
        for (std::size_t j = 0; j < i; j++) 
        { 
        // rows we will operate on 
        arma::mat MAT = convert(); 
        arma::rowvec row1 = MAT.row(i);   // get the row of arma matrix 
        arma::rowvec row2 = MAT.row(j); 
   
        // compute the average using std::tranform from the STL 
        std::vector<double> avg(row1.n_elem); 
        std::transform(row1.begin(), row1.end(), // input range 1 
            row2.begin(),    // input range 2 
            avg.begin(),    // output range 
            average);     // function to apply 
   
        // calculate divergences 
        double d1 = kl_divergence(row1.begin(), row1.end(), avg.begin()); 
        double d2 = kl_divergence(row2.begin(), row2.end(), avg.begin()); 
   
        // write to output matrix 
        tmp_rmat(i,j) = sqrt(.5 * (d1 + d2)); 
        } 
       } 
       } 
   }; 
   
   // [[Rcpp::export]] 
   NumericMatrix rcpp_parallel_js_distance_modify(const Rcpp::NumericMatrix& matrix_input, int N_cores) 
   { 
       // allocate the matrix we will return 
       NumericMatrix matrix_output(matrix_input.nrow(), matrix_input.nrow()); 
       std::size_t row_size = matrix_input.nrow(); 
       std::size_t col_size = matrix_input.ncol(); 
   
       // create the worker 
       JsDistance jsDistance(matrix_input, matrix_output, row_size, col_size); 
   
       // call it with parallelFor 
       parallelFor(0, matrix_input.nrow(), jsDistance, matrix_input.nrow()/N_cores);   // parallelFor with grain size setting 
   
       return matrix_output; 
   } 
   
   // Example compare: 
   n_row = 1E3; 
   n_col = 1E2; 
   m = matrix(runif(n_row*n_col), nrow = n_row, ncol = n_col); 
   m = m/rowSums(m); 
   
   res <- benchmark(rcpp_parallel_js_distance(m, 6), 
        rcpp_parallel_js_distance_orignal(m, 6), 
        order="relative") 
   res[,1:4]; 
   
   #test         #elapsed #relative 
   rcpp_parallel_js_distance_orignal(m, 6) 128.069 1.000 
   rcpp_parallel_js_distance(m, 6)   129.210 1.009 
   

Come possiamo vedere, il tipo di dati all'interno di operator sarà C++ arma, e ora possiamo tranquillamente e velocemente chiamare il nostro funzione definita dall'utente utilizzando direttamente l'oggetto anziché solo i puntatori, che non può essere generico o facilmente progettato.

Ora, questa struttura parallelFor condividerà la stessa fonte di dati senza una copia extra nel calcolo parallelo, e quindi potremo apportare qualche lieve cambiamento per il backtest utilizzando l'idea menzionata nella domanda precedente.