A função Rcpp trava
Estou usando o R.3.0.1 em conjunto com o RStudio 0.97.551 em um PC com Windows7 de 64 bits e comecei a terceirizar uma função para C / C ++ usando o Rcpp. A função é compilada, mas a avaliação dentro de uma função R produz um erro de tempo de execução. Não consigo descobrir por que e como corrigir isso.
DetalhesAbaixo está o meu arquivo cpp ... digamos que seja chamado "vector.cpp"
#include <Rcpp.h>
using namespace Rcpp;
// [[Rcpp::export]]
NumericVector l5(double z, double k, double s, double K, double theta, double x, double h, NumericVector m){
int n = m.size();
NumericVector a(n);
NumericVector bu(n);
NumericVector b(n);
NumericVector c(n);
for(int i=0; i<n+1; i++){
a[i] = pow(z,m[i]) * (pow((x*pow(h,m[i])/K), theta) - 1) * (K/theta);
for (int j=0; j<i; j++){
bu[i] += pow(z,j) * (1 - z) * fmax(((pow((s/K), theta) - 1) * (K/theta)), ((pow((x*pow(h, j)/K), theta) - 1) * (K/theta)));
}
b[i] = k *bu[i];
c[i] = k * pow(z, m[i]) * fmax(((pow((s/K), theta) - 1) * (K/theta)), ((pow((x*pow(h, m[i])/K), theta) - 1) * (K/theta)));
}
return wrap(a-b-c);
}
que eu compilar em R (ou RStudio) usando o comando
sourceCpp(<path to file>/vector.cpp)
Compila - até agora tudo bem. No entanto, quando uso a função l5 em outras funções R, muitas vezes leva R a travar (tanto no RStudio quanto na GUI simples R). De fato, também a avaliação em si não é mais estável. Para reproduzir isso, p. tente avaliar l6 várias vezes
l6 <- function(zs, ks, ss, Ks, thetas, xs, hs){
z=zs
k=ks
s=ss
K=Ks
theta=thetas
x=xs
h=hs
m=0:30
res <- l5(z, k, s, K, theta,x, h, m)
return(res)
}
e corra
l6(0.9, 0.1, 67, 40, 0.5, 44, 1.06)
Especificamente, produz o seguinte erro de tempo de execução
This application has requested Runtime to terminate it in an unusual way.
Então, o que há de errado com minha função?
SoluçãoComo Dirk sugeriu abaixo, há um erro elementar no loop for, em que i é executado de 0 a n e, portanto, possui n + 1 elementos, mas eu inicializei apenas vetores de comprimento n. Para evitar esse erro, agora implementei a função usando iteradores
#include <Rcpp.h>
using namespace Rcpp;
// [[Rcpp::export]]
NumericVector l5(double z, double k, double s, double K, double theta, double x, double h, NumericVector m){
int n = m.size();
NumericVector a(n);
NumericVector bu(n);
NumericVector b(n);
NumericVector c(n);
for(NumericVector::iterator i = m.begin(); i != m.end(); ++i){
a[*i] = pow(z, m[*i]) * (pow((x*pow(h, m[*i])/K), theta) - 1) * (K/theta);
for(int j=0; j<*i; j++){
bu[*i] += pow(z,j) * (1 - z) * fmax(((pow((s/K), theta) - 1) * (K/theta)), ((pow((x*pow(h, j)/K), theta) - 1) * (K/theta)));
}
b[*i] = k *bu[*i];
c[*i] = k * pow(z, m[*i]) * fmax(((pow((s/K), theta) - 1) * (K/theta)), ((pow((x*pow(h, m[*i])/K), theta) - 1) * (K/theta)));
}
return wrap(a-b-c);
}
Muito obrigado novamente!