Bielajew A.F. Fundamentals of the Monte Carlo method for neutral and charged particle transport

Подождите немного. Документ загружается.

1.2. SOME HISTORY 11

1.2 Some history

The usual ﬁrst reference to the Monte Carlo method is usually that of Comte de Buf-

fon [dB77] who proposed a Monte Carlo-like method to evaluate the probability of tossing a

needle onto a ruled sheet. This reference goes back to 1777, well before the contemplation of

automatic calculating machines. Buﬀon calculated that a needle of length L tossed randomly

on a plane ruled with parallel lines of distance d apart where d>Lwould have a probability

p =

πd

. (1.1)

A computer simulation of 50 needles (where d/L =

) on a ﬁnite grid of 5 lines is shown in

Figure 1.9.

The Buffon needle simulation

Figure 1.9: A computer simulation of the Buﬀon needle problem.

Later, Laplace [Lap86] suggested that this procedure could be employed to determine the

value of π, albeit slowly. Several other historical uses of Monte Carlo predating computers

12 CHAPTER 1. WHAT IS THE MONTE CARLO METHOD?

are cited by Kalos and Whitlock [KW86]. The modern Monte Carlo age was ushered in by

von Neumann and Ulam during the initial development of thermonuclear weapons

.Ulam

and von Neumann coined the phrase “Monte Carlo” and were pioneers in the development

of the Monte Carlo technique and its realizations on digital computers

The two books by Richard Rhodes, The making of the Atomic Bomb and Dark Sun are excellent historical

sources for this period.

Despite its chequered history, thermonuclear weapons have never been deployed in conﬂict. So Monte

Carlo calculations have not been employed in a destructive way. In contrast, Monte Carlo calculations are

employed for beneﬁcial purposes, such as in the prediction of dose in cancer treatments [Bie94] and should

be credited with saving lives.

Bibliography

[BHNR94] A. F. Bielajew, H. Hirayama, W. R. Nelson, and D. W. O. Rogers. His-

tory, overview and recent improvements of EGS4. National Research Council

of Canada Report PIRS-0436, 1994.

[Bie94] A. F. Bielajew. Monte Carlo Modeling in External Electron-Beam Radiother-

apy — Why Leave it to Chance? In “Proceedings of the XI’th Conference on

the Use of Computers in Radiotherapy” (Medical Physics Publishing, Madison,

Wisconsin), pages 2 – 5, 1994.

[Bie95] A. F. Bielajew. EGS4 timing benchmark results: Why Monte Carlo is a viable

option for radiotherapy treatment planning. In “Proceedings of the International

Conference on Mathematics and Computations, Reactor Physics, and Environ-

mental Analyses” (American Nuclear Society Press, La Grange Park, Illinois,

U.S.A.), pages 831 – 837, 1995.

[BR92] A. F. Bielajew and D. W. O. Rogers. A standard timing benchmark for EGS4

Monte Carlo calculations. Medical Physics, 19:303 – 304, 1992.

[BW91] A. F. Bielajew and P. E. Weibe. EGS-Windows - A Graphical Interface to EGS.

NRCC Report: PIRS-0274, 1991.

