Why do two players see different results from the same inputs?

Because one of them called Math.sin. Across 200,000 identical inputs, Node 23 and Python 3.14 on the same Apple M3 returned different bits for sin on 4.40% of them, exp on 9.61%, and tan on 41.16% — each by one to three ULP. Dropping a single sin term into an otherwise identical particle simulation put the two clients 1 world unit apart after 464 steps and 133 units apart after 600 — 7.7 and 10 seconds of play at 60 Hz. The same simulation using only + - * / and Math.sqrt ran 5,000 steps in both languages and produced the identical 32-bit state hash, 1210764c.

Hardware: Apple M3, 8 cores, 16 GB, macOS 26.4.1 (build 25E253). Node v23.5.0 (V8 12.9.202.28-node.12) and Python 3.14.6. Everything here was measured on one machine and one CPU architecture — ARM64. I could not test x86, ARM32, a different libm, Unity's Mono against IL2CPP, or a GPU. What that limitation means is spelled out near the end, because it changes which conclusions transfer.

The short answer

• IEEE 754 requires + - * / and sqrt to be correctly rounded, and they behaved that way: 200,000 inputs each, zero bit differences between Node 23 and Python 3.14. Transcendental functions carry no such requirement, and sin differed on 4.40% of inputs, pow on 10.17%, tan on 41.16%.
• A one-ULP difference in a single particle's starting x became a full world unit of position error in 413 steps — 6.9 seconds — and 100 units in 482. Divergence is exponential, not gradual; there is no "small desync".
• Reversing the order in which forces are accumulated, with no other change, broke the state hash on step 1 and reached one world unit by step 408. Same code, same machine, same inputs.
• 16.16 fixed point in Int32 cost 6.66x — 21.10 µs against 140.54 µs per step in Node — and produced byte-identical results in Node and Python, hash 6467a425 in both.
• Hashing the state every step is affordable, but not with your own hash function. Node's native SHA-256 over the 6,400-byte state took 3.07 µs against a 21.21 µs step (14.5%). A hand-written FNV-1a in JavaScript took 26.92 µs — more than the simulation step it was checking.

What is being measured

Two hundred particles in a 512x512 world. Each step tests all 19,900 pairs, applies a linear repulsion inside a radius of 24, adds gravity, integrates at dt = 1/60 and reflects off the walls. The state is 800 float64 values — x, y, vx, vy — hashed with FNV-1a over the raw bytes, so one bit changing anywhere changes the hash. Starting positions come from mulberry32(12345).

The physics is deliberately restricted to + - * / and Math.sqrt. That restriction is the experiment: those five operations are the ones IEEE 754-2019 requires to be correctly rounded. Everything else in Math is unregulated.

Is float math deterministic on one machine?

Yes, and more robustly than I expected. Eight separate Node processes, 5,000 steps each, all produced hash 1210764c and x[0] = 429.534394289. Then I forced V8 down completely different code paths:

A 53x spread in execution speed and not one bit of difference in the result. The interpreter, the baseline compiler and the optimising compiler all agree, because for these five operations the standard leaves them nothing to disagree about. This is the useful half of the story: within one binary on one machine, a float simulation is reproducible, which is why your replays work locally and your desync only shows up in a real match.

Does the same float code match across languages?

It did. I ported the simulation to Python 3.14 line for line — same mulberry32, same loop order, same struct-packed hash — and ran 5,000 steps:

node sim.js   -> float64 steps=5000 hash=1210764c ms=105.5   x0=429.534394289
python3 sim.py -> py float64 steps=5000 hash=1210764c ms=5904.1 x0=429.534394289

Two unrelated runtimes, 99.5 million pair tests, byte-identical state. That is what "correctly rounded" buys you, and it is worth stating plainly because much of the writing on this topic implies float arithmetic is inherently unreliable. The basic operations are not. The library around them is.

Do Math.sin results match across languages?

No. Same machine, same 200,000 input doubles (verified bit-identical on both sides before comparing):

V8 computes Math.sin with its own polynomial kernel; CPython calls the platform libm. Neither is wrong — IEEE 754 only recommends correct rounding for these — and both are within a few ULP of the true value. That is irrelevant. A lockstep simulation does not need accuracy, it needs agreement, and one ULP of disagreement is a desync.

tan at 41% is the worst case, because it amplifies argument-reduction error near its poles — and it is exactly the kind of function that ends up in a camera or aiming calculation.

How fast does one wrong bit become a different game?

Three divergence experiments, each comparing two simulations step by step and recording the largest position difference across all 200 particles.

The growth is the part worth internalising. In the 1-ULP run the maximum position delta was 5.68e-14 at step 10, 9.02e-12 at step 100, 106 units at step 500 and 474 units at step 20,000 — the width of the world. A colliding n-body simulation is chaotic, so the error doubles every 6 to 12 steps once the first contact resolves differently. There is no window in which a desync is minor and will settle. You have about seven seconds between the first wrong bit and two players watching different games — which is also seven seconds in which to catch it.

Does the order you add things up matter?

Yes — but the first version of this test said it did not, and the reason is worth more than the result. Summing one million random values scaled to ±500, forward against backward against a four-way chunked reduction, gave exactly the same double every time. Kahan summation agreed too.

