乱数発生のサンプルプログラム(種の設定方法)

Fortran90で乱数を使用するサンプルプログラムです。乱数の種(random seed)の設定方法が少しわかりづらかったので、微妙に条件を変えて3通りやってみました。

ちょっと変に思ったのは、seedが配列である必要があること。random_seed関数でnrandにサイズを代入していますが、中身を見てみると1なので、seedという一次元配列は要素が一つしかallocateされていません。何か理由があるのでしょうか。

3つのseedの設定はそれぞれ、デフォルト値、システム時間、任意の数、で設定しています。システム時間を使用した場合は呼び出される度に乱数が変わることが分かります。数値計算では再現性が重要になることもあるので、任意に乱数の種を決めたい時もありますよね。そういうときはシステム時間と同じくらい巨大な数にしないと、乱数のばらつきが悪くなるようです。デフォルト値では乱数の種(seed変数)は0になっているので、そのまま使用することはできません。ここでは任意のランダムシードを1234567891234657という数にしていますが、ちょうどいい乱数を出すための特別な数があると聞いたことが…。それはまた今度調べよう。


rand.f90

program rand
  implicit none

  integer(kind=8), allocatable, dimension(:) :: seed
  integer(kind=8) :: nrand
  integer(kind=8), parameter :: SEED_SELF = 1234567891234567
  integer(kind=8) :: clock
  integer(kind=8) :: i
  integer(kind=8), parameter :: NCHECK = 10
  real(kind=8) :: x

  ! initialize
  call random_seed(size=nrand)
  allocate(seed(nrand))
  call system_clock(count=clock)

  ! show parameter
  write(*,*) "## parameters"
  write(*,'(a,I20)') "    clock       = ", clock
  write(*,'(a,I20)') "    SEED_SELF   = ", SEED_SELF
  write(*,'(a,I20)') "    nrand       = ", nrand
  do i = 1, size(seed)
    write(*,'(a,I5,a,I20)') "    seed(", i, ") = ", seed(i)
  enddo
  write(*,'(a,I20)') "    NCHECK      = ", NCHECK

  ! default random seed
  write(*,*) "## default random seed"
  call random_seed(put=seed)

  do i = 1, NCHECK
    call random_number(x)
    write(*, '(I5,a,E15.5)') i, " : ", x
    !print *, i, " : ", x
  enddo

  ! system clock random seed
  write(*,*) "## system clock random seed"
  seed = clock
  call random_seed(put=seed)

  do i = 1, NCHECK
    call random_number(x)
    write(*, '(I5,a,E15.5)') i, " : ", x
    !print *, i, " : ", x
  enddo

  ! personalized random seed
  write(*,*) "## personalized random seed"
  seed = SEED_SELF
  call random_seed(put=seed)

  do i = 1, NCHECK
    call random_number(x)
    write(*, '(I5,a,E15.5)') i, " : ", x
    !print *, i, " : ", x
  enddo

end program rand

実行方法

ifort -o rand rand.f90
./rand

出力結果(例)

 ## parameters
    clock       =     1366181897413064
    SEED_SELF   =     1234567891234567
    nrand       =                    1
    seed(    1) =                    0
    NCHECK      =                   10
 ## default random seed
    1 :     0.39209E-06
    2 :     0.25480E-01
    3 :     0.35252E+00
    4 :     0.66691E+00
    5 :     0.96306E+00
    6 :     0.83829E+00
    7 :     0.33536E+00
    8 :     0.91533E+00
    9 :     0.79586E+00
   10 :     0.83269E+00
 ## system clock random seed
    1 :     0.10319E+00
    2 :     0.56227E+00
    3 :     0.74591E+00
    4 :     0.58723E+00
    5 :     0.33721E+00
    6 :     0.83088E+00
    7 :     0.54183E+00
    8 :     0.43535E-01
    9 :     0.16427E+00
   10 :     0.83179E+00
 ## personalized random seed
    1 :     0.22664E+00
    2 :     0.84194E+00
    3 :     0.12193E+00
    4 :     0.62328E+00
    5 :     0.66848E+00
    6 :     0.25764E+00
    7 :     0.68095E+00
    8 :     0.78163E+00
    9 :     0.17263E+00
   10 :     0.82049E+00