code wiki / (root) / nx_pendulum.nx

nx_pendulum.nx source

↩ module page · 131 lines · 6238 B

1// nx_pendulum.nx -- Simulates a nonlinear pendulum using integer arithmetic and symplectic integration for accurate physical behavior. 2const O_MAGIC_3141593: i64 = 3141593 3const O_MAGIC_6283185: i64 = 6283185 4const O_MAGIC_1570796: i64 = 1570796 5const O_MAGIC_1047198: i64 = 1047198 6const O_MAGIC_4096: i64 = 4096 7// nx_pendulum.nx -- BROWSER-RUNNABLE real pendulum simulator (base-relative, integer/no-float, no JS logic). SAME 8// code runs native (base = mmap -> PNG) and in the browser (base = 0 = sovereign WASM memory). Integrates the TRUE 9// nonlinear pendulum theta'' = -(g/L) sin(theta) with a symplectic (semi-implicit) integer step and our own Taylor 10// sine. The real, experimentally-established pendulum laws emerge: (a) Galileo's T proportional to sqrt(L), (b) 11// small-amplitude ISOCHRONISM (period ~ independent of amplitude), (c) large-amplitude period GROWTH (~+18% at 90 12// deg) -- the nonlinear effect a LINEARIZED fake (-k theta) cannot reproduce. Restoring strength GLN encodes g/L 13// (GLN=980 => g=9.8, L=1; period ~2007 steps at dt=1ms). nx_wasm auto-exports every func + the linear memory. 14// license_tier: ORIGINAL 15const W: i64 = 256 16const H: i64 = 256 17const O_FB: i64 = 0 18const O_TH: i64 = 524288 // theta (micro-rad, x1e6) 19const O_OM: i64 = 524296 // omega (micro-rad / s) 20const O_GLN: i64 = 524304 // restoring numerator = 980/L (g/L * dt * 1e5, g=9.8 dt=1ms) 21const O_AMP: i64 = 524312 // initial amplitude (micro-rad) 22const O_LIN: i64 = 524320 // 0 = true sin pendulum, 1 = linear approximation (for the liar-kill) 23const O_STEP: i64 = 524328 24const O_THEAD: i64 = 524336 25const O_TCNT: i64 = 524344 26const O_TR: i64 = 524352 // bob trail (display px,py) 27const TRAILN: i64 = 120 28const SC: i64 = 1000000 // micro-unit scale 29const DT_DIV: i64 = 1000 // theta += omega/1000 (dt = 1ms) 30const GL_DIV: i64 = 100000 // omega -= GLN*sin(theta)/1e5 31const CX: i64 = 128 32const PIVY: i64 = 40 33const LPX: i64 = 155 // rod length in pixels 34const STEPS_PER_TICK: i64 = 16 35 36func rgb(r: i64, g: i64, b: i64) -> i64 { return (r & 255) | ((g & 255) << 8) | ((b & 255) << 16) } 37func sin_fixed(Xin: i64) -> i64 { 38 let PI: i64=O_MAGIC_3141593; let TWO_PI: i64=O_MAGIC_6283185; let HALF: i64=O_MAGIC_1570796 39 var sign: i64=1; var X: i64=Xin 40 if X<0 { X=0-X; sign=0-1 } 41 X = X % TWO_PI 42 if X > PI { X = TWO_PI - X; sign = 0-sign } 43 if X > HALF { X = PI - X } 44 let x2: i64 = X*X/SC 45 var term: i64=X; var sum: i64=X; var k: i64=1 46 while k<=8 { term = (0-term)*x2/SC/((2*k)*(2*k+1)); sum=sum+term; k=k+1 } 47 return sum*sign 48} 49func cos_fixed(X: i64) -> i64 { return sin_fixed(X+O_MAGIC_1570796) } 50 51func putpx(base: i64, px: i64, py: i64, c: i64) -> i64 { 52 if px<0 { return 0 } if px>=W { return 0 } if py<0 { return 0 } if py>=H { return 0 } 53 let fb: *i64 = (base + O_FB) as *i64; fb[py*W+px]=c; return 0 54} 55func disc(base: i64, px: i64, py: i64, r: i64, c: i64) -> i64 { 56 var dy: i64=0-r 57 while dy<=r { var dx: i64=0-r 58 while dx<=r { if dx*dx+dy*dy <= r*r { putpx(base, px+dx, py+dy, c) } dx=dx+1 } 59 dy=dy+1 } 60 return 0 61} 62func line(base: i64, x0: i64, y0: i64, x1: i64, y1: i64, c: i64) -> i64 { 63 var dx: i64=x1-x0; if dx<0 { dx=0-dx } 64 var dy: i64=y1-y0; if dy<0 { dy=0-dy } 65 var n: i64=dx; if dy>n { n=dy } 66 if n==0 { putpx(base,x0,y0,c); return 0 } 67 var i: i64=0 68 while i<=n { putpx(base, x0+(x1-x0)*i/n, y0+(y1-y0)*i/n, c); i=i+1 } 69 return 0 70} 71func 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 } 72 73func init_impl(base: i64) -> i64 { 74 let TH: *i64=(base+O_TH) as *i64; let OM: *i64=(base+O_OM) as *i64; let GLN: *i64=(base+O_GLN) as *i64 75 let AMP: *i64=(base+O_AMP) as *i64; let LIN: *i64=(base+O_LIN) as *i64 76 GLN[0]=980; AMP[0]=O_MAGIC_1047198; LIN[0]=0 // L=1m, amplitude 60deg, true pendulum 77 TH[0]=AMP[0]; OM[0]=0 78 let stp: *i64=(base+O_STEP) as *i64; stp[0]=0 79 let th2: *i64=(base+O_THEAD) as *i64; th2[0]=0 80 let tc: *i64=(base+O_TCNT) as *i64; tc[0]=0 81 return 0 82} 83// symplectic (semi-implicit Euler) step: update omega using theta_n, then theta using omega_{n+1} 84func step_one_impl(base: i64) -> i64 { 85 let TH: *i64=(base+O_TH) as *i64; let OM: *i64=(base+O_OM) as *i64; let GLN: *i64=(base+O_GLN) as *i64; let LIN: *i64=(base+O_LIN) as *i64 86 var th: i64=TH[0]; var om: i64=OM[0] 87 var rest: i64=0 88 if LIN[0]==1 { rest = GLN[0]*th/GL_DIV } else { rest = GLN[0]*sin_fixed(th)/GL_DIV } 89 om = om - rest 90 th = th + om/DT_DIV 91 TH[0]=th; OM[0]=om 92 return 0 93} 94func record_trail(base: i64) -> i64 { 95 let TH: *i64=(base+O_TH) as *i64; let TR: *i64=(base+O_TR) as *i64 96 let th: i64=TH[0] 97 let bx: i64=CX + LPX*sin_fixed(th)/SC; let by: i64=PIVY + LPX*cos_fixed(th)/SC 98 let h0: *i64=(base+O_THEAD) as *i64; let tc: *i64=(base+O_TCNT) as *i64 99 let h: i64=h0[0]; TR[h*2]=bx; TR[h*2+1]=by 100 h0[0]=(h+1)%TRAILN; if tc[0]<TRAILN { tc[0]=tc[0]+1 } 101 return 0 102} 103func tick_impl(base: i64, cmd: i64) -> i64 { 104 if cmd==3 { init_impl(base); return 0 } 105 var s: i64=0; while s<STEPS_PER_TICK { step_one_impl(base); s=s+1 } 106 record_trail(base) 107 let stp: *i64=(base+O_STEP) as *i64; stp[0]=stp[0]+1 108 return 0 109} 110func render_impl(base: i64) -> i64 { 111 let TH: *i64=(base+O_TH) as *i64; let TR: *i64=(base+O_TR) as *i64; let tc: *i64=(base+O_TCNT) as *i64 112 clear_fb(base, rgb(10,12,20)) 113 let th: i64=TH[0] 114 let bx: i64=CX + LPX*sin_fixed(th)/SC; let by: i64=PIVY + LPX*cos_fixed(th)/SC 115 // bob trail 116 var t: i64=0 117 while t<tc[0] { putpx(base, TR[t*2], TR[t*2+1], rgb(60,90,120)); t=t+1 } 118 disc(base, CX, PIVY, 3, rgb(180,180,190)) // pivot 119 line(base, CX, PIVY, bx, by, rgb(150,160,180)) // rod 120 disc(base, bx, by, 7, rgb(90,200,235)) // bob 121 return 0 122} 123 124// ---- wasm interface ---- 125func ww() -> i64 { return W } 126func hh() -> i64 { return H } 127func fb_off() -> i64 { return O_FB } 128func init() -> i64 { return init_impl(0) } 129func tick(cmd: i64) -> i64 { return tick_impl(0, cmd) } 130func render() -> i64 { return render_impl(0) } 131func mem_bytes() -> i64 { return O_TR + TRAILN*2*8 + O_MAGIC_4096 }