常见伪随机数常生

#include "stdlib.h"
#include "stdio.h"
#include "math.h"

double uniform(double a,double b,long int* seed);
double gauss(double mean,double sigma,long int *seed);
double exponent(double beta,long int *seed);
double laplace(double beta,long int* seed);
double rayleigh(double sigma,long int *seed);
double weibull(double a,double b,long int*seed);
int bn(double p,long int*seed);
int bin(int n,double p,long int*seed);
int poisson(double lambda,long int *seed);



void main()
{
double a,b,x,mean;
int i,j;
long int s;

a=4;
b=0.7;
s=13579;
mean=0;
for(i=0;i<10;i++)
{
for(j=0;j<5;j++)
{
x=poisson(a,&s);
mean+=x;
printf("%-13.7f",x);
}
printf("/n");
}
mean/=50;
printf("平均值为:%-13.7f/n",mean);
}




/*******************************************************************
* 求[a,b]上的均匀分布
* 输入: a--双精度实型变量，给出区间的下限
*       b--双精度实型变量，给出区间的上限
*    seed--长整型指针变量，*seed为随机数的种子  
********************************************************************/
double uniform(double a,double b,long int*seed)
{
double t;
*seed=2045*(*seed)+1;
*seed=*seed-(*seed/1048576)*1048576;
t=(*seed)/1048576.0;
t=a+(b-a)*t;

return(t);
}

/*******************************************************************
* 正态分布
* 输入: mean--双精度实型变量，正态分布的均值
*      sigma--双精度实型变量，正态分布的均方差
*       seed--长整型指针变量，*seed为随机数的种子  
********************************************************************/
double gauss(double mean,double sigma,long int*seed)
{
int i;
double x,y;
for(x=0,i=0;i<12;i++)
x+=uniform(0.0,1.0,seed);
x=x-6.0;
y=mean+x*sigma;

return(y);
}


/*******************************************************************
* 指数分布
* 输入: beta--指数分布均值
*       seed--种子
*******************************************************************/
double exponent(double beta,long int *seed)
{
double u,x;
u=uniform(0.0,1.0,seed);
x=-beta*log(u);

return(x);
}

/*******************************************************************
* 拉普拉斯随机分布
* beta--拉普拉斯分布的参数
* *seed--随机数种子
*******************************************************************/
double laplace(double beta,long int* seed)
{
double u1,u2,x;

u1=uniform(0.,1.,seed);
u2=uniform(0.,1.,seed);
if(u1<=0.5)
x=-beta*log(1.-u2);
else
x=beta*log(u2);

return(x);
}



/********************************************************************
* 瑞利分布
*
********************************************************************/
double rayleigh(double sigma,long int *seed)
{
double u,x;
u=uniform(0.,1.,seed);
x=-2.0*log(u);
x=sigma*sqrt(x);
return(x);
}

/************************************************************************/
/* 韦伯分布                                                                     */
/************************************************************************/
double weibull(double a,double b,long int*seed)
{
double u,x;

u=uniform(0.0,1.0,seed);
u=-log(u);
x=b*pow(u,1.0/a);

return(x);
}

/************************************************************************/
/* 贝努利分布                                                           */
/************************************************************************/
int bn(double p,long int*seed)
{
int x;
double u;
u=uniform(0.0,1.0,seed);
x=(u<=p)?1:0;
return(x);
}

/************************************************************************/
/* 二项式分布                                                           */
/************************************************************************/
int bin(int n,double p,long int*seed)
{
int i,x;
for(x=0,i=0;i<n;i++)
x+=bn(p,seed);
return(x);
}

/************************************************************************/
/* 泊松分布                                                             */
/************************************************************************/
int poisson(double lambda,long int *seed)
{
int i,x;
double a,b,u;

a=exp(-lambda);
i=0;
b=1.0;
do {
u=uniform(0.0,1.0,seed);
b*=u;
i++;
} while(b>=a);
x=i-1;
return(x);
}

 

/************************************************************************
* 离散傅立叶变换与反变换
* 输入: x--要变换的数据的实部
* y--要变换的数据的虚部
*       a--变换结果的实部
*       b--变换结果的虚部
*       n--数据长度
*    sign--sign=1时，计算离散傅立叶正变换;sign=-1时;计算离散傅立叶反变换
************************************************************************/
void dft(double x[],double y[],double a[],double b[],int n,int sign)
{
int i,k;
double c,d,q,w,s;

q=6.28318530718/n;
for(k=0;k<n;k++)
{
w=k*q;
a[k]=b[k]=0.0;
for(i=0;i<n;i++)
{
d=i*w;
c=cos(d);
s=sin(d)*sign;
a[k]+=c*x+s*y;
b[k]+=c*y-s*x;
}
}
if(sign==-1)
{
c=1.0/n;
for(k=0;k<n;k++)
{
a[k]=c*a[k];
b[k]=c*b[k];
}
}
}