Random Number Generators 
Author Message
 Random Number Generators

I have been trying to use the 'rand()' random number generator,
but it doesn't work adequately.  Now, I am trying to find other
random number generators that might be of use.  I have found
through the man pages that 'random()', 'ran()', and 'lcrans'
are available.  It seems that 'lcrans' would be best  for my
purposes, but I am programming in fortran and this is a C Library
Function.  Any ideas on how I can access this function -- or
any other function in the C Library -- from a FORTRAN program?

Also, does anyone have any ideas for random number generators that
might be of use to me?  I would prefer a pseudo-random number
generator...

                            Thanks,
                              Eva S. Simmons

--
***THE SIMMONS FACTOR ****EVA SABRINA SIMMONS******UT-AUSTIN GRAD STUDENT*****

         WWW Personal Page:   http://www.*-*-*.com/ ~weevey
****************** WATCH IT, OR IT MIGHT ATTACK!! ****************************



Fri, 29 Aug 1997 09:11:43 GMT  
 Random Number Generators

Quote:
>I have been trying to use the 'rand()' random number generator,
>but it doesn't work adequately.  Now, I am trying to find other
>random number generators that might be of use.  I have found
>through the man pages that 'random()', 'ran()', and 'lcrans'
>are available.  It seems that 'lcrans' would be best  for my
>purposes, but I am programming in FORTRAN and this is a C Library
>Function.  Any ideas on how I can access this function -- or
>any other function in the C Library -- from a FORTRAN program?

By reading the documentation of your compilers.  Since you didn't mention
the system you're using, this is all we can provide as help.

On VMS, you simply call your function and specify the C library at link
time.

On Unix, you write a wrapper function in C and don't specify any library
at link time (the C library will be searched by default).

Quote:

>Also, does anyone have any ideas for random number generators that
>might be of use to me?  I would prefer a pseudo-random number
>generator...

Without specialized hardware (e.g. sources and detectors of radiations)
all the generators are pseudo-random.  One of the volumes of Knuth's
"The Art of Computer Programming" is discussing the problem in great
detail.  Check your library.

Dan
--
Dan Pop
CERN, CN Division

Mail:  CERN - PPE, Bat. 31 R-004, CH-1211 Geneve 23, Switzerland



Fri, 29 Aug 1997 19:39:34 GMT  
 Random Number Generators

Quote:

> I have been trying to use the 'rand()' random number generator,
> but it doesn't work adequately.  Now, I am trying to find other
> random number generators that might be of use.  I have found
> through the man pages that 'random()', 'ran()', and 'lcrans'
> are available.  It seems that 'lcrans' would be best  for my
> purposes, but I am programming in FORTRAN and this is a C Library
> Function.  Any ideas on how I can access this function -- or
> any other function in the C Library -- from a FORTRAN program?

> Also, does anyone have any ideas for random number generators that
> might be of use to me?  I would prefer a pseudo-random number
> generator...

>                        Thanks,
>                          Eva S. Simmons

Here's an earlier Marsaglia generator that is easier
to implement.

-bb

======================================================
C This random number generator originally appeared in "Toward a Universal
C Random Number Generator" by George Marsaglia and Arif Zaman.
C Florida State University Report: FSU-SCRI-87-50 (1987)
C
C It was later modified by F. James and published in "A Review of Pseudo-
C random Number Generators"
C
C THIS IS THE BEST KNOWN RANDOM NUMBER GENERATOR AVAILABLE.
C       (However, a newly discovered technique can yield
C         a period of 10^600. But that is still in the development stage.)
C
C It passes ALL of the tests for random number generators and has a period
C   of 2^144, is completely portable (gives bit identical results on all
C   machines with at least 24-bit mantissas in the floating point
C   representation).
C
C The algorithm is a combination of a Fibonacci sequence (with lags of 97
C   and 33, and operation "subtraction plus one, modulo one") and an
C   "arithmetic sequence" (using subtraction).
C
C On a Vax 11/780, this random number generator can produce a number in
C    13 microseconds.
C========================================================================

      PROGRAM TstRAN
      REAL temp(100)
      INTEGER IJ, KL, len
C These are the seeds needed to produce the test case results
      IJ = 1802
      KL = 9373

C Do the initialization
      call rmarin(ij,kl)

