Physics & Mathematics
A guided walkthrough of every equation and physical concept that powers the simulator — from special relativity to Monte Carlo integration.
Monte Carlo Methods
Monte Carlo methods are a class of computational algorithms that rely on repeated random sampling to obtain numerical results. Named after the Monte Carlo Casino in Monaco, these techniques were formalised during the Manhattan Project by Stanisław Ulam and John von Neumann in the late 1940s.
In high-energy physics, Monte Carlo simulations are indispensable. They are used to:
- Generate pseudo-data that mimics detector output (event generators like Pythia, Herwig)
- Integrate high-dimensional phase-space integrals that are analytically intractable
- Propagate uncertainties by sampling nuisance parameters (systematic error estimation)
The fundamental idea: estimate an expectation value by averaging over N random samples:
The statistical uncertainty decreases as 1/√N regardless of dimensionality — a key advantage over quadrature methods in high dimensions:
Special Relativity & Kinematics
At the energies we simulate (GeV scale), particles travel at speeds close to c. Newtonian mechanics breaks down and we must use Einstein's special relativity. We work in natural units where ℏ = c = 1, so energy, momentum, and mass all have units of GeV.
The cornerstone is the relativistic energy–momentum relation, which replaces E = ½mv²:
From this we extract the 3-momentum magnitude:
We decompose the momentum into Cartesian components using spherical coordinates (θ = polar, φ = azimuthal):
Why natural units?
Setting ℏ = c = 1 removes factors of 3×10⁸ from every equation. Energies, momenta, and masses all share the same unit (GeV), which simplifies both algebra and code. Every formula in this simulator uses natural units.
Invariant Mass & Mandelstam s
When two particles collide, the quantity that determines what new particles can be created is the centre-of-mass energy √s. It is a Lorentz invariant — the same in every reference frame — and is defined through the Mandelstam variable s:
Expanding the squared 3-momentum sum with the opening angle θ between the two beams:
Head-on (θ = π)
cos θ = −1 → the cross-term is maximised. This is the collider configuration (LHC, RHIC) that yields the highest √s for given beam energies.
Fixed-target (θ ≈ 0)
cos θ ≈ 1 → the cross-term is minimised. Much of the beam energy goes into boosting the CM system rather than creating particles.
Centre-of-Mass Frame & Final States
In the centre-of-mass (CM) frame the total 3-momentum is zero. This makes the kinematics of the outgoing particles particularly simple. For a 2→2 process with outgoing masses m₁ and m₂, energy–momentum conservation uniquely determines the final-state energies:
When both outgoing particles are massless (e.g. γγ) this simplifies to an equal split:
Lab-frame conservation in the simulator
Our simulator stores lab-frame quantities. To ensure exact energy conservation we split the lab-frame total energy E_A + E_B between the two outgoing particles with a stochastic asymmetry, rather than splitting √s. This guarantees E₁ + E₂ = E_A + E_B to machine precision (~10⁻¹⁷ GeV).
Collider Geometry
The simulator models a symmetric collider: particle A travels along the +z axis and particle B along −z. The crossing half-angle θ parameterises how far each beam deviates from perfect head-on alignment:
The azimuthal angles φ_A and φ_B are drawn independently and uniformly in [0, 2π), breaking cylindrical symmetry when θ ≠ 0. For θ = 0 the beams are perfectly head-on and √s is maximised.
Simulation Algorithm
Putting everything together, here is the complete event-generation algorithm executed for each of the N collision events. In the implementation every step is fully vectorised with NumPy, so all N events are processed simultaneously.
Sample beam energies
Draw E_A and E_B independently from a uniform distribution over the configured energy range.
Sample collision angle (isotropic)
Use inverse-transform sampling on cos θ to respect the solid-angle element dΩ = sin θ dθ dφ.
Compute 3-momenta
Use the energy–momentum relation and spherical-to-Cartesian conversion.
Compute invariant mass
Evaluate √s from the total four-momentum of the system.
Final-state energies
Split the lab-frame total energy with a small stochastic asymmetry to model the random boost direction.
Record event
Store (event_id, E_A, E_B, θ, p_x, p_y, p_z, E₁, E₂, √s) in the dataset.