150 likes | 374 Views
Monte Carlo Methods. When you can’t do the math, simulate the process with random numbers Numerical integration to get areas/volumes Particle physics Interactions (e.g., biological) in which the level of response is variable “Monte Carlo” as in throwing dice to determine outcomes.
E N D
Monte Carlo Methods When you can’t do the math, simulate the process with random numbers • Numerical integration to get areas/volumes • Particle physics • Interactions (e.g., biological) in which the level of response is variable “Monte Carlo” as in throwing dice to determine outcomes
Random Numbers We don’t really want random numbers, we want “pseudorandom” numbers, because we want reproducible results in experiments. We want numbers that behave the way that random numbers would, but that are produced in a deterministic way. What do we mean by “behave the way that random numbers would”?
Random Numbers (2) What do we mean by “behave the way that random numbers would”? If the numbers run 0.0 to 1.0, then 1/n of them should fall into any interval of length 1/n Similarly for 2 dimensions, 3 dimensions, etc. Shouldn’t be able to predict the next from the previous, or the next from the two previous, or the next from the three previous, etc. But we don’t want to require that the numbers be exactly evenly balanced either… NIST has a set of standard statistical tests
Random Numbers (3) The RNG should generate the same numbers regardless of the machine used Different sequences of random numbers should come from different “seed” values The RNG should use minimal CPU and memory
Random Numbers (4) RNGs usually produce integers that are uniformlydistributed in the range of 1 to some maximum value M We want a maximal period, i.e., M distinct integers before a repetition We convert to reals, uniformly distributed in the range 0.0 to 1.0, by dividing by M (with appropriate mode conversion) We convert to reals with some chosen probabilistic distribution by using the random numbers in the distribution function
Congruential Random Number Generators Consider powers of 2 modulo 13, say Power N Power N 1 2 7 24=11 2 4 8 22=9 3 8 9 18=5 4 16=3 10 10 5 6 11 20=7 6 12 12 14=1 2**12 = 1 modulo 13
Congruential RNGs (2) (Fermat’s Little) Theorem: If p is a prime number, and a is any integer, then a**(p-1) = 1 modulo p So the 2**12 = 1 mod 13 is not a coincidence Observation: The sequence of values assumed as we power up is more or less random among the values from 1 through p-1 This would be called a multiplicativecongruential random number generator
Congruential RNGs (3) Choose a prime p Choose a multiplier a Choose an addend b Compute xn+1 = (a*xn + b) modulo p If done correctly, then we get all the numbers 1 through p-1 before we repeat
Congruential RNGs (4) Subprogram for your assignment For a long time it was highly useful that 2**31 – 1 was a prime The “usual” 2s complement arithmetic for 32-bit words worked out nicely, so the whole multiplicative congruential stuff required no extra effort
Lagged Fibonacci Generators Fill an array of values Xj for j = 1, …, N (Note this means N seeds, not just one) Compute Xj = Xj-m + Xj-n, for m,n < N “Fibonacci” because the series looks sort of like a Fibonacci series Fn+2 = Fn+1 + Fn “Lagged” from looking at the subscripts
Parallel RNGs For a very large Monte Carlo computation, we might need lots of random numbers We don’t want to have to compute these sequentially (very slow) Lagged Fibonacci RNGs allow at least N at a time to be computed in parallel Parallel RNGs are still an active research topic.
Difference and Differential Equations The derivative of a function is the rate of change of that function A differentialequation is an equation involving a derivative If one is driving at constant velocity, then the velocity is the rate of change of position dx/dt = v Starting at position zero, and travelling at velocity v for t time units, one travels v*t distance
Euler’s Method dy/dt = f(y,t) (the differential equation) y(0) = y0 (the initial condition) Approximate dy/dt by Δy/Δt = (yi+1 – yi)/ Δt = f(yi,ti) so yi+1 = yi + Δt * f(yi,ti) This is called Euler’s method. If we can compute the function f(y,t), then we can step forward a solution of the differential equation from time 0 to time whatever. In essence, approximate the function by trapezoids.