[Rd] R-gui sessions end when executing C-code
William Dunlap
wdunlap at tibco.com
Fri Feb 2 21:28:15 CET 2018
SEXP eta = PROTECT(allocVector(REALSXP,H_c)); n_prot++;
double *eta_c; eta_c = REAL(eta);
for (i=0;i<N_c;i++)
{
start_c = i * H_c;
for (j=0; j<H_c; j++)
{
eta_c[start_c + j] = ...
It looks you expect to be able to write into the N_c*H_c element
of eta, but you only allocated H_c elements for it.
Bill Dunlap
TIBCO Software
wdunlap tibco.com
On Fri, Feb 2, 2018 at 6:22 AM, Jesper Hybel Pedersen <jeshyb at dtu.dk> wrote:
> Hi
>
> I'm trying to develop some C code to find the fixpoint of a contraction
> mapping, the code compiles and gives the right results when executed in R.
> However R-gui session is frequently terminated. I'm suspecting some access
> violation error due to the exception code 0xc0000005
> In the error report windows 10 gives me.
>
> It is the first time I'm writing any C-code so I'm guessing I have done
> something really stupid. I have been trying to debug the code for a couple
> of days now,
> But I simply can't figure out what generates the problem. Could it be
> something particular to my windows set up and security stuff?
>
>
> I'm in the process of reading Writing R Extensions and Hadley Wickham's
> Advanced R but might have missed something.
>
> The windows error report:
>
> Faulting application name: Rgui.exe, version: 3.33.6774.0, time stamp:
> 0x58bd6d26
> Faulting module name: R.dll, version: 3.33.6774.0, time stamp: 0x58bd6d0b
> Exception code: 0xc0000005
> Fault offset: 0x000000000010b273
> Faulting process id: 0x1d14
> Faulting application start time: 0x01d39aede45c96e9
> Faulting application path: C:\Program Files\R\R-3.3.3\bin\x64\Rgui.exe
> Faulting module path: C:\Program Files\R\R-3.3.3\bin\x64\R.dll
> Report Id: c78d7c52-72c5-40f3-a3cc-927323d2af07
> Faulting package full name:
> Faulting package-relative application ID:
>
>
> ####### How I call the C-function in R
>
> dyn.load("C://users//jeshyb//desktop//myC//lce_fixpoint_cc.dll")
>
>
> N = 10
> H = 3
> v <- rnorm(N*H)
> mu <- 0.1
> psi <- matrix(c(1,0,0.5,0.5,0,1),nrow=2)
> K <- dim(psi)[1]
> p <- rep(1/H,N*H)
> error <- 1e-10
>
>
> f<-function(p,v,mu,psi,N,H,K)
> {
>
> .Call("lce_fixpoint_cc",p, v, mu, psi, as.integer(N), as.integer(H),
> as.integer(K),error)
> }
>
>
> for (i in 1:100)
> {
> v <- rnorm(N*H)
> p <- rep(1/H,N*H)
>
> a<-f(p,v,mu,psi,N,H,K)
> }
>
>
> a<-f(p,v,mu,psi,N,H,K)
> plot(a)
>
>
>
> ######## The C- function
>
>
>
> #include <R.h>
> #include <Rinternals.h>
>
>
> SEXP lce_fixpoint_cc(SEXP q, SEXP v, SEXP mu, SEXP psi, SEXP N,SEXP H,
> SEXP K, SEXP err)
> {
>
> int n_prot = 0;
> /* Make ready integer and double constants */
> PROTECT(N); n_prot++;
> PROTECT(H); n_prot++;
> PROTECT(K); n_prot++;
> int N_c = asInteger(N);
> int H_c = asInteger(H);
> int K_c = asInteger(K);
>
> double mu_c = asReal(mu);
> double mu2_c = 1.0 - mu_c;
> double error_c = asReal(err);
> double lowest_double = 1e-15;
> double tmp_c;
> double denom;
> double error_temp;
> double error_i_c;
>
>
> /* Make ready vector froms input */
> PROTECT(q); n_prot++;
> PROTECT(v); n_prot++;
> PROTECT(psi); n_prot++;
> double *v_c; v_c = REAL(v);
> double *psi_c; psi_c = REAL(psi);
>
> /* Initialize new vectors not given as input */
> SEXP q_copy = PROTECT(duplicate(q)); n_prot++;
> double *q_c; q_c = REAL(q_copy);
>
> SEXP q_new = PROTECT(allocVector(REALSXP,length(q)));
> n_prot++;
> double *q_new_c; q_new_c = REAL(q_new);
>
> SEXP eta = PROTECT(allocVector(REALSXP,H_c));
> n_prot++;
> double *eta_c; eta_c = REAL(eta);
>
> SEXP exp_eta = PROTECT(allocVector(REALSXP,H_c));
> n_prot++;
> double *exp_eta_c; exp_eta_c = REAL(exp_eta);
>
> SEXP psi_ln_psiq =
> PROTECT(allocVector(REALSXP,H_c)); n_prot++;
> double *psi_ln_psiq_c; psi_ln_psiq_c =
> REAL(psi_ln_psiq);
>
> int not_converged;
> int maxIter = 10000;
> int iter;
> int start_c;
>
> /* loop indeces */
> int i;
> int j;
> int k;
>
> /* loop over observational units to find choice
> probabilities for i=1,...,N */
> for (i=0;i<N_c;i++)
> {
>
> start_c = i * H_c;
> not_converged = 1;
> iter = 0;
>
> while(iter < maxIter
> && not_converged)
> {
>
> /* make v_ij
> + (1-mu)*ln(q_ij) */
>
> for (j=0; j<H_c; j++)
>
> {
>
> eta_c[start_c + j] = v_c[start_c + j] +
> mu2_c * log(q_c[start_c + j]);
>
> psi_ln_psiq_c[j] = 0.0;
>
> }
>
>
> /* Make psi_ln_psiq_c vector for individual i */
>
> for (k=0;k<K_c;k++)
>
> {
>
> tmp_c = 0.0;
>
>
> /* Calculate row k of psi %*% q */
>
> for (j=0;j<H_c;j++)
>
> {
>
> tmp_c +=
> psi_c[k + j*K_c] * q_c[start_c +j];
>
> }
>
>
> tmp_c = mu2_c
> * log(tmp_c);
>
>
> for (j=0;j<H_c;j++)
>
> {
>
>
> psi_ln_psiq_c[j] += psi_c[k + j*K_c] * tmp_c;
>
> }
>
> }
>
>
> denom = 0.0;
>
> for (j=0;j<H_c;j++)
>
> {
>
> exp_eta_c[start_c + j] = exp(
> eta_c[start_c + j] - psi_ln_psiq_c[j] ) + lowest_double;
>
> denom += exp_eta_c[start_c + j];
>
> }
>
>
> error_i_c = 0.0;
>
> error_temp = 0.0;
>
>
> /* calculate error and update choice prob */
>
> for (j=0;j<H_c;j++)
>
> {
>
> q_new_c[start_c + j] = exp_eta_c[start_c
> + j]/denom;
>
> error_temp = fabs(q_new_c[start_c + j] -
> q_c[start_c + j]);
>
> if (error_temp>error_i_c)
>
> {
>
>
> error_i_c = error_temp;
>
> }
>
> q_c[start_c + j] = q_new_c[start_c + j];
>
> }
>
>
> not_converged = error_i_c > error_c;
>
> iter++;
> } /* End while loop
> for individual i to solve for q_i */
>
>
> } /* End loop over individuals */
>
>
> UNPROTECT(n_prot);
> return(q_new);
> }
>
>
>
> Best regards
> Jesper Hybel
>
>
> [[alternative HTML version deleted]]
>
> ______________________________________________
> R-devel at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-devel
>
[[alternative HTML version deleted]]
More information about the R-devel
mailing list