Random numbers and Monte Carlo

1. Random numbers and Monte Carlo Inmetro, 8/7/2003 João R. T. de Mello Neto IF - UFRJ

2. References • The Nature of Mathematical Modeling, N. Gershenfelder, Cambridge, 1999; • Numerical Recipes in C, Second Edition, W.H Press et al., Cambridge, 1992; • Statistical Data Analysis, G. Cowan, Oxford, 1998 • Computational Physics, Dean Karlen (online), 1998

3. Random O acaso não existe. car sticker

4. Random numbers Important task: generate random variables from known probability distributions. random numbers: produced by the computer in a strictly deterministic way – pseudorandom (→ Monte Carlo Method) random number generators: linear congruential generators: c=3, a=5, m=16 ex: 0,2,3,13,4,7,6,1,8,11,10,5,12,15,14,9,0,2….

5. Generators Arguments from number theory: good values for a, c and m. • longest possible period • individual elements within a period should • follow each other “randomly” Good generators: Ex. 1 RANDU a = 65539, m=231, c=0

6. Generators

7. Generators

8. Generators

9. Generators “Random numbers fall mainly in the planes” Marsaglia, Proc. Acad. Sci. 61, 25, 1968 What is seen in RANDU is present in any multiplicative congruential generator. In 32 bit machines: maximum number of hyperplanes in the space of d-dimensions is: d=3 2953 d=4 566 d=6 120 d=10 41 RANDU has much less than the maximum!

10. Generators Ex. 2 Minimal standard generator (Num. Recp. ran0) RAN1 and RAN2, given in the first edition, are “at best mediocre” a = 75 =16807, m=231-1

11. Generators cernlib (mathlib V115, F.James ) A portable high-quality random generator for lattice field theory simulation, M. Lüscher, Comp. Phys. Comm. 79, 100, 1994 Ex. 3 RANLUX period ≥ 10165 Functional definition: an algorithm that generates uniform numbers in acceptable if it is not rejected by a set of tests. Lots of literature about random generators testing. See a set of programs (Die Hard) http://stat.fsu.edu/~geo/diehard.html • Recommendations: • Do some simple tests. • Check the results with other generator

12. Inverse transform p(x) uniform probability distribution x: uniform deviate Generate y according to g(y). 1 G(y) uniform deviate in analytically or numerically. x g(y) 0 y out

13. F(x) 1 0.5 x 0 1 2 3 Inverse transform Ex. 4 discrete case f(x=0)=0.3 f(x=1)=0.2 f(x=2)=0.5 u 0.0 ≤ u ≤ 0.3 x=0 0.3 ≤ u ≤ 0.5 x=1 0.5 ≤ u ≤ 1.0 x=2 for 2000 deviates: x=0 0.2880 X=1 0.2075 X=2 0.5045

14. Inverse transform • amount of time until a specific event occurs; • time between independent events (time until a part fails); Exponential

15. ymax f(x) 0 xmin xmax Acceptance-rejection method if (x,y) under the curve, accept the point, else, discard.

16. Acceptance-rejection method

17. Importance sampling g(x): random numbers are easy to generate g(x) • generate random x according to g(x) • generate u uniformily between 0 and g(x); • if u < f(x), accept x, if not, reject. f(x)

18. Multivariate Normal Very useful Starts with a d-dimensional vector of standard normal Random numbers; Transformed to the desired distribution using: dX1 vector of standard normal random numbers dX1 vector representing the mean covariance matrix desired dXd matrix such that (Cholesky factorization of the cov. matrix)

19. Multivariate Normal

20. Random walk 10 simulations 1K points each