C Generate 20000 random numbers
      len = 100
      do 10 i = 1, 200
         call RANMAR(temp, len)
10    continue

C If the random number generator is working properly, the next six random
C numbers should be:
C           6533892.0  14220222.0  7275067.0
C           6172232.0  8354498.0   10633180.0

      len = 6
      call RANMAR(temp, len)

      write(6,20) (4096.0*4096.0*temp(I), I=1,6)
20    format (3f12.1)
      end

      subroutine RMARIN(IJ,KL)
C This is the initialization routine for the random number generator RANMAR()
C NOTE: The seed variables can have values between:    0 <= IJ <= 31328
C                                                      0 <= KL <= 30081
C The random number sequences created by these two seeds are of sufficient
C length to complete an entire calculation with. For example, if sveral
C different groups are working on different parts of the same calculation,
C each group could be assigned its own IJ seed. This would leave each group
C with 30000 choices for the second seed. That is to say, this random
C number generator can create 900 million different subsequences -- with
C each subsequence having a length of approximately 10^30.
C
C Use IJ = 1802 & KL = 9373 to test the random number generator. The
C subroutine RANMAR should be used to generate 20000 random numbers.
C Then display the next six random numbers generated multiplied by 4096*4096
C If the random number generator is working properly, the random numbers
C should be:
C           6533892.0  14220222.0  7275067.0
C           6172232.0  8354498.0   10633180.0

      real U(97), C, CD, CM
      integer I97, J97
      logical TEST
      data TEST /.FALSE./
      common /raset1/ U, C, CD, CM, I97, J97, TEST

      if( IJ .lt. 0  .or.  IJ .gt. 31328  .or.
     *    KL .lt. 0  .or.  KL .gt. 30081 ) then
          print '(A)', ' The first random number seed must have a value
     *between 0 and 31328'
          print '(A)',' The second seed must have a value between 0 and        
     *30081'
            stop
      endif

      i = mod(IJ/177, 177) + 2
      j = mod(IJ    , 177) + 2
      k = mod(KL/169, 178) + 1
      l = mod(KL,     169)

      do 2 ii = 1, 97
         s = 0.0
         t = 0.5
         do 3 jj = 1, 24
            m = mod(mod(i*j, 179)*k, 179)
            i = j
            j = k
            k = m
            l = mod(53*l+1, 169)
            if (mod(l*m, 64) .ge. 32) then
               s = s + t
            endif
            t = 0.5 * t
3        continue
         U(ii) = s
2     continue

      C = 362436.0 / 16777216.0
      CD = 7654321.0 / 16777216.0
      CM = 16777213.0 /16777216.0

      I97 = 97
      J97 = 33

      TEST = .TRUE.
      return
      end

      subroutine ranmar(RVEC, LEN)
C This is the random number generator proposed by George Marsaglia in
C Florida State University Report: FSU-SCRI-87-50
C It was slightly modified by F. James to produce an array of pseudorandom
C numbers.
      REAL RVEC(*)
      real U(97), C, CD, CM
      integer I97, J97
      logical TEST
      common /raset1/ U, C, CD, CM, I97, J97, TEST

      integer ivec

      if( .NOT. TEST ) then
         print '(A)',' Call the init routine (RMARIN) before calling RAN        
     *MAR'  
         stop
      endif

      do 100 ivec = 1, LEN
         uni = U(I97) - U(J97)
         if( uni .lt. 0.0 ) uni = uni + 1.0
         U(I97) = uni
         I97 = I97 - 1
         if(I97 .eq. 0) I97 = 97
         J97 = J97 - 1
         if(J97 .eq. 0) J97 = 97
         C = C - CD
         if( C .lt. 0.0 ) C = C + CM
         uni = uni - C
         if( uni .lt. 0.0 ) uni = uni + 1.0
         RVEC(ivec) = uni
100   continue
      return
      end

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



Fri, 29 Aug 1997 22:14:44 GMT  
 Random Number Generators

Quote:

> I have been trying to use the 'rand()' random number generator,
> but it doesn't work adequately.  Now, I am trying to find other
> random number generators that might be of use.  I have found
> through the man pages that 'random()', 'ran()', and 'lcrans'
> are available.  It seems that 'lcrans' would be best  for my
> purposes, but I am programming in FORTRAN and this is a C Library
> Function.  Any ideas on how I can access this function -- or
> any other function in the C Library -- from a FORTRAN program?

