gaussrand2.c

生成均匀随机数列算法实现和生成高斯随即序列的两种算法实现

/*        高斯分布随机序列的产生和验证                    */

/*  产生含有N个值的高斯分布随机序列               
    求出产生的随机序列的均值和方差
    对所求的自相关函数值来判断是否符合白噪声模型*/

#include "stdio.h"
#include "math.h"
#include "gauss2.c"

#define N 1500                /*产生的高斯分布随机序列含有N个值*/
#define M 50                 /*产生的高斯分布随机序列的自相关函数的自变量取值范围，即(0,49)*/

main()
{
 int i,j;
 double X[N],Rx[M];  /*设置数组X[N]存放每次产生的随机序列，设置数组Rx[M]存放自相关函数*/
 long int s;                  /*随机种子*/
 double mean,dev;             /*均值mean,方差dev*/
 double gauss2(double,double,long int *);   /*产生高斯分布的伪随机数函数*/
 int num1,num2,num3;

 for(i=0;i<M;i++)
 {
  Rx[i]=0;
 }

 scanf("%ld",&s);
  for(i=0;i<(N/5);i++)
   {for(j=0;j<5;j++)
      {X[i*5+j]=gauss2(0.0,1.0,&s)*sqrt(0.5)+1;
        printf("%13.7f",X[i*5+j]);
      }
   printf("\n");
   }
  printf("\n");

  mean=0;
  for (i=0;i<N;i++)  mean=mean+X[i];    /*计算随机序列的均值*/              
  mean=mean/N;
  printf("mean=%13.7f    ",mean);
  dev=0;
  for (i=0;i<N;i++)  dev=dev+(X[i]-mean)*(X[i]-mean);    /*计算随机序列的方差*/
  dev=dev/N;
  printf("dev=%13.7f\n\n",dev);

  num1=0;                         /*分别计算产生的随机数在区间(mean-sqrt(dev),mean+sqrt(dev))，*/
  num2=0;                         /*(mean-2*sqrt(dev),mean=2*sqrt(dev))，                      */
  num3=0;                         /*(mean-3*sqrt(dev),mean+3*sqrt(dev))中的概率                */
  for(i=0;i<N;i++)  if((X[i]>=(mean-sqrt(dev)))&&(X[i]<=(mean+sqrt(dev))))  num1++;  
  printf("num1=%d  ",num1);
  printf("the percent of (u-p,u+p) is: %13.7f\n",(double)num1/(double)N);
  for(i=0;i<N;i++)  if((X[i]>=(mean-2*sqrt(dev)))&&(X[i]<=(mean+2*sqrt(dev))))  num2++;
  printf("num2=%d  ",num2);
  printf("the percent of (u-2p,u+2p) is: %13.7f\n",(double)num2/(double)N);
  for(i=0;i<N;i++)  if((X[i]>=(mean-3*sqrt(dev)))&&(X[i]<=(mean+3*sqrt(dev))))  num3++;
  printf("num3=%d  ",num3);
  printf("the percent of (u-3p,u+3p) is: %13.7f\n\n",(double)num3/(double)N);

  for(i=0;i<M;i++)              /*计算自相关函数*/
  {for(j=0;j<N-i;j++)  Rx[i]=Rx[i]+X[j]*X[j+i];                
   Rx[i]=Rx[i]/(N-i); 
   printf("%13.7f  \n",Rx[i]);
  }
}