rgamma function

[Chandler Zuo]

Hi,

Has anyone encountered the problem of rgamma function in C? The
following simplified program always dies for me, and I wonder if anyone
can tell me the reason.

#include <Rmath.h>
#include <time.h>
#include <Rinternals.h>

SEXP generateGamma ()
{
     srand(time(NULL));
     return (rgamma(5000,1));
}

Has anyone encountered a similar problem before? Is there another way of
generating Gamma random variable in C?

P.S. I have no problem compiling and loading this function in R.

Thanks for suggestions in advance!

--Chandler

______________________________________________
[hidden email] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[Peter Dalgaard-2]

On Jul 14, 2012, at 04:55 , Chandler Zuo wrote:

> Hi,
>
> Has anyone encountered the problem of rgamma function in C? The following simplified program always dies for me, and I wonder if anyone can tell me the reason.
>
> #include <Rmath.h>
> #include <time.h>
> #include <Rinternals.h>
>
> SEXP generateGamma ()
> {
>    srand(time(NULL));
>    return (rgamma(5000,1));
> }
>
> Has anyone encountered a similar problem before? Is there another way of generating Gamma random variable in C?
>
> P.S. I have no problem compiling and loading this function in R.

It doesn't even give off a warning??

The prototype in Rmath.h is

double  rgamma(double, double);

and you should be returning an SEXP. As soon as something tries to interpret the double value as a pointer -- Poof!

Notice that rgamma in C is not the same function as the R counterpart, in particular it isn't vectorized, so only generates one random number at a time. The long and the short of it is that you need to read up on sections 5.9 and 5.10 of Writing R Extensions.

>
> Thanks for suggestions in advance!
>
> --Chandler
>
> ______________________________________________
> [hidden email] mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.

--
Peter Dalgaard, Professor,
Center for Statistics, Copenhagen Business School
Solbjerg Plads 3, 2000 Frederiksberg, Denmark
Phone: [phone]
Email: [hidden email]  Priv: [hidden email]

______________________________________________
[hidden email] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[Thomas Lumley-2]

On Fri, Jul 13, 2012 at 7:55 PM, Chandler Zuo <[hidden email]> wrote:

> Hi,
>
> Has anyone encountered the problem of rgamma function in C? The following
> simplified program always dies for me, and I wonder if anyone can tell me
> the reason.
>
> #include <Rmath.h>
> #include <time.h>
> #include <Rinternals.h>
>
> SEXP generateGamma ()
> {
>     srand(time(NULL));
>     return (rgamma(5000,1));
> }

rgamma doesn't return an SEXP, it returns a double.  Also, the srand()
call is pointless.

> Has anyone encountered a similar problem before? Is there another way of
> generating Gamma random variable in C?
>
> P.S. I have no problem compiling and loading this function in R.

Strange. You should get compiler warnings that the return type is
incompatible.  I get

foo.c: In function ‘generateGamma’:
foo.c:7: warning: implicit declaration of function ‘srand’
foo.c:8: error: incompatible types in return

I thought the ANSI standard actually *required* a diagnostic for the
incompatible return types.

   -thomas

--
Thomas Lumley
Professor of Biostatistics
University of Auckland

______________________________________________
[hidden email] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[Chandler Zuo]

Sorry that I posted the wrong syntax. My initial program is very long
and I tried to post the section where I have been narrowed to locate the
problem.

In this sample I am simulating a Gamma with size parameter 5000 and
scale parameter 1.

I plugged in many breaks in my initial program and what I found was that
most of the time the program stops when encountering a call to rgamma()
function. It just freezes without popping out any error.

#include<Rmath.h>
#include<time.h>
#include<Rinternals.h>

SEXP generateGamma ()
{
     SEXP a;
     PROTECT(a=allocVector(REALSXP,1));
     srand(time(NULL));
     REAL(a)[0]=rgamma(5000,1);
     UNPROTECT(1);
     return (a);
}


Thanks for your suggestions!


On 07/15/12 07:07, peter dalgaard wrote:

> On Jul 14, 2012, at 04:55 , Chandler Zuo wrote:
>
>> Hi,
>>
>> Has anyone encountered the problem of rgamma function in C? The following simplified program always dies for me, and I wonder if anyone can tell me the reason.
>>
>> #include<Rmath.h>
>> #include<time.h>
>> #include<Rinternals.h>
>>
>> SEXP generateGamma ()
>> {
>>     srand(time(NULL));
>>     return (rgamma(5000,1));
>> }
>>
>> Has anyone encountered a similar problem before? Is there another way of generating Gamma random variable in C?
>>
>> P.S. I have no problem compiling and loading this function in R.
> It doesn't even give off a warning??
>
> The prototype in Rmath.h is
>
> double  rgamma(double, double);
>
> and you should be returning an SEXP. As soon as something tries to interpret the double value as a pointer -- Poof!
>
> Notice that rgamma in C is not the same function as the R counterpart, in particular it isn't vectorized, so only generates one random number at a time. The long and the short of it is that you need to read up on sections 5.9 and 5.10 of Writing R Extensions.
>
>
>
>
>> Thanks for suggestions in advance!
>>
>> --Chandler
>>
>> ______________________________________________
>> [hidden email] mailing list
>> https://stat.ethz.ch/mailman/listinfo/r-help
>> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
>> and provide commented, minimal, self-contained, reproducible code.

______________________________________________
[hidden email] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[Thomas Lumley-2]

On Mon, Jul 16, 2012 at 12:49 PM, Chandler Zuo <[hidden email]> wrote:

> Sorry that I posted the wrong syntax. My initial program is very long and I
> tried to post the section where I have been narrowed to locate the problem.
>
> In this sample I am simulating a Gamma with size parameter 5000 and scale
> parameter 1.
>
> I plugged in many breaks in my initial program and what I found was that
> most of the time the program stops when encountering a call to rgamma()
> function. It just freezes without popping out any error.
>
> #include<Rmath.h>
> #include<time.h>
> #include<Rinternals.h>
>
> SEXP generateGamma ()
> {
>     SEXP a;
>     PROTECT(a=allocVector(REALSXP,1));
>     srand(time(NULL));
>     REAL(a)[0]=rgamma(5000,1);
>     UNPROTECT(1);
>     return (a);
> }
>

srand() isn't relevant to the R random number generator, and you
haven't included the header file that defines it (this should lead to
a warning about implicit declaration).

More importantly, you haven't included the R equivalents GetRNGstate()
and PutRNGstate().  You want something like

#include "R.h"
#include<Rmath.h>
#include<Rinternals.h>

SEXP generateGamma ()
{
    SEXP a;
    PROTECT(a=allocVector(REALSXP,1));
    GetRNGstate();
    REAL(a)[0]=rgamma(5000,1);
    PutRNGstate();
    UNPROTECT(1);
    return (a);
}

    - thomas

--
Thomas Lumley
Professor of Biostatistics
University of Auckland

______________________________________________
[hidden email] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.