nx_gassim.nx source
↩ module page · 133 lines · 6922 B
1// nx_gassim.nx -- Simulates elastic collisions and motion of particles in a 2D box to demonstrate kinetic theory and thermalization.
2const O_MAGIC_1600: i64 = 1600
3const O_MAGIC_3200: i64 = 3200
4const O_MAGIC_2000: i64 = 2000
5const O_MAGIC_3600: i64 = 3600
6const O_MAGIC_4096: i64 = 4096
7// nx_gassim.nx -- BROWSER-RUNNABLE kinetic-theory gas / molecular-dynamics simulator (base-relative, integer/no
8// float, no JS logic). The SAME code runs native (base = mmap -> PNG verify) and in the browser (base = 0 =
9// sovereign WASM linear memory). This is the computational-CHEMISTRY pillar: N=20 equal-mass hard disks in a box,
10// moving ballistically with EXACT elastic wall bounces (sign flip) and pairwise ELASTIC collisions (equal-mass
11// normal-component exchange -- v_i -= (dv.dp/|dp|^2) dp, which conserves total momentum EXACTLY per collision and
12// kinetic energy up to integer rounding). All particles start at the SAME speed (200); collisions thermalize them
13// so a velocity DISTRIBUTION EMERGES from the delta -- the onset of Maxwell-Boltzmann, the foundation of kinetic
14// theory / statistical mechanics. Particles are coloured by speed (blue slow -> red fast) so you watch the gas
15// thermalize. nx_wasm auto-exports every func + the linear memory. license_tier: ORIGINAL
16const W: i64 = 256
17const H: i64 = 256
18const O_FB: i64 = 0 // framebuffer W*H i64 packed-RGB
19const N: i64 = 20 // particles
20const O_PX: i64 = 524288 // px[N] (sim coords, box 0..BW)
21const O_PY: i64 = 524448
22const O_VX: i64 = 524608
23const O_VY: i64 = 524768
24const O_COL: i64 = 524928 // [0] cumulative collision count
25const O_STEP: i64 = 524936 // [0] step counter
26const BW: i64 = 16000 // sim box width/height
27const RAD: i64 = 240 // particle radius (sim units)
28const TWO_R2: i64 = 230400 // (2*RAD)^2
29const STEPS_PER_TICK: i64 = 2
30
31func rgb(r: i64, g: i64, b: i64) -> i64 { return (r & 255) | ((g & 255) << 8) | ((b & 255) << 16) }
32func isqrt(N0: i64) -> i64 { if N0<2 { return N0 } var x: i64=N0; var y: i64=(x+1)/2; while y<x { x=y; y=(x + N0/x)/2 } return x }
33
34func putpx(base: i64, px: i64, py: i64, c: i64) -> i64 {
35 if px<0 { return 0 } if px>=W { return 0 } if py<0 { return 0 } if py>=H { return 0 }
36 let fb: *i64 = (base + O_FB) as *i64; fb[py*W+px]=c; return 0
37}
38func disc(base: i64, px: i64, py: i64, r: i64, c: i64) -> i64 {
39 var dy: i64=0-r
40 while dy<=r { var dx: i64=0-r
41 while dx<=r { if dx*dx+dy*dy <= r*r { putpx(base, px+dx, py+dy, c) } dx=dx+1 }
42 dy=dy+1 }
43 return 0
44}
45func clear_fb(base: i64, c: i64) -> i64 { let fb: *i64=(base+O_FB) as *i64; var i: i64=0; while i<W*H { fb[i]=c; i=i+1 } return 0 }
46// round-to-nearest integer division (den>0). Truncation-toward-zero in the collision factor biases every collision
47// weak -> systematic energy LOSS; rounding to nearest removes the bias so energy does a small bounded walk instead.
48func rdiv(num: i64, den: i64) -> i64 { if num>=0 { return (num + den/2)/den } return 0 - ((0-num + den/2)/den) }
49
50func energy(base: i64) -> i64 {
51 let VX: *i64=(base+O_VX) as *i64; let VY: *i64=(base+O_VY) as *i64
52 var e: i64=0; var i: i64=0; while i<N { e=e+VX[i]*VX[i]+VY[i]*VY[i]; i=i+1 } return e
53}
54
55func init_impl(base: i64) -> i64 {
56 let PX: *i64=(base+O_PX) as *i64; let PY: *i64=(base+O_PY) as *i64; let VX: *i64=(base+O_VX) as *i64; let VY: *i64=(base+O_VY) as *i64
57 // positions: 5x4 grid, spacing 3200 >> 2*RAD (no initial overlap)
58 var gx: i64=0
59 while gx<5 { var gy: i64=0
60 while gy<4 { let p: i64=gx*4+gy; PX[p]=O_MAGIC_1600+gx*O_MAGIC_3200; PY[p]=O_MAGIC_2000+gy*O_MAGIC_3600; gy=gy+1 }
61 gx=gx+1 }
62 // velocities: 20 distinct directions, ALL magnitude exactly 200 (Pythagorean: 200,160-120,120-160,56-192,192-56)
63 VX[0]=200; VY[0]=0; VX[1]=0-200; VY[1]=0; VX[2]=0; VY[2]=200; VX[3]=0; VY[3]=0-200
64 VX[4]=160; VY[4]=120; VX[5]=160; VY[5]=0-120; VX[6]=0-160; VY[6]=120; VX[7]=0-160; VY[7]=0-120
65 VX[8]=120; VY[8]=160; VX[9]=120; VY[9]=0-160; VX[10]=0-120; VY[10]=160; VX[11]=0-120; VY[11]=0-160
66 VX[12]=56; VY[12]=192; VX[13]=56; VY[13]=0-192; VX[14]=0-56; VY[14]=192; VX[15]=0-56; VY[15]=0-192
67 VX[16]=192; VY[16]=56; VX[17]=192; VY[17]=0-56; VX[18]=0-192; VY[18]=56; VX[19]=0-192; VY[19]=0-56
68 let col: *i64=(base+O_COL) as *i64; col[0]=0
69 let stp: *i64=(base+O_STEP) as *i64; stp[0]=0
70 return 0
71}
72
73func step_one_impl(base: i64) -> i64 {
74 let PX: *i64=(base+O_PX) as *i64; let PY: *i64=(base+O_PY) as *i64; let VX: *i64=(base+O_VX) as *i64; let VY: *i64=(base+O_VY) as *i64
75 let col: *i64=(base+O_COL) as *i64
76 // 1) ballistic move + EXACT elastic wall bounce (sign flip conserves energy exactly)
77 var i: i64=0
78 while i<N {
79 PX[i]=PX[i]+VX[i]; PY[i]=PY[i]+VY[i]
80 if PX[i] < RAD { PX[i]=RAD; VX[i]=0-VX[i] }
81 if PX[i] > BW-RAD { PX[i]=BW-RAD; VX[i]=0-VX[i] }
82 if PY[i] < RAD { PY[i]=RAD; VY[i]=0-VY[i] }
83 if PY[i] > BW-RAD { PY[i]=BW-RAD; VY[i]=0-VY[i] }
84 i=i+1
85 }
86 // 2) pairwise elastic collisions (equal mass, normal-component exchange)
87 i=0
88 while i<N {
89 var j: i64=i+1
90 while j<N {
91 let dpx: i64=PX[i]-PX[j]; let dpy: i64=PY[i]-PY[j]; let d2: i64=dpx*dpx+dpy*dpy
92 if d2 < TWO_R2 { if d2 > 0 {
93 let dvx: i64=VX[i]-VX[j]; let dvy: i64=VY[i]-VY[j]; let dot: i64=dvx*dpx+dvy*dpy
94 if dot < 0 { // approaching -> collide
95 let fx: i64=rdiv(dot*dpx, d2); let fy: i64=rdiv(dot*dpy, d2)
96 VX[i]=VX[i]-fx; VY[i]=VY[i]-fy; VX[j]=VX[j]+fx; VY[j]=VY[j]+fy
97 col[0]=col[0]+1
98 }
99 } }
100 j=j+1
101 }
102 i=i+1
103 }
104 return 0
105}
106func tick_impl(base: i64, cmd: i64) -> i64 {
107 if cmd==3 { init_impl(base); return 0 }
108 var s: i64=0; while s<STEPS_PER_TICK { step_one_impl(base); s=s+1 }
109 let stp: *i64=(base+O_STEP) as *i64; stp[0]=stp[0]+1
110 return 0
111}
112func render_impl(base: i64) -> i64 {
113 let PX: *i64=(base+O_PX) as *i64; let PY: *i64=(base+O_PY) as *i64; let VX: *i64=(base+O_VX) as *i64; let VY: *i64=(base+O_VY) as *i64
114 clear_fb(base, rgb(12,12,20))
115 var i: i64=0
116 while i<N {
117 let dx: i64=PX[i]*W/BW; let dy: i64=PY[i]*H/BW
118 let spd: i64=isqrt(VX[i]*VX[i]+VY[i]*VY[i])
119 var rr: i64=spd*255/360; if rr>255 { rr=255 } // blue(slow) -> red(fast)
120 disc(base, dx, dy, 3, rgb(rr, 80, 255-rr))
121 i=i+1
122 }
123 return 0
124}
125
126// ---- wasm interface ----
127func ww() -> i64 { return W }
128func hh() -> i64 { return H }
129func fb_off() -> i64 { return O_FB }
130func init() -> i64 { return init_impl(0) }
131func tick(cmd: i64) -> i64 { return tick_impl(0, cmd) }
132func render() -> i64 { return render_impl(0) }
133func mem_bytes() -> i64 { return O_STEP + O_MAGIC_4096 }