mercredi 4 mars 2020

How to randomly generate integers following a certain distribution (power-law, normal, poisson...) in c

I'm pretty much ignorant when it comes to probabilities and statistics, and I'm trying to create random integer generators following certain distributions in c.
From what I gathered from wikipedia and questions already answered in here, I wrote the following functions that create an array of uint64_t following those distributions.

#define EXP 2.718281828459045



uint64_t *uniform_dist_ints(uint64_t nb_vals,uint64_t upper_bound) {
    uint64_t *arr=malloc(nb_vals*sizeof(uint64_t ));
    for(uint64_t i=0;i<nb_vals;i++)
        arr[i]=1+arc4random_uniform(upper_bound-1);//pas de 0
    return arr;
}

//from http://mathworld.wolfram.com/RandomNumber.html
uint64_t *power_dist_ints(uint64_t nb_vals,uint64_t lower_bound,uint64_t upper_bound,uint64_t p) {
    uint64_t *arr=malloc(nb_vals*sizeof(uint64_t));
    uint64_t x0=lower_bound,
        x1=upper_bound,y=0;
    for(uint64_t i=0;i<nb_vals;i++){
        y=(1+arc4random_uniform(upper_bound-1));
        arr[i]=(uint64_t)pow(((pow(x1,p+1)-pow(x0,p+1))*y+pow(x0,p+1)),1./(p+1));
    }
    return arr;
}

/*
algorithm poisson random number (Knuth):
    init:
        Let L ← e−λ, k ← 0 and p ← 1.
    do:
        k ← k + 1.
        Generate uniform random number u in [0,1] and let p ← p × u.
    while p > L.
    return k − 1.
*/

//from wikipedia https://en.wikipedia.org/wiki/Poisson_distribution#Generating_Poisson-distributed_random_variables
uint64_t *poisson_dist_ints_v1(uint64_t nb_vals,double lambda){
    //Knuth's method
    double p,l,u;
    uint64_t *res=malloc(sizeof(uint64_t)*nb_vals);
    uint64_t i,k;

    for(i=0;i<nb_vals;i++){
        l=pow(EXP,-lambda);
        k=0;
        p=1.;
        do{
            k++;
            u=((double)arc4random() / UINT32_MAX);
            p=p*u;
        }while(p>l);
        res[i]=k-1;
    }
    return res;
}

/*wikipedia
algorithm Poisson generator based upon the inversion by sequential search:[54]:505
    init:
        Let x ← 0, p ← e−λ, s ← p.
        Generate uniform random number u in [0,1].
    while u > s do:
        x ← x + 1.
        p ← p × λ / x.
        s ← s + p.
    return x.
*/

//also from https://en.wikipedia.org/wiki/Poisson_distribution#Generating_Poisson-distributed_random_variables
uint64_t *poisson_dist_ints_v2(uint64_t nb_vals,double lambda){
    double p,s,u;
    uint64_t *res=malloc(sizeof(uint64_t)*nb_vals);
    uint64_t i,x;

    for(i=0;i<nb_vals;i++){
        x=0;
        p=pow(EXP,-lambda);
        s=p;
        u=((double)arc4random() / UINT32_MAX);
        while(u>s){
            x++;
            p=p*lambda/x;
            s+=p;
        }
        res[i]=x;
    }
    return res;  
}

My questions are, what do you think about these functions first? and where can I find similar algorithms but to generate integers following a normal distribution and a power-law distribution (other than the one I already used)? Thanks.
P.S: I don't want to use a predefined library or something like that, I'm looking for some straight forward algorithms that work, I'm also not looking to read a +2000 page book for it, these functions will only be used to test another program, I don't want to spend too much time on them.




Aucun commentaire:

Enregistrer un commentaire