Bei dieser Methode werden zwei statistisch voneinander unabhängige gaußverteilte Zufallsgrößen x und y aus den beiden zwischen 0 und 1 gleichverteilten und statistisch unabhängigen Zufallsgrößen u und υ durch nichtlineare Transformation erzeugt:
Das Verfahren nach Box und Muller – abgekürzt BM – kann wie folgt charakterisiert werden:
Der theoretische Hintergrund für die Gültigkeit obiger Generierungsvorschriften basiert auf den Gesetzmäßigkeiten für zweidimensionale Zufallsgrößen (siehe Kapitel 4.1).
Die obigen Gleichungen liefern sukszessive zwei Gaußsche Zufallswerte x und y ohne statistische Bindungen untereinander. Diese Tatsache kann man zur Verkürzung der Simulationszeit nutzen, indem man bei jedem Funktionsaufruf ein Tupel (x, y) von Gaußwerten generiert.
Ein Vergleich der Rechenzeiten zeigt, dass – bei jeweils bestmöglicher Implementierung – das BM-Verfahren der Additionsmethode mit I = 12 etwa um den Faktor 3 überlegen ist.
Der Wertebereich ist hier weniger begrenzt als bei der Additionsmethode, so dass bei einer Simulation auch kleine Wahrscheinlichkeiten besser wiedergegeben werden. Aber auch mit dem BM-Verfahren lassen sich keine beliebig kleinen Fehlerwahrscheinlichkeiten simulieren.
Beispiel: Bei einem 32 Bit-Rechner ist die kleinste darstellbare Floatzahl ca. 2–31 ≈ 0.466 · 10–9. Der Maximalwert für den Wurzelausdruck in der Generierungsvorschrift des BM-Verfahrens kann somit nicht größer als ca. 6.55 werden und ist äußerst unwahrscheinlich. Da sowohl die Cosinus- als auch die Sinusfunktion betragsmäßig auf 1 beschränkt ist, wäre das gleichzeitig der maximal mögliche Wert für die Zufallsgrößen x und y (mx = my = 0 und σx = σy = 1 vorausgesetzt).
Eine in[ES96]dokumentierte Simulation über eine Milliarde Abtastwerte hat aber gezeigt, dass das BM-Verfahren nur bis zu Fehlerwahrscheinlichkeiten von 10–5 die Q-Funktion sehr gut approximiert, dann aber der Kurvenverlauf steil abbricht. Der maximal auftretende Wert war bei der in [ES96] beschriebenen Simulation nicht 6.55, sondern aufgrund der aktuellen Zufallsgrößen u und υ nur etwa 4.6, womit sich der schlagartige Abfall ab etwa 10–5 erklären lässt. Bei 64 Bit-Rechenoperationen kann dieses Verfahren natürlich noch deutlich verbessert werden.