home *** CD-ROM | disk | FTP | other *** search
/ The Fred Fish Collection 1.5 / ffcollection-1-5-1992-11.iso / ff_progs / libs / matlab.lzh / MATLAB / MATLAB.LZH / Source / MatLab / URAND.FOR < prev    next >
Encoding:
Text File  |  1991-04-13  |  1.7 KB  |  58 lines
  1.       DOUBLE PRECISION FUNCTION URAND (IY)
  2.       IMPLICIT NONE
  3. C
  4. C UNIFORM RANDOM NUMBER GENERATOR BASED ON THEORY AND SUGGESTIONS GIVEN IN
  5. C  D.E. KNUTH (1969), VOL 2.  THE INTEGER IY SHOULD BE INITIALIZED TO AN
  6. C  ARBITRARY INTEGER PRIOR TO THE FIRST CALL TO URAND.  THE CALLING PROGRAM
  7. C  SHOULD NOT ALTER THE VALUE OF IY BETWEEN SUBSEQUENT CALLS TO URAND.
  8. C  VALUES OF URAND WILL BE RETURNED IN THE INTERVAL (0,1).
  9. C
  10.       INTEGER IY
  11. C
  12.       INTEGER IA, IC, ITWO, M2, M, MIC
  13.       DOUBLE PRECISION HALFM, S
  14. C
  15.       DOUBLE PRECISION DATAN, DSQRT, DFLOAT
  16. C
  17.       DATA M2 / 0 /, ITWO / 2 /
  18. C
  19. C
  20.       IF (M2.NE.0) GO TO 20
  21. C
  22. C ***      IF FIRST ENTRY, COMPUTE MACHINE INTEGER WORD LENGTH
  23.       M = 1
  24. 10    CONTINUE
  25.       M2 = M
  26.       M = ITWO*M2
  27.       IF (M.GT.M2) GO TO 10
  28.       HALFM = M2
  29. C
  30. C ***      COMPUTE MULTIPLIER AND INCREMENT FOR LINEAR CONGRUENTIAL METHOD
  31.       IA = 8*IDINT (HALFM*DATAN (1.D0)/8.D0)+5
  32.       IC = 2*IDINT (HALFM*(0.5D0-DSQRT (3.D0)/6.D0))+1
  33.       MIC = (M2-IC)+M2
  34. C
  35. C ***      S IS THE SCALE FACTOR FOR CONVERTING TO FLOATING POINT
  36.       S = 0.5D0/HALFM
  37. C
  38. C ***      COMPUTE NEXT RANDOM NUMBER
  39. 20    CONTINUE
  40.       IY = IY*IA
  41. C
  42. C ***      THE FOLLOWING STATEMENT IS FOR COMPUTERS WHICH DO NOT ALLOW
  43. C           INTEGER OVERFLOW ON ADDITION
  44.       IF (IY.GT.MIC) IY = (IY-M2)-M2
  45.       IY = IY+IC
  46. C
  47. C ***      THE FOLLOWING STATEMENT IS FOR COMPUTERS WHERE THE WORD LENGTH
  48. C           FOR ADDITION IS GREATER THAN FOR MULTIPLICATION
  49.       IF (IY/2.GT.M2) IY = (IY-M2)-M2
  50. C
  51. C ***      THE FOLLOWING STATEMENT IS FOR COMPUTERS WHERE INTEGER OVERFLOW
  52. C           AFFECTS THE SIGN BIT
  53.       IF (IY.LT.0) IY = (IY+M2)+M2
  54.       URAND = DFLOAT (IY)*S
  55. C
  56.       RETURN
  57.       END
  58.