#define _DEFAULT_SOURCE
#define  _ISOC99_SOURCE
#include <stdlib.h>
#include <math.h>
/*****************************************************************************/
/* 
The alea() function returns a pseudo-random double in the range [0,1[ (0. inclusive to 1. exclusive).
It uses random() function (rather than rand() function) because it is supposed to be more efficient.
*/
double alea(){
#ifdef __osf__
  double rand_max=1024.*1024.*1024.*2.;
#elif defined __linux__
  double rand_max=(double)RAND_MAX+1.;
#else 
  double rand_max=(double)RAND_MAX+1.;
#endif

  return ((double)random()/(double)(rand_max));

}
/*****************************************************************************/
/* 
The aleai() function returns a pseudo-random double in the range [0,1] (0. inclusive to 1. inclusive).
It uses random() function (rather than rand() function) because it is supposed to be more efficient.
*/
double aleai(){
#ifdef __osf__
  double rand_max=1024.*1024.*1024.*2.;
#elif defined __linux__
  double rand_max=RAND_MAX;
#else 
  double rand_max=RAND_MAX;
#endif

  return ((double)random()/(double)(rand_max));

}
/*****************************************************************************/
/* 
   initialization of pseudo-random process used by alea() function
 */
int salea(unsigned seed){

#ifdef __osf__
  return srandom(seed);
#elif defined __linux__
  srandom(seed);
  return 0;
#else
  srandom(seed);
  return 0;
#endif
}

/***********************************************************************************/
/* Pseudo-random generator of Poisson distribution.
   Lambda is the mean of the random process.
   The uniform pseudo-random process random() is supposed to have been initialized with
   salea(seed) or srandom(seed).
 */
int random_poisson(double lambda){
  int n;
  double t;
  long int u;
  
  n = -1;
  t=0.;
  do{

    do{
      u = random();
    } while(u==0);
    /* doing so, 0 is exluded */

    t -=log( u / (double)(RAND_MAX));
    n++;

  } while(t<lambda);
  
  return n;
  
}