diff options
Diffstat (limited to 'src/intersim2/.svn/pristine/e2')
| -rw-r--r-- | src/intersim2/.svn/pristine/e2/e216ca41a2789704029b85d462988590d43813f3.svn-base | 108 |
1 files changed, 0 insertions, 108 deletions
diff --git a/src/intersim2/.svn/pristine/e2/e216ca41a2789704029b85d462988590d43813f3.svn-base b/src/intersim2/.svn/pristine/e2/e216ca41a2789704029b85d462988590d43813f3.svn-base deleted file mode 100644 index c27740f..0000000 --- a/src/intersim2/.svn/pristine/e2/e216ca41a2789704029b85d462988590d43813f3.svn-base +++ /dev/null @@ -1,108 +0,0 @@ -/* This program by D E Knuth is in the public domain and freely copyable. - * It is explained in Seminumerical Algorithms, 3rd edition, Section 3.6 - * (or in the errata to the 2nd edition --- see - * http://www-cs-faculty.stanford.edu/~knuth/taocp.html - * in the changes to Volume 2 on pages 171 and following). */ - -/* N.B. The MODIFICATIONS introduced in the 9th printing (2002) are - included here; there's no backwards compatibility with the original. */ - -/* This version also adopts Brendan McKay's suggestion to - accommodate naive users who forget to call ran_start(seed). */ - -/* If you find any bugs, please report them immediately to - * (and you will be rewarded if the bug is genuine). Thanks! */ - -/************ see the book for explanations and caveats! *******************/ -/************ in particular, you need two's complement arithmetic **********/ - -#define KK 100 /* the long lag */ -#define LL 37 /* the short lag */ -#define MM (1L<<30) /* the modulus */ -#define mod_diff(x,y) (((x)-(y))&(MM-1)) /* subtraction mod MM */ - -long ran_x[KK]; /* the generator state */ - -#ifdef __STDC__ -void ran_array(long aa[],int n) -#else -void ran_array(aa,n) /* put n new random numbers in aa */ - long *aa; /* destination */ - int n; /* array length (must be at least KK) */ -#endif -{ - register int i,j; - for (j=0;j<KK;j++) aa[j]=ran_x[j]; - for (;j<n;j++) aa[j]=mod_diff(aa[j-KK],aa[j-LL]); - for (i=0;i<LL;i++,j++) ran_x[i]=mod_diff(aa[j-KK],aa[j-LL]); - for (;i<KK;i++,j++) ran_x[i]=mod_diff(aa[j-KK],ran_x[i-LL]); -} - -/* the following routines are from exercise 3.6--15 */ -/* after calling ran_start, get new randoms by, e.g., "x=ran_arr_next()" */ - -#define QUALITY 1009 /* recommended quality level for high-res use */ -long ran_arr_buf[QUALITY]; -long ran_arr_dummy=-1, ran_arr_started=-1; -long *ran_arr_ptr=&ran_arr_dummy; /* the next random number, or -1 */ - -#define TT 70 /* guaranteed separation between streams */ -#define is_odd(x) ((x)&1) /* units bit of x */ - -#ifdef __STDC__ -void ran_start(long seed) -#else -void ran_start(seed) /* do this before using ran_array */ - long seed; /* selector for different streams */ -#endif -{ - register int t,j; - long x[KK+KK-1]; /* the preparation buffer */ - register long ss=(seed+2)&(MM-2); - for (j=0;j<KK;j++) { - x[j]=ss; /* bootstrap the buffer */ - ss<<=1; if (ss>=MM) ss-=MM-2; /* cyclic shift 29 bits */ - } - x[1]++; /* make x[1] (and only x[1]) odd */ - for (ss=seed&(MM-1),t=TT-1; t; ) { - for (j=KK-1;j>0;j--) x[j+j]=x[j], x[j+j-1]=0; /* "square" */ - for (j=KK+KK-2;j>=KK;j--) - x[j-(KK-LL)]=mod_diff(x[j-(KK-LL)],x[j]), - x[j-KK]=mod_diff(x[j-KK],x[j]); - if (is_odd(ss)) { /* "multiply by z" */ - for (j=KK;j>0;j--) x[j]=x[j-1]; - x[0]=x[KK]; /* shift the buffer cyclically */ - x[LL]=mod_diff(x[LL],x[KK]); - } - if (ss) ss>>=1; else t--; - } - for (j=0;j<LL;j++) ran_x[j+KK-LL]=x[j]; - for (;j<KK;j++) ran_x[j-LL]=x[j]; - for (j=0;j<10;j++) ran_array(x,KK+KK-1); /* warm things up */ - ran_arr_ptr=&ran_arr_started; -} - -#define ran_arr_next() (*ran_arr_ptr>=0? *ran_arr_ptr++: ran_arr_cycle()) -long ran_arr_cycle() -{ - if (ran_arr_ptr==&ran_arr_dummy) - ran_start(314159L); /* the user forgot to initialize */ - ran_array(ran_arr_buf,QUALITY); - ran_arr_buf[KK]=-1; - ran_arr_ptr=ran_arr_buf+1; - return ran_arr_buf[0]; -} - -#include <stdio.h> -int main() -{ - register int m; long a[2009]; - ran_start(310952L); - for (m=0;m<=2009;m++) ran_array(a,1009); - printf("%ld\n", a[0]); /* 995235265 */ - ran_start(310952L); - for (m=0;m<=1009;m++) ran_array(a,2009); - printf("%ld\n", a[0]); /* 995235265 */ - return 0; -} |
