mercredi 15 juin 2016

Rounding of double precision to single precision: Forcing an upper bound

I'm using a Mersenne Twister implementation which provides me numbers with double precision.

http://ift.tt/1ZQZ1zr (implementation in Fortran 77 by Tsuyoshi Tada, I'm using genrand_real2)

However, my application needs, in order to avoid warnings while multiplying numbers with different precisions, a single precision random number. So, I wrote a small function to convert between the two data types:

    function genrand_real()

    real   genrand_real
    real*8 genrand_real2

    genrand_real = real(genrand_real2())

    return
    end

I'm using real and real*8 to be consistent with the code I'm working on. It works perfectly most of the time (besides de fact that I'm not sure about how fast real() is), however it changes the upper bound of my RNG, since the conversion changes the [0,1) to [0,1]. I've never thought about that until I've had problems with it.

My question is, how can I ensure the upper bound in an efficient way, or even how could I write a function similar to genrand_real2 (the original one) that provides me single precision reals. My guess is I only need to replace the divisor 4294967296.d0 but I don't know by which number

  function genrand_real2()

  double precision genrand_real2,r
  integer genrand_int32
  r=dble(genrand_int32())
  if(r.lt.0.d0)r=r+2.d0**32
  genrand_real2=r/4294967296.d0

  return
  end




Aucun commentaire:

Enregistrer un commentaire