The test data was the problem. mulberry32 returns k / 2^32, so (rnd() - 0.5) * 1000 is always an integer multiple of 2^-29, and every partial sum stayed under 2^53 — every addition was exact, so no ordering could change anything. Sorting the array first did change the answer, by 2.96e-7, because ascending order pushes the running total to about 1.25e8 and out of that exact range. Regenerating with full 53-bit mantissas produced the expected result:

Then the version that actually bites. Four worker_threads summing quarters of a 4-million-element SharedArrayBuffer, with the main thread adding the partial sums in the order the messages arrive, repeated 30 times: 15 distinct arrival orders and two distinct sums — 18 runs giving one value, 12 the other, both about 3e-9 from the serial result. Same binary, same data, same machine, two answers. If your simulation fans work out to a job system and reduces in completion order, it is already nondeterministic.

What else desyncs that is not floating point at all?

• Hash-set iteration. In four separate Python processes, iterating a set of eight unit names gave four different orders — PYTHONHASHSEED is randomised per process. The equivalent dict gave insertion order every time; a set of small ints gave sorted order every time. In JavaScript, Map and Set are insertion-ordered by spec, but removing and re-adding an element moves it to the end — exactly what a kill-and-respawn path does.
• float32 accumulation. 100,000 additions of 0.1: float64 landed on 10000.000000018848, Math.fround at every step on 9998.556640625 — an error of 1.44, not 1.9e-8. Storing simulation state as float32 to save bandwidth puts you at a precision where errors show up with no cross-machine effect at all.
• The PRNG, usually not. mulberry32(12345) produced identical first six outputs to 17 significant digits in both languages, because it is integer arithmetic converted to a double at the end. Seeded PRNGs desync when the number of calls differs between clients, not because the algorithm drifts.

What does fixed point cost?

16.16 fixed point in Int32 — one part in 65,536 — with multiplication split into 16-bit halves to stay inside the exact range of a double, truncating division, and an integer square root corrected so it never depends on how Math.sqrt rounds.

6.66x slower, and identical in both languages. In C# or C++ the ratio would be smaller — much of the JavaScript penalty is |0 coercions and the lack of a 64-bit intermediate — but it is not free anywhere, because fixed-point multiply and divide are several integer operations each.

Note what fixed point does not buy: accuracy. The fixed simulation's trajectory differs completely from the float one — x[0] = 282.41 against 429.53 after 5,000 steps. It is not a better approximation of the physics. It is a reproducible one, which is the only property lockstep actually needs.

How often can you afford to check?

Compare a state hash between clients and you catch the desync on the frame it happens. The cost, for the 6,400-byte state:

The hand-rolled byte-at-a-time hash costs more than the physics it is checking. Node's native SHA-256 is 8.8x faster than my JavaScript FNV-1a on the same buffer, so the hash you run every frame should be the one you did not write.

Quantising positions before hashing, so tiny differences do not raise false alarms, is the usual advice. Its cost in detection latency, on the 1-ULP run:

Quantising to a whole world unit — absurdly coarse — delayed detection by only 367 steps, because once divergence starts it covers nine orders of magnitude in about 230 steps. Quantise as coarsely as you like; chaos finds you within seven seconds either way. The choice costs six seconds, not six minutes.

What this could not measure

One machine, one architecture. Every cross-implementation result above is Node against Python on ARM64 — the easy case: same CPU, same rounding mode, same 64-bit doubles. I could not test x86 against ARM, 32-bit x87 with its 80-bit intermediates, FMA contraction, -ffast-math, Unity Mono against IL2CPP against Burst, or any GPU. That is where the folklore says it gets worse, and I have no numbers for it. What these measurements do establish is a lower bound: the divergence is already present between two runtimes on one chip, and it comes from the math library, not the arithmetic.

For the simulation side of the same problem, see at how many objects naive collision detection breaks and how to path 500 units without dropping frames; for streaming state through a hash without holding it all in memory, see parsing huge JSON.

Check it yourself

Two files, no dependencies. Each runs the simulation with sqrt only, then again with one Math.sin term added. Run both and compare.

// det.js  — node det.js
const N = 200, W = 512, H = 512, R = 24, K = 60, G = 30, DT = 1/60, D = 0.999;
function rng(a){return()=>{a=(a+0x6D2B79F5)|0;let t=a;t=Math.imul(t^(t>>>15),t|1);
  t^=t+Math.imul(t^(t>>>7),t|61);return((t^(t>>>14))>>>0)/4294967296;};}
function fnv(b){let h=0x811c9dc5;for(let i=0;i<b.length;i++){h^=b[i];
  h=Math.imul(h,0x01000193)>>>0;}return h>>>0;}
function init(){const r=rng(12345),s={x:new Float64Array(N),y:new Float64Array(N),
  vx:new Float64Array(N),vy:new Float64Array(N)};
  for(let i=0;i<N;i++){s.x[i]=r()*W;s.y[i]=r()*H;s.vx[i]=(r()-0.5)*80;s.vy[i]=(r()-0.5)*80;}
  return s;}
