JOVANA
Explore Library Glossary Getting Started Three Levels Fields How it works Mission
Join the mission
All guides

Physics by Dice: Monte Carlo and Molecular Dynamics

Some problems have too many dimensions to grid up — a gas of 10^{23} particles, a hundred-dimensional integral. The escape is to compute with deliberate randomness. Meet the Monte Carlo method, its dimension-blind 1/\sqrt{N} error, the Metropolis trick for thermal systems, and molecular dynamics — physics simulated one force at a time.

Why randomness beats a grid

Grid methods die in high dimensions. Integrating over the positions of N gas particles is a $3N-dimensional integral; even a coarse 10 points per axis gives 10^{3N}$ evaluations — hopeless once N exceeds a handful. This is the curse of dimensionality. The Monte Carlo method sidesteps it: estimate an integral by averaging the integrand over randomly chosen sample points.

\int_V f\,dV \;\approx\; \frac{V}{N}\sum_{i=1}^{N} f(\mathbf{x}_i)

The integral is the volume times the average of f over N random points — a statistical estimate of a deterministic quantity.

The magic is in the error. By the central limit theorem the statistical error of the estimate falls as 1/\sqrt{N}independent of the number of dimensions. A grid's error scales like N^{-k/d}, getting worse as the dimension d grows; Monte Carlo's N^{-1/2} does not care about d at all. Above a handful of dimensions, the dice win.

\sigma_{\text{est}} \;\propto\; \frac{\sigma_f}{\sqrt{N}}

To halve the error, quadruple the samples — slow, but dimension-blind convergence that no grid can match in high dimensions.

Throwing darts at pi

The cleanest illustration: estimate \pi by throwing darts. Scatter random points in the unit square and count the fraction that land inside the quarter circle of radius 1. That fraction approaches the area ratio \pi/4.

\pi \;\approx\; 4\,\frac{N_{\text{inside}}}{N_{\text{total}}}

The quarter disk's area over the square's is \pi/4; counting hits estimates that ratio, and four times it estimates \pi.

Add samples and watch the running estimate of \pi converge, with the scatter shrinking like 1/\sqrt{N} — a hands-on Monte Carlo integration you drive yourself.

Watch two things in the widget. The estimate wanders toward \pi but never sits exactly on it; and the wobble narrows only slowly — to gain one more decimal digit you need roughly 100\times more darts. Monte Carlo is robust and general but not fast; its virtue is that it keeps working where nothing else can.

Sampling what matters: importance sampling and Metropolis

Uniform darts waste effort where f is tiny. Importance sampling instead draws points preferentially where the integrand is large and reweights them — slashing the variance for the same N. In statistical mechanics we need averages weighted by the Boltzmann factor e^{-E/k_BT}, where almost every configuration has negligible weight; uniform sampling is hopeless.

The Metropolis algorithm solves this by building a random walk through configurations that visits each in proportion to its Boltzmann weight. Propose a small change; then accept or reject it by one simple rule.

P_{\text{accept}} = \min\!\left(1,\; e^{-\Delta E / k_B T}\right)

Always accept a move that lowers the energy; accept an uphill move only with probability e^{-\Delta E/k_BT}. The walk settles into the canonical distribution.

This single rule lets you compute thermal averages for systems like the Ising model of magnetism, where you can watch the magnetization vanish at the critical temperature. Metropolis Monte Carlo is how much of statistical physics is actually done — it samples the ensemble behind the partition function without ever computing the astronomically large sum Z directly.

Molecular dynamics: physics one force at a time

The complementary approach is molecular dynamics (MD): don't sample randomly — actually integrate Newton's equations for every particle, forward in time, and let the system explore its states dynamically. Here Guide 2's Verlet integrator earns its keep: for thousands to billions of atoms, each step computes the force on every atom from its neighbours and Verlet-steps them all.

\mathbf{a}_i = \frac{1}{m_i}\sum_{j \neq i} \mathbf{F}_{ij}

The force on each atom is the sum over all pair interactions; this O(N^2) force sum (tamed to nearly O(N) with cutoffs and neighbour lists) is the computational bottleneck.

By the equipartition theorem, the average kinetic energy per degree of freedom fixes the temperature, so an MD box is a genuine thermodynamic system you can heat, cool, and measure. MD gives you dynamics — transport, diffusion, real time-dependence — while Monte Carlo gives you equilibrium averages efficiently. Together they simulate everything from protein folding to the melting of a crystal.