Hello, 

This stands out 

double U; 
  for(int i=0; i<n; i++){
        for(int j=0; j<n; j++){
                U = U + A_kl(i,j) * B_kl(i,j)/n/n ;
        }
  } 
  U = sqrt(U);
  
U is not initialized. so all bets are on. Anything can happen. Same holds for 
aa and bb. 


Also, in terms of how you use attributes, you can at your advantage replace 
things like this: 

//[[Rcpp::export]]
double mydcov(SEXP X, SEXP Y) {
   NumericVector x = as<NumericVector>(X);
   NumericVector y = as<NumericVector>(Y);
   …
}

by things like this: 

//[[Rcpp::export]]
double mydcov(NumericVector x, NumericVector y){
   …
}

Romain

Le 19 janv. 2014 à 14:09, 晔张 <zhang...@gmail.com> a écrit :

> Hello, everyone. 
> I guess the probrem maybe occur because of sourceCpp(), so I have made 
> another test. I make a simple package using 
> RcppArmadillo.package.skeletion(). 
> But the same problem appear.
> My code is on Github:
> https://github.com/ZhangYet/RcppTest
> Is there anything related to my compiler or OS? I use 64-bit CentOS and g++ 
> 4.4.5
> 
> 
> On Sat, Jan 18, 2014 at 2:01 AM, 张晔 <zhang...@gmail.com> wrote:
> Thanks, Dirk.
> I have upload my code to Github: https://github.com/ZhangYet/RcppTest/
> And I made a mistake:
> 
> 
> | res1 <- mydcov(x,y)/sqrt(mydcov(x,x)*mydcov(y,y))
> |
> | a <- mydcov(x,y)
> | b <- mydcov(x,x)
> | c <- mydcov(y,y)
> | res2 <- a/sqrt(b*c)
> |
> | and res2 is different from res1.(res2 is the right answer.)
> 
> I checked my test result. res2 is not the right answer, but it's different 
> from res.
> 
> 于 2014/1/17 21:36, Dirk Eddelbuettel 写道:
> 
> On 17 January 2014 at 20:51, 晔张 wrote:
> | Dear all,
> |
> | I have a strange question such that I can't conclude it in the title.
> | I have write an function:
> | //[[Rcpp::export]]
> | double mydcov(NumericVector x, NumericVector y)
> 
> It would help us to see the function to have any guess at what may be wrong.
>   | This function is computing  ditance covariance. This function is OK. It 
> give
> | the right answer.
> | I use this function to compute distance correltion.
> |
> | //[[Rcpp::export]]
> | double mydcor(NumericVector x, NumericVector y)
> | {
> |     return mydcov(x,y)/sqrt(mydcov(x,x)*mydcov(y,y));
> | }
> |
> | Then, I compile the source code using sourceCpp.
> | And then mydcor(x,y) give a wrong answer.
> |
> | Here comes the strange part:
> | In R, I run code as follow
> |
> | res1 <- mydcov(x,y)/sqrt(mydcov(x,x)*mydcov(y,y))
> |
> | a <- mydcov(x,y)
> | b <- mydcov(x,x)
> | c <- mydcov(y,y)
> | res2 <- a/sqrt(b*c)
> |
> | and res2 is different from res1.(res2 is the right answer.)
> |
> | And I have another test, rewriting the source code for mydcor():
> | version1:
> | //[[Rcpp::export]]
> | double mydcor(NumericVector x, NumericVector y)
> | {
> |     double a = mydcov(x,y)
> |     double b = mydcov(x,x);
> |     double c = mydcov(y,y);
> |     return (a/sqrt(b*c));
> | }
> |
> | version2:
> | //[[Rcpp::export]]
> | double mydcor(NumericVector x, NumericVector y)
> | {
> |     double a = mydcov(x,y)
> |     Rcout<<a<<"\n";
> |     double b = mydcov(x,x);
> |     Rcout<<b<<"\n";
> |     double c = mydcov(y,y);
> |     Rcout<<c<<"\n"
> |     return (a/sqrt(b*c));
> | }
> | This time, version2 make a right answer.
> | It's too weird. I havn't meet such situation before.
> 
> I am not sure we are in a position to help you here.
> 
> Regards, Dirk
> 
> |
> | Happy Lunar New Year (the most important festival in China), and Thanks in
> | advance.
> |
> | Best regards
> | ZhangYet
> |
> | ----------------------------------------------------------------------
> | _______________________________________________
> | Rcpp-devel mailing list
> | Rcpp-devel@lists.r-forge.r-project.org
> | https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/rcpp-devel
> 
> 
> _______________________________________________
> Rcpp-devel mailing list
> Rcpp-devel@lists.r-forge.r-project.org
> https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/rcpp-devel

_______________________________________________
Rcpp-devel mailing list
Rcpp-devel@lists.r-forge.r-project.org
https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/rcpp-devel

Reply via email to