转自:http://blog.csdn.net/chl033/article/details/3936584
#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);
}
自己补充:最新的项目需要产生0-100万内的为随机数,之前用rand()函数产生的伪随机数的范围为0-0x7FFF,无法满足大随机数的产生条件,
故需要另外选择大随机数的产生的方法。
以下是关于该随机函数的解释
seed 是指针,内容可以被函数uniform修改 *seed=2045*(*seed)+1; 修改了seed的内容。 *seed=*seed-(*seed/1048576)*1048576;是取*seed的低20位(2进制)计算机作移位处理1048576是0x100000 t=(*seed)/1048576.0;中*seed 被转化成float并且小于1048576所以t=0.xxxx 就是所谓的随机数
/*******************************************************************
* 正态分布
* 输入: 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];
}
}
}