home *** CD-ROM | disk | FTP | other *** search
/ OS/2 Shareware BBS: 22 gnu / 22-gnu.zip / db02_src.zip / rand.cc < prev    next >
C/C++ Source or Header  |  1993-11-05  |  3KB  |  170 lines
  1. //
  2. // Source     : rand.cc
  3. // Author     : Darren Platt
  4. // Version    : 1.1
  5. // Date Begun : 10.5.91
  6. // Last Rev   : 92/05/18
  7. // Id String  : @(#)rand.cc    1.1
  8. // 
  9. // Modification log
  10. //-----------------
  11. // 0.0  tjp        Tom was here.
  12. // 1.0    daz        Expanded to use a class
  13. // 1.1    daz        Added uniformInt and put in SCCS
  14. //
  15. // Portable random number generator.
  16. // Modified lahmer, generates integers in range min_long to max_long
  17. // 1.6.91 Expanded to include a random number class.
  18.  
  19. static char *WhatString = "@(#)rand.cc    1.1";
  20.  
  21. #include "rand.h"
  22. #include <assert.h>
  23. #include <math.h>
  24. #ifndef MALLOC_H_MISSING
  25. #include <malloc.h>
  26. #endif
  27.  
  28. unsigned long     global_table[128];
  29. unsigned long     *rand_table = global_table;
  30. unsigned long     lahmer_seed;
  31. long             taus_seed;
  32.  
  33. void seed_lahmer(long seed)
  34. {
  35.     lahmer_seed = seed;
  36. }
  37. unsigned long lahmer(void)
  38. {
  39.     return lahmer_seed=(69069l*lahmer_seed+1);
  40. }
  41. void fill_table(void)
  42. {
  43.     for(int i=0;i<128;i++)
  44.         rand_table[i] = lahmer();
  45. }
  46. void seed_taus(long seed)
  47. {
  48.     taus_seed = seed;
  49. }
  50. long taus(void)
  51. {
  52.     const long k_value = 123456;
  53.  
  54.     if (taus_seed >=0) return (taus_seed = (taus_seed<<1)%1000003);
  55.     else return taus_seed = ((taus_seed<<1)^k_value)%1000003;
  56. }
  57. void seed_random(unsigned long seed)
  58. {
  59.     seed_lahmer(seed);
  60.     seed_taus((signed long)seed);
  61.     fill_table();
  62. }
  63.  
  64. unsigned long my_rand(void)
  65. {
  66.     int choice = taus()%128;
  67.     if (choice<0) choice = -choice;
  68.     unsigned long ret_value = rand_table[choice];
  69.     rand_table[choice] = lahmer();
  70.     return ret_value;
  71. }
  72. //
  73. // Method definitions for the randcl class
  74. //
  75.  
  76. inline double randcl::quick_unit(void)
  77. {
  78.     return (my_rand()+0.5)/(double)max_ulong;
  79. }
  80. inline void randcl::ungrab(void)
  81. {
  82.     taus_seed = ::taus_seed;
  83.     lahmer_seed = ::lahmer_seed;
  84.     rand_table = ::rand_table;
  85. }
  86. inline void randcl::grab(void)
  87. {
  88.     ::taus_seed = taus_seed;
  89.     ::lahmer_seed = lahmer_seed;
  90.     ::rand_table = rand_table;
  91. }
  92.  
  93. void randcl::seed(unsigned long seed)
  94. {
  95.     grab();
  96.     seed_random(seed);
  97.     ungrab();
  98. }
  99. randcl::randcl(void)
  100. {
  101.     assert(rand_table = new unsigned long[128]);
  102.     seed(1);
  103.     phase = 0;
  104. }
  105. randcl::randcl(unsigned long seed)
  106. {
  107.     assert(rand_table = new unsigned long[128]);
  108.     randcl::seed(seed);
  109.     phase=0;
  110. }
  111. double randcl::neg_exp(double mean)
  112. {
  113.     return -mean*log(unit());
  114. }
  115. double randcl::normal(double mean,double sd)
  116. {
  117.     if (phase) {    // Already have one stored up.
  118.         phase = 0;
  119.         return (sqratio * q * sd)+mean;
  120.     }
  121.     double    p,v;
  122.     grab();
  123.     do {
  124.         p = quick_unit()*2-1; q = quick_unit()*2-1;
  125.         v = p*p + q*q;
  126.     } while(v > 1.0 || v <0.25);
  127.     sqratio = sqrt(-2*log(quick_unit()) / v);
  128.     ungrab();
  129.     phase = 1;
  130.     return (sd * sqratio * p)+mean;
  131. }
  132. long randcl::uniformInt(long lower,long upper)
  133. {
  134.     return (ulong()%(upper-lower+1))+lower;
  135. }
  136.  
  137. double randcl::uniform(double lower,double upper)
  138. {
  139.     return (upper-lower)*unit()+lower;
  140. }
  141. randcl::~randcl()
  142. {
  143.     delete rand_table;
  144. }
  145. double randcl::unit(void)
  146. {
  147.     grab();
  148.     unsigned long    ul = my_rand();
  149.     ungrab();
  150.     return (ul+0.5)/(double)max_ulong;
  151. }
  152. unsigned long randcl::ulong(void)
  153. {
  154.     grab();
  155.     unsigned long ul = my_rand();
  156.     ungrab();
  157.     return ul;
  158. }
  159.  
  160. #ifdef RAND_TEST
  161. #include <iostream.h>
  162. main()
  163. {
  164.     randcl Myrand(21);
  165.     for(int i=0;i<20;i++) {
  166.         cout << Myrand.uniformInt(-1,1) << endl;
  167.     }
  168. }
  169. #endif
  170.