[Rcpp-devel] Best way to create random normal matrix

[Julian Smith]

Thanks for the feedback. This should be easy enough to port over to
RcppEigen, right?

On Tue, Dec 4, 2012 at 6:17 PM, Dirk Eddelbuettel <edd at debian.org> wrote:

>
> Julian,
>
> Your prorgam was actually fine as such, but you "used" it wrong -- what
> comes
> back from cxxfunction() is already an R function, see below.
>
> If you use the R RNGs (which is easy from Rcpp) your results are
> reproducible
> with R.  That is all documented so I won't repeat it here.
>
> Dirk
>
> R> library(inline)
> R>
> R> src <- '
> +     Rcpp::NumericMatrix Xr(Xs);
> +     int q = Rcpp::as<int>(ys);
> +     int n = Xr.nrow(), k = Xr.ncol();
> +     arma::mat X(Xr.begin(), n, k, false);
> +     arma::mat G, Y, B;
> +     G = arma::randn(n,q);
> +     Y = X*G;
> +     arma::mat Q, R;
> +     arma::qr(Q,R,Y);
> +     Q = Q(arma::span::all,arma::span(0,q-1));
> +     B = arma::trans(Q)*X;
> +     arma::vec d;
> +     arma::mat u, v;
> +     arma::svd(u, d, v, B );
> +     arma::mat U = Q*u;
> +     return Rcpp::List::create(Rcpp::Named("d")=d,
> +                               Rcpp::Named("U")=U,
> +                               Rcpp::Named("V")=v);'
> R>
> R> rsvd <- cxxfunction(signature(Xs="numeric", ys="integer"), body=src,
> plugin= "RcppArmadillo")
> R> rsvd(matrix(1:4,2,2), 2)
> $d
>          [,1]
> [1,] 5.464986
> [2,] 0.365966
>
> $U
>          [,1]      [,2]
> [1,] 0.576048 -0.817416
> [2,] 0.817416  0.576048
>
> $V
>          [,1]      [,2]
> [1,] 0.404554  0.914514
> [2,] 0.914514 -0.404554
>
> R>
>
>
> --
> Dirk Eddelbuettel | edd at debian.org | http://dirk.eddelbuettel.com
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.r-forge.r-project.org/pipermail/rcpp-devel/attachments/20121204/11822b99/attachment.html>