> Also, does anyone have any ideas for random number generators that
> might be of use to me?  I would prefer a pseudo-random number
> generator...

>                        Thanks,
>                          Eva S. Simmons

Here's a supposedly superior generator that was posted a while back.
I haven't tried it yet.

-bb

========================================================
Hello everyone!

Well, here's what you've all been waiting for.  Again, would
someone please crosspost this to the relevant groups (numerical
analysis and so on)?

I am including, separated by a row of 72 stars, my response
(received today) from Fred James from CERN, who wrote the code,
and my original email from James from about a year ago.  What
you should cut out and compile starts right after the  row of
plus signs in this original email.  This consists of the
SUBROUTINE RANLUX, up until the END statement, and then the
PROGRAM RLXTST, which you can use to make sure everything is
working correctly on you computer, which includes the output
which should be produced in the comment lines at the end.

All I did was delete James's opening remarks, separate the
subroutine and the program, and compile them separately.
Except for the C's at the beginning of the lines, of course,
the test output on my computer (VAX 3100) was identical to the
test output in the RLXTST comments.  All in all, less than 2
minutes work.  Which speaks a lot for coding in absolutely
standard FORTRAN, as well as for the VAX compiler.  (While most
compilers, if not all, which offer extensions of course support
`only' standard FORTRAN, some, such as IBM AIX on a RISC 6000,
have to be given options in order _NOT_ to support extensions
which are INCOMPATIBLE with the standard (like folding to lower
case).  If you start a compiler with the defaults, then my view
is that it is OK to support extensions, but it should at least
compile standard FORTRAN without having to be explicitly told that
this is standard FORTRAN or otherwise produce results incompatible
with the standard.)

The Martin Luescher mentioned is a professor at DESY/University
of Hamburg who does theoretical elementary particle physics and
who wrote an article (DESY 93--133, September 1993, ISSN 0418-9833)
entitled `A Portable High-Quality Random Number Generator for
Lattice Field Theory Simulations' which I accidentally ran into
in the Hamburg Observatory library.  Luescher directed me to Fred
James, and since then we've used the generator quite a bit, tested
it extensively just to be safe and have noticed no problems.  For
most purposes, one could probably use the highest level with no
problems as far as CPU time goes.

Have fun!!

Phillip Helbig               Tel. .............. +49 40 7252 4110

Gojenbergsweg 112            Fax ............... +49 40 7252 4198
D-21029 Hamburg              Telex ............... 217884 hamst d

************************************************************************

Hello,
     Thank you for your message about RANLUX.  I am happy to hear of
satisfied users.  Yes, I know a little about some of the various rng's
that are being proposed these days, in particular by the Numerical
Recipes people.  I have told Bill Press about RANLUX and he has the
references. The area of rng's was always one of the weakest points in
Num. Recip., and I guess that will continue to be true, although their
algorithms will improve of course.  They do a very good job in so many
areas, but you can't really expect them to be best in everything.

   Press seems to be working closely with Marsaglia, which I also did at
one time, until I discovered Martin Luescher and Pierre L'Ecuyer who have
a more rigorous and more convincing approach to the problem.  You may
know that Marsaglia once said: "Random numbers are like sex: ... Even
when they are bad they are still pretty good."  However, we must remember
that the actual RANLUX algorithm owes much to several people including
Marsaglia:

  1. Marsaglia invented the "subtract-with-borrow" algorithm which has
many practical advantages (fast and portable because all computations are
done in portable floating-point, very long period, and reasonably easy to
initialize, restart, etc.), but is not sufficiently random.
  2. L'Ecuyer and Tuzaka (and perhaps others) recognized and proved that
Marsaglia's algorithm is in fact equivalent to a linear congruential
generator (which is bad) with a very big multiplier (which is good).
This gave a first basis for understanding the randomness properties and
in particular the flaws.
  3. Luescher, using Kolmogorov's chaos theory, showed rigorously how to
improve the algorithm by skipping, and how much you have to skip to make
it lose all memory of its past.  He also can calculate the period
exactly, with or without skipping.
  4. I wrote the Fortran version of the generator.

   Feel free to propagate RANLUX, but make sure to give proper credit
always to Luescher and the two articles in Computer Physics Communications
where it all is published.  Thanks again,   Fred.

************************************************************************

I guess that theoretical cosmology is easier than experimental cosmology
(less time-consuming) but this random number generator should be good enough
for both.  The references are Lu:scher, Computer Physics Communications
79 (1994) 100, and James, CPC 79 (1994) 111.

+++++++++++++++++++++++++++++++++++
      SUBROUTINE RANLUX(RVEC,LENV)
C         Subtract-and-borrow random number generator proposed by
C         Marsaglia and Zaman, implemented by F. James with the name
C         RCARRY in 1991, and later improved by Martin Luescher
C         in 1993 to produce "Luxury Pseudorandom Numbers".
C     Fortran 77 coded by F. James, 1993
C
C       references:
C  M. Luscher, Computer Physics Communications  79 (1994) 100
C  F. James, Computer Physics Communications 79 (1994) 111
C
C   LUXURY LEVELS.
C   ------ ------      The available luxury levels are:
C
C  level 0  (p=24): equivalent to the original RCARRY of Marsaglia
C           and Zaman, very long period, but fails many tests.
C  level 1  (p=48): considerable improvement in quality over level 0,
C           now passes the gap test, but still fails spectral test.
C  level 2  (p=97): passes all known tests, but theoretically still
C           defective.
C  level 3  (p=223): DEFAULT VALUE.  Any theoretically possible
C           correlations have very small chance of being observed.
C  level 4  (p=389): highest possible luxury, all 24 bits chaotic.
C
C!!! ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
C!!!  Calling sequences for RANLUX:                                  ++
C!!!      CALL RANLUX (RVEC, LEN)   returns a vector RVEC of LEN     ++
C!!!                   32-bit random floating point numbers between  ++
C!!!                   zero (not included) and one (also not incl.). ++
C!!!      CALL RLUXGO(LUX,INT,K1,K2) initializes the generator from  ++
C!!!               one 32-bit integer INT and sets Luxury Level LUX  ++
C!!!               which is integer between zero and MAXLEV, or if   ++
C!!!               LUX .GT. 24, it sets p=LUX directly.  K1 and K2   ++
C!!!               should be set to zero unless restarting at a break++
C!!!               point given by output of RLUXAT (see RLUXAT).     ++
C!!!      CALL RLUXAT(LUX,INT,K1,K2) gets the values of four integers++
C!!!               which can be used to restart the RANLUX generator ++
C!!!               at the current point by calling RLUXGO.  K1 and K2++
C!!!               specify how many numbers were generated since the ++
C!!!               initialization with LUX and INT.  The restarting  ++
C!!!               skips over  K1+K2*E9   numbers, so it can be long.++
C!!!   A more efficient but less convenient way of restarting is by: ++
C!!!      CALL RLUXIN(ISVEC)    restarts the generator from vector   ++
C!!!                   ISVEC of 25 32-bit integers (see RLUXUT)      ++
C!!!      CALL RLUXUT(ISVEC)    outputs the current values of the 25 ++
C!!!                 32-bit integer seeds, to be used for restarting ++
C!!!      ISVEC must be dimensioned 25 in the calling program        ++
C!!! ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
      DIMENSION RVEC(LENV)
      DIMENSION SEEDS(24), ISEEDS(24), ISDEXT(25)
      PARAMETER (MAXLEV=4, LXDFLT=3)
      DIMENSION NDSKIP(0:MAXLEV)
      DIMENSION NEXT(24)
      PARAMETER (TWOP12=4096., IGIGA=1000000000,JSDFLT=314159265)
      PARAMETER (ITWO24=2**24, ICONS=2147483563)
      SAVE NOTYET, I24, J24, CARRY, SEEDS, TWOM24, TWOM12, LUXLEV
      SAVE NSKIP, NDSKIP, IN24, NEXT, KOUNT, MKOUNT, INSEED
      INTEGER LUXLEV
      LOGICAL NOTYET
      DATA NOTYET, LUXLEV, IN24, KOUNT, MKOUNT /.TRUE., LXDFLT, 0,0,0/
      DATA I24,J24,CARRY/24,10,0./
C                               default
C  Luxury Level   0     1     2   *3*    4
      DATA NDSKIP/0,   24,   73,  199,  365 /
Corresponds to p=24    48    97   223   389
C     time factor 1     2     3     6    10   on slow workstation
C                 1    1.5    2     3     5   on fast mainframe
C
C  NOTYET is .TRUE. if no initialization has been performed yet.
C              Default Initialization by Multiplicative Congruential
      IF (NOTYET) THEN
         NOTYET = .FALSE.
         JSEED = JSDFLT
         INSEED = JSEED
         WRITE(6,'(A,I12)') '
...

read more »



Fri, 29 Aug 1997 22:05:34 GMT  
 Random Number Generators

Quote:
>I have been trying to use the 'rand()' random number generator,
>but it doesn't work adequately.  Now, I am trying to find other
>random number generators that might be of use.  I have found
>through the man pages that 'random()', 'ran()', and 'lcrans'
>are available.  It seems that 'lcrans' would be best  for my
>purposes, but I am programming in FORTRAN and this is a C Library
>Function.  Any ideas on how I can access this function -- or
>any other function in the C Library -- from a FORTRAN program?

>Also, does anyone have any ideas for random number generators that
>might be of use to me?  I would prefer a pseudo-random number
>generator...

>                        Thanks,
>                          Eva S. Simmons

The book "Numerical Recipes in Fortran" by William H. Press et al.
contains an entire chapter on random number generators.  This
will probably have all that you need.

Erik

- Show quoted text -

Quote:
>--
>***THE SIMMONS FACTOR ****EVA SABRINA SIMMONS******UT-AUSTIN GRAD STUDENT*****

>     WWW Personal Page:  http://ccwf.cc.utexas.edu/~weevey
>****************** WATCH IT, OR IT MIGHT ATTACK!! ****************************



Sat, 30 Aug 1997 08:07:58 GMT  
 Random Number Generators

[snip]
w;Also, does anyone have any ideas for random number generators that
w;might be of use to me?  I would prefer a pseudo-random number
w;generator...

Eva,

Find a copy of Numerical Recipes - it should be available in a good
college bookstore or engineering library. It gives some reasonable
algorithms including one based on the Data Encryption Standard.

Bob
---
~ CMPQwk #1.42-20~ UNREGISTERED EVALUATION COPY

Lousiana - Darwin Takes A Holiday!
'[1;35;40m-=> Delphi Internet Jet SST v2.009 - (C) PBE



Sun, 31 Aug 1997 11:13:46 GMT  
 Random Number Generators
Thanks for the helpful advice about random number generators.  I've
decided to implement the RAN4 routine from the Numerical Recipes book
for FORTRAN.  So far, it has worked just fine...:)  

