I've got C code which I'm calling from R:
.C("giveProb",as.double(2),as.double(2),as.double(c(0,1,0,1,1,0,1,0)))
I'd like to call it a bunch (hundreds of thousands or millions) of times (with differing third arguments), and it works fine when I put it in a for loop for around 100 times, but anything above that crashes R.
I've got a feeling that it's a memory problem from using R_alloc. I've got six arrays allocated in C, e.g.:
newCoefArray = (double *)R_alloc(1,curSize * sizeof(double));
But according to the R Manual "Writing R Extensions":
This memory is taken from the heap, and released at the end of the .C, .Call or .External call.
which I took to mean that the memory would be freed during each iteration of the loop. But in the next sentence:
Users can also manage it, by noting the current position with a call to vmaxget and clearing memory allocated subsequently by a call to vmaxset. This is only recommended for experts.
As someone who is nowhere near an expert, I was hoping someone here could help. The C code, in its entirety is below.
#include <R.h>
#include <Rmath.h>
#include <stdio.h>
void giveProb(double *k, double *q,double *order){
double curSize;
double tmpSize;
double *newCoefArray;
double *oldCoefArray;
double *newAArray;
double *oldAArray;
double *newBArray;
double *oldBArray;
int position=0;
long int factorial(int n){
if(n==0||n==1){
return(1);
}
int tmp=1,i=1;
while(i<=n){
tmp=tmp*i;
i++;
}
return(tmp);
}
void expander(double a, double b,double c,double d,double coeff){
double leadingTerm=beta(a,b);
int bb=b;
double index[bb], sumLeaders[bb];
for(int i=0;i<bb;i++){
index[i]=a+i;
sumLeaders[i]=factorial(a+b-1)/(factorial(index[i])*factorial(a+b-1-index[i]));
newCoefArray[i+position]=coeff*leadingTerm*sumLeaders[i];
newAArray[i+position]=index[i]+c+1;
newBArray[i+position]=a+b+d-index[i];
}
position=position+bb;
curSize=position;
}
void separator(double e, double f){
double a, b, coeff;
for(int i=0; i<tmpSize; i++){
coeff=oldCoefArray[i];
a=oldAArray[i];
b=oldBArray[i];
expander(a,b,e,f,coeff);
}
}
void condenser(){
tmpSize=0;
for(int i=1; i<curSize; i++){
for(int j=0; j<i; j++){
if(newAArray[j]==newAArray[i]){
newCoefArray[j]=newCoefArray[j]+newCoefArray[i];
newCoefArray[i]=0;
}
}
}
for(int i=0; i<curSize; i++){
tmpSize=tmpSize+(newCoefArray[i]!=0);
}
oldCoefArray =(double *) R_alloc(1,tmpSize * sizeof(double));
oldAArray = (double *)R_alloc(1,tmpSize * sizeof(double));
oldBArray = (double *)R_alloc(1,tmpSize * sizeof(double));
for(int i=0; i<tmpSize; i++){
oldCoefArray[i]=newCoefArray[i];
oldAArray[i]=newAArray[i];
oldBArray[i]=newBArray[i];
}
curSize=tmpSize;
}
long double coefficient=1;
for(int i=0;i<*k;i++){
coefficient=coefficient*factorial(*k)/(factorial(i)*factorial(*k-i-1));
}
for(int i=0;i<*q;i++){
coefficient=coefficient*factorial(*q)/(factorial(i)*factorial(*q-i-1));
}
double numObs=*k+*q;
double out=0;
curSize=order[1]+1;
newCoefArray = (double *)R_alloc(1,curSize * sizeof(double));
newAArray = (double *)R_alloc(1,curSize * sizeof(double));
newBArray = (double *)R_alloc(1,curSize * sizeof(double));
expander(order[0]+1,order[1]+1,order[2],order[3],coefficient);
oldCoefArray = (double *)R_alloc(1,curSize * sizeof(double));
oldAArray = (double *)R_alloc(1,curSize * sizeof(double));
oldBArray = (double *)R_alloc(1,curSize * sizeof(double));
for(int i=0;i<curSize; i++){
oldCoefArray[i]=newCoefArray[i];
oldAArray[i]=newAArray[i];
oldBArray[i]=newBArray[i];
}
for(int i=4;i<2*numObs;i+=2){
position=0;
tmpSize=curSize;
separator(order[i],order[i+1]);
condenser();
}
position=0;
for(int i=0;i<curSize;i++){
out=out+newCoefArray[i]*beta(newAArray[i],newBArray[i]);
}
*k=out;
}
UPDATE: Using the suggestion in the comment below, I get the following (which sort of confirms what I had previously thought, right?):
R -d valgrind -f test_script.R ==11131== Memcheck, a memory error detector ==11131== Copyright (C) 2002-2010, and GNU GPL'd, by Julian Seward et al. ==11131== Using Valgrind-3.6.0 and LibVEX; rerun with -h for copyright info ==11131== Command: /usr/lib64/R/bin/exec/R -f test_script.R ==11131==
R version 2.15.1 (2012-06-22) -- "Roasted Marshmallows"
Copyright (C) 2012 The R Foundation for Statistical Computing
ISBN 3-900051-07-0
Platform: x86_64-redhat-linux-gnu (64-bit)
R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.
==11131== Conditional jump or move depends on uninitialised value(s)
==11131== at 0x3A6A685F80: __GI___strcasecmp_l (in /lib64/libc-2.12.so)
==11131== by 0x3A6A61FF24: __gconv_open (in /lib64/libc-2.12.so)
==11131== by 0x3A6A62D3B7: _nl_find_msg (in /lib64/libc-2.12.so)
==11131== by 0x3A6A62DB83: __dcigettext (in /lib64/libc-2.12.so)
==11131== by 0x3A6C3BD2DF: ??? (in /usr/lib64/R/lib/libR.so)
==11131== by 0x3A6C313FC8: setup_Rmainloop (in /usr/lib64/R/lib/libR.so)
==11131== by 0x3A6C315278: Rf_mainloop (in /usr/lib64/R/lib/libR.so)
==11131== by 0x40084A: main (in /usr/lib64/R/bin/exec/R)
==11131==
==11131== Use of uninitialised value of size 8
==11131== at 0x3A6A6863A4: __GI___strcasecmp_l (in /lib64/libc-2.12.so)
==11131== by 0x3A6A61FF24: __gconv_open (in /lib64/libc-2.12.so)
==11131== by 0x3A6A62D3B7: _nl_find_msg (in /lib64/libc-2.12.so)
==11131== by 0x3A6A62DB83: __dcigettext (in /lib64/libc-2.12.so)
==11131== by 0x3A6C3BD2DF: ??? (in /usr/lib64/R/lib/libR.so)
==11131== by 0x3A6C313FC8: setup_Rmainloop (in /usr/lib64/R/lib/libR.so)
==11131== by 0x3A6C315278: Rf_mainloop (in /usr/lib64/R/lib/libR.so)
==11131== by 0x40084A: main (in /usr/lib64/R/bin/exec/R)
==11131==
==11131== Use of uninitialised value of size 8
==11131== at 0x3A6A6863A8: __GI___strcasecmp_l (in /lib64/libc-2.12.so)
==11131== by 0x3A6A61FF24: __gconv_open (in /lib64/libc-2.12.so)
==11131== by 0x3A6A62D3B7: _nl_find_msg (in /lib64/libc-2.12.so)
==11131== by 0x3A6A62DB83: __dcigettext (in /lib64/libc-2.12.so)
==11131== by 0x3A6C3BD2DF: ??? (in /usr/lib64/R/lib/libR.so)
==11131== by 0x3A6C313FC8: setup_Rmainloop (in /usr/lib64/R/lib/libR.so)
==11131== by 0x3A6C315278: Rf_mainloop (in /usr/lib64/R/lib/libR.so)
==11131== by 0x40084A: main (in /usr/lib64/R/bin/exec/R)
==11131==
Natural language support but running in an English locale
R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.
Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.
[Previously saved workspace restored]
> dyn.load("SchWolfenew.so")
> for(i in 1:1000){
+ .C("giveProb",as.double(2),as.double(2),as.double(c(0,1,0,1,1,0,1,0)))
+ }
==29371== Invalid read of size 1
==29371== at 0x3A6C31A2E9: ??? (in /usr/lib64/R/lib/libR.so)
==29371== by 0x3A6C31CBC9: Rf_cons (in /usr/lib64/R/lib/libR.so)
==29371== by 0x3A6C2D2B74: ??? (in /usr/lib64/R/lib/libR.so)
==29371== by 0x3A6C2DD121: Rf_eval (in /usr/lib64/R/lib/libR.so)
==29371== by 0x3A6C2DF830: Rf_applyClosure (in /usr/lib64/R/lib/libR.so)
==29371== by 0x3A6C32B828: Rf_usemethod (in /usr/lib64/R/lib/libR.so)
==29371== by 0x3A6C32BAE7: ??? (in /usr/lib64/R/lib/libR.so)
==29371== by 0x3A6C2D290B: ??? (in /usr/lib64/R/lib/libR.so)
==29371== by 0x3A6C2DD121: Rf_eval (in /usr/lib64/R/lib/libR.so)
==29371== by 0x3A6C2DF830: Rf_applyClosure (in /usr/lib64/R/lib/libR.so)
==29371== by 0x3A6C2DD3F7: Rf_eval (in /usr/lib64/R/lib/libR.so)
==29371== by 0x3A6C2DEF4F: ??? (in /usr/lib64/R/lib/libR.so)
==29371== Address 0x3ff0000000000003 is not stack'd, malloc'd or (recently) free'd
==29371==
*** caught segfault ***
address (nil), cause 'unknown'
==29371== Invalid read of size 1
==29371== at 0x3A6C31AF0B: ??? (in /usr/lib64/R/lib/libR.so)
==29371== by 0x3A6C31CBC9: Rf_cons (in /usr/lib64/R/lib/libR.so)
==29371== by 0x3A6C31CC71: Rf_allocList (in /usr/lib64/R/lib/libR.so)
==29371== by 0x3A6C2C8CD4: R_GetTraceback (in /usr/lib64/R/lib/libR.so)
==29371== by 0x3A6C313472: ??? (in /usr/lib64/R/lib/libR.so)
==29371== by 0x3A6B20F4FF: ??? (in /lib64/libpthread-2.12.so)
==29371== by 0x3A6C31A2E8: ??? (in /usr/lib64/R/lib/libR.so)
==29371== by 0x3A6C31CBC9: Rf_cons (in /usr/lib64/R/lib/libR.so)
==29371== by 0x3A6C2D2B74: ??? (in /usr/lib64/R/lib/libR.so)
==29371== by 0x3A6C2DD121: Rf_eval (in /usr/lib64/R/lib/libR.so)
==29371== by 0x3A6C2DF830: Rf_applyClosure (in /usr/lib64/R/lib/libR.so)
==29371== by 0x3A6C32B828: Rf_usemethod (in /usr/lib64/R/lib/libR.so)
==29371== Address 0x4020000000000003 is not stack'd, malloc'd or (recently) free'd
==29371==
==29371==
==29371== Process terminating with default action of signal 11 (SIGSEGV)
==29371== General Protection Fault
==29371== at 0x3A6C31AF0B: ??? (in /usr/lib64/R/lib/libR.so)
==29371== by 0x3A6C31CBC9: Rf_cons (in /usr/lib64/R/lib/libR.so)
==29371== by 0x3A6C31CC71: Rf_allocList (in /usr/lib64/R/lib/libR.so)
==29371== by 0x3A6C2C8CD4: R_GetTraceback (in /usr/lib64/R/lib/libR.so)
==29371== by 0x3A6C313472: ??? (in /usr/lib64/R/lib/libR.so)
==29371== by 0x3A6B20F4FF: ??? (in /lib64/libpthread-2.12.so)
==29371== by 0x3A6C31A2E8: ??? (in /usr/lib64/R/lib/libR.so)
==29371== by 0x3A6C31CBC9: Rf_cons (in /usr/lib64/R/lib/libR.so)
==29371== by 0x3A6C2D2B74: ??? (in /usr/lib64/R/lib/libR.so)
==29371== by 0x3A6C2DD121: Rf_eval (in /usr/lib64/R/lib/libR.so)
==29371== by 0x3A6C2DF830: Rf_applyClosure (in /usr/lib64/R/lib/libR.so)
==29371== by 0x3A6C32B828: Rf_usemethod (in /usr/lib64/R/lib/libR.so)
==29371==
==29371== HEAP SUMMARY:
==29371== in use at exit: 29,307,724 bytes in 12,896 blocks
==29371== total heap usage: 28,845 allocs, 15,949 frees, 48,495,252 bytes allocated
==29371==
==29371== LEAK SUMMARY:
==29371== definitely lost: 0 bytes in 0 blocks
==29371== indirectly lost: 0 bytes in 0 blocks
==29371== possibly lost: 0 bytes in 0 blocks
==29371== still reachable: 29,307,724 bytes in 12,896 blocks
==29371== suppressed: 0 bytes in 0 blocks
==29371== Rerun with --leak-check=full to see details of leaked memory
==29371==
==29371== For counts of detected and suppressed errors, rerun with: -v
==29371== Use --track-origins=yes to see where uninitialised values come from
==29371== ERROR SUMMARY: 5 errors from 5 contexts (suppressed: 21 from 9)
Segmentation fault (core dumped)
The realloc() function changes the size of a previously reserved storage block. The ptr argument points to the beginning of the block. The size argument gives the new size of the block, in bytes. The contents of the block are unchanged up to the shorter of the new and old sizes.
The function realloc is used to resize the memory block which is allocated by malloc or calloc before. Here is the syntax of realloc in C language, void *realloc(void *pointer, size_t size) Here, pointer − The pointer which is pointing the previously allocated memory block by malloc or calloc.
malloc is not required, you can use realloc only. malloc(n) is equivalent to realloc(NULL, n) . However, it is often clearer to use malloc instead of special semantics of realloc . It's not a matter of what works, but not confusing people reading the code.
Use of realloc function Using realloc function, we can resize the memory area which is already created by malloc or calloc. It's is also declared in stdlib. h library.
I put your C code in a file memory.c, and ran R CMD SHLIB memory.c. Here's my test script
dyn.load("/tmp/memory.so")
set.seed(123L)
while (TRUE)
.C("giveProb",as.double(2),as.double(2), sample(c(0, 1), 8, TRUE))
And R -d valgrind -f test_script.R says
> dyn.load("/tmp/memory.so")
> set.seed(123L)
> while (TRUE)
+ .C("giveProb",as.double(2),as.double(2), sample(c(0, 1), 8, TRUE))
==3461== Invalid write of size 8
==3461== at 0xBF29D78: expander.4631 (memory.c:35)
==3461== by 0xBF29EBB: separator.4643 (memory.c:49)
==3461== by 0xBF29AAB: giveProb (memory.c:108)
==3461== by 0x4EEAF49: do_dotCode (dotcode.c:1689)
==3461== by 0x4F1F4E2: Rf_eval (eval.c:493)
==3461== by 0x4F218F6: do_for (eval.c:1310)
==3461== by 0x4F1F2E8: Rf_eval (eval.c:467)
==3461== by 0x4F6B5FB: Rf_ReplIteration (main.c:256)
==3461== by 0x4F6B7B2: R_ReplConsole (main.c:305)
==3461== by 0x4F6D022: run_Rmainloop (main.c:987)
==3461== by 0x4F6D037: Rf_mainloop (main.c:994)
==3461== by 0x400845: main (Rmain.c:32)
(plus a lot more). memory.c line 35 is
newCoefArray[i+position]=coeff*leadingTerm*sumLeaders[i];
so i + position is larger than the allocated array.
If you love us? You can donate to us via Paypal or buy me a coffee so we can maintain and grow! Thank you!
Donate Us With