nx_research_orbit_gate.nx source

1// nx_research_orbit_gate.nx -- the capstone of the dynamics thread: a 2-D gravitational orbit (inverse-square 2// central force, a = -GM r / |r|^3) integrated in sovereign fixed-point (using our isqrt for |r|). THE canonical 3// reason symplectic integrators exist: over many orbits the SYMPLECTIC leapfrog keeps the orbit BOUNDED and energy 4// conserved, while naive forward-EULER pumps energy in and the orbit SPIRALS OUTWARD. 5// Fixed-point precision trick: store velocity scaled by VS=2D, so the leapfrog half-kick v += a*dt/2 = a/(2D) 6// becomes v_scaled += a EXACTLY (no truncation -- the dominant error source). Position step then divides by 2D*D. 7// Verifier = energy conservation + orbit-radius bound; liar-kill = Euler's plausible orbit fails the bound the 8// symplectic orbit holds. No float, deterministic. GREEN iff 6/6. license_tier: ORIGINAL 9import "nx_syscalls.nx" 10 11func g_w(s: *u8) -> i64 { var n: i64=0; while s[n]!=(0 as u8){n=n+1} sys_write(1,s,n); return 0 } 12func g_n(v: i64) -> i64 { var m: i64=v; if m<0{g_w("-");m=0-m} let t:*u8=sys_mmap(24); var k:i64=0; if m==0{t[0]=48 as u8;k=1}; while m>0{t[k]=(48+(m%10)) as u8;m=m/10;k=k+1}; var i:i64=0; let o:*u8=sys_mmap(24); while i<k{o[i]=t[k-1-i];i=i+1}; sys_write(1,o,k); return 0 } 13func g_row(id: *u8, ok: i64, pass: *i64) -> i64 { g_w(" "); g_w(id); g_w(": "); if ok==1 { g_w("OK\n"); pass[0]=pass[0]+1 } else { g_w("FAIL\n") } return 0 } 14func iabs(x: i64) -> i64 { if x<0 { return 0-x } return x } 15func isqrt(N: i64) -> i64 { if N<2 { return N } var x: i64=N; var y: i64=(x+1)/2; while y<x { x=y; y=(x + N/x)/2 } return x } 16 17func accel(x: i64, y: i64, GM: i64, out: *i64) -> i64 { 18 let r2: i64 = x*x + y*y 19 let r: i64 = isqrt(r2) 20 let r3: i64 = r2*r 21 out[0] = 0 - GM*x/r3 22 out[1] = 0 - GM*y/r3 23 return 0 24} 25// 2x specific energy using scaled velocities (v_real = vs/VS): (vsx^2+vsy^2)/VS^2 - 2GM/r. 26func orbit_E2(x: i64, y: i64, vsx: i64, vsy: i64, GM: i64, VS: i64) -> i64 { let r: i64=isqrt(x*x+y*y); let vrx: i64=vsx/VS; let vry: i64=vsy/VS; return vrx*vrx+vry*vry - 2*GM/r } 27 28// leapfrog with velocity scaled by VS=2D (exact half-kick). out[0]=max|E-E0|, out[1]=min r, out[2]=max r. 29func sim_lf_orbit(x0: i64, y0: i64, vx0: i64, vy0: i64, GM: i64, D: i64, N: i64, out: *i64) -> i64 { 30 let VS: i64=2*D; let DD2: i64=2*D*D 31 var x: i64=x0; var y: i64=y0; var vsx: i64=vx0*VS; var vsy: i64=vy0*VS 32 let E0: i64=orbit_E2(x,y,vsx,vsy,GM,VS) 33 let a: *i64=sys_mmap(8*4) as *i64 34 var maxdev: i64=0; var minr: i64=isqrt(x*x+y*y); var maxr: i64=minr 35 var i: i64=0 36 while i<N { 37 accel(x,y,GM,a) 38 let vhx: i64=vsx + a[0]; let vhy: i64=vsy + a[1] // exact: a*VS/(2D)=a 39 x = x + vhx/DD2; y = y + vhy/DD2 40 accel(x,y,GM,a) 41 vsx = vhx + a[0]; vsy = vhy + a[1] 42 let E: i64=orbit_E2(x,y,vsx,vsy,GM,VS); let dd: i64=iabs(E-E0); if dd>maxdev { maxdev=dd } 43 let r: i64=isqrt(x*x+y*y); if r<minr { minr=r } if r>maxr { maxr=r } 44 i=i+1 45 } 46 out[0]=maxdev; out[1]=minr; out[2]=maxr 47 return 0 48} 49// forward-Euler (same VS scaling). out[0]=max|E-E0|, out[1]=max r. 50func sim_eu_orbit(x0: i64, y0: i64, vx0: i64, vy0: i64, GM: i64, D: i64, N: i64, out: *i64) -> i64 { 51 let VS: i64=2*D; let DD2: i64=2*D*D 52 var x: i64=x0; var y: i64=y0; var vsx: i64=vx0*VS; var vsy: i64=vy0*VS 53 let E0: i64=orbit_E2(x,y,vsx,vsy,GM,VS) 54 let a: *i64=sys_mmap(8*4) as *i64 55 var maxdev: i64=0; var maxr: i64=isqrt(x*x+y*y) 56 var i: i64=0 57 while i<N { 58 accel(x,y,GM,a) 59 let nx: i64=x + vsx/DD2; let ny: i64=y + vsy/DD2 60 let nvsx: i64=vsx + 2*a[0]; let nvsy: i64=vsy + 2*a[1] // euler: a*VS/D = 2a 61 x=nx; y=ny; vsx=nvsx; vsy=nvsy 62 let E: i64=orbit_E2(x,y,vsx,vsy,GM,VS); let dd: i64=iabs(E-E0); if dd>maxdev { maxdev=dd } 63 let r: i64=isqrt(x*x+y*y); if r>maxr { maxr=r } 64 i=i+1 65 } 66 out[0]=maxdev; out[1]=maxr 67 return 0 68} 69 70func main() -> i64 { 71 let pass: *i64 = sys_mmap(8) as *i64; pass[0]=0 72 g_w("=== NX-RESEARCH-ORBIT GATE (gravitational orbit: symplectic stays bounded, Euler spirals out) ===\n") 73 let SC: i64=100000; let D: i64=64; let N: i64=40000 74 let GM: i64=15625000000000 // (SC/8)^2 * SC -> circular speed SC/8 at r=SC 75 let x0: i64=SC; let y0: i64=0; let vx0: i64=0; let vy0: i64=12500 76 let VS: i64=2*D 77 let E0: i64=orbit_E2(x0,y0,vx0*VS,vy0*VS,GM,VS); let aE0: i64=iabs(E0) 78 79 let lp: *i64=sys_mmap(8*4) as *i64; sim_lf_orbit(x0,y0,vx0,vy0,GM,D,N,lp) 80 let eu: *i64=sys_mmap(8*4) as *i64; sim_eu_orbit(x0,y0,vx0,vy0,GM,D,N,eu) 81 let lp_pm: i64=lp[0]*1000/aE0; let eu_pm: i64=eu[0]*1000/aE0 82 83 g_w(" steps="); g_n(N); g_w(" (~12 orbits) r0="); g_n(SC); g_w(" E0(x2)="); g_n(E0); g_w("\n") 84 g_w(" leapfrog: Edev_permille="); g_n(lp_pm); g_w(" rmin/SC*1000="); g_n(lp[1]/(SC/1000)); g_w(" rmax/SC*1000="); g_n(lp[2]/(SC/1000)); g_w("\n") 85 g_w(" euler: Edev_permille="); g_n(eu_pm); g_w(" rmax/SC*1000="); g_n(eu[1]/(SC/1000)); g_w("\n") 86 87 g_row("SETUP: bound orbit initialized (E<0), isqrt gives r0=SC" as *u8, ((E0<0) as i64)*((isqrt(x0*x0+y0*y0)==SC) as i64), pass) 88 g_row("BOUND (leapfrog): orbit radius stays bounded over ~12 orbits (rmax < 1.3*SC, rmin > 0.7*SC)" as *u8, ((lp[2] < SC+SC*3/10) as i64)*((lp[1] > SC*7/10) as i64), pass) 89 g_row("CONSERVATION (leapfrog): orbital energy stays bounded (deviation < 5%)" as *u8, (lp_pm < 50) as i64, pass) 90 g_row("SPIRAL (euler): forward-Euler pumps energy and the radius grows outward (rmax > 1.2*SC vs leapfrog ~1.0*SC)" as *u8, (eu[1] > SC+SC/5) as i64, pass) 91 g_row("DISCRIMINATION: leapfrog energy deviation far below Euler's (>=3x better)" as *u8, (eu_pm > lp_pm*3) as i64, pass) 92 var liar: i64=0; if lp_pm < 50 { if eu_pm >= 50 { liar=1 } } 93 g_row("LIAR-KILL: Euler's spiraling orbit FAILS the energy bound the symplectic orbit holds" as *u8, liar, pass) 94 95 g_w("RESEARCH-ORBIT-GATE rows=6 pass="); g_n(pass[0]) 96 if pass[0]==6 { g_w(" verdict=GREEN\n"); sys_exit(0); return 0 } 97 g_w(" verdict=RED\n"); sys_exit(1); return 1 98}