home *** CD-ROM | disk | FTP | other *** search
/ GEMini Atari / GEMini_Atari_CD-ROM_Walnut_Creek_December_1993.iso / files / math / ols / ran1.c < prev    next >
Encoding:
C/C++ Source or Header  |  1993-07-28  |  1.2 KB  |  65 lines
  1. #include <stdio.h>
  2. #include <string.h>
  3. #include "utils.h"
  4.  
  5. #define M1 259200
  6. #define IA1 7141
  7. #define IC1 54773
  8. #define RM1 (1.0/M1)
  9. #define M2 134456
  10. #define IA2 8121
  11. #define IC2 28411
  12. #define RM2 (1.0/M2)
  13. #define M3 243000
  14. #define IA3 4561
  15. #define IC3 51349
  16.  
  17. double 
  18. ran1 (idum)
  19.      int *idum;
  20. {
  21.   static long ix1, ix2, ix3;
  22.   static double r[98];
  23.   double temp;
  24.   static int iff = 0;
  25.   int j;
  26.  
  27.   if (*idum < 0 || iff == 0)
  28.     {
  29.       iff = 1;
  30.       ix1 = (IC1 - (*idum)) % M1;
  31.       ix1 = (IA1 * ix1 + IC1) % M1;
  32.       ix2 = ix1 % M2;
  33.       ix1 = (IA1 * ix1 + IC1) % M1;
  34.       ix3 = ix1 % M3;
  35.       for (j = 1; j <= 97; j++)
  36.     {
  37.       ix1 = (IA1 * ix1 + IC1) % M1;
  38.       ix2 = (IA2 * ix2 + IC2) % M2;
  39.       r[j] = (ix1 + ix2 * RM2) * RM1;
  40.     }
  41.       *idum = 1;
  42.     }
  43.   ix1 = (IA1 * ix1 + IC1) % M1;
  44.   ix2 = (IA2 * ix2 + IC2) % M2;
  45.   ix3 = (IA3 * ix3 + IC3) % M3;
  46.   j = 1 + ((97 * ix3) / M3);
  47.   if (j > 97 || j < 1)
  48.     FatalError ("RAN1: This cannot happen.");
  49.   temp = r[j];
  50.   r[j] = (ix1 + ix2 * RM2) * RM1;
  51.   return temp;
  52. }
  53.  
  54. #undef M1
  55. #undef IA1
  56. #undef IC1
  57. #undef RM1
  58. #undef M2
  59. #undef IA2
  60. #undef IC2
  61. #undef RM2
  62. #undef M3
  63. #undef IA3
  64. #undef IC3
  65.