CODEWITHSUNDEEP
c
Damien
Damien

Reputation: 43

Generating random numbers of exponential distribution

#include <stdio.h>
#include <stdlib.h>
#include <math.h>

double ran_expo(lambda){
    double u;

    u = rand() / (RAND_MAX + 1.0);

    return -log(1- u) / lambda;
}

I am using this (I do not put other parts of the codes here) to generate random numbers of an exponential distribution lambda = 0.05. (Probability density function being lambda * exp(-lambda * x)).

However, I always get very small numbers, such as 0.000041, or something like 1.#INF00 (what is this?).

In fact, for exponential distribution with lambda = 0.05, numbers generated should be quite large, namely most are not far from 30. My results are very weird. In addition, the degree of precision is not satisfactory, only to 10^(-6). I tried long double, but it is still like this.

I am using DEV C++ under windows.

Upvotes: 5

Views: 10733

Answers (2)

Weather Vane
Weather Vane

Reputation: 34575

This is because you did not declare the type for lambda, having corrected that the results are in the range you seek. If undeclared, old compilers will assume it to be int.

#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <time.h>

double ran_expo(double lambda){
    double u;
    u = rand() / (RAND_MAX + 1.0);
    return -log(1- u) / lambda;
}

int main(void)
{
    int i;
    srand((unsigned)time(NULL));
    for (i=0; i<20; i++)
        printf("%f\n", ran_expo(0.05));
    return 0;
}

Program output:

0.025040
16.582459
4.296027
33.079902
17.589123
13.073084
8.624299
45.254803
34.611211
27.454302
3.825699
39.168172
24.790600
14.411160
7.247698
0.301951
1.917010
9.065004
3.187146
3.627885

Upvotes: 5

djhaskin987
djhaskin987

Reputation: 10077

I threw your code verbatim into code chef, it looks like it works fine. Compiler is gcc c++ 4.9.2 . This code prints out 36.6751:

#include <iostream>
#include <string>
#include <cstdlib>
#include <cmath>

using namespace std;

double ran_expo(double lambda){
    double u;

    u = rand() / (RAND_MAX + 1.0);

    return -log(1- u) / lambda;
}


int main() {
    cout << ran_expo(0.05);
    return 0;
}

Works the same way in C, too (again, gcc 4.9.2):

#include <stdlib.h>
#include <math.h>

double ran_expo(double lambda){
    double u;

    u = rand() / (RAND_MAX + 1.0);

    return -log(1- u) / lambda;
}


int main() {
    printf("%f\n", ran_expo(0.05));
    return 0;
}

Must be your compiler/libraries?

Upvotes: 0

Related Questions