Random_number

[mohammad_b_]

Hi

I use ifort compiler. in my code i need random integer number between 1<= i <=64.

i used random_number(x) and it is between 0<x<1. then for my purpose i use this transformation

i = floor(64*x) + 1

but when i use this transformation sometimes i = 65. after debug my code i find that when x = 0.9999

this mistake is happened.I guess the problem is due to rounding error. how to avoid from this mistake?

[Laura_S_3]

That seems like very unexpected behavior. I find this very concerning since I depend on floor working properly in my programs.

Can you provide a small reproducer? Also, what version of compiler?

[mecej4]

Indeed, you need to provide a reproducer. The only way for floor(64*x) to equal 64 is for x to be equal to 1, and RANDOM_NUMBER is supposed to never return 1; call random_number(x) is supposed to return x such that 0 <= x < 1. Here is a test program that does not use RANDOM_NUMBER that you can use to test your compiler:

program tstflr
implicit none
real x,dx
integer i,ix
!
dx=16.0*epsilon(x)
do
   x = 1.0 -dx
   if(x == 1.0)stop
   if(floor(64*x) >= 64)then
      write(*,*)'When x = ',x,' floor(64*x) exceeds 63'
      stop
   end if
   write(*,'(1x,2ES15.8)')x,dx
   dx=0.5*dx
end do
end program

The output that I get from it:

  9.99998093E-01 1.90734863E-06
  9.99999046E-01 9.53674316E-07
  9.99999523E-01 4.76837158E-07
  9.99999762E-01 2.38418579E-07
  9.99999881E-01 1.19209290E-07
  9.99999940E-01 5.96046448E-08

This shows that x went through several values > 0.9999 before it became indistinguishable from 1.0 in single precision.

[mohammad_b_]

Laura S. wrote:

That seems like very unexpected behavior. I find this very concerning since I depend on floor working properly in my programs.

Can you provide a small reproducer? Also, what version of compiler?

my compiler version is 15.0.0      and  20140723

this is my code and if this mistake don't occur then in line 56 n1 must be printed as 0

but n1 printed as 2

program simple_ising

  implicit none
  integer, parameter :: n = 64
  real, parameter :: t = 2.269
  integer :: s( 1:n , 1:n )
  real :: x, y
  real :: prob
  real :: finish, start
  integer :: i, j, ix, iy, r, l, u, d, mag, m1, mn, delta, en, n1

  open(unit=1, file='a3')
  
  call cpu_time(start)
  
  mag = 0
  en = 0
!   intial state is random
  call random_seed
  do i = 1, n, 1

    do j = 1, n, 1

        s(i,j) = -1

        call random_number(x)

        if ( x .gt. 0.5 ) then
          s(i,j) = 1
        end if

    end do

  end do

  mag = sum(s)
  print*,mag/real(n*n)
    n1=0

  do i = 1, n*n, 1

    do j = 1, n*n, 1

      call random_number(x)
      ix = floor( n*x + 1.0 )

      call random_number(y)
      iy = floor( n*y + 1.0) 

      r = ix + 1 ; l = ix -1 ; u = iy + 1 ; d = iy - 1

      if ( ix + 1 .gt. n ) r = 1
      if ( ix - 1 .lt. 1 ) l = n
      if ( iy + 1 .gt. n ) u = 1
      if ( iy - 1 .lt. 1 ) d = n
      if ( ix .gt. n .or. iy .gt. n) n1=n1+1
      mn = s( r, iy) + s( l, iy) + s( ix, u) + s( ix, d )

      delta = mn * s( ix, iy )
      !delta=0
      prob = exp((-2*delta)/t)
      !prob=1.1
      call random_number(x)

      if ( x .lt. prob ) then
         s(ix,iy) = -s(ix,iy)
         en = en + 2*delta
         mag = mag + 2*s(ix,iy)
      end if

    end do
    write(1,*) i, mag/real(n*n), en/real(n*n)
    !pause

  end do
  call cpu_time(finish)
  print*,n1, (finish-start)

end program simple_ising

[mecej4]

On Line-48, instead of ix = floor( n*x + 1.0 ), use ix = floor( n*x) + 1. Similarly for iy on Line-48. The following test program illustrates the difference and the importance of truncation immediately after multiplying and before adding 1:

program tflr
implicit none
integer ix1,ix2
real x
!
x=nearest(1.0,-1.0)
ix1=floor(64*x+1)
ix2=floor(64*x)+1
print *,ix1,ix2
end program

Note that the FLOOR function does not in general have the property

     floor(A + B) = floor(A) + floor(B)

and, therefore, you have to write your program accordingly.

[mohammad_b_]

mecej4 wrote:

On Line-48, instead of ix = floor( n*x + 1.0 ), use ix = floor( n*x) + 1. Similarly for iy on Line-48. The following test program illustrates the difference and the importance of truncation immediately after multiplying and before adding 1:

program tflr
implicit none
integer ix1,ix2
real x
!
x=nearest(1.0,-1.0)
ix1=floor(64*x+1)
ix2=floor(64*x)+1
print *,ix1,ix2
end program

Note that the FLOOR function does not in general have the property

     floor(A + B) = floor(A) + floor(B)

and, therefore, you have to write your program accordingly.

yes. but in this case this property is true

[mecej4]

mohammad b. wrote:

Quote:

mecej4 wrote:

On Line-48, instead of ix = floor( n*x + 1.0 ), use ix = floor( n*x) + 1. Similarly for iy on Line-48. Note that the FLOOR function does not in general have the property

     floor(A + B) = floor(A) + floor(B)

and, therefore, you have to write your program accordingly.

yes. but in this case this property is true

No, it is not; if you can advance logical arguments to support that claim, I'll revise my opinion. On the other hand, making the changes (to replace floor(x*64.0+1.0) by floor(x*64)+1, and similarly with y) to your program in #4 and running shows that my explanation is correct (of how ix or iy exceeded the intended range). The result displayed for n1 is 0 after the changes, but you may see various non-zero values in your original program. If you wish to track down this issue more fully, you may set the random number seed to a predetermined value at the beginning of the program.