fortran random number seed [fortran] by 바죠

포트란 언어에서 랜덤 넘버를 사용하는 방법이다. 기본적으로 씨드를 먼저 설정해 주고 랜덤넘버를 사용하게 된다. 씨드를 사용자가 설정하지 않는 경우, 디폴트 씨드가 사용되게 된다.

가장 간단한 방법은 call random_seed()를 제일 먼저 불러서 초기화를 시킨다음에 call random_number(tmp) 처럼 사용하는것이다.
이 때, call random_seed()는 가장 초기에 한 번만 부르면된다.

http://incredible.egloos.com/7378157
http://fortranwiki.org/fortran/show/init_seed
http://incredible.egloos.com/4790690
http://incredible.egloos.com/4494190
http://incredible.egloos.com/7378157


Sobol sequence:

http://incredible.egloos.com/7378157



cat seed.f90
   program main
   implicit none
!  integer, parameter :: wp = selected_real_kind(15,307)
   integer, parameter :: ndiscard = 100

   integer :: nstate_size, i
   integer, allocatable, dimension(:) :: istate
!  real(wp) :: ran, oldran
   real*8 :: ran, oldran

   call random_seed(size=nstate_size)
   allocate(istate(nstate_size))
   istate = 20180815
   istate = 20180814
   call random_seed(put=istate)
   deallocate(istate)
!  ran = 0.5_wp
   ran = 0.5d0
   do i=1,ndiscard
      oldran = ran
      call random_number(ran)
   enddo
   end program main



단순하게 초기화하는 방법:
call random_seed()



subroutine init_seed()
  integer :: n, ival(8), v(3), i
  integer, allocatable :: seed(:)
  call date_and_time(values=ival)
  v(1) = ival(8) + 2048*ival(7)
  v(2) = ival(6) + 64*ival(5)     ! value(4) isn't really 'random'
  v(3) = ival(3) + 32*ival(2) + 32*8*ival(1)
  call random_seed(size=n)
  allocate(seed(n))
  call random_seed()   ! Give the seed an implementation-dependent kick
  call random_seed(get=seed)
  do i=1, n
     seed(i) = seed(i) + v(mod(i-1, 3) + 1)
  enddo
  call random_seed(put=seed)
  deallocate(seed)
end subroutine
random number에 대한 이미지 검색결과
--------------------------------------------------------------------------------------------------------------------

!234567890
      program test
      implicit none
      real*8 tmp
      integer i

      call init_seed()

      do i=1,10
      call random_number(tmp)
      write(6,*) tmp
      enddo

      end program test
      subroutine init_seed()
      implicit none
      integer :: n, ival(8), v(3), i
      integer, allocatable :: seed(:)
        call date_and_time(values=ival)
        v(1) = ival(8) + 2048*ival(7)
        v(2) = ival(6) + 64*ival(5)     ! value(4) isn't really 'random'
        v(3) = ival(3) + 32*ival(2) + 32*8*ival(1)
        call random_seed(size=n)
        allocate(seed(n))
        call random_seed()   ! Give the seed an implementation-dependent kick
        call random_seed(get=seed)
        do i=1, n
           seed(i) = seed(i) + v(mod(i-1, 3) + 1)
        enddo
        call random_seed(put=seed)
        deallocate(seed)
      end subroutine


 ./a.out
  0.546180230297762    
  0.696882323005700    
  0.768761350468132    
  0.845658298526402    
  0.882858783958011    
  0.933070317521213    
  0.813091367535648    
  0.639849523728344    
  0.852785819902453    
  0.279603133334902    
[ihlee@master source1]$ ./a.out
  0.323252338672265    
  0.903215112990365    
  0.503255362984121    
  0.804420441564049    
  3.073751675555892E-002
  3.593891768474497E-002
  0.890067805375792    
  0.538413460722725    
  0.354780752750282    
  0.178050204708365    
[ihlee@master source1]$ ./a.out
  0.105318159773966    
  0.924913085353378    
  0.263212870514529    
  0.659453458643305    
  0.561370073685634    
  0.653400796716597    
  0.902523513284742    
  0.704490006380550    
  0.103280202848286    
  0.708269233444185    

