r语言 - Rcpp 继续为一项看似简单的任务运行

Question

我已经想了一整天，仍然无法弄清楚为什么会发生这种情况。我的目标很简单：STEP1，生成一个函数 S（h，p）;STEP2，通过梯形规则对S（h，p）相对于p进行数值积分，得到一个新的函数SS（h）。我编写了代码并通过sourceCpp获取源代码，它成功地在R中生成了两个函数S（h，p）和SS（h）。但是当我尝试通过计算 SS（1） 来测试它时，R 只是继续运行并且从未给出结果，这很奇怪，因为计算量不是那么大。知道为什么会这样吗？

我的代码在这里：

#include <Rcpp.h>
using namespace Rcpp;
//generate the first function that gives S(h,p)           
// [[Rcpp::export]]
double S(double h, double p){
  double out=2*(h+p+h*p);
  return out;
}
//generate the second function that gives the numerically integreation of S(h,p) w.r.t p
//[[Rcpp::export]]
double SS(double h){
  double out1=0;
  double sum=0;
  for (int i=0;i<1;i=i+0.01){
    sum=sum+S(h,i);
  }
  out1=0.01/2*(2*sum-S(h,0)-S(h,1));
  return out1;
}

Answer 1

问题是你把i当作不是这个语句中的int：

for (int i=0;i<1;i=i+0.01){
    sum=sum+S(h,i);
}

每次迭代后，您尝试将0.01添加到整数中，该整数当然会立即截断为 0，这意味着i始终等于零，并且您有一个无限循环。一个突出问题的最小示例，有几个可能的解决方案：

#include <Rcpp.h>
// [[Rcpp::export]]
void bad_loop() {
    for (int i = 0; i < 1; i += 0.01) {
        std::printf("i = %dn", i);
        Rcpp::checkUserInterrupt();
    }
}
// [[Rcpp::export]]
void good_loop() {
    for (int i = 0; i < 100; i++) {
        std::printf("i = %dn", i);
        Rcpp::checkUserInterrupt();
    }
}
// [[Rcpp::export]]
void good_loop2() {
    for (double j = 0.0; j < 1.0; j += 0.01) {
        std::printf("j = %.2fn", j);
        Rcpp::checkUserInterrupt();
    }
}

---

第一种选择（good_loop）是适当地缩放你的步长——从0到99乘1循环的迭代次数与从0到0.99乘0.01的循环次数相同。此外，您可以只使用double而不是int，如good_loop2。无论如何，这里的主要要点是，在选择C++变量类型时需要更加小心。与R不同，当你声明i是一个int时，它将被视为一个int，而不是一个浮点数。

Answer 2

正如@nrussell非常专业地指出的那样，当持有的类型是double时，将i视为int存在问题。发布此答案的目的是强调需要避免使用 double 或 float 作为循环增量器。为了提高可读性，我选择将其作为答案而不是评论发布。

请注意，由于精度问题，循环增量不应作为double或float给出。 例如，由于i = 0.981111111等原因，很难获得i = .99...

相反，我会选择将循环处理为int并尽快将其转换为double/float，例如

for (int i=0; i < 100; i++){
    // Make sure to use double division 
    // (e.g. either numerator or denominator is a floating / double) 
    sum += S(h, i/100.0);
}

进一步说明：

RCPP犰狳和C++师问题

使用浮点数/双精度作为循环变量