poisson.cpp

C语言

#include "stdio.h"
#include "math.h"
#define UL unsigned long
#define znew  (z=36969*(z&65535)+(z>>16))
#define wnew  (w=18000*(w&65535)+(w>>16))
#define MWC   ((znew<<16)+wnew)
#define SHR3  (jsr^=(jsr<<17),jsr^=(jsr>>13),jsr^=(jsr<<5))
#define CONG  (jcong=69069*jcong+1234567)
#define FIB   ((b=a+b),(a=b-a))
#define KISS  ((MWC^CONG)+SHR3)
#define LFIB4 (c++,t[c]=t[c]+t[UC(c+58)]+t[UC(c+119)]+t[UC(c+178)])
#define SWB   (c++,bro=(x<y),t[c]=(x=t[UC(c+58)])-(y=t[UC(c+19)]+bro))
#define UNI   (KISS*2.328306e-10)
#define VNI   ((long) KISS)*4.656613e-10
#define UC    (unsigned char)/*a cast operation*/
#define pi     3.14159
#define lamd   0.2
/*  Global static variables: */
 static UL z=362436069, w=521288629, jsr=123456789, jcong=380116160;
 static UL a=224466889,b=7584631,t[256];
 static UL x=0,y=0,bro; static unsigned char c=0;

/* Random seeds must be used to reset z,w,jsr,jcong and
the table t[256]  Here is an example procedure, using KISS: */
FILE  *f;
 void settable(UL i1,UL i2,UL i3,UL i4,UL i5,UL i6)
 { int i; z=i1;w=i2,jsr=i3; jcong=i4,a=i5,b=i6;
 for(i=0;i<256;i++)  t[i]=SHR3;        }
void main()
 {
	f=fopen("out.txt","w");
	 int k;
     double a[1000];
	 settable(12345,65435,34221,12345,9983651,95746118);
	 for(k=0;k<1000;k++)
	 {
		 a[k]=UNI;
	 }

	double y[1000];
	for(k=0;k<1000;k++)
		y[k]=-log(a[k])/lamd;
	double sum=0.000;
	int	j=0;
	int i=1;
    do
	{
      sum=sum+y[i-1];
	  if(sum<1)  j++;
	  else if(sum>=1)
	  {
		  printf("%d\n ",j);
		  fprintf(f,"%d\n ",j);
		  j=1;sum=y[i-1];
	  }
	  i++;
	}while(i<=1000);
	printf("\n");
	fclose(f);

}