Enclosed is a copy of this code:

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

c    **********************************************************************
      function ran4(idum)
c    **********************************************************************
c     PURPOSE:  Returns a uniform random deviate in the range 0.0 to 1.0
c      generated by pseudo-DES (DES-like) hasing of the 64-bit word
c      (idums,idum), where idums was set by a previous call with negative
c      idum.  Also increments idum.  Routine can be used to generate a
c      random sequence by successive calls, leaving idum unaltered between
c      calls; or it can randomly access the nth deviate in a sequence by
c      calling with idum=n.  Different sequences are initialized by calls
c      with differing negative values of idum.
c    **********************************************************************
      integer idum
      real ran4
      integer idums,irword,itemp,jflmsk,jflone,lword
      real ftemp
      equivalence (itemp,ftemp)
      save idums,jflone,jflmsk
      data idums/0/, jflone/Z'3F800000'/, jflmsk /Z'007FFFFF'/
c    **********************************************************************
      if (idum.lt.0) then
         idums=-idum
         idum=1
      endif
      irword=idum
      lword=idums
      call psdes(lword,irword)
      itemp=ior(jflone,iand(jflmsk,irword))
      ran4=ftemp-1.0
      idum=idum+1
c    **********************************************************************
      return
      end
c    **********************************************************************

