home *** CD-ROM | disk | FTP | other *** search
/ BMUG Revelations / BMUG Revelations.toast / Programming / Programming Utilities / Randoms / Randoms.p < prev    next >
Encoding:
Text File  |  1991-09-21  |  6.6 KB  |  250 lines  |  [TEXT/PJMM]
  1. unit Randoms;
  2.  
  3.  
  4. {This code given to me courtesy of my Operating Systems instructor, Gerald B. Blanton.}
  5. {'Liberated' from the MS-dos world on 9/20/91 by David W. Bock}
  6. {If you use these routines, I'd like to hear about it!  Drop me E-Mail at}
  7. {David Bock or IC Dave on America Online, BOCKD@ITHACA on Bitnet, or snail mail at:}
  8. {Fuzzy Navel Software}
  9. {PO Box 862}
  10. {Great Falls, VA 22066}
  11. {}
  12. {Thanks!!!  (And I'd appreciate any credit you could give me or my instructor in your}
  13. {docs or 'About...' dialog.)}
  14. {}
  15. {PROGRAMMERS!   If you have any good sample code, RELEASE IT!  That's what I'm}
  16. {doing...  I'd like to create an atmosphere where mac programmers help each other}
  17. {out.  I'm not asking you to give away any proprietary secrets, but if you have a clever}
  18. {little routine or a better mouse trap, Release it...  I'm interested in creating a}
  19. {P/D Library of sample code snippets.  If you have something you'd like to ad, send it to me}
  20. {or tell me about it.  You can reach me at any of the addresses above. - Thanks!}
  21. {                                                                                                               -db}
  22.  
  23.  
  24. {-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-}
  25. {the random number routines - uses random number generator from CACM.}
  26. {Includes the user distribution routines (uniform distribution, exponential}
  27. {distribution and normal distribution.}
  28. {  The Random Number (rn) routine uses a byte argument to select one of}
  29. {8 possible seeds from the ran array.  All user distribution routines use this}
  30. {same convention.  rn returns a real value between 0 and 1.}
  31. {-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-}
  32. interface
  33.  
  34.     const
  35.         NUMRANDOMS = 8;
  36.     var
  37.         ran: array[1..NUMRANDOMS] of longint;
  38.         norm: array[1..79] of record
  39.                 z, cp: real;
  40.             end;
  41.  
  42.  
  43.  
  44.     procedure InitRandoms;
  45. {Call this routine before using any of the three functions below.}
  46. {data structures are set up and globals are initialized.}
  47.  
  48.     function Uniform (low, hi, rnIndex: integer): Integer;
  49.  
  50. {a standard Random Number Generator.  when passed integers for low and high,}
  51. {a number between low and high will be returned (with psedo-equal probability)}
  52. { rnIndex is a number from 1 to 8 and is used as a seed. (it actually indexes }
  53. {an array of seeds below.}
  54.  
  55.     function Exponent (mean: real; rnIndex: integer): Integer;
  56.  
  57. {a Random Number Generator that passes back an integer.  The probability of an}
  58. {integer coming back is on the exponential curve with the mean passed in 'mean'.}
  59. {(see the sample program... this is a hard one to explain.)  rnIndex is used as above.}
  60.  
  61.     function Normal (mean, stdDev, rnIndex: integer): Integer;
  62.  
  63. {a Random number generator that passes back an integer.  Passed a mean and a standard}
  64. {deviation, the probability of a certain integer coming back is drawn by a bell curve}
  65. {around the mean. Standard deviation controls the 'width' of the bell.  (again, see the}
  66. {sample program...)  enIndex is used as above.}
  67.  
  68. implementation
  69.  
  70.     function rn (ranNum: Byte): Real;
  71.         const
  72.             a = 16807;
  73.             m = 2147483647;
  74.             q = 127773;
  75.             r = 3826;
  76.         var
  77.             lo, hi, test: Longint;
  78.  
  79.     begin
  80.         hi := ran[ranNum] div q;
  81.         lo := ran[ranNum] mod q;
  82.         test := a * lo - r * hi;
  83.         if test > 0 then
  84.             ran[ranNum] := test
  85.         else
  86.             ran[ranNum] := test + m;
  87.         rn := ran[ranNum] / m;
  88.     end;
  89.  
  90. {-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-}
  91.     function Uniform (low, hi, rnIndex: integer): Integer;
  92.  
  93.  
  94.     begin
  95.         Uniform := trunc(low + (hi - low + 1) * rn(rnIndex));
  96.     end;
  97.  
  98. {-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-}
  99.     function Exponent (mean: real; rnIndex: integer): Integer;
  100.  
  101.     begin
  102.         Exponent := trunc(mean * (-ln(1 - rn(rnIndex))));
  103.     end;
  104.  
  105. {-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-}
  106.     procedure InitNorm;
  107.         var
  108.             i: integer;
  109.  
  110.     begin
  111.         norm[40].z := 0.0;
  112.         norm[40].cp := 0.5;
  113.         norm[41].z := 0.1;
  114.         norm[41].cp := 0.53983;
  115.         norm[42].z := 0.2;
  116.         norm[42].cp := 0.57926;
  117.         norm[43].z := 0.3;
  118.         norm[43].cp := 0.61791;
  119.         norm[44].z := 0.4;
  120.         norm[44].cp := 0.65542;
  121.         norm[45].z := 0.5;
  122.         norm[45].cp := 0.69146;
  123.         norm[46].z := 0.6;
  124.         norm[46].cp := 0.72575;
  125.         norm[47].z := 0.7;
  126.         norm[47].cp := 0.75804;
  127.         norm[48].z := 0.8;
  128.         norm[48].cp := 0.78814;
  129.         norm[49].z := 0.9;
  130.         norm[49].cp := 0.81594;
  131.  
  132.         norm[50].z := 1.0;
  133.         norm[50].cp := 0.84134;
  134.         norm[51].z := 1.1;
  135.         norm[51].cp := 0.86433;
  136.         norm[52].z := 1.2;
  137.         norm[52].cp := 0.88493;
  138.         norm[53].z := 1.3;
  139.         norm[53].cp := 0.90320;
  140.         norm[54].z := 1.4;
  141.         norm[54].cp := 0.91924;
  142.         norm[55].z := 1.5;
  143.         norm[55].cp := 0.93319;
  144.         norm[56].z := 1.6;
  145.         norm[56].cp := 0.94520;
  146.         norm[57].z := 1.7;
  147.         norm[57].cp := 0.95543;
  148.         norm[58].z := 1.8;
  149.         norm[58].cp := 0.96407;
  150.         norm[59].z := 1.9;
  151.         norm[59].cp := 0.97128;
  152.  
  153.         norm[60].z := 2.0;
  154.         norm[60].cp := 0.97725;
  155.         norm[61].z := 2.1;
  156.         norm[61].cp := 0.98214;
  157.         norm[62].z := 2.2;
  158.         norm[62].cp := 0.98610;
  159.         norm[63].z := 2.3;
  160.         norm[63].cp := 0.98928;
  161.         norm[64].z := 2.4;
  162.         norm[64].cp := 0.99180;
  163.         norm[65].z := 2.5;
  164.         norm[65].cp := 0.99379;
  165.         norm[66].z := 2.6;
  166.         norm[66].cp := 0.99534;
  167.         norm[67].z := 2.7;
  168.         norm[67].cp := 0.99653;
  169.         norm[68].z := 2.8;
  170.         norm[68].cp := 0.99744;
  171.         norm[69].z := 2.9;
  172.         norm[69].cp := 0.99813;
  173.  
  174.         norm[70].z := 3.0;
  175.         norm[70].cp := 0.99865;
  176.         norm[71].z := 3.1;
  177.         norm[71].cp := 0.99903;
  178.         norm[72].z := 3.2;
  179.         norm[72].cp := 0.99931;
  180.         norm[73].z := 3.3;
  181.         norm[73].cp := 0.99952;
  182.         norm[74].z := 3.4;
  183.         norm[74].cp := 0.99966;
  184.         norm[75].z := 3.5;
  185.         norm[75].cp := 0.99977;
  186.         norm[76].z := 3.6;
  187.         norm[76].cp := 0.99984;
  188.         norm[77].z := 3.7;
  189.         norm[77].cp := 0.99989;
  190.         norm[78].z := 3.8;
  191.         norm[78].cp := 0.99993;
  192.         norm[79].z := 3.9;
  193.         norm[79].cp := 0.99995;
  194.  
  195.         for i := 1 to 39 do
  196.             begin
  197.                 norm[i].z := -norm[80 - i].z;
  198.                 norm[i].cp := 1.0 - norm[80 - i].cp;
  199.             end;
  200.     end;
  201.  
  202. {-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-}
  203.     function Normal (mean, stdDev, rnIndex: integer): Integer;
  204.  
  205.  
  206.         function GetZ (rnIndex: integer): real;
  207.             var
  208.                 lo, hi: integer;
  209.                 rancp: Real;
  210.  
  211.         begin
  212.             rancp := rn(rnIndex);
  213.             if rancp < norm[1].cp then
  214.                 GetZ := -4.0
  215.             else if rancp > norm[79].cp then
  216.                 GetZ := 4.0
  217.             else
  218.                 begin
  219.                     lo := 1;
  220.                     hi := 79;
  221.                     while hi - lo > 1 do
  222.                         if rancp < norm[(hi + lo) div 2].cp then
  223.                             hi := (hi + lo) div 2
  224.                         else
  225.                             lo := (hi + lo) div 2;
  226.                     GetZ := norm[lo].z;
  227.                 end;
  228.         end;
  229.  
  230.     begin
  231.         Normal := trunc(GetZ(rnIndex) * stdDev + mean);
  232.     end;
  233.  
  234. {-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-}
  235.  
  236.  
  237.     procedure InitRandoms;
  238.     begin
  239.         InitNorm;
  240.         ran[1] := 37584381;
  241.         ran[2] := 1909996635;
  242.         ran[3] := 1964463183;
  243.         ran[4] := 1235671459;
  244.         ran[5] := 1480745561;
  245.         ran[6] := 442596621;
  246.         ran[7] := 340029185;
  247.         ran[8] := 2030226625;
  248.     end;
  249.  
  250. end.