random_numbers.F90

!***********************************************************************
!*                   GNU Lesser General Public License
!*
!* This file is part of the GFDL Flexible Modeling System (FMS).
!*
!* FMS is free software: you can redistribute it and/or modify it under
!* the terms of the GNU Lesser General Public License as published by
!* the Free Software Foundation, either version 3 of the License, or (at
!* your option) any later version.
!*
!* FMS is distributed in the hope that it will be useful, but WITHOUT
!* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
!* FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
!* for more details.
!*
!* You should have received a copy of the GNU Lesser General Public
!* License along with FMS.  If not, see <http://www.gnu.org/licenses/>.
!***********************************************************************
 
!> @defgroup random_numbers_mod random_numbers_mod
!> @ingroup random_numbers
!> @brief Generic module to wrap random number generators.
!!
!> The module defines a type that identifies the particular stream of random
!!  numbers, and has procedures for initializing it and getting real numbers
!!  in the range 0 to 1.
!!  This version uses the Mersenne Twister to generate random numbers on [0, 1].
 
module random_numbers_mod
 
  use mersennetwister_mod, only: randomnumbersequence, & ! The random number engine.
                                 new_randomnumbersequence, getrandomreal
  use time_manager_mod, only: time_type, get_date
  use platform_mod, only: r4_kind, r8_kind
 
  implicit none
  private
 
  !> @brief Type to hold a stream of randomly generated numbers
  !> @ingroup random_numbers_mod
  type randomnumberstream
    type(randomNumberSequence) :: theNumbers
  end type randomnumberstream
 
  !> Returns scalar, 1 or 2 D random real numbers
  !!
  !> @param stream @ref randomNumberStream to generate from
  !> @param[out] number output number(s)
  !> @ingroup random_numbers_mod
  interface getrandomnumbers
    module procedure :: get_random_number_0d_r4, get_random_number_0d_r8
    module procedure :: get_random_number_1d_r4, get_random_number_1d_r8
    module procedure :: get_random_number_2d_r4, get_random_number_2d_r8
  end interface getrandomnumbers
 
  !> Initializes stream for generating random numbers.
  !> @ingroup random_numbers_mod
  interface initializerandomnumberstream
    module procedure initializerandomnumberstream_s, initializerandomnumberstream_v
  end interface initializerandomnumberstream
 
  public :: randomnumberstream,                             &
            initializerandomnumberstream, getrandomnumbers, &
            constructseed
 
!> @addtogroup random_numbers_mod
!> @{
 
contains
 
  !> Initialization
  function initializerandomnumberstream_s(seed) result(new)
    integer, intent( in)     :: seed
    type(randomnumberstream) :: new
 
    new%theNumbers = new_randomnumbersequence(seed)
 
  end function initializerandomnumberstream_s
 
  function initializerandomnumberstream_v(seed) result(new)
    integer, dimension(:), intent( in) :: seed
    type(randomnumberstream)           :: new
 
    new%theNumbers = new_randomnumbersequence(seed)
  end function initializerandomnumberstream_v
 
  !> Constructs a unique seed from grid cell index and model date/time
  !!   The perm is supplied we generate a different seed by
  !!   circularly shifting the bits of the seed - this is useful
  !!   if we want to create more than one seed for a given
  !!   column and model date/time.
  !!   Note that abs(perm) must be <= the number of bits used
  !!   to represent the default integer (likely 32)
  function constructseed(i, j, time, perm) result(seed)
    integer,           intent( in)  :: i, j
    type(time_type),   intent( in) :: time
    integer, optional, intent( in) :: perm
    integer, dimension(8) :: seed
 
    ! Local variables
    integer :: year, month, day, hour, minute, second
 
 
    call get_date(time, year, month, day, hour, minute, second)
    seed = (/ i, j, year, month, day, hour, minute, second /)
    if(present(perm)) seed = ishftc(seed, perm)
  end function constructseed
 
#include "random_numbers_r4.fh"
#include "random_numbers_r8.fh"
 
end module random_numbers_mod
 
!> @}
! close documentation grouping