@@ -3,111 +3,36 @@ suppressMessages(library(Rcpp))
33suppressMessages(library(inline ))
44suppressMessages(library(rbenchmark ))
55
6- # # NOTE: This is the old way to compile Rcpp code inline.
7- # # The code here has left as a historical artifact and tribute to the old way.
8- # # Please use the code under the "new" inline compilation section.
9-
10- rcppGamma_old <- cxxfunction(signature(xs = " numeric" ), plugin = " Rcpp" , body = '
11- NumericVector x(xs);
12- int n = x.size();
13-
14- // Initialize Random number generator
15- RNGScope scope;
16-
17- const double y = 1.234;
18- for (int i=0; i<n; i++) {
19- x[i] = ::Rf_rgamma(3.0, 1.0/(y*y+4));
20- }
21-
22- // Return to R
23- return x;
24- ' )
25-
26-
27- gslGamma_old <- cxxfunction(signature(xs = " numeric" ), plugin = " RcppGSL" ,
28- include = ' #include <gsl/gsl_rng.h>
29- #include <gsl/gsl_randist.h>' ,
30- body = '
31- NumericVector x(xs);
32- int n = x.size();
33-
34- gsl_rng *r = gsl_rng_alloc(gsl_rng_mt19937);
35- const double y = 1.234;
36- for (int i=0; i<n; i++) {
37- x[i] = gsl_ran_gamma(r,3.0,1.0/(y*y+4));
38- }
39- gsl_rng_free(r);
40-
41- // Return to R
42- return x;
43- ' )
44-
45-
46- rcppNormal_old <- cxxfunction(signature(xs = " numeric" ), plugin = " Rcpp" , body = '
47- NumericVector x(xs);
48- int n = x.size();
49-
50- // Initialize Random number generator
51- RNGScope scope;
52-
53- const double y = 1.234;
54- for (int i=0; i<n; i++) {
55- x[i] = ::Rf_rnorm(1.0/(y+1),1.0/sqrt(2*y+2));
56- }
57-
58- // Return to R
59- return x;
60- ' )
61-
62-
63- gslNormal_old <- cxxfunction(signature(xs = " numeric" ), plugin = " RcppGSL" ,
64- include = ' #include <gsl/gsl_rng.h>
65- #include <gsl/gsl_randist.h>' ,
66- body = '
67- NumericVector x(xs);
68- int n = x.size();
69-
70- gsl_rng *r = gsl_rng_alloc(gsl_rng_mt19937);
71- const double y = 1.234;
72- for (int i=0; i<n; i++) {
73- x[i] = 1.0/(y+1)+gsl_ran_gaussian(r,1.0/sqrt(2*y+2));
74- }
75- gsl_rng_free(r);
76-
77- // Return to R
78- return x;
79- ' )
80-
81-
826# # NOTE: Within this section, the new way to compile Rcpp code inline has been
83- # # written. Please use the code next as a template for your own project.
7+ # # written. Please use the code next as a template for your own project, and
8+ # # do NOT use the old code below
849
8510cppFunction('
8611NumericVector rcppGamma(NumericVector x){
8712 int n = x.size();
88-
13+
8914 const double y = 1.234;
9015 for (int i=0; i<n; i++) {
9116 x[i] = R::rgamma(3.0, 1.0/(y*y+4));
9217 }
93-
18+
9419 // Return to R
9520 return x;
9621}' )
9722
98- # # This approach is a bit sloppy. Generally, you will want to use
23+ # # This approach is a bit sloppy. Generally, you will want to use
9924# # sourceCpp() if there are additional includes that are required.
10025cppFunction('
10126NumericVector gslGamma(NumericVector x){
10227 int n = x.size();
103-
28+
10429 gsl_rng *r = gsl_rng_alloc(gsl_rng_mt19937);
10530 const double y = 1.234;
10631 for (int i=0; i<n; i++) {
10732 x[i] = gsl_ran_gamma(r,3.0,1.0/(y*y+4));
10833 }
10934 gsl_rng_free(r);
110-
35+
11136 // Return to R
11237 return x;
11338}' , includes = ' #include <gsl/gsl_rng.h>
@@ -118,18 +43,18 @@ NumericVector gslGamma(NumericVector x){
11843cppFunction('
11944NumericVector rcppNormal(NumericVector x){
12045 int n = x.size();
121-
46+
12247 const double y = 1.234;
12348 for (int i=0; i<n; i++) {
12449 x[i] = R::rnorm(1.0/(y+1),1.0/sqrt(2*y+2));
12550 }
126-
51+
12752 // Return to R
12853 return x;
12954}' )
13055
13156
132- # # Here we demonstrate the use of sourceCpp() to show the continuity
57+ # # Here we demonstrate the use of sourceCpp() to show the continuity
13358# # of the code artifact.
13459
13560sourceCpp(code = '
@@ -144,14 +69,14 @@ using namespace Rcpp;
14469// [[Rcpp::export]]
14570NumericVector gslNormal(NumericVector x){
14671 int n = x.size();
147-
72+
14873 gsl_rng *r = gsl_rng_alloc(gsl_rng_mt19937);
14974 const double y = 1.234;
15075 for (int i=0; i<n; i++) {
15176 x[i] = 1.0/(y+1)+gsl_ran_gaussian(r,1.0/sqrt(2*y+2));
15277 }
15378 gsl_rng_free(r);
154-
79+
15580 // Return to R
15681 return x;
15782}' )
@@ -167,3 +92,85 @@ res <- benchmark(rcppGamma(x),
16792print(res )
16893
16994
95+
96+
97+
98+ # #
99+ # #
100+ # # Old code below. Do not use in new projects, it is here solely for comparison
101+ # #
102+ # #
103+
104+
105+ # # NOTE: This is the old way to compile Rcpp code inline.
106+ # # The code here has left as a historical artifact and tribute to the old way.
107+ # # Please use the code under the "new" inline compilation section.
108+
109+ rcppGamma_old <- cxxfunction(signature(xs = " numeric" ), plugin = " Rcpp" , body = '
110+ NumericVector x(xs);
111+ int n = x.size();
112+
113+ RNGScope scope; // Initialize Random number generator. Not needed when Attributes used.
114+
115+ const double y = 1.234;
116+ for (int i=0; i<n; i++) {
117+ x[i] = ::Rf_rgamma(3.0, 1.0/(y*y+4));
118+ }
119+
120+ // Return to R
121+ return x;
122+ ' )
123+
124+
125+ gslGamma_old <- cxxfunction(signature(xs = " numeric" ), plugin = " RcppGSL" ,
126+ include = ' #include <gsl/gsl_rng.h>
127+ #include <gsl/gsl_randist.h>' ,
128+ body = '
129+ NumericVector x(xs);
130+ int n = x.size();
131+
132+ gsl_rng *r = gsl_rng_alloc(gsl_rng_mt19937);
133+ const double y = 1.234;
134+ for (int i=0; i<n; i++) {
135+ x[i] = gsl_ran_gamma(r,3.0,1.0/(y*y+4));
136+ }
137+ gsl_rng_free(r);
138+
139+ // Return to R
140+ return x;
141+ ' )
142+
143+
144+ rcppNormal_old <- cxxfunction(signature(xs = " numeric" ), plugin = " Rcpp" , body = '
145+ NumericVector x(xs);
146+ int n = x.size();
147+
148+ RNGScope scope; // Initialize Random number generator. Not needed when Attributes used.
149+
150+ const double y = 1.234;
151+ for (int i=0; i<n; i++) {
152+ x[i] = ::Rf_rnorm(1.0/(y+1),1.0/sqrt(2*y+2));
153+ }
154+
155+ // Return to R
156+ return x;
157+ ' )
158+
159+
160+ gslNormal_old <- cxxfunction(signature(xs = " numeric" ), plugin = " RcppGSL" ,
161+ include = ' #include <gsl/gsl_rng.h>
162+ #include <gsl/gsl_randist.h>' ,
163+ body = '
164+ NumericVector x(xs);
165+ int n = x.size();
166+
167+ gsl_rng *r = gsl_rng_alloc(gsl_rng_mt19937);
168+ const double y = 1.234;
169+ for (int i=0; i<n; i++) {
170+ x[i] = 1.0/(y+1)+gsl_ran_gaussian(r,1.0/sqrt(2*y+2));
171+ }
172+ gsl_rng_free(r);
173+
174+ // Return to R
175+ return x;
176+ ' )
0 commit comments