WebGPU Monte-Carlo simulation: Ising Model

Using WebGPU, we simulate an Ising Model for magnetic spins, with no external field. The critical parameter is beta = J / kT, where J is the coupling constant between spins (J>0: ferromagnetism, J<0: anti-ferromagnetism).
Since the cubic lattice can be decomposed into two independent sub-lattices (depending on the parity of (x + y + z)), both cases are symmetric and share the same critical temperature Tc=4.5115...
The anti-ferromagnetic case is visually interesting with the alternating planes of +1 and -1 spins.
All the Monte-Carlo simulation is done as a WebGPU compute-shader (flipping the spins, drawing a random number and rejecting the move or not, etc.) which run several loops internally before a pass of visualisation is done.
Just for fun, we scan the whole range of temperature [Tmin, Tmax], performing a warm-up steps, followed by a proper sampling of the magnetic field and internal energy.
The resulting magnetic field |M| and specific heat Cv are plotted on the bottom-left. Instantaneous magnetic field is plotted on the bottom-right.
Of course, a proper scientific Monte-Carlo experiment should let run the system at a single temperature value for millions of thermalisation steps before any measurement.
An often overlooked part of the Monte-Carlo simulation is the Random Number Generator, which needs to be high-quality, given the millions of draws involved, which can easily exhaust the period of the RNG if not careful.
Here, each spins carries its own RNG state as a triplet of u32 which are used for three Tausworthe steps with different periods, and then combined into a final 24bit result.
This is how you do to stay massively parallel!