[Rcpp-devel] Pointer Problem ?

bbonit at tin.it bbonit at tin.it
Fri Jun 22 21:54:39 CEST 2012


i have another problem maby cause by my inexperience   in the following code in Rccp inline  in the middle of the double cicle  (just after  ...
 for (j=0; j<n; j++) {....    in the abiut middle of code). maby the code is imprecise but i need to undertand (if it is possible) why this troule aris
there is an assegnation   th1=th0  but i have seen that it cames in mathematical sense : th0= th1   is the same of th1 =th0  but i dont want it becaus  for (j=0; j<n; j++) {
e this make for all iteratin exp(f1-f0) = 0 (except the number 1) because th0 came always the dame th1    like that th0 is the same of th1 ..
and so vanish the metropolis ratio....this not happen in R "classic"    i think it is a problem about pointer  but im at the startin study c++ .Thank you and sorry for noise.....
The folloing code

require(inline)
code <- '
RNGScope scope;
NumericVector start_= as<NumericVector>(start);
int n=Rf_length(start_);
int N = as<int>(nsim);
NumericVector tune = as<NumericVector>(scale);
NumericMatrix U(N,n);
NumericVector init (n);
NumericVector th0 (n);
NumericVector th1 (n);
NumericMatrix store(N,n);
NumericVector acount (n);
NumericVector f0;
NumericVector f1;
Function ltd(fun);
int i,j;
init = start_;
bool u;
for( i=0; i<n; i++ ){
th0[i]   = init[i];
                    }
f0=ltd(th0);  

for (i=0; i<N; i++) {
  for (j=0; j<n; j++) {
     th1=th0 ;
     th1[j] = Rf_rnorm(th0[j],tune[j]);
     f1=ltd(th1); 
      u= Rf_runif(0.0,1.0) <= exp(f1[0] - f0[0]);
         if(u) {th0[j] = th1[j];
                   f0 = f1;
                   acount[j]++;
                   }
       store(i,j)=th0[j];
                   }
}
return  List::create (Named("sim") = store, Named("acount") = acount,th0,th1,U);
'

AM<- cxxfunction(signature(fun="function",start="numeric",nsim="numeric",scale="numeric"),
code,  plugin="Rcpp")
 
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.r-forge.r-project.org/pipermail/rcpp-devel/attachments/20120622/0d960ed1/attachment.html>


More information about the Rcpp-devel mailing list