Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

Can Random Number Generator of Fortran 90 be trusted for Monte Carlo Integration?

I have written a short monte carlo integration algorithm to calculate an integral in Fortran 90. I once compared the result obtained by solving the integral with respect to some parameter using the intrinsic random number generator with the random number generator method ran1 presented in Numerical Recipes for Fortran90 Volume 2.

Running the same algorithm twice, once calling the intrinsic random_seed(), then always call random_number() and once calling the ran1() method provided in the Numerical Recipe book I obtain as result in principal the same shape but the intrinsic result is a continuous curve in contrast to the ran1 result. In both cases I call the function with random parameters 10,000 times for a parameter value q, add it and then go on to the next q value and call the function 10,000 times etc.

A comparative image of the result can be found here: http://i.imgur.com/gZVFdOP.png

If I increase the number of calls both curves converge. But I was wondering: why does the intrinsic random number generator generate this smoothness? Is it still generally advised to use it or are there are other more advised RNG? I suppose the continuous result is a result of the "less" randomness of the intrinsic number generator.

(I left out the source code as I don't think that there is a lot of input from it. If somebody cares I can hand it in later.)

like image 545
DomiD Avatar asked Jul 08 '16 09:07

DomiD


People also ask

How are random numbers used in Monte Carlo simulation?

Monte Carlo simulation is one of the main applications involving the use of random number generators. It is also one of the best methods of testing the randomness properties of such generators, by comparing results of simulations using different generators with each other, or with analytic results.

How do I generate a random number in Fortran 90?

Fortran 90 introduced a pseudo-random number generator (PRNG) to the language standard. The routine random_number() returns a single number or an array of numbers from the uniform distribution over the range 0 ≤ x < 1.

How do I generate a random number in Fortran 95?

Description: RAND(FLAG) returns a pseudo-random number from a uniform distribution between 0 and 1. If FLAG is 0, the next number in the current sequence is returned; if FLAG is 1, the generator is restarted by CALL SRAND(0) ; if FLAG has any other value, it is used as a new seed with SRAND .

How do I use random seed in Fortran?

The Fortran random number generator has a "seed" value which it uses to start its sequence from. You can change this value using using the RANDOM_SEED intrinsic subroutine. The way this is done in a completely portable way (i.e. that can be run on any machine) is complicated.


2 Answers

There are NO guarantees about the quality of the pseudo random generator in standard Fortran. If you care about some particular quality of implementation for cryptography or science sensitive to random numbers (Monte-Carlo), you should use some library which you have control about.

You can study the manual of your compiler to find out what it says about the random number generator, but every compiler can implement a completely different algorithm to generate random numbers.

Numerical Recipes is actually not well received by some people in the numerical mathematics community http://www.uwyo.edu/buerkle/misc/wnotnr.html

This site is not for software recommendation, but this article (link given by roygvib in a comment): https://arxiv.org/abs/1005.4117 is a good review with examples of bad and good algorithms, methods how to test them, how to generate arbitrary number distributions and examples of calls of two example libraries in C (one of them can be called from Fortran as well).

Personally I use this https://bitbucket.org/LadaF/elmm/src/master/src/rng_par_zig.f90 parallel PRNG, but I didn't test the quality, I personally just need speed. But this is not a software recommendation site.

like image 199
Vladimir F Героям слава Avatar answered Oct 13 '22 23:10

Vladimir F Героям слава


The particular random number generator used depends on the compiler. Perhaps documented, perhaps not. Perhaps subject to change. For high quality work, I would use library / source code from elsewhere. A webpage with Fortran implementations of the Mersenne Twister: http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/VERSIONS/FORTRAN/fortran.html. I have used http://theo.phys.sci.hiroshima-u.ac.jp/~ishikawa/PRNG/mt_stream_en.html and verified against the original implementation. This version is useful for multi-threaded programs.

like image 26
M. S. B. Avatar answered Oct 14 '22 01:10

M. S. B.