function step(s,t,useSin){const{x,y,vx,vy}=s,ax=new Float64Array(N),ay=new Float64Array(N);
  for(let i=0;i<N;i++)for(let j=i+1;j<N;j++){const dx=x[j]-x[i],dy=y[j]-y[i],d2=dx*dx+dy*dy;
    if(d2<R*R&&d2>1e-9){const d=Math.sqrt(d2),f=(K*(R-d))/d,fx=dx*f,fy=dy*f;
      ax[i]-=fx;ay[i]-=fy;ax[j]+=fx;ay[j]+=fy;}}
  for(let i=0;i<N;i++){if(useSin)ax[i]+=20*Math.sin(x[i]*0.01+t*DT);
    vx[i]=(vx[i]+ax[i]*DT)*D;vy[i]=(vy[i]+(ay[i]+G)*DT)*D;x[i]+=vx[i]*DT;y[i]+=vy[i]*DT;
    if(x[i]<0){x[i]=-x[i];vx[i]=-vx[i];}if(x[i]>W){x[i]=2*W-x[i];vx[i]=-vx[i];}
    if(y[i]<0){y[i]=-y[i];vy[i]=-vy[i];}if(y[i]>H){y[i]=2*H-y[i];vy[i]=-vy[i];}}}
function hash(s){const b=new Uint8Array(4*N*8),v=new Float64Array(b.buffer);
  v.set(s.x,0);v.set(s.y,N);v.set(s.vx,2*N);v.set(s.vy,3*N);return fnv(b);}
for(const useSin of[false,true]){const s=init(),n=useSin?700:2000;
  for(let k=1;k<=n;k++)step(s,k,useSin);
  console.log(`${useSin?'sin ':'sqrt'} steps=${n} hash=${hash(s).toString(16)
    .padStart(8,'0')} x0=${s.x[0].toFixed(9)}`);}

# det.py  — python3 det.py
import math, struct
N,W,H,R,K,G,DT,D = 200,512.0,512.0,24.0,60.0,30.0,1.0/60.0,0.999
M32 = 0xFFFFFFFF
def rng(a):
    st=[a & M32]
    def n():
        st[0]=(st[0]+0x6D2B79F5)&M32; t=st[0]
        t=((t^(t>>15))*(t|1))&M32
        t=(t^((t+(((t^(t>>7))*(t|61))&M32))&M32))&M32
        return ((t^(t>>14))&M32)/4294967296.0
    return n
def fnv(b):
    h=0x811c9dc5
    for c in b: h=((h^c)*0x01000193)&M32
    return h
def init():
    r=rng(12345); x=[0.0]*N; y=[0.0]*N; vx=[0.0]*N; vy=[0.0]*N
    for i in range(N):
        x[i]=r()*W; y[i]=r()*H; vx[i]=(r()-0.5)*80.0; vy[i]=(r()-0.5)*80.0
    return [x,y,vx,vy]
def step(s,t,use_sin):
    x,y,vx,vy=s; ax=[0.0]*N; ay=[0.0]*N
    for i in range(N):
        for j in range(i+1,N):
            dx=x[j]-x[i]; dy=y[j]-y[i]; d2=dx*dx+dy*dy
            if d2<R*R and d2>1e-9:
                d=math.sqrt(d2); f=(K*(R-d))/d; fx=dx*f; fy=dy*f
                ax[i]-=fx; ay[i]-=fy; ax[j]+=fx; ay[j]+=fy
    for i in range(N):
        if use_sin: ax[i]+=20*math.sin(x[i]*0.01+t*DT)
        vx[i]=(vx[i]+ax[i]*DT)*D; vy[i]=(vy[i]+(ay[i]+G)*DT)*D
        x[i]+=vx[i]*DT; y[i]+=vy[i]*DT
        if x[i]<0: x[i]=-x[i]; vx[i]=-vx[i]
        if x[i]>W: x[i]=2*W-x[i]; vx[i]=-vx[i]
        if y[i]<0: y[i]=-y[i]; vy[i]=-vy[i]
        if y[i]>H: y[i]=2*H-y[i]; vy[i]=-vy[i]
def hsh(s):
    return fnv(struct.pack('<%dd'%(4*N), *(s[0]+s[1]+s[2]+s[3])))
for use_sin in (False,True):
    s=init(); n=700 if use_sin else 2000
    for k in range(1,n+1): step(s,k,use_sin)
    print('%s steps=%d hash=%08x x0=%.9f' % ('sin ' if use_sin else 'sqrt',
                                             n, hsh(s), s[0][0]))

On the machine above, in 0.11 s and 3.17 s respectively:

$ node det.js
sqrt steps=2000 hash=8878d217 x0=408.511548854
sin  steps=700 hash=65a26d59 x0=377.522038117

$ python3 det.py
sqrt steps=2000 hash=8878d217 x0=408.511548854
sin  steps=700 hash=6596e6d1 x0=429.015297437

Same hash on the sqrt line. Different hash on the sin line, and the two runtimes disagree about where particle zero is by 51.5 world units after 700 steps. If your netcode has one transcendental call in the simulation path, that is your bug.