mtprng_rand64 Function

public function mtprng_rand64(state) result(r)

NAME

mtprng_rand64(3f) - [M_random:MERSENNE TWISTER] Obtain the next 64-bit integer in the pseudo-random sequence
(LICENSE:CUSTOM OPEN)

SYNOPSIS

function mtprng_rand64(state) result(r)
type(mtprng_state), intent(inout) :: state
integer(INT64) :: r

DESCRIPTION

Obtain the next 64-bit integer in the pseudo-random sequence in the range 0 to 2^32-1.
Note that the range is considerably below the value of HUGE(0_int64).

OPTIONS

state  generator state initialized by mtprng_init(3f) or mtprng_init_array(3f)

RETURNS

r      next pseudo-random value in the range 0 to 2^32-1

EXAMPLE

Sample program:

program demo_mtprng_rand64
use M_random, only : mtprng_state, mtprng_init, mtprng_rand64
use, intrinsic :: iso_fortran_env, only : int32, int64
implicit none
integer(INT32) :: seed
type(mtprng_state) :: state
  GET_SEED: block
  integer :: count
  integer :: count_rate
     call system_clock(count, count_rate)
  seed = count
  endblock GET_SEED
  call mtprng_init(seed, state)
  write(*,*) mtprng_rand64(state)
end program demo_mtprng_rand64

!r = ieor(r,ishft(r,-11))

Arguments

Type IntentOptional Attributes Name
type(mtprng_state), intent(inout) :: state

Return Value integer(kind=INT64)


Source Code

function mtprng_rand64(state) result(r)

! ident_10="@(#) M_random mtprng_rand64(3f) Obtain the next 64-bit integer in the pseudo-random sequence"

! arguments
type(mtprng_state), intent(inout) :: state
!return type
integer(INT64) :: r

   ! internal constants
   integer(INT64), dimension(0:1), parameter :: mag01 = (/ 0_INT64, -1727483681_INT64 /)

   ! Period parameters
   integer(INT64), parameter :: UPPER_MASK =  2147483648_INT64
   integer(INT64), parameter :: LOWER_MASK =  2147483647_INT64

   ! Tempering parameters
   integer(INT64), parameter :: TEMPERING_B = -1658038656_INT64
   integer(INT64), parameter :: TEMPERING_C =  -272236544_INT64

   ! Note: variable names match those in original example
   integer(INT32) :: kk

   ! Generate N words at a time
   if (state%mti >= N) then
      ! The value -1 acts as a flag saying that the seed has not been set.
      if (state%mti == -1) call mtprng_init(4357_INT32,state)

      ! Fill the mt array
      do kk = 0, N - M - 1
         r = ior(iand(state%mt(kk),UPPER_MASK),iand(state%mt(kk+1),LOWER_MASK))
         state%mt(kk) = ieor(ieor(state%mt(kk + M),ishft(r,-1_INT64)),mag01(iand(r,1_INT64)))
      enddo

      do kk = N - M, N - 2
         r = ior(iand(state%mt(kk),UPPER_MASK),iand(state%mt(kk+1),LOWER_MASK))
         state%mt(kk) = ieor(ieor(state%mt(kk + (M - N)),ishft(r,-1_INT64)),mag01(iand(r,1_INT64)))
      enddo

      r = ior(iand(state%mt(N-1),UPPER_MASK),iand(state%mt(0),LOWER_MASK))
      state%mt(N-1) = ieor(ieor(state%mt(M-1),ishft(r,-1)),mag01(iand(r,1_INT64)))

      ! Start using the array from first element
      state%mti = 0
   endif

   ! Here is where we actually calculate the number with a series of
   !   transformations
   r = state%mt(state%mti)
   state%mti = state%mti + 1

 !-------------------------
 !*!r = ieor(r,ishft(r,-11))
   r = ieor(r,ishft(iand(4294967295_INT64,r),-11)) ! Added a 32-bit mask to first r shift
 !-------------------------

   r = iand(4294967295_INT64,ieor(r,iand(ishft(r, 7),TEMPERING_B)))
   r = iand(4294967295_INT64,ieor(r,iand(ishft(r,15),TEMPERING_C)))
   r = ieor(r,ishft(r,-18))

end function mtprng_rand64