Pseudo-Random Numbers and Seeds in Simulation

Deterministic numbers that behave randomly

A pseudo-random number generator (PRNG) is a small deterministic machine: a fixed-size internal state and a function that updates the state and emits a number. Start it from the same state and it emits the same sequence, always. What makes it useful is that the sequence passes statistical tests for randomness — uniform frequencies, no detectable correlation between successive values, no short cycles — well enough that a simulation cannot tell the difference.

The seed is the value used to fill the initial state. Two consequences follow, and simulation work depends on both:

• Reproducibility. The same seed gives the same run. A surprising result can be re-run, debugged and shown to someone else.
• Replication. Different seeds give independent-looking runs. Estimating how much a result varies means running several seeds, not one.

Why not Math.random()?

JavaScript’s Math.random() is a reasonable generator (V8, the engine in Chrome and Node, uses xorshift128+), but it cannot be seeded. Every page load produces a different stream, so a run cannot be repeated and a browser result cannot be compared with a Python one. Python’s random module can be seeded, but it uses the Mersenne Twister, a different algorithm, so the same seed means different numbers in the two languages.

OpenSim.zip therefore ships its own generator, written twice with identical arithmetic. The interactive pages, every simulation.js and every simulation.py in the model packages use it, and the build checks that all three agree for every model.

The generator: xoshiro128**

xoshiro128** (Blackman and Vigna) keeps four 32-bit words of state and updates them with shifts, rotations and exclusive-ors, then scrambles one word with two multiplications and a rotation to produce the output. It has a period of 2¹²⁸ − 1, passes the standard statistical test batteries, and needs only 32-bit integer arithmetic — which JavaScript does exactly through Math.imul and >>> 0.

The four state words are filled from the seed with a SplitMix-style hash, so nearby seeds (1, 2, 3…) still start from unrelated states.

const rotl = (x, k) => ((x << k) | (x >>> (32 - k))) >>> 0;

const next = () => {
  const r = Math.imul(rotl(Math.imul(b, 5) >>> 0, 7), 9) >>> 0;
  const t = (b << 9) >>> 0;
  c = (c ^ a) >>> 0;
  d = (d ^ b) >>> 0;
  b = (b ^ c) >>> 0;
  a = (a ^ d) >>> 0;
  c = (c ^ t) >>> 0;
  d = rotl(d, 11);
  return r;
};
const uniform = () => next() / 4294967296;   // [0, 1), 32-bit resolution

M32 = 0xFFFFFFFF

def _rotl(x, k):
    return ((x << k) | (x >> (32 - k))) & M32

def next(self):
    r = (_rotl((self.b * 5) & M32, 7) * 9) & M32
    t = (self.b << 9) & M32
    self.c ^= self.a
    self.d ^= self.b
    self.b ^= self.c
    self.a ^= self.d
    self.c ^= t
    self.d = _rotl(self.d, 11)
    return r

Python integers are unbounded, so every product and shift is masked back to 32 bits; JavaScript gets the same effect from Math.imul and the unsigned shift. With seed 42 both produce 2837322924, 544945897, 479756282 as their first three outputs.

From 32-bit integers to the distributions models need

• Uniform on [0, 1): divide by 2³². The smallest step is about 2.3 × 10⁻¹⁰, far finer than any model here resolves.
• Integer in [0, n): ⌊u · n⌋. For n that is not a power of two some integers are hit once more than others out of 2³² draws — a relative bias below n/2³², negligible for dice and grid positions, not acceptable for cryptography.
• Exponential with rate λ: inversion, −ln(1 − u)/λ. Using 1 − u keeps the logarithm away from zero.
• Standard normal: the Box–Muller transform, √(−2 ln u₁) · cos(2πu₂). The site draws one normal per call and discards the second value Box–Muller could produce; that wastes a little speed but keeps the call sequence trivial to mirror in another language.

Pitfalls that bias results

• Treating one seed as the answer. A single run is one sample. Report results across seeds, or report a standard error — see how accurate a Monte Carlo estimate is.
• Choosing the seed after seeing the result. Fix seeds before running experiments. Hunting for a seed that "works" is a way of fabricating a result.
• Sharing one stream between parts of a model carelessly. If adding a feature changes how many numbers an early component consumes, every later draw shifts and the whole run changes. That is expected, not a bug — but compare model variants with several seeds, not one.
• Seeding from the clock in parallel jobs. Jobs started in the same second get the same seed and produce identical "independent" runs.
• Using a simulation PRNG for security. xoshiro, xorshift and the Mersenne Twister are predictable from their outputs. Use crypto.getRandomValues or Python’s secrets module for anything an adversary might want to guess.

When bit-for-bit agreement is not guaranteed

The integer stream is identical everywhere. Floating-point results built on it are identical when a model uses only +, −, ×, ÷ and square roots in the same order, because IEEE 754 requires those to be correctly rounded. Functions such as exp, log, sin and cos are not required to be correctly rounded, and the JavaScript and Python runtimes can differ in the last bit. For most models that difference never changes a decision; for chaotic ones like the double pendulum it eventually would, and the package documents how that model handles it.