lrand.c 1.2 KB

1234567891011121314151617181920212223242526272829303132333435363738394041424344454647484950515253545556575859606162636465666768697071727374757677787980818283848586
  1. /* lrand.c
  2. *
  3. * Pseudorandom number generator
  4. *
  5. *
  6. *
  7. * SYNOPSIS:
  8. *
  9. * long y, drand();
  10. *
  11. * drand( &y );
  12. *
  13. *
  14. *
  15. * DESCRIPTION:
  16. *
  17. * Yields a long integer random number.
  18. *
  19. * The three-generator congruential algorithm by Brian
  20. * Wichmann and David Hill (BYTE magazine, March, 1987,
  21. * pp 127-8) is used. The period, given by them, is
  22. * 6953607871644.
  23. *
  24. *
  25. */
  26. #include <math.h>
  27. /* Three-generator random number algorithm
  28. * of Brian Wichmann and David Hill
  29. * BYTE magazine, March, 1987 pp 127-8
  30. *
  31. * The period, given by them, is (p-1)(q-1)(r-1)/4 = 6.95e12.
  32. */
  33. static int sx = 1;
  34. static int sy = 10000;
  35. static int sz = 3000;
  36. /* This function implements the three
  37. * congruential generators.
  38. */
  39. long lrand()
  40. {
  41. int r, s;
  42. unsigned long ans;
  43. /*
  44. if( arg )
  45. {
  46. sx = 1;
  47. sy = 10000;
  48. sz = 3000;
  49. }
  50. */
  51. /* sx = sx * 171 mod 30269 */
  52. r = sx/177;
  53. s = sx - 177 * r;
  54. sx = 171 * s - 2 * r;
  55. if( sx < 0 )
  56. sx += 30269;
  57. /* sy = sy * 172 mod 30307 */
  58. r = sy/176;
  59. s = sy - 176 * r;
  60. sy = 172 * s - 35 * r;
  61. if( sy < 0 )
  62. sy += 30307;
  63. /* sz = 170 * sz mod 30323 */
  64. r = sz/178;
  65. s = sz - 178 * r;
  66. sz = 170 * s - 63 * r;
  67. if( sz < 0 )
  68. sz += 30323;
  69. ans = sx * sy * sz;
  70. return(ans);
  71. }