Annotation of capa/capa51/pProj/com.c, revision 1.1
1.1 ! albertel 1: #include "ranlib.h"
! 2: #include <stdio.h>
! 3: #include <stdlib.h>
! 4: void advnst(long k)
! 5: /*
! 6: **********************************************************************
! 7: void advnst(long k)
! 8: ADV-a-N-ce ST-ate
! 9: Advances the state of the current generator by 2^K values and
! 10: resets the initial seed to that value.
! 11: This is a transcription from Pascal to Fortran of routine
! 12: Advance_State from the paper
! 13: L'Ecuyer, P. and Cote, S. "Implementing a Random Number Package
! 14: with Splitting Facilities." ACM Transactions on Mathematical
! 15: Software, 17:98-111 (1991)
! 16: Arguments
! 17: k -> The generator is advanced by2^K values
! 18: **********************************************************************
! 19: */
! 20: {
! 21: #define numg 32L
! 22: extern void gsrgs(long getset,long *qvalue);
! 23: extern void gscgn(long getset,long *g);
! 24: extern long Xm1,Xm2,Xa1,Xa2,Xcg1[],Xcg2[];
! 25: static long g,i,ib1,ib2;
! 26: static long qrgnin;
! 27: /*
! 28: Abort unless random number generator initialized
! 29: */
! 30: gsrgs(0L,&qrgnin);
! 31: if(qrgnin) goto S10;
! 32: puts(" ADVNST called before random generator initialized - ABORT");
! 33: exit(1);
! 34: S10:
! 35: gscgn(0L,&g);
! 36: ib1 = Xa1;
! 37: ib2 = Xa2;
! 38: for(i=1; i<=k; i++) {
! 39: ib1 = mltmod(ib1,ib1,Xm1);
! 40: ib2 = mltmod(ib2,ib2,Xm2);
! 41: }
! 42: setsd(mltmod(ib1,*(Xcg1+g-1),Xm1),mltmod(ib2,*(Xcg2+g-1),Xm2));
! 43: /*
! 44: NOW, IB1 = A1**K AND IB2 = A2**K
! 45: */
! 46: #undef numg
! 47: }
! 48: void getsd(long *iseed1,long *iseed2)
! 49: /*
! 50: **********************************************************************
! 51: void getsd(long *iseed1,long *iseed2)
! 52: GET SeeD
! 53: Returns the value of two integer seeds of the current generator
! 54: This is a transcription from Pascal to Fortran of routine
! 55: Get_State from the paper
! 56: L'Ecuyer, P. and Cote, S. "Implementing a Random Number Package
! 57: with Splitting Facilities." ACM Transactions on Mathematical
! 58: Software, 17:98-111 (1991)
! 59: Arguments
! 60: iseed1 <- First integer seed of generator G
! 61: iseed2 <- Second integer seed of generator G
! 62: **********************************************************************
! 63: */
! 64: {
! 65: #define numg 32L
! 66: extern void gsrgs(long getset,long *qvalue);
! 67: extern void gscgn(long getset,long *g);
! 68: extern long Xcg1[],Xcg2[];
! 69: static long g;
! 70: static long qrgnin;
! 71: /*
! 72: Abort unless random number generator initialized
! 73: */
! 74: gsrgs(0L,&qrgnin);
! 75: if(qrgnin) goto S10;
! 76: printf(
! 77: " GETSD called before random number generator initialized -- abort!\n");
! 78: exit(0);
! 79: S10:
! 80: gscgn(0L,&g);
! 81: *iseed1 = *(Xcg1+g-1);
! 82: *iseed2 = *(Xcg2+g-1);
! 83: #undef numg
! 84: }
! 85: long ignlgi(void)
! 86: /*
! 87: **********************************************************************
! 88: long ignlgi(void)
! 89: GeNerate LarGe Integer
! 90: Returns a random integer following a uniform distribution over
! 91: (1, 2147483562) using the current generator.
! 92: This is a transcription from Pascal to Fortran of routine
! 93: Random from the paper
! 94: L'Ecuyer, P. and Cote, S. "Implementing a Random Number Package
! 95: with Splitting Facilities." ACM Transactions on Mathematical
! 96: Software, 17:98-111 (1991)
! 97: **********************************************************************
! 98: */
! 99: {
! 100: #define numg 32L
! 101: extern void gsrgs(long getset,long *qvalue);
! 102: extern void gssst(long getset,long *qset);
! 103: extern void gscgn(long getset,long *g);
! 104: extern void inrgcm(void);
! 105: extern long Xm1,Xm2,Xa1,Xa2,Xcg1[],Xcg2[];
! 106: extern long Xqanti[];
! 107: static long ignlgi,curntg,k,s1,s2,z;
! 108: static long qqssd,qrgnin;
! 109: /*
! 110: IF THE RANDOM NUMBER PACKAGE HAS NOT BEEN INITIALIZED YET, DO SO.
! 111: IT CAN BE INITIALIZED IN ONE OF TWO WAYS : 1) THE FIRST CALL TO
! 112: THIS ROUTINE 2) A CALL TO SETALL.
! 113: */
! 114: gsrgs(0L,&qrgnin);
! 115: if(!qrgnin) inrgcm();
! 116: gssst(0,&qqssd);
! 117: if(!qqssd) setall(1234567890L,123456789L);
! 118: /*
! 119: Get Current Generator
! 120: */
! 121: gscgn(0L,&curntg);
! 122: s1 = *(Xcg1+curntg-1);
! 123: s2 = *(Xcg2+curntg-1);
! 124: k = s1/53668L;
! 125: s1 = Xa1*(s1-k*53668L)-k*12211;
! 126: if(s1 < 0) s1 += Xm1;
! 127: k = s2/52774L;
! 128: s2 = Xa2*(s2-k*52774L)-k*3791;
! 129: if(s2 < 0) s2 += Xm2;
! 130: *(Xcg1+curntg-1) = s1;
! 131: *(Xcg2+curntg-1) = s2;
! 132: z = s1-s2;
! 133: if(z < 1) z += (Xm1-1);
! 134: if(*(Xqanti+curntg-1)) z = Xm1-z;
! 135: ignlgi = z;
! 136: return ignlgi;
! 137: #undef numg
! 138: }
! 139: void initgn(long isdtyp)
! 140: /*
! 141: **********************************************************************
! 142: void initgn(long isdtyp)
! 143: INIT-ialize current G-e-N-erator
! 144: Reinitializes the state of the current generator
! 145: This is a transcription from Pascal to Fortran of routine
! 146: Init_Generator from the paper
! 147: L'Ecuyer, P. and Cote, S. "Implementing a Random Number Package
! 148: with Splitting Facilities." ACM Transactions on Mathematical
! 149: Software, 17:98-111 (1991)
! 150: Arguments
! 151: isdtyp -> The state to which the generator is to be set
! 152: isdtyp = -1 => sets the seeds to their initial value
! 153: isdtyp = 0 => sets the seeds to the first value of
! 154: the current block
! 155: isdtyp = 1 => sets the seeds to the first value of
! 156: the next block
! 157: **********************************************************************
! 158: */
! 159: {
! 160: #define numg 32L
! 161: extern void gsrgs(long getset,long *qvalue);
! 162: extern void gscgn(long getset,long *g);
! 163: extern long Xm1,Xm2,Xa1w,Xa2w,Xig1[],Xig2[],Xlg1[],Xlg2[],Xcg1[],Xcg2[];
! 164: static long g;
! 165: static long qrgnin;
! 166: /*
! 167: Abort unless random number generator initialized
! 168: */
! 169: gsrgs(0L,&qrgnin);
! 170: if(qrgnin) goto S10;
! 171: printf(
! 172: " INITGN called before random number generator initialized -- abort!\n");
! 173: exit(1);
! 174: S10:
! 175: gscgn(0L,&g);
! 176: if(-1 != isdtyp) goto S20;
! 177: *(Xlg1+g-1) = *(Xig1+g-1);
! 178: *(Xlg2+g-1) = *(Xig2+g-1);
! 179: goto S50;
! 180: S20:
! 181: if(0 != isdtyp) goto S30;
! 182: goto S50;
! 183: S30:
! 184: /*
! 185: do nothing
! 186: */
! 187: if(1 != isdtyp) goto S40;
! 188: *(Xlg1+g-1) = mltmod(Xa1w,*(Xlg1+g-1),Xm1);
! 189: *(Xlg2+g-1) = mltmod(Xa2w,*(Xlg2+g-1),Xm2);
! 190: goto S50;
! 191: S40:
! 192: printf("isdtyp not in range in INITGN");
! 193: exit(1);
! 194: S50:
! 195: *(Xcg1+g-1) = *(Xlg1+g-1);
! 196: *(Xcg2+g-1) = *(Xlg2+g-1);
! 197: #undef numg
! 198: }
! 199: void inrgcm(void)
! 200: /*
! 201: **********************************************************************
! 202: void inrgcm(void)
! 203: INitialize Random number Generator CoMmon
! 204: Function
! 205: Initializes common area for random number generator. This saves
! 206: the nuisance of a BLOCK DATA routine and the difficulty of
! 207: assuring that the routine is loaded with the other routines.
! 208: **********************************************************************
! 209: */
! 210: {
! 211: #define numg 32L
! 212: extern void gsrgs(long getset,long *qvalue);
! 213: extern long Xm1,Xm2,Xa1,Xa2,Xa1w,Xa2w,Xa1vw,Xa2vw;
! 214: extern long Xqanti[];
! 215: static long T1;
! 216: static long i;
! 217: /*
! 218: V=20; W=30;
! 219: A1W = MOD(A1**(2**W),M1) A2W = MOD(A2**(2**W),M2)
! 220: A1VW = MOD(A1**(2**(V+W)),M1) A2VW = MOD(A2**(2**(V+W)),M2)
! 221: If V or W is changed A1W, A2W, A1VW, and A2VW need to be recomputed.
! 222: An efficient way to precompute a**(2*j) MOD m is to start with
! 223: a and square it j times modulo m using the function MLTMOD.
! 224: */
! 225: Xm1 = 2147483563L;
! 226: Xm2 = 2147483399L;
! 227: Xa1 = 40014L;
! 228: Xa2 = 40692L;
! 229: Xa1w = 1033780774L;
! 230: Xa2w = 1494757890L;
! 231: Xa1vw = 2082007225L;
! 232: Xa2vw = 784306273L;
! 233: for(i=0; i<numg; i++) *(Xqanti+i) = 0;
! 234: T1 = 1;
! 235: /*
! 236: Tell the world that common has been initialized
! 237: */
! 238: gsrgs(1L,&T1);
! 239: #undef numg
! 240: }
! 241: void setall(long iseed1,long iseed2)
! 242: /*
! 243: **********************************************************************
! 244: void setall(long iseed1,long iseed2)
! 245: SET ALL random number generators
! 246: Sets the initial seed of generator 1 to ISEED1 and ISEED2. The
! 247: initial seeds of the other generators are set accordingly, and
! 248: all generators states are set to these seeds.
! 249: This is a transcription from Pascal to Fortran of routine
! 250: Set_Initial_Seed from the paper
! 251: L'Ecuyer, P. and Cote, S. "Implementing a Random Number Package
! 252: with Splitting Facilities." ACM Transactions on Mathematical
! 253: Software, 17:98-111 (1991)
! 254: Arguments
! 255: iseed1 -> First of two integer seeds
! 256: iseed2 -> Second of two integer seeds
! 257: **********************************************************************
! 258: */
! 259: {
! 260: #define numg 32L
! 261: extern void gsrgs(long getset,long *qvalue);
! 262: extern void gssst(long getset,long *qset);
! 263: extern void gscgn(long getset,long *g);
! 264: extern long Xm1,Xm2,Xa1vw,Xa2vw,Xig1[],Xig2[];
! 265: static long T1;
! 266: static long g,ocgn;
! 267: static long qrgnin;
! 268: T1 = 1;
! 269: /*
! 270: TELL IGNLGI, THE ACTUAL NUMBER GENERATOR, THAT THIS ROUTINE
! 271: HAS BEEN CALLED.
! 272: */
! 273: gssst(1,&T1);
! 274: gscgn(0L,&ocgn);
! 275: /*
! 276: Initialize Common Block if Necessary
! 277: */
! 278: gsrgs(0L,&qrgnin);
! 279: if(!qrgnin) inrgcm();
! 280: *Xig1 = iseed1;
! 281: *Xig2 = iseed2;
! 282: initgn(-1L);
! 283: for(g=2; g<=numg; g++) {
! 284: *(Xig1+g-1) = mltmod(Xa1vw,*(Xig1+g-2),Xm1);
! 285: *(Xig2+g-1) = mltmod(Xa2vw,*(Xig2+g-2),Xm2);
! 286: gscgn(1L,&g);
! 287: initgn(-1L);
! 288: }
! 289: gscgn(1L,&ocgn);
! 290: #undef numg
! 291: }
! 292: void setant(long qvalue)
! 293: /*
! 294: **********************************************************************
! 295: void setant(long qvalue)
! 296: SET ANTithetic
! 297: Sets whether the current generator produces antithetic values. If
! 298: X is the value normally returned from a uniform [0,1] random
! 299: number generator then 1 - X is the antithetic value. If X is the
! 300: value normally returned from a uniform [0,N] random number
! 301: generator then N - 1 - X is the antithetic value.
! 302: All generators are initialized to NOT generate antithetic values.
! 303: This is a transcription from Pascal to Fortran of routine
! 304: Set_Antithetic from the paper
! 305: L'Ecuyer, P. and Cote, S. "Implementing a Random Number Package
! 306: with Splitting Facilities." ACM Transactions on Mathematical
! 307: Software, 17:98-111 (1991)
! 308: Arguments
! 309: qvalue -> nonzero if generator G is to generating antithetic
! 310: values, otherwise zero
! 311: **********************************************************************
! 312: */
! 313: {
! 314: #define numg 32L
! 315: extern void gsrgs(long getset,long *qvalue);
! 316: extern void gscgn(long getset,long *g);
! 317: extern long Xqanti[];
! 318: static long g;
! 319: static long qrgnin;
! 320: /*
! 321: Abort unless random number generator initialized
! 322: */
! 323: gsrgs(0L,&qrgnin);
! 324: if(qrgnin) goto S10;
! 325: printf(
! 326: " SETANT called before random number generator initialized -- abort!\n");
! 327: exit(1);
! 328: S10:
! 329: gscgn(0L,&g);
! 330: Xqanti[g-1] = qvalue;
! 331: #undef numg
! 332: }
! 333: void setsd(long iseed1,long iseed2)
! 334: /*
! 335: **********************************************************************
! 336: void setsd(long iseed1,long iseed2)
! 337: SET S-ee-D of current generator
! 338: Resets the initial seed of the current generator to ISEED1 and
! 339: ISEED2. The seeds of the other generators remain unchanged.
! 340: This is a transcription from Pascal to Fortran of routine
! 341: Set_Seed from the paper
! 342: L'Ecuyer, P. and Cote, S. "Implementing a Random Number Package
! 343: with Splitting Facilities." ACM Transactions on Mathematical
! 344: Software, 17:98-111 (1991)
! 345: Arguments
! 346: iseed1 -> First integer seed
! 347: iseed2 -> Second integer seed
! 348: **********************************************************************
! 349: */
! 350: {
! 351: #define numg 32L
! 352: extern void gsrgs(long getset,long *qvalue);
! 353: extern void gscgn(long getset,long *g);
! 354: extern long Xig1[],Xig2[];
! 355: static long g;
! 356: static long qrgnin;
! 357: /*
! 358: Abort unless random number generator initialized
! 359: */
! 360: gsrgs(0L,&qrgnin);
! 361: if(qrgnin) goto S10;
! 362: printf(
! 363: " SETSD called before random number generator initialized -- abort!\n");
! 364: exit(1);
! 365: S10:
! 366: gscgn(0L,&g);
! 367: *(Xig1+g-1) = iseed1;
! 368: *(Xig2+g-1) = iseed2;
! 369: initgn(-1L);
! 370: #undef numg
! 371: }
! 372: long Xm1,Xm2,Xa1,Xa2,Xcg1[32],Xcg2[32],Xa1w,Xa2w,Xig1[32],Xig2[32],Xlg1[32],
! 373: Xlg2[32],Xa1vw,Xa2vw;
! 374: long Xqanti[32];
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>
![[BACK]](/icons/back.gif){kind=icon}
Return to com.c CVS log ![[TXT]](/icons/text.gif){kind=icon}
![[DIR]](/icons/dir.gif){kind=icon}
Up to [LON-CAPA] / capa / capa51 / pProj