c    **********************************************************************
      subroutine psdes(lword,irword)
c    **********************************************************************
c     PURPOSE:  "Pseudo-DES" hashing of the 64-bit word (lword,irword).
c      Both 32-bit arguments are turned hashed on all bits.  NOTE:  This
c      routine assumes that arbitrary 32-bit integers can be added without
c      overflow.  To accomplish this, you may need to compile with a
c      special directive (e.g., /check=nooverflow for VMS).  In other
c      languages, such as C, one can instead type the integers as "unsigned".
c    **********************************************************************
      integer i,ia,ib,iswap,itmph,itmpl,c1(4),c2(4),irword,lword,NITER
      parameter (NITER=4)
      save c1,c2
      data c1 /Z'BAA96887',Z'1E17D32C',Z'03BCDC3C',Z'0F33D1B2'/,
     +  c2 /Z'4B0F3B58',Z'E874F0C3',Z'6955C5A6',Z'55A7CA46'/
c    **********************************************************************
      do 500 i=1,NITER
         iswap=irword
         ia=ieor(irword,c1(i))
         itmpl=iand(ia,65535)
         itmph=iand(ishft(ia,-16),65535)
         ib=itmpl**2+not(itmph**2)
         ia=ior(ishft(ib,i6),iand(ishft(ib,-16),65535))
         irword=ieor(lword,ieor(c2(i),ia)+itmpl*itmph)
         lword=iswap
