On 10 December 2014 at 14:54, Maxime To wrote: | Ok, thanks, I'll try to build on it. In the example I tried to isolate the problem, but in my real program I have lot of other matrix step using armadillo, that's why I put it in that way... I'd like to avoid armadillo, but it makes matrix calculus messier... | | About the memory, it does not seem to work, when I run a long loop the program just crash because of full memory.
The Writing R Extensions manual has some material on how to use tools like valgrind. Dirk | | -----Message d'origine----- | De : "Romain François" <[email protected]> | Envoyé : 10/12/2014 14:32 | À : "Maxime To" <[email protected]> | Cc : "Dirk Eddelbuettel" <[email protected]>; "[email protected]" <[email protected]> | Objet : Re: [Rcpp-devel] Rcpp Parallel and Rcpp Armadillo | | Some pointers. | | | When you use an arma::mat passed by value in an Rcpp::export, this means copying all of the data of the underlying R object into armadillo. I’d suggest you use a reference to const to avoid that, i.e. | | | mat contrib1(const mat& X1) { … } | | | Then in pQnorm, you do: | | | NumericMatrix x_q = Rcpp::as<Rcpp::NumericMatrix>(wrap(xx_q)); | | | That is yet again, copying all of the data from the arma::mat into an Rcpp matrix. | | | You then return a arma::mat, which data is copied implicitly as the return of contrib1. | | | I’d suggest you do all this without armadillo, which you don’t really use except for inducing a lot of extra copies of data. | | | To anwser your last question, R uses a garbage collector, so the memory is not automatically reclaimed as soon as it is no longer needed. | | | Hope this helps. | | | Romain | | | Le 10 déc. 2014 à 15:01, Maxime To <[email protected]> a écrit : | | | Hi, | | I changed the function as indicated by Dirk and I modify the functions and the program does work now. | However, I am still puzzled by the memory use of the program. when I run a loop of my function in R as in the code below, it seems that the program does not free the memory used in the previous iterations... which is annoying when I need to optimize on my final object. | | So I was wondering whether it was a question of declaration of object in my code? | | ------------------------------------------------------------------------------------------------------------------ | | sourceCpp("Rcpp/test.cpp") # | qwe = matrix(runif(10000), nrow = 100) | a = contrib1(qwe) | b = qnorm(qwe) | a - b | | for (i in 1:20000) a = contrib1(qwe) | ---------------------------------------------------------- | // test.cpp | | #include <RcppArmadillo.h> | #include <cmath> | #include <algorithm> | #include <RcppParallel.h> | #include <boost/math/distributions/inverse_gaussian.hpp> | | using namespace Rcpp; | using namespace arma; | using namespace std; | using namespace RcppParallel; | | // [[Rcpp::depends(RcppArmadillo, RcppParallel, BH)]] | | double qnorm_f(const double& x_q) { | boost::math::normal s; | return boost::math::quantile(s, x_q); | }; | | | | struct Qnorm : public Worker | { | // source matrix | const RMatrix<double> input_q; | | // destination matrix | RMatrix<double> output_q; | | // initialize with source and destination | Qnorm(const NumericMatrix input_q, NumericMatrix output_q) | : input_q(input_q), output_q(output_q) {} | | // take the Pnorm of the range of elements requested | void operator()(std::size_t begin, std::size_t end) { | std::transform(input_q.begin() + begin, | input_q.begin() + end, | output_q.begin() + begin, | ::qnorm_f); | } | }; | | mat pQnorm(mat xx_q) { | | NumericMatrix x_q = Rcpp::as<Rcpp::NumericMatrix>(wrap(xx_q)); | | // allocate the output matrix | const NumericMatrix output_q(x_q.nrow(), x_q.ncol()); | | // Pnorm functor (pass input and output matrices) | Qnorm qnorm_temp(x_q, output_q); | | // call parallelFor to do the work | parallelFor(0, x_q.length(), qnorm_temp); | | // return the output matrix | mat outmat_q(output_q.begin(), output_q.nrow(),output_q.ncol()); | return outmat_q; | | } | | // [[Rcpp::export]] | mat contrib1(mat X1) { | | mat test = pQnorm(X1); | mat results = test; | | return results; | } | | ---------------------------------------------------------- | | | > Date: Tue, 9 Dec 2014 09:07:10 -0600 | > To: [email protected] | > CC: [email protected]; [email protected] | > Subject: Re: [Rcpp-devel] Rcpp Parallel and Rcpp Armadillo | > From: [email protected] | > | > | > On 9 December 2014 at 09:46, Qiang Kou wrote: | > | What do you mean by "doesn't work" ? Compiling error or the result is not | > | right? | > | | > | I just tried the code, and it seems the code can compile and work. | > | > I am generally very careful about calling back to anything related to R from | > functions to be parallelized. So for | > | > inline double f(double x) { return ::Rf_pnorm5(x, 0.0, 1.0, 1, 0); } | > | > I think going with an equivalent pnorm() function from Boost / Bh may be better. | > | > But I am shooting from my hip here as I have not had time to look at this, | > having been out way too late at a nice concert :) | > | > Dirk | > | > -- | > http://dirk.eddelbuettel.com | @eddelbuettel | [email protected] | | _______________________________________________ | Rcpp-devel mailing list | [email protected] | https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/rcpp-devel | _______________________________________________ | Rcpp-devel mailing list | [email protected] | https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/rcpp-devel -- http://dirk.eddelbuettel.com | @eddelbuettel | [email protected] _______________________________________________ Rcpp-devel mailing list [email protected] https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/rcpp-devel
