如果添加print语句,相同的Rcpp函数返回不同的输出

Same Rcpp function returns different output if print statement added

本文关键字:函数 返回 输出 Rcpp 添加 语句 如果 print      更新时间:2023-10-16

我使用Rcpp编写的c++函数根据代码中是否有RcoutRprintf语句给出不同的输出。下面的代码返回1,这是正确的值,对于具有打印语句H_sigma_1()的函数。但是,对于没有print语句的H_sigma_2()函数,该函数返回2。我在Ubuntu 16.04.1和CentOS 6.8上进行了测试。不过,我没能在Windows 10上重现这个bug。因此,这似乎是一个Linux问题。

Rcpp代码
library(Rcpp)
cppFunction ( 
  "double H_sigma_1(IntegerVector sigma, NumericMatrix J, NumericVector h)
  {
    double first_sum, second_sum = 0;
    int n = sigma.size();
    for(int i = 0; i < n; i++)
    {
      for(int j = 0; j < n; j++)
      {
        // skip inside loop if i >= j to stop double counting
        if(i >= j) {continue;}
        first_sum += J(i, j) * sigma[i] * sigma[j];
        Rcout << first_sum << std::endl;
      }
      second_sum += h[i] * sigma[i];
    }
    return(-1.0 * first_sum - second_sum);
  }"
)
cppFunction (
  "double H_sigma_2(IntegerVector sigma, NumericMatrix J, NumericVector h)
  {
    double first_sum, second_sum = 0;
    int n = sigma.size();
    for(int i = 0; i < n; i++)
    {
      for(int j = 0; j < n; j++)
      {
        // skip inside loop if i >= j to stop double counting
        if(i >= j) {continue;}
        first_sum += J(i, j) * sigma[i] * sigma[j];
        // Rcout << first_sum << std::endl;
      }
      second_sum += h[i] * sigma[i];
    }
    return(-1.0 * first_sum - second_sum);
  }"
)

示例调用

设置:

n = 2
params = rep(1, n) 
h = rep(params[1], n)
J = toeplitz(c(0, params[2], rep(0, n - 2)))
测试1:

H_sigma_1(c(-1, -1), J, h)
输出:

1
[1] 1
测试2

H_sigma_2(c(-1, -1), J, h)

[1] 2

您遇到的问题是在使用求和之前没有显式声明它们的起始值。

double first_sum, second_sum = 0;

不一样
double first_sum = 0, second_sum = 0;

参见编译器警告标志:

file11d36f564b15.cpp:17:5: warning: variable 'first_sum' is uninitialized when used here [-Wuninitialized]
    first_sum += J(i, j) * sigma[i] * sigma[j];
    ^~~~~~~~~
file11d36f564b15.cpp:8:21: note: initialize the variable 'first_sum' to silence this warning
    double first_sum, second_sum = 0;

立即使用:

first_sum += J(i, j) * sigma[i] * sigma[j];

不设置first_sum = ...;


此外,另一个问题是:

second_sum = 0;

用整数值初始化双精度对象。虽然这个问题的范围很小,但要纠正这个问题,只需使用0.0而不是0

second_sum = 0.0;

这也适用于first_sum


修复以上问题的代码:

library(Rcpp)
cppFunction ( 
    "double H_sigma_1(IntegerVector sigma, NumericMatrix J, NumericVector h)
    {
    double first_sum = 0.0, second_sum = 0.0;
    int n = sigma.size();
    for(int i = 0; i < n; i++) {
      for(int j = 0; j < n; j++) {
        // skip inside loop if i >= j to stop double counting
        if(i >= j) {continue;}
        first_sum += J(i, j) * sigma[i] * sigma[j];
        Rcout << first_sum << std::endl;
      }
      second_sum += h[i] * sigma[i];
    }
    return(-1.0 * first_sum - second_sum);
    }"
)
cppFunction ( 
    "double H_sigma_2(IntegerVector sigma, NumericMatrix J, NumericVector h)
    {
    double first_sum = 0.0, second_sum = 0.0;
    int n = sigma.size();
    for(int i = 0; i < n; i++) {
      for(int j = 0; j < n; j++) {
        // skip inside loop if i >= j to stop double counting
        if(i >= j) {continue;}
        first_sum += J(i, j) * sigma[i] * sigma[j];
      }
      second_sum += h[i] * sigma[i];
    }
    return(-1.0 * first_sum - second_sum);
    }"
)

测试:

n = 2
params = rep(1, n) 
h = rep(params[1], n)
J = toeplitz(c(0, params[2], rep(0, n - 2)))
H_sigma_1(c(-1, -1), J, h)
输出:

1
[1] 1
测试2:

H_sigma_2(c(-1, -1), J, h)
输出:

[1] 1