500   continue
c    **********************************************************************
      return
      end
c    **********************************************************************

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

                        Thanks again for the help!:)
                                Eva S. Simmons

--
***THE SIMMONS FACTOR ****EVA SABRINA SIMMONS******UT-AUSTIN GRAD STUDENT*****

         WWW Personal Page:  http://ccwf.cc.utexas.edu/~weevey
****************** WATCH IT, OR IT MIGHT ATTACK!! ****************************



Mon, 01 Sep 1997 06:16:41 GMT  
 Random Number Generators
On 15 Mar 1995 16:16:41 -0600 Eva Sabrona Simmons pontificated:

Quote:
> Thanks for the helpful advice about random number generators.  I've
> decided to implement the RAN4 routine from the Numerical Recipes book
> for FORTRAN.  So far, it has worked just fine...:)  
> Enclosed is a copy of this code:

Oops! Posting Numerical Recipes code is technically a violation of their copyright;
they (Flannery, Press, et. al.) do not permit free distribution of the code in any
form, especially electronically.  The book is quite clear about this; you are expected
to pay royalties for each routine you distribute.  Besides, better random number routines
have been posted to this newsgroup, and/or are available from NETLIB.

There is also some work by that venerable sage Knuth, who studied the problem in
the same detailed way as he did with sorting.

--
Rick Venable                  =====\     |=|    "Eschew Obfuscation"
FDA/CBER Biophysics Lab       |____/     |=|
Bethesda, MD  U.S.A.          |   \  /   |=|  / Not an official statement \



Tue, 02 Sep 1997 12:53:33 GMT  
 Random Number Generators
"Numerical Recipes" codes are in many cases good pedagogically, but not
the bestimplementations of (semi-)numerical computations to be used for
for production/research work.In general, copying code from textbooks,
even if diskettes are available, should be done with caution.


Tue, 02 Sep 1997 18:18:42 GMT  
 Random Number Generators

Quote:
>"Numerical Recipes" codes are in many cases good pedagogically, but not
>the bestimplementations of (semi-)numerical computations to be used for
>for production/research work.In general, copying code from textbooks,
>even if diskettes are available, should be done with caution.

This is particularly true when the book in question is "Numerical Recipes".
Those who use or plan to use this book, should have a look at:

http://www.lysator.liu.se/c/numerical-recipes.html

Dan
--
Dan Pop
CERN, CN Division

Mail:  CERN - PPE, Bat. 31 R-004, CH-1211 Geneve 23, Switzerland



Wed, 03 Sep 1997 06:03:18 GMT  
 
 [ 10 post ] 

 Relevant Pages 

1. Random Number Generator to produce SAME random number from 12:00am-11:59pm

2. Pentium III processor number and random number generator

3. J random-number generator: what is used?

4. Random number generators

5. Random number generators

6. random number generator?

7. Random number generator for Smalltalk/GemStone

8. random number generator

9. Problem Solved: Need A Random Number Generator

10. Seeding the Random Number generator

11. Random Number Generator

12. Wanted: random number generator

 

 
Powered by phpBB® Forum Software