[dB77] G. Comte de Buﬀon. Essai d’arithm´etique morale, volume 4. Suppl´ement `a

l’Histoire Naturelle, 1777.

[KW86] M. H. Kalos and P. A. Whitlock. Monte Carlo methods, Volume I: Basics. John

Wiley and Sons, New York, 1986.

[Lap86] P. S. Laplace. Theorie analytique des probabilit´es, Livre 2. In Oeuvres compl´etes

de Laplace, volume 7, Part 2, pages 365 – 366. L’acad´emie des Sciences, Paris,

1886.

[MRR95] M. S. MacPherson, C. K. Ross, and D. W. O. Rogers. A technique for accurate

measurement of electron stopping powers. Med. Phys. (abs), 22:950, 1995.

[MRR96] M. S. MacPherson, C. K. Ross, and D. W. O. Rogers. Measured electron stopping

powers for elemental absorbers. Med. Phys. (abs), 23:797, 1996.

14 BIBLIOGRAPHY

[NBRH94] W. R. Nelson, A. F. Bielajew, D. W. O. Rogers, and H. Hirayama. EGS4 in ’94:

A decade of enhancements. Stanford Linear Accelerator Report SLAC-PUB-6625

(Stanford, Calif), 1994.

[NHR85] W. R. Nelson, H. Hirayama, and D. W. O. Rogers. The EGS4 Code System.

Report SLAC–265, Stanford Linear Accelerator Center, Stanford, Calif, 1985.

Chapter 2

Elementary probability theory

Fundamental to understanding the operation of a Monte Carlo process and interpreting the

results of Monte Carlo calculations, is some understanding of elementary probability theory.

In this chapter we introduce some elementary probability theory to facilitate the discussions

in the later chapters.

2.1 Continuous random variables

2.1.1 One-dimensional probability distributions

A probability distribution function, p(x) is a measure of the likelihood of observing x.For

example, x could be the position at which a photon interacts via the Compton interaction.

If p(x

)=2p(x

) then what this actually means is that an observation of x in a diﬀerential

interval dx is twice as likely to be observed than in an interval the same size about x

This really only has meaning in the limit that dx goes to zero, or if dx is ﬁnite, when p(x)

varies slowly over the interval in the vicinity of x

and x

. A more concrete example for

p(x)=exp(−x) is shown in Figure 2.1.

In general, p(x) has some special properties that distinguishes it from other functions:

• p(x) ≥ 0 since negative probabilities have no meaning (at least none that we will

mention herein),

• p(x) is normalized in the following fashion,

max

min

dxp(x)=1, (2.1)

16 CHAPTER 2. ELEMENTARY PROBABILITY THEORY

0.0 1.0 2.0 3.0 4.0 5.0

0.0

0.2

0.4

0.6

0.8

1.0

p(x)

Probability distribution function

p(x) = exp(−x)

(log 2, 0.5)

(log 4, 0.25)

Figure 2.1: The probability function p(x)=exp(−x). p(log 2) = 2p(log 4)

•−∞<x

min

max

< +∞,thatis,x

min

and x

max

can be any real number so long as

min

is less than x

max

These are the only restrictions on p(x). Note that the statement of normalization above

implies that p(x) is integrable over its range of deﬁnition. Other than that, it can be

discontinuous and even inﬁnite! (The above deﬁnition in terms of intervals preserves our

interpretation

Certain probability functions can be characterized in terms of their moments,

i =

max

min

dxx

p(x) , (2.2)

however, the existence of these moments is not guaranteed nor even necessary. Assuming

that hxi and hx

i exist, the variance associated with the probability function may be deﬁned

An example of an inﬁnite probability would be p(x)=aδ(x − x

)+(1− a)δ(x − x

), where δ() is the

Dirac delta function.

2.1. CONTINUOUS RANDOM VARIABLES 17

to be:

var{x} = hx

i−hxi

. (2.3)

The variance is a measure of the width of the distribution of x.var{x} is zero for the Dirac

delta function and greater than zero for all other probability distribution functions, even

combinations of delta functions.

An example of a probability distribution function that has no moments is the Cauchy or

Lorentz distribution function

deﬁned over the inﬁnite range:

p(x)=

+ x

−∞<x<∞. (2.4)

The mean value of the Cauchy/Lorentz distribution can be deﬁned to be zero if the “principle

value” of the ﬁrst moment is obtained in the following way:

hxi = lim

a−→ ∞

−a

+ x

, (2.5)

but the second moment and hence the variance can not be deﬁned in any fashion.

Another interesting probability distribution function of great importance in the simulation

of electron transport is the screened Rutherford or Wentzel distribution function:

p(µ)=

a(2 + a)

(1 − µ + a)

; − 1 ≤ µ ≤ 1, (2.6)

where µ is the cosine of the scattering angle, cos Θ.

The conventional small-angle form of the screened Rutherford or Wentzel distribution func-

tion is:

p(Θ) = 4a

(Θ

+2a)

;0≤ Θ < ∞, (2.7)

Itsﬁrstmomentexists,hΘi = π

a/2 but its second moment is inﬁnite! This strange

behavior, the nonexistence of an angular variance is responsible for the theory of electron

transport being so problematic. Of course, one could restrict the range of integration to

physical angles, 0 ≤ Θ ≤ π, but the problems persist.

2.1.2 Two-dimensional probability distributions

Consideration of two and higher-dimensional probability distributions follows from a gen-

eralization of one-dimensional distributions with the added features of correlation between

observables and conditional probabilities.

More pathological examples could be invented but this distribution arises from an interesting

application—the intrinsic probability distribution of the energy of a quantum from an excited atomic state

of ﬁnite lifetime.

18 CHAPTER 2. ELEMENTARY PROBABILITY THEORY

Consider a two-dimensional (or joint) probability function p(x,y). A tangible example is the

distribution in energy and angle of a photon undergoing an inelastic collision with an atom

Another example is the two-dimensional probability distribution presented in Figure 2.2.

The meaning of two-dimensional probability distributions is as follows: Hold one variable,

Figure 2.2: A two-dimensional probability function f(x, y)

say x ﬁxed, and the resulting distribution is a probability distribution function in the other

variable, y. It as if you cut through the two-dimensional probability distribution at a given

point in x and then displayed the “cross cuts”. Several examples of these cross cuts are

shown in Figure 2.3.

The notions of normalization and moments follow directly. Thus,

i =

max

min

max

min

dxdyx

p(x, y) , (2.8)

with the normalization condition hx

i = 1, the only “moment” that must be deﬁned.

“Inelastic” in this deﬁnition relates to the energy of the photon changing, not the energy of the atom.

2.1. CONTINUOUS RANDOM VARIABLES 19

Figure 2.3: Cross cuts of the two-dimensional probability function given in Figure 2.2.

Higher order moments may or may not exist. If they exist, we deﬁne the covariance:

cov{x, y} = hxyi−hxihyi , (2.9)

which can be positive or negative. Note that cov{x, x} =var{x}.

The covariance is a measure of the independence of observing x or y.Ifx and y are inde-

pendent random variables, then p(x, y)=p

(x)p

(y)andcov{x, y} = 0. A related function

is the correlation coeﬃcient:

ρ{x, y} =

cov{x, y}

var{x}var{y}

, (2.10)

where −1 ≤ ρ{x, y}≤1.

Some interesting relations involving the variances and the covariances may be found. For

example,

var{x ± y} =var{x} +var{y}±2cov{x, y} , (2.11)

20 CHAPTER 2. ELEMENTARY PROBABILITY THEORY

or simply var{x} +var{y} if x and y are independent.

The marginal probabilities are deﬁned by integrating out the other variables.

m(x)=

max

min

dyp(x, y); m(y)=

max

min

dxp(x, y) . (2.12)

Note that the marginal probability distributions are properly normalized. For the example

of joint energy and angular distributions, one marginal probability distribution relates to the

distribution in energy irrespective of angle and the other refers to the distribution in angle

irrespective of energy. This the joint probability distribution function may be written:

p(x, y)=m(x)p(y|x) , (2.13)

where the conditional probability is deﬁned by

p(y|x)=

p(x, y)

m(x)

. (2.14)

The interpretation of the conditional probability is that given x, what is the probability

that y occurs. The appearance of m(x) in the denominator guaranteed the normalization of

p(y|x).

2.1.3 Cumulative probability distributions

Associated with each one-dimensional probability distribution function is its cumulative

probability distribution function

c(x)=

min

p(x

) . (2.15)

Cumulative probability distribution functions have the following properties which follow

directly from its deﬁnition and the properties of probability distribution functions:

• p(x)andc(x) are related by a derivative:

p(x)=

dc(x)

, (2.16)

• c(x) is zero at the beginning of its range of deﬁnition

c(x

min

)=0, (2.17)

• and unity at the end of its range of deﬁnition

c(x

max

)=1, (2.18)