5.5 Monte Carlo integration

Source: Ideas from Computational Physics (Landau, Páez, Bordeianu), §5.5. Rewritten in original wording. Continues 5.4 Gaussian quadrature. Labs: 5.5b Monte Carlo integration labs.

Pond area by stone throwing

Imagine a pond inside a rectangular fence of known area $A_{\mathrm{box}}$. Throw $N$ stones uniformly at random into the box; count how many land in the water ($N_{\mathrm{pond}}$). Walk the hit-fraction ratio, the pond-area estimator, and $\pi\simeq 4 N_{\mathrm{hit}}/N$ — pick the matching equation at each step.

No antiderivative, no mesh — only a uniform sampler and a membership test (“inside the pond?”). That is the Monte Carlo idea for integration: replace geometry by counting.

$\pi$ from a unit circle in a square

Put the unit disk $x^{2}+y^{2}\le 1$ inside the square $[-1,1]^{2}$ of side $2$ (area $A_{\mathrm{box}}=4$). The guide’s last two steps are the $\pi/4$ hit fraction and the factor-of-four estimator. With the chapter LCG and $N=10^{4}$ throws you typically land near $3.14$ — close enough to see the method, not a replacement for $\mathrm{atan}(1)$ in production. The statistical error scales like $1/\sqrt{N}$: quadrupling the sample count halves the typical fluctuation.

Why Monte Carlo for integrals?

In one dimension, trapezoid, Simpson, and Gauss usually beat crude Monte Carlo for smooth $f$. Monte Carlo shines when:

• the domain is awkward (irregular pond, complicated acceptance region);
• the dimension $D$ is large — a grid with $M$ points per axis needs $M^{D}$ evaluations, which explodes;
• you only need a rough answer with a transparent error bar $\sim 1/\sqrt{N}$.

The next lesson (§5.6) makes the mean-value form precise and shows the $D$-independence of the $1/\sqrt{N}$ law.

Error investigation: trap / Simpson / Gauss / MC

A standard calibration integral is

$$I=\int_0^{1} e^{-t}\,dt=1-e^{-1}.$$

Run trapezoid, Simpson, Gauss–Legendre, and a mean-value Monte Carlo estimate at several $N$, then plot $\log_{10}|I_{\mathrm{approx}}-I|$ versus $\log_{10} N$. On that log–log plane a clean power law

$$\varepsilon\sim C\,N^{\alpha}$$

appears as a straight segment with slope $\alpha$. Deterministic rules for smooth $f$ often show steep negative slopes (Simpson $\alpha\approx -4$ in $h$, Gauss even steeper in $N$ for analytic integrands) until round-off flattens the floor — algorithms stop improving when floating-point noise dominates (the Fig. 5.5 idea). Monte Carlo’s slope is near $\alpha=-1/2$; it loses in 1D on this smooth problem, but that shallow slope does not get worse when $D$ grows.

Use the explorers below (and §5.6’s mean-value chart) to compare regimes rather than memorizing a single “best” rule.

Quiz

Monte Carlo integration

1 / 3

Stone-throwing estimates A_pond by…

Pick an answer to see if you are right.

Next: 5.5b Monte Carlo integration labs.