#include <sys/time.h>
#include <math.h>
#include <gsl/gsl_randist.h>


int main(void)
{
  int nMu = 30;
  int nSim = 1000;
  double minLogMu=8;
  double maxLogMu=11;
  unsigned long iSim,iMu;
  double x,mu;
  double xMean;
  
  struct timeval tv;  

  //Set Random number generator
  gsl_rng* rng;
  rng = gsl_rng_alloc(gsl_rng_default);
  gettimeofday(&tv,NULL);
  gsl_rng_set(rng,tv.tv_usec);


  for(iMu=0;iMu<nMu;iMu++)
  {
    mu=pow(10.,minLogMu + (maxLogMu - minLogMu) * iMu / (nMu - 1.));
    xMean = 0;
    
    for(iSim=0;iSim<nSim;iSim++)
    {
      xMean += gsl_ran_poisson(rng,mu);
    }
    xMean /= (double)nSim;
    printf(" mu = %f    Mean(x) = %f   Ratio=%f\n",mu,xMean,xMean/mu);
    
  }
}