!! uniform : float function (seeds : in out array [1..2] of integer)
!! Version 1.1

!! Purpose:
!! Generate a uniform random number using a combined generator.  The
!! combined generator uses three 16-bit integer values (and updates them
!! each time a new number is generated) stored in global variables.
!!
!! This generator has a period approximately equal to 8.125_44e12.
!!
!! This generator will never return either 0.0 or 1.0, as long as the
!! FLOAT (or REAL) variables have at least 23 bits for the fraction
!! part (commonly, but erroneously, called the mantissa).

!!######################################################################
!!#                                                                    #
!!# The procedure InitU initializes the values of the global variables #
!!# and must be called before Unifrm is called.                        #
!!#                                                                    #
!!######################################################################

!! Source:
!!     Pierre L'Ecuyer, Efficient and Portable Combined Random Number
!!     Generators, Communications of the ACM, Vol. 31, No. 6 (June 1988),
!!     742-749, 774.

!! Author(s):  James J. Fullerton
!! Creation Date:  22 JUN 88
!! Date Last Approved:  22 JUN 88
!! Revisions (Date, Name, and summary of revision(s)):

!! Pierre L'Ecuyer's Pascal code:
!!
!! s1, s2, s3 : integer {# 16-bit integers #}
!!
!! . . .
!!
!! function Uniform : real
!!   ;var z, k : integer  {# 16-bit integers #}
!! ;begin
!!    k := s1 div 206
!!   ;s1 := 157 * (s1 - k * 206) - k * 21
!!   ;if s1 < 0 then s1 := s1 + 32_363
!!
!!   ;k := s2 div 217
!!   ;s2 := 146 * (s2 - k * 217) - k * 45
!!   ;if s2 < 0 then s2 := s2 + 31_727
!!
!!   ;k := s3 div 222
!!   ;s3 := 142 * (s3 - k * 222) - k * 133
!!   ;if s3 < 0 then s3 := s3 + 31_657
!!
!!   ;z := s1 - s2
!!   ;if z > 706 then z := z - 32_362
!!   ;z := z + s3
!!   ;if z < 1 then z := z + 32_362
!!
!!   ;Uniform := z * 3.089_9e-5
!!
!! end { Uniform }

      real*4 function Unifrm (Seeds)

      implicit character*255 (a-z)

      !!----------------------------------------------------------------------
      !! Parameters.
      !!----------------------------------------------------------------------

      !! IN OUT:

      integer*2 Seeds (1:3)	!! Seed values used (and updated) each time
				!! a new random variate is computed.

      !!----------------------------------------------------------------------
      !! Local variables.
      !!----------------------------------------------------------------------

      integer*2 z
      integer*2 k

!! begin

      k= Seeds(1) / 206
      Seeds(1)= 157 * (Seeds(1) - k * 206) - k * 21
      if (Seeds(1) .lt. 0) Seeds(1)= Seeds(1) + 32363

      k= Seeds(2) / 217
      Seeds(2)= 146 * (Seeds(2) - k * 217) - k * 45
      if (Seeds(2) .lt. 0) Seeds(2)= Seeds(2) + 31727

      k= Seeds(3) / 222
      Seeds(3)= 142 * (Seeds(3) - k * 222) - k * 133
      if (Seeds(3) .lt. 0) Seeds(3)= Seeds(3) + 31657

      z= Seeds(1) - Seeds(2)
      if (z .gt. 706) z= z - 32362
      z= z + Seeds(3)
      if (z .lt. 1) z= z + 32362

      Unifrm= z * 3.0899e-5

      return

      end !! Unifrm