mercredi 20 juillet 2016

overflow in implementation of linear congruential random number generator using fortran

The famous linear congruential random number generator also known as minimal standard use formula

x(i+1)=16807*x(i) mod (2^31-1)

I want to implement this using fortran.

However, as pointed out by "Numerical Recipes", directly implement the formula with default Integer type (32bit) will cause 16807*x(i) to overflow.

So the book recommend Schrage’s algorithm is based on an approximate factorization of m. This method can still implemented with default integer type.

However, I am wondering fortran actually has Integer(8) type whose range is -9,223,372,036,854,775,808 to 9,223,372,036,854,775,807 which is much bigger than 16807*x(i) could be.

but the book even said the following sentence

It is not possible to implement equations (7.1.2) and (7.1.3) directly in a high-level language, since the product of a and m − 1 exceeds the maximum value for a 32-bit integer.

So why can't we just use Integer(8) type to implement the formula directly?




Aucun commentaire:

Enregistrer un commentaire