code wiki / _hdl_build / nx_research_dynamics_gate.nx
nx_research_dynamics_gate.nx source
↩ module page · 94 lines · 5604 B
1// nx_research_dynamics_gate.nx -- TIME-STEPPING dynamics simulation, sovereign + fixed-point (no float), with a
2// conservation-law verifier (operator NEXT: "time-stepping dynamics sims"). Simulates a harmonic oscillator
3// (a = -x, unit mass/stiffness) two ways and lets the verifier tell them apart by ENERGY conservation:
4// LEAPFROG (symplectic): v_half=v - x/(2D); x'=x + v_half/D; v'=v_half - x'/(2D) -- energy stays BOUNDED.
5// FORWARD-EULER (naive): x'=x + v/D; v'=v - x/D -- energy DRIFTS upward.
6// The well-known physics: a symplectic integrator conserves energy to bounded error over arbitrarily many steps;
7// Euler injects energy and the orbit spirals out. The verifier measures max energy deviation; LIAR-KILL: Euler
8// produces a perfectly plausible trajectory but FAILS the conservation bar leapfrog passes -- a drifting
9// integrator cannot be rubber-stamped as energy-conserving. Fixed-point integers, deterministic. GREEN iff 6/6.
10// license_tier: ORIGINAL
11import "nx_syscalls.nx"
12
13func 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 }
14func 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 }
15func 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 }
16func iabs(x: i64) -> i64 { if x<0 { return 0-x } return x }
17
18// run leapfrog; out[0]=max energy deviation, out[1]=max|x|, out[2]=min x.
19func sim_leapfrog(x0: i64, v0: i64, D: i64, N: i64, out: *i64) -> i64 {
20 var x: i64=x0; var v: i64=v0
21 let E0: i64=x0*x0+v0*v0
22 var maxdev: i64=0; var maxx: i64=iabs(x0); var minx: i64=x0
23 var i: i64=0
24 while i<N {
25 let vh: i64 = v - x/(2*D)
26 x = x + vh/D
27 v = vh - x/(2*D)
28 let E: i64 = x*x + v*v
29 let d: i64 = iabs(E-E0); if d>maxdev { maxdev=d }
30 if iabs(x)>maxx { maxx=iabs(x) }
31 if x<minx { minx=x }
32 i=i+1
33 }
34 out[0]=maxdev; out[1]=maxx; out[2]=minx
35 return 0
36}
37// run forward-Euler; out[0]=max energy deviation, out[1]=max|x|.
38func sim_euler(x0: i64, v0: i64, D: i64, N: i64, out: *i64) -> i64 {
39 var x: i64=x0; var v: i64=v0
40 let E0: i64=x0*x0+v0*v0
41 var maxdev: i64=0; var maxx: i64=iabs(x0)
42 var i: i64=0
43 while i<N {
44 let xn: i64 = x + v/D
45 let vn: i64 = v - x/D
46 x=xn; v=vn
47 let E: i64 = x*x + v*v
48 let d: i64 = iabs(E-E0); if d>maxdev { maxdev=d }
49 if iabs(x)>maxx { maxx=iabs(x) }
50 i=i+1
51 }
52 out[0]=maxdev; out[1]=maxx
53 return 0
54}
55// permille of E0: dev*1000/E0, via integer compare helper returns dev*1000 and caller compares to E0*thr.
56func permille(dev: i64, E0: i64) -> i64 { return dev*1000/E0 } // dev*1000 < 2^63 for our magnitudes
57
58func main() -> i64 {
59 let pass: *i64 = sys_mmap(8) as *i64; pass[0]=0
60 g_w("=== NX-RESEARCH-DYNAMICS GATE (time-stepping sim: symplectic leapfrog conserves energy, Euler drifts) ===\n")
61
62 let SC: i64=1000000; let D: i64=64; let N: i64=2000
63 let x0: i64=SC; let v0: i64=0; let E0: i64=x0*x0+v0*v0
64
65 let lp: *i64=sys_mmap(8*8) as *i64; sim_leapfrog(x0,v0,D,N,lp)
66 let eu: *i64=sys_mmap(8*8) as *i64; sim_euler(x0,v0,D,N,eu)
67
68 let lp_dev: i64=lp[0]; let lp_maxx: i64=lp[1]; let lp_minx: i64=lp[2]
69 let eu_dev: i64=eu[0]; let eu_maxx: i64=eu[1]
70 let lp_pm: i64=permille(lp_dev,E0); let eu_pm: i64=permille(eu_dev,E0)
71
72 g_w(" steps="); g_n(N); g_w(" dt=1/"); g_n(D); g_w(" E0="); g_n(E0); g_w("\n")
73 g_w(" leapfrog: maxEdev_permille="); g_n(lp_pm); g_w(" max|x|/SC*1000="); g_n(lp_maxx/(SC/1000)); g_w(" minx/SC*1000="); g_n(lp_minx/(SC/1000)); g_w("\n")
74 g_w(" euler: maxEdev_permille="); g_n(eu_pm); g_w(" max|x|/SC*1000="); g_n(eu_maxx/(SC/1000)); g_w("\n")
75
76 // row1: leapfrog actually oscillates (swings to the negative side ~ -SC)
77 g_row("SIM: leapfrog oscillator swings through equilibrium (x reaches the negative half-amplitude)" as *u8, (lp_minx < (0-SC/2)) as i64, pass)
78 // row2: leapfrog conserves energy (bounded deviation < 50 permille)
79 g_row("CONSERVATION: symplectic leapfrog keeps energy bounded (max deviation < 5%)" as *u8, (lp_pm < 50) as i64, pass)
80 // row3: euler drifts (deviation > 150 permille)
81 g_row("DRIFT-DETECT: the verifier catches forward-Euler injecting energy (deviation > 15%)" as *u8, (eu_pm > 150) as i64, pass)
82 // row4: discrimination -- euler dev clearly larger than leapfrog (>=5x)
83 g_row("DISCRIMINATION: leapfrog energy deviation is far smaller than Euler's (symplectic measurably better)" as *u8, (eu_pm > lp_pm*5) as i64, pass)
84 // row5: amplitude -- leapfrog amplitude bounded ~SC, euler amplitude grows beyond it
85 var amp: i64=0; if lp_maxx < (SC + SC/5) { if eu_maxx > (SC + SC/8) { amp=1 } }
86 g_row("AMPLITUDE: leapfrog orbit stays bounded (~A) while Euler's spirals outward (amplitude grows)" as *u8, amp, pass)
87 // row6: liar-kill -- Euler fails the SAME conservation bar leapfrog passes
88 var liar: i64=0; if lp_pm < 50 { if eu_pm >= 50 { liar=1 } }
89 g_row("LIAR-KILL: Euler's plausible trajectory FAILS the energy-conservation bar leapfrog passes" as *u8, liar, pass)
90
91 g_w("RESEARCH-DYNAMICS-GATE rows=6 pass="); g_n(pass[0])
92 if pass[0]==6 { g_w(" verdict=GREEN\n"); sys_exit(0); return 0 }
93 g_w(" verdict=RED\n"); sys_exit(1); return 1
94}