--------------------------------------------------------------------------------------------------------------------

          subroutine init_random_seed()
            implicit none
            integer, allocatable :: seed(:)
            integer :: i, n, un, istat, dt(8), pid, t(2), s
            integer(8) :: count, tms
         
            call random_seed(size = n)
            allocate(seed(n))
            ! First try if the OS provides a random number generator
            open(newunit=un, file="/dev/urandom", access="stream", &
                 form="unformatted", action="read", status="old", iostat=istat)
            if (istat == 0) then
               read(un) seed
               close(un)
            else
               ! Fallback to XOR:ing the current time and pid. The PID is
               ! useful in case one launches multiple instances of the same
               ! program in parallel.
               call system_clock(count)
               if (count /= 0) then
                  t = transfer(count, t)
               else
                  call date_and_time(values=dt)
                  tms = (dt(1) - 1970) * 365_8 * 24 * 60 * 60 * 1000 &
                       + dt(2) * 31_8 * 24 * 60 * 60 * 1000 &
                       + dt(3) * 24 * 60 * 60 * 60 * 1000 &
                       + dt(5) * 60 * 60 * 1000 &
                       + dt(6) * 60 * 1000 + dt(7) * 1000 &
                       + dt(8)
                  t = transfer(tms, t)
               end if
               s = ieor(t(1), t(2))
               pid = getpid() + 1099279 ! Add a prime
               s = ieor(s, pid)
               if (n >= 3) then
                  seed(1) = t(1) + 36269
                  seed(2) = t(2) + 72551
                  seed(3) = pid
                  if (n > 3) then
                     seed(4:) = s + 37 * (/ (i, i = 0, n - 4) /)
                  end if
               else
                  seed = s + 37 * (/ (i, i = 0, n - 1 ) /)
               end if
            end if
            call random_seed(put=seed)
          end subroutine init_random_seed


http://gcc.gnu.org/onlinedocs/gfortran/RANDOM_005fSEED.html#RANDOM_005fSEED

subroutine init_random_seed()
integer :: i, n, clock
integer, dimension(:), allocatable :: seed

call random_seed(size = n)
allocate(seed(n))

call system_clock(count=clock)

seed = clock + 37 * (/ (i - 1, i = 1, n) /)
call random_seed(put = seed)

deallocate(seed)
end subroutine init_random_see

program test_nint
real(4) x4
real(8) x8
x4 = 1.234E0_4
x8 = 4.321_8
print *, nint(x4), idnint(x8)
end program test_nint


program test_floor
real :: x = 63.29
real :: y = -63.59
print *, floor(x) ! returns 63
print *, floor(y) ! returns -64
end program test_floor

program test_ceiling
real :: x = 63.29
real :: y = -63.59
print *, ceiling(x) ! returns 64
print *, ceiling(y) ! returns -63
end program test_ceiling


http://gcc.gnu.org/onlinedocs/gfortran/CEILING.html#CEILING
 program test_nint            
real(4) x4
real(8) x8


x4 = 1.234E0_4
x8 = 4.321_8
print *, nint(x4), idnint(x8)
end program test_nint



program test_random_number           
REAL :: r(5,5)         
  
CALL init_random_seed()         ! see example of RANDOM_SEED           
CALL RANDOM_NUMBER(r)     
    
end program
------------------------------------------
 CALL Random_Seed
 CALL RANDOM_NUMBER(rRAND)
------------------------------------------
  subroutine random_init()

    implicit none

    integer :: i, n
    integer :: clock


    i = 0
    call random_seed(size=n)
    allocate(iSeed(n))
    call system_clock(count=clock)
    iSeed = clock + 37 * (/ (i - 1, i = 1, n) /)
    call random_seed(put=iSeed)


  end subroutine random_init
  subroutine random_final()

    implicit none

    deallocate(iSeed)

  end subroutine random_final

-------------------------------------------------------------------------------------------------------------------

1,2,3,4,5, 6 중 하나의 정수를 무작위로 추출하는 방법, 균등 확률:

call random_number(tmp)
i= tmp*6 +1
http://incredible.egloos.com/4790690
http://incredible.egloos.com/4494190
http://incredible.egloos.com/7378157


------------------------------------------------------------------------------------------------------------------
Dice
------------------------------------------------------------------------------------------------------------------
integer  i
real*8  tmp

call random_number(tmp)
i=tmp*10+1

i는 1, 2, 3,..., 10 중의 한 숫자가 됨.
동등한 확률로 분포함.

------------------------------------------------------------------------------------------------------------------
------------------------------------------------------------------------------------------------------------------


핑백

덧글

댓글 입력 영역

최근 포토로그



MathJax