diff --git a/.trinity/seals/analog_AnalogSarComparator.json b/.trinity/seals/analog_AnalogSarComparator.json new file mode 100644 index 0000000000..62ec02f8a1 --- /dev/null +++ b/.trinity/seals/analog_AnalogSarComparator.json @@ -0,0 +1,37 @@ +{ + "built_by": "t27c-bootstrap@0.5.2+9b556b7c3", + "closure": [ + { + "spec": "specs/analog/linear.t27", + "spec_hash": "sha256:3ddef7be183aba583ed285334ee389b48a65b41335e98914514235db7c4e5d62" + }, + { + "spec": "specs/analog/mos_table.t27", + "spec_hash": "sha256:168cd93374cec92ece02ddab1be51f49308d9efd3288ac0427f912d4a7521782" + } + ], + "config": "gen=c,rust,verilog,zig;tests=zig", + "gen_hash_c": "sha256:e6f3677d7f54078eabbdd3e3786cd6028cb0445eb4cf8f05b313a3d4ff618225", + "gen_hash_rust": "sha256:0b7f0fe869197f808776e207bceee2075fd591d253ced3f2ac9f6c6c6cd62f97", + "gen_hash_verilog": "sha256:bab3512235c19a5adbe9a9149e589999fc93bbac49484ff62864ed7dae37e140", + "gen_hash_zig": "sha256:ac367deed0470ab29568a4caa691a067ad8d341b9f6e5d356be531a026506510", + "module": "AnalogSarComparator", + "ring": 12, + "seal_schema": 2, + "sealed_at": "2026-10-18T18:44:35Z", + "sealed_by": "t27c-bootstrap@0.5.2", + "spec_hash": "sha256:2df42ed0afe41add0b43762dfa05c8aa38ecbc022bea8f7bbd72f55917bcb806", + "spec_path": "specs/analog/sar_comparator.t27", + "tests": { + "failed": 0, + "failing": [], + "forced": false, + "passed": 20, + "total": 20, + "vacuous": 0 + }, + "toolchain": { + "test_runner": "sha256:8512b80da533b79289110fca4a65f3336abd7a059ca4b72ea595d5e33ff8ad10", + "zig": "0.16.0" + } +} \ No newline at end of file diff --git a/specs/analog/decks/sar_cmp_m05.cir b/specs/analog/decks/sar_cmp_m05.cir new file mode 100644 index 0000000000..5bb8f87ce4 --- /dev/null +++ b/specs/analog/decks/sar_cmp_m05.cir @@ -0,0 +1,63 @@ +SAR ADC dynamic comparator, m05: VDP - VDN = -0.5 mV, backward Euler on a fixed 0.1 ps grid, MOSFET LEVEL=1 +* specs/analog/decks/sar_cmp_m05.cir -- test data for specs/analog/sar_comparator.t27 (#10583). A plain netlist, no .control block. +* Circuit: the dynamic_comparator subcircuit of the uberWAVE golden SAR deck (wp2-parser/golden/top.spice; GPL-2.0, read at +* run time and never vendored) -- thirteen transistors, M1..M13 here for its XM1..XM13, with the W, L and m of that deck -- +* plus the inverter x1 that makes n_clkc from clkc (Mip, Min). Square law: MOSFET LEVEL=1 (Shichman-Hodges), KP 300 uA/V^2 (n) +* and 100 uA/V^2 (p), VTO +-0.35 V, LAMBDA 0, GAMMA 0. beta = KP W / L = 666.7 uA/V^2 per unit m for both types (n: 1u / 0.45u, +* p: 3u / 0.45u). The model sets no TOX, CGSO, CGDO or CJ, so the only capacitance is the explicit one: 4 fF on ap, an, dp, dn; +* 6 fF on net1; 2 fF on net2, net3 and nclk. This is the SQUARE LAW, not PSP103. +* Inputs: VDD 1.5 V, common mode 0.75 V, VDP and VDN as written below (units of 2^-24 V: 4194 per side, so the difference +* is 8388 units = 0.49996 mV). The clock is a B source (it declares no breakpoint): 0 for the first 10 grid points, then a +* ramp of 16 steps of 1.5/16 V, then 1.5 V. nint(time/dt) picks the grid index whatever the offset of ngspice's opening ramp. +* The circuit rests until the clock moves: ap = an = net2 = net3 = nclk = 1.5 V, dp = dn = 0, net1 = max(VDP, VDN) - 0.35 V +* (the higher input device sits at threshold, so no current flows; 0.35 V is taken as 5872026 units of 2^-24 V, the spec's own +* rounding of it, so that the deck and the spec start from the same word). +* Fixed step: gear with maxord 1 is backward Euler, trtol 1e7 stops the truncation test from shrinking a step; reltol, vntol, +* abstol and gmin are tightened so Newton ends far below the 1e-6 V the spec compares. Step j of the run (j = 1..140) ends at the +* grid point with nint(time/dt) = j + 10. Probe sources: mx = max(dp, dn), diff = |dp - dn| (a B source senses without loading). +* Measured (printed by the .meas lines): tck, the clock through 0.75 V; tfirst, the first output through 0.75 V; tres, |dp - dn| +* through 0.75 V; dpmax and dnmax, the highest value each output reaches (the loser's peak); dpend and dnend, the outputs at the +* last point. Decision time = tfirst - tck; resolution time = tres - tck. +* Run: ngspice -b sar_cmp_m05.cir (do not add -r: ngspice refuses .meas together with a raw file) +.model nmos1 nmos level=1 vto=0.35 kp=300u lambda=0 gamma=0 +.model pmos1 pmos level=1 vto=-0.35 kp=100u lambda=0 +Vdd vdd 0 1.5 +Vdp vdp 0 0.749750018119812 +Vdn vdn 0 0.750249981880188 +Bclk clkc 0 V = 1.5*min(max(nint(time/1e-13)-10,0),16)/16 +M1 ap clkc vdd vdd pmos1 w=3u l=0.45u m=9 +M2 an clkc vdd vdd pmos1 w=3u l=0.45u m=9 +M3 ap vdn net1 0 nmos1 w=1u l=0.45u m=36 +M4 an vdp net1 0 nmos1 w=1u l=0.45u m=36 +M5 net1 clkc 0 0 nmos1 w=1u l=0.45u m=18 +M6 net2 dp vdd vdd pmos1 w=3u l=0.45u m=30 +M7 net3 dn vdd vdd pmos1 w=3u l=0.45u m=30 +M10 dn dp 0 0 nmos1 w=1u l=0.45u m=12 +M11 dp dn 0 0 nmos1 w=1u l=0.45u m=12 +M12 dn nclk 0 0 nmos1 w=1u l=0.45u m=6 +M13 dp nclk 0 0 nmos1 w=1u l=0.45u m=6 +M8 dp an net3 vdd pmos1 w=3u l=0.45u m=30 +M9 dn ap net2 vdd pmos1 w=3u l=0.45u m=30 +Mip nclk clkc vdd vdd pmos1 w=3u l=0.45u m=1 +Min nclk clkc 0 0 nmos1 w=1u l=0.45u m=1 +Cap ap 0 4f +Can an 0 4f +Cn1 net1 0 6f +Cn2 net2 0 2f +Cn3 net3 0 2f +Cdp dp 0 4f +Cdn dn 0 4f +Cnc nclk 0 2f +Bmx mx 0 V = max(v(dp),v(dn)) +Bdiff diff 0 V = abs(v(dp)-v(dn)) +.ic v(ap)=1.5 v(an)=1.5 v(net1)=0.4002499580383301 v(net2)=1.5 v(net3)=1.5 v(dp)=0 v(dn)=0 v(nclk)=1.5 +.options method=gear maxord=1 trtol=1e7 reltol=1e-10 vntol=1e-13 abstol=1e-16 gmin=1e-15 itl4=200 +.tran 0.1p 15.5p 0 0.1p uic +.meas tran tck when v(clkc)=0.75 rise=1 +.meas tran tfirst when v(mx)=0.75 rise=1 +.meas tran tres when v(diff)=0.75 rise=1 +.meas tran dpmax max v(dp) +.meas tran dnmax max v(dn) +.meas tran dpend find v(dp) at=15.5p +.meas tran dnend find v(dn) at=15.5p +.end diff --git a/specs/analog/decks/sar_cmp_m2.cir b/specs/analog/decks/sar_cmp_m2.cir new file mode 100644 index 0000000000..41744b1012 --- /dev/null +++ b/specs/analog/decks/sar_cmp_m2.cir @@ -0,0 +1,63 @@ +SAR ADC dynamic comparator, m2: VDP - VDN = -2 mV, backward Euler on a fixed 0.1 ps grid, MOSFET LEVEL=1 +* specs/analog/decks/sar_cmp_m2.cir -- test data for specs/analog/sar_comparator.t27 (#10583). A plain netlist, no .control block. +* Circuit: the dynamic_comparator subcircuit of the uberWAVE golden SAR deck (wp2-parser/golden/top.spice; GPL-2.0, read at +* run time and never vendored) -- thirteen transistors, M1..M13 here for its XM1..XM13, with the W, L and m of that deck -- +* plus the inverter x1 that makes n_clkc from clkc (Mip, Min). Square law: MOSFET LEVEL=1 (Shichman-Hodges), KP 300 uA/V^2 (n) +* and 100 uA/V^2 (p), VTO +-0.35 V, LAMBDA 0, GAMMA 0. beta = KP W / L = 666.7 uA/V^2 per unit m for both types (n: 1u / 0.45u, +* p: 3u / 0.45u). The model sets no TOX, CGSO, CGDO or CJ, so the only capacitance is the explicit one: 4 fF on ap, an, dp, dn; +* 6 fF on net1; 2 fF on net2, net3 and nclk. This is the SQUARE LAW, not PSP103. +* Inputs: VDD 1.5 V, common mode 0.75 V, VDP and VDN as written below (units of 2^-24 V: 16777 per side, so the difference +* is 33554 units = 1.99997 mV). The clock is a B source (it declares no breakpoint): 0 for the first 10 grid points, then a +* ramp of 16 steps of 1.5/16 V, then 1.5 V. nint(time/dt) picks the grid index whatever the offset of ngspice's opening ramp. +* The circuit rests until the clock moves: ap = an = net2 = net3 = nclk = 1.5 V, dp = dn = 0, net1 = max(VDP, VDN) - 0.35 V +* (the higher input device sits at threshold, so no current flows; 0.35 V is taken as 5872026 units of 2^-24 V, the spec's own +* rounding of it, so that the deck and the spec start from the same word). +* Fixed step: gear with maxord 1 is backward Euler, trtol 1e7 stops the truncation test from shrinking a step; reltol, vntol, +* abstol and gmin are tightened so Newton ends far below the 1e-6 V the spec compares. Step j of the run (j = 1..140) ends at the +* grid point with nint(time/dt) = j + 10. Probe sources: mx = max(dp, dn), diff = |dp - dn| (a B source senses without loading). +* Measured (printed by the .meas lines): tck, the clock through 0.75 V; tfirst, the first output through 0.75 V; tres, |dp - dn| +* through 0.75 V; dpmax and dnmax, the highest value each output reaches (the loser's peak); dpend and dnend, the outputs at the +* last point. Decision time = tfirst - tck; resolution time = tres - tck. +* Run: ngspice -b sar_cmp_m2.cir (do not add -r: ngspice refuses .meas together with a raw file) +.model nmos1 nmos level=1 vto=0.35 kp=300u lambda=0 gamma=0 +.model pmos1 pmos level=1 vto=-0.35 kp=100u lambda=0 +Vdd vdd 0 1.5 +Vdp vdp 0 0.7490000128746033 +Vdn vdn 0 0.7509999871253967 +Bclk clkc 0 V = 1.5*min(max(nint(time/1e-13)-10,0),16)/16 +M1 ap clkc vdd vdd pmos1 w=3u l=0.45u m=9 +M2 an clkc vdd vdd pmos1 w=3u l=0.45u m=9 +M3 ap vdn net1 0 nmos1 w=1u l=0.45u m=36 +M4 an vdp net1 0 nmos1 w=1u l=0.45u m=36 +M5 net1 clkc 0 0 nmos1 w=1u l=0.45u m=18 +M6 net2 dp vdd vdd pmos1 w=3u l=0.45u m=30 +M7 net3 dn vdd vdd pmos1 w=3u l=0.45u m=30 +M10 dn dp 0 0 nmos1 w=1u l=0.45u m=12 +M11 dp dn 0 0 nmos1 w=1u l=0.45u m=12 +M12 dn nclk 0 0 nmos1 w=1u l=0.45u m=6 +M13 dp nclk 0 0 nmos1 w=1u l=0.45u m=6 +M8 dp an net3 vdd pmos1 w=3u l=0.45u m=30 +M9 dn ap net2 vdd pmos1 w=3u l=0.45u m=30 +Mip nclk clkc vdd vdd pmos1 w=3u l=0.45u m=1 +Min nclk clkc 0 0 nmos1 w=1u l=0.45u m=1 +Cap ap 0 4f +Can an 0 4f +Cn1 net1 0 6f +Cn2 net2 0 2f +Cn3 net3 0 2f +Cdp dp 0 4f +Cdn dn 0 4f +Cnc nclk 0 2f +Bmx mx 0 V = max(v(dp),v(dn)) +Bdiff diff 0 V = abs(v(dp)-v(dn)) +.ic v(ap)=1.5 v(an)=1.5 v(net1)=0.4009999632835388 v(net2)=1.5 v(net3)=1.5 v(dp)=0 v(dn)=0 v(nclk)=1.5 +.options method=gear maxord=1 trtol=1e7 reltol=1e-10 vntol=1e-13 abstol=1e-16 gmin=1e-15 itl4=200 +.tran 0.1p 15.5p 0 0.1p uic +.meas tran tck when v(clkc)=0.75 rise=1 +.meas tran tfirst when v(mx)=0.75 rise=1 +.meas tran tres when v(diff)=0.75 rise=1 +.meas tran dpmax max v(dp) +.meas tran dnmax max v(dn) +.meas tran dpend find v(dp) at=15.5p +.meas tran dnend find v(dn) at=15.5p +.end diff --git a/specs/analog/decks/sar_cmp_p05.cir b/specs/analog/decks/sar_cmp_p05.cir new file mode 100644 index 0000000000..15fe15d062 --- /dev/null +++ b/specs/analog/decks/sar_cmp_p05.cir @@ -0,0 +1,64 @@ +SAR ADC dynamic comparator, p05: VDP - VDN = +0.5 mV, backward Euler on a fixed 0.1 ps grid, MOSFET LEVEL=1 +* specs/analog/decks/sar_cmp_p05.cir -- test data for specs/analog/sar_comparator.t27 (#10583). A plain netlist, no .control block. +* Circuit: the dynamic_comparator subcircuit of the uberWAVE golden SAR deck (wp2-parser/golden/top.spice; GPL-2.0, read at +* run time and never vendored) -- thirteen transistors, M1..M13 here for its XM1..XM13, with the W, L and m of that deck -- +* plus the inverter x1 that makes n_clkc from clkc (Mip, Min). Square law: MOSFET LEVEL=1 (Shichman-Hodges), KP 300 uA/V^2 (n) +* and 100 uA/V^2 (p), VTO +-0.35 V, LAMBDA 0, GAMMA 0. beta = KP W / L = 666.7 uA/V^2 per unit m for both types (n: 1u / 0.45u, +* p: 3u / 0.45u). The model sets no TOX, CGSO, CGDO or CJ, so the only capacitance is the explicit one: 4 fF on ap, an, dp, dn; +* 6 fF on net1; 2 fF on net2, net3 and nclk. This is the SQUARE LAW, not PSP103. +* Inputs: VDD 1.5 V, common mode 0.75 V, VDP and VDN as written below (units of 2^-24 V: 4194 per side, so the difference +* is 8388 units = 0.49996 mV). The clock is a B source (it declares no breakpoint): 0 for the first 10 grid points, then a +* ramp of 16 steps of 1.5/16 V, then 1.5 V. nint(time/dt) picks the grid index whatever the offset of ngspice's opening ramp. +* The circuit rests until the clock moves: ap = an = net2 = net3 = nclk = 1.5 V, dp = dn = 0, net1 = max(VDP, VDN) - 0.35 V +* (the higher input device sits at threshold, so no current flows; 0.35 V is taken as 5872026 units of 2^-24 V, the spec's own +* rounding of it, so that the deck and the spec start from the same word). +* Fixed step: gear with maxord 1 is backward Euler, trtol 1e7 stops the truncation test from shrinking a step; reltol, vntol, +* abstol and gmin are tightened so Newton ends far below the 1e-6 V the spec compares. Step j of the run (j = 1..140) ends at the +* grid point with nint(time/dt) = j + 10. Probe sources: mx = max(dp, dn), diff = |dp - dn| (a B source senses without loading). +* Measured (printed by the .meas lines): tck, the clock through 0.75 V; tfirst, the first output through 0.75 V; tres, |dp - dn| +* through 0.75 V; dpmax and dnmax, the highest value each output reaches (the loser's peak); dpend and dnend, the outputs at the +* last point. Decision time = tfirst - tck; resolution time = tres - tck. +* Run: ngspice -b sar_cmp_p05.cir (do not add -r: ngspice refuses .meas together with a raw file) +.model nmos1 nmos level=1 vto=0.35 kp=300u lambda=0 gamma=0 +.model pmos1 pmos level=1 vto=-0.35 kp=100u lambda=0 +Vdd vdd 0 1.5 +Vdp vdp 0 0.750249981880188 +Vdn vdn 0 0.749750018119812 +Bclk clkc 0 V = 1.5*min(max(nint(time/1e-13)-10,0),16)/16 +M1 ap clkc vdd vdd pmos1 w=3u l=0.45u m=9 +M2 an clkc vdd vdd pmos1 w=3u l=0.45u m=9 +M3 ap vdn net1 0 nmos1 w=1u l=0.45u m=36 +M4 an vdp net1 0 nmos1 w=1u l=0.45u m=36 +M5 net1 clkc 0 0 nmos1 w=1u l=0.45u m=18 +M6 net2 dp vdd vdd pmos1 w=3u l=0.45u m=30 +M7 net3 dn vdd vdd pmos1 w=3u l=0.45u m=30 +M10 dn dp 0 0 nmos1 w=1u l=0.45u m=12 +M11 dp dn 0 0 nmos1 w=1u l=0.45u m=12 +M12 dn nclk 0 0 nmos1 w=1u l=0.45u m=6 +M13 dp nclk 0 0 nmos1 w=1u l=0.45u m=6 +M8 dp an net3 vdd pmos1 w=3u l=0.45u m=30 +M9 dn ap net2 vdd pmos1 w=3u l=0.45u m=30 +Mip nclk clkc vdd vdd pmos1 w=3u l=0.45u m=1 +Min nclk clkc 0 0 nmos1 w=1u l=0.45u m=1 +Cap ap 0 4f +Can an 0 4f +Cn1 net1 0 6f +Cn2 net2 0 2f +Cn3 net3 0 2f +Cdp dp 0 4f +Cdn dn 0 4f +Cnc nclk 0 2f +Bmx mx 0 V = max(v(dp),v(dn)) +Bdiff diff 0 V = abs(v(dp)-v(dn)) +.ic v(ap)=1.5 v(an)=1.5 v(net1)=0.4002499580383301 v(net2)=1.5 v(net3)=1.5 v(dp)=0 v(dn)=0 v(nclk)=1.5 +.options method=gear maxord=1 trtol=1e7 reltol=1e-10 vntol=1e-13 abstol=1e-16 gmin=1e-15 itl4=200 +.tran 0.1p 15.5p 0 0.1p uic +.meas tran tck when v(clkc)=0.75 rise=1 +.meas tran tfirst when v(mx)=0.75 rise=1 +.meas tran tres when v(diff)=0.75 rise=1 +.meas tran dpmax max v(dp) +.meas tran dnmax max v(dn) +.meas tran dpend find v(dp) at=15.5p +.meas tran dnend find v(dn) at=15.5p +.print tran v(ap) v(an) v(net1) v(net2) v(net3) v(dp) v(dn) v(nclk) +.end diff --git a/specs/analog/decks/sar_cmp_p2.cir b/specs/analog/decks/sar_cmp_p2.cir new file mode 100644 index 0000000000..8ed6f5f39f --- /dev/null +++ b/specs/analog/decks/sar_cmp_p2.cir @@ -0,0 +1,64 @@ +SAR ADC dynamic comparator, p2: VDP - VDN = +2 mV, backward Euler on a fixed 0.1 ps grid, MOSFET LEVEL=1 +* specs/analog/decks/sar_cmp_p2.cir -- test data for specs/analog/sar_comparator.t27 (#10583). A plain netlist, no .control block. +* Circuit: the dynamic_comparator subcircuit of the uberWAVE golden SAR deck (wp2-parser/golden/top.spice; GPL-2.0, read at +* run time and never vendored) -- thirteen transistors, M1..M13 here for its XM1..XM13, with the W, L and m of that deck -- +* plus the inverter x1 that makes n_clkc from clkc (Mip, Min). Square law: MOSFET LEVEL=1 (Shichman-Hodges), KP 300 uA/V^2 (n) +* and 100 uA/V^2 (p), VTO +-0.35 V, LAMBDA 0, GAMMA 0. beta = KP W / L = 666.7 uA/V^2 per unit m for both types (n: 1u / 0.45u, +* p: 3u / 0.45u). The model sets no TOX, CGSO, CGDO or CJ, so the only capacitance is the explicit one: 4 fF on ap, an, dp, dn; +* 6 fF on net1; 2 fF on net2, net3 and nclk. This is the SQUARE LAW, not PSP103. +* Inputs: VDD 1.5 V, common mode 0.75 V, VDP and VDN as written below (units of 2^-24 V: 16777 per side, so the difference +* is 33554 units = 1.99997 mV). The clock is a B source (it declares no breakpoint): 0 for the first 10 grid points, then a +* ramp of 16 steps of 1.5/16 V, then 1.5 V. nint(time/dt) picks the grid index whatever the offset of ngspice's opening ramp. +* The circuit rests until the clock moves: ap = an = net2 = net3 = nclk = 1.5 V, dp = dn = 0, net1 = max(VDP, VDN) - 0.35 V +* (the higher input device sits at threshold, so no current flows; 0.35 V is taken as 5872026 units of 2^-24 V, the spec's own +* rounding of it, so that the deck and the spec start from the same word). +* Fixed step: gear with maxord 1 is backward Euler, trtol 1e7 stops the truncation test from shrinking a step; reltol, vntol, +* abstol and gmin are tightened so Newton ends far below the 1e-6 V the spec compares. Step j of the run (j = 1..140) ends at the +* grid point with nint(time/dt) = j + 10. Probe sources: mx = max(dp, dn), diff = |dp - dn| (a B source senses without loading). +* Measured (printed by the .meas lines): tck, the clock through 0.75 V; tfirst, the first output through 0.75 V; tres, |dp - dn| +* through 0.75 V; dpmax and dnmax, the highest value each output reaches (the loser's peak); dpend and dnend, the outputs at the +* last point. Decision time = tfirst - tck; resolution time = tres - tck. +* Run: ngspice -b sar_cmp_p2.cir (do not add -r: ngspice refuses .meas together with a raw file) +.model nmos1 nmos level=1 vto=0.35 kp=300u lambda=0 gamma=0 +.model pmos1 pmos level=1 vto=-0.35 kp=100u lambda=0 +Vdd vdd 0 1.5 +Vdp vdp 0 0.7509999871253967 +Vdn vdn 0 0.7490000128746033 +Bclk clkc 0 V = 1.5*min(max(nint(time/1e-13)-10,0),16)/16 +M1 ap clkc vdd vdd pmos1 w=3u l=0.45u m=9 +M2 an clkc vdd vdd pmos1 w=3u l=0.45u m=9 +M3 ap vdn net1 0 nmos1 w=1u l=0.45u m=36 +M4 an vdp net1 0 nmos1 w=1u l=0.45u m=36 +M5 net1 clkc 0 0 nmos1 w=1u l=0.45u m=18 +M6 net2 dp vdd vdd pmos1 w=3u l=0.45u m=30 +M7 net3 dn vdd vdd pmos1 w=3u l=0.45u m=30 +M10 dn dp 0 0 nmos1 w=1u l=0.45u m=12 +M11 dp dn 0 0 nmos1 w=1u l=0.45u m=12 +M12 dn nclk 0 0 nmos1 w=1u l=0.45u m=6 +M13 dp nclk 0 0 nmos1 w=1u l=0.45u m=6 +M8 dp an net3 vdd pmos1 w=3u l=0.45u m=30 +M9 dn ap net2 vdd pmos1 w=3u l=0.45u m=30 +Mip nclk clkc vdd vdd pmos1 w=3u l=0.45u m=1 +Min nclk clkc 0 0 nmos1 w=1u l=0.45u m=1 +Cap ap 0 4f +Can an 0 4f +Cn1 net1 0 6f +Cn2 net2 0 2f +Cn3 net3 0 2f +Cdp dp 0 4f +Cdn dn 0 4f +Cnc nclk 0 2f +Bmx mx 0 V = max(v(dp),v(dn)) +Bdiff diff 0 V = abs(v(dp)-v(dn)) +.ic v(ap)=1.5 v(an)=1.5 v(net1)=0.4009999632835388 v(net2)=1.5 v(net3)=1.5 v(dp)=0 v(dn)=0 v(nclk)=1.5 +.options method=gear maxord=1 trtol=1e7 reltol=1e-10 vntol=1e-13 abstol=1e-16 gmin=1e-15 itl4=200 +.tran 0.1p 15.5p 0 0.1p uic +.meas tran tck when v(clkc)=0.75 rise=1 +.meas tran tfirst when v(mx)=0.75 rise=1 +.meas tran tres when v(diff)=0.75 rise=1 +.meas tran dpmax max v(dp) +.meas tran dnmax max v(dn) +.meas tran dpend find v(dp) at=15.5p +.meas tran dnend find v(dn) at=15.5p +.print tran v(ap) v(an) v(net1) v(net2) v(net3) v(dp) v(dn) v(nclk) +.end diff --git a/specs/analog/sar_comparator.t27 b/specs/analog/sar_comparator.t27 new file mode 100644 index 0000000000..27bc593e96 --- /dev/null +++ b/specs/analog/sar_comparator.t27 @@ -0,0 +1,2368 @@ +// SPDX-License-Identifier: Apache-2.0 +; specs/analog/sar_comparator.t27 -- the SAR ADC dynamic comparator: multi-node Newton, square law, against ngspice LEVEL=1 +; M4b of the SAR ADC milestone (#10524, epic #10486, issue #10583). M3a (specs/analog/mos_table.t27) solved ONE nonlinear node; the +; comparator of the golden SAR deck is a regenerative latch with eight internal nodes, so this file adds the multi-node Newton +; step. THE DEVICE LEVEL IS THE SQUARE LAW (cpphdl's MOS law, SPICE MOSFET LEVEL=1), NOT PSP103: IHP SG13G2 tables replace the +; law in M3b (#10516) and change nothing in the solver. Everything hand-written is this file and the four decks in +; specs/analog/decks/sar_cmp_*.cir (plain SPICE netlists, data). +; - +; THE CIRCUIT. The subckt dynamic_comparator of uberWAVE's golden SAR deck (wp2-parser/golden/top.spice, GPL-2.0: read at run time, +; never vendored; this file holds the W, L and m of its thirteen transistors, which are facts of the circuit, and none of its text). +; A StrongARM-style latch. The clock clkc is low in reset: the PMOS XM1, XM2 (m 9) hold the precharge nodes ap, an at Vdd, the +; NMOS XM12, XM13 (m 6, gate n_clkc, the inverter x1's output) hold the outputs dp, dn at 0, and the tail NMOS XM5 (m 18) is off. +; When clkc rises the tail conducts and the input pair XM3 (gate VDN) and XM4 (gate VDP), both m = 36, discharge ap and an at rates +; set by the inputs; the PMOS XM8 (gate an) and XM9 (gate ap), m = 30, through XM7 and XM6 (gates dn, dp, m = 30), pull dp and dn +; up; the cross-coupled NMOS XM10, XM11 (m 12) and the PMOS pair regenerate the difference until one output is at Vdd and the other +; at 0. Larger VDP discharges an faster, so XM8 conducts first and dp wins: the decision is "VDP > VDN" when dp ends high. +; Eight internal nodes: ap, an, net1 (the tail node), net2, net3 (the sources of XM9, XM8), dp, dn and n_clkc. Fifteen devices +; (the thirteen and the inverter's PMOS and NMOS, m = 1). The golden top adds 4 fF on each of dp, dn, an inverter on each (x5, x6) +; and a NAND of the two inverted outputs, "Valid": it is high once the FIRST output has passed an inverter's switching point, +; which is here taken as Vdd / 2. That event, the first output through Vdd / 2, is the DECISION TIME below. +; - +; THE DEVICES. cpphdl's law (analog/include/solver.hpp), as in M3a: ov = Vgs - Vth, no current for ov <= 0; Vds < ov: I = beta +; (ov Vds - Vds^2 / 2); else I = beta ov^2 / 2; Vth = 0.35 V; the PMOS is the mirror; the channel is symmetric (the higher end +; acts as drain for an NMOS and as source for a PMOS). beta = KP W / L m with KP = 300 uA/V^2 (n), 100 uA/V^2 (p), W = 1 um (n), +; 3 um (p), L = 0.45 um as in the golden deck: 666.7 uA/V^2 per unit of m for BOTH types, so one table serves both. LAMBDA = 0, no +; body effect. Capacitance is explicit and nothing else: 4 fF on ap, an, dp, dn (the deck's own 4 fF on dp, dn), 6 fF on net1, 2 fF +; on net2, net3 and n_clkc. Supply 1.5 V, common mode 0.75 V. The ngspice decks carry exactly these numbers (MOSFET LEVEL=1: KP, +; VTO = +-0.35, LAMBDA = 0, GAMMA = 0; no TOX, CGSO or CJ, so no other capacitance). +; - +; THE DEVICE IS A TABLE, AND THE TABLE IS M3A'S. `use analog::mos_table` gives the reader unchanged: id_of, gm_cell, gds_cell, a +; 25 x 25 table read by exact-integer bilinear interpolation, gm and gds the slopes of the same interpolant. Its grid is 24 cells +; of 2^20 units; here the unit is 2^-24 V, so the 24 cells span 1.5 V and a cell is 62.5 mV (M3a's 0.05 V cells span its 1.2 V; the +; reader only sees 24 cells of 2^20). Vdd = 25165824 units = 1.5 * 2^24, and so 0.75 V = 12582912. A table entry is T(x, y) = the +; unit device's current at Vgs = x, Vds = y times dt / C_unit, in the same unit, with dt = 0.1 ps and a charge unit of 0.125 fF: +; dt / C_unit = 800 V/A, the largest entry (1.5 V, 1.5 V) is 5916765 < 2^23, which is M3a's rule for the u64 interpolation. A +; device of multiplier m reads the table once and multiplies by m. THE PRICE: 0.35 V is 5.6 cells, so the threshold lies INSIDE +; cell 5 (0.3125 - 0.375 V) and the table has a soft knee. The chord from 0 at node 5 to 2796 at node 6 reads 783 units at 0.33 V +; where the law is exactly 0 (a test). The input pair, sitting at threshold at rest, therefore conducts a little: see THE REST. +; - +; THE STEP. Backward Euler on all eight nodes at once. Node n, c_n the capacitance in charge units (32, 32, 48, 16, 16, 32, 32, +; 16), v0 the state before the step: G_n(v) = c_n (v_n - v0_n) + (channel currents leaving n) - (channel currents entering n) +; = 0, every current m T(x, y) in the same units, J = dG/dv with the device slopes m gm_cell and m gds_cell (a pure number is +; the cell value over 2^20). Each channel is stamped from its higher end to its lower end, read as x = Vg - Vlow (NMOS) or +; Vhigh - Vg (PMOS) and y = Vhigh - Vlow; the partials of a device are (gds, gm, -gm - gds) on (high, gate, low) for an NMOS and +; (gm + gds, -gm, -gds) for a PMOS; row high takes them, row low takes them negated. K = 3 FIXED Newton iterations from v = v0: +; J x = -G by Gaussian elimination with PARTIAL PIVOTING on the dense 8 x 9 matrix, v <- v + x clamped to [0, Vdd]. The residual +; after the last iteration is computed at every step as the Newton move still owed, max over nodes of round(|G_n| 2^20 / J_nn), +; in units (59.6 nV); the bound is SC_RES_BOUND = 2 units (0.12 uV), met at 31 of the 840 steps of the six runs and never exceeded. +; - +; THE NUMBERS ARE UNSIGNED. The registers J+, J-, G+, G- hold the positive and negative terms as u64 sums (as linear.t27 does for +; its state): a current adds to one register by the direction it flows. The matrix and the solution are sign-magnitude u64, bit 63 +; the direction (sc_pack, sc_sadd); a multiplier is round(|a_rc| 2^20 / |pivot|) <= 2^20; back substitution forms +; (b 2^20 - sum A_ij x_j) / A_ii in two u64 registers. No signed type is used (`t27b verilog` refuses one) and no float is used +; from sc_clk_at to sc_sign. The widest word any stage of the solve forms is 4.95e14, under 2^49 (a test), of 2^63. +; - +; THE RUN. The rest state is: ap = an = net2 = net3 = n_clkc = Vdd, dp = dn = 0, net1 = max(VDP, VDN) - 0.35 V (5872026 units), +; clock at 0. Step j = 1..140 (0.1 ps each) has the clock at Vdd min(j, 16) / 16, a rise of 1.6 ps with its 50 % point exactly at +; step 8. The inputs are held at 0.75 V +- h: h = 16777 units (the differential 2 h = 33554 units = 2 mV less 0.4 uV, "+-2 mV"), +; 4194 (0.5 mV), 1049 (0.125 mV) and 262 (31 uV). The run has no idle steps before the clock: the decks idle 10 grid points at +; rest, which changes nothing there. S, an order-sensitive fold of all 1120 samples, is pinned for each run. +; MEASURED, all by linear interpolation between samples: the sign (dp above dn at step 140); the DECISION TIME, the first output +; through Vdd / 2, counted from the clock's 50 % point; the RESOLUTION TIME, |dp - dn| through Vdd / 2 (the point at which the +; pair has separated by half the supply); the residual. Times carry 12 fractional bits of a step. +; - +; THE REST. At rest every node is in balance except three: the input devices sit at the threshold, inside the knee, so their +; chords conduct. ap and an sag by 91 and 96 uV and net1 rises by 145 uV in the first 0.1 ps (a test), and net1 keeps rising, +; slowly, until it is 0.3125 V under the input. The real circuit and the decks rest exactly (the law is 0 at threshold). This is +; the table's knee and the run includes it: it is part of the table error below, and the PSP103 tables of M3b will have their own +; subthreshold tail. +; - +; THE RESULT, K = 3, fixed point, against ngspice MOSFET LEVEL=1 on the same grid (the .meas lines of the decks, 6 digits): +; decision time (first output) resolution time (|dp - dn| = Vdd / 2) sign +; +-2 mV 5.489 ps ngspice 5.496 ps -0.14 % 6.963 ps ngspice 6.931 ps +0.47 % right +; +-0.5 mV 5.558 ps ngspice 5.569 ps -0.20 % 7.734 ps ngspice 7.694 ps +0.53 % right +; The mirror inputs (-2 mV, -0.5 mV) give the same times and the opposite sign, in ngspice and here. The issue asks for 5 %; the +; worst of the four numbers is 0.53 %, a factor 9 inside it. Two more inputs, 0.125 mV and 31 uV, decide the right way at 8.501 ps +; and 9.267 ps (resolution), the decision time moving only from 5.558 to 5.584 ps: the first output is the common-mode rise of +; both, and only the resolution depends on the input, as the logarithm of it (THE METASTABILITY). +; - +; THE THREE ERRORS, kept apart (W10): +; arithmetic the fixed-point run against a REAL run of the same iteration (the same table words, K = 3, the same clamps, no +; rounding, a real Gauss): under 3 units (0.2 uV) up to step 50, then amplified by the regeneration, about e-fold per 5.5 +; steps, to 68.05 units (4.1 uV) at 2 mV (bound 70, from both sides) and 250.8 units (15 uV) at 0.5 mV (bound 255); +; the decision and resolution times differ by under 3.4e-4 steps (the time is truncated to 1/4096). +; table the real table run against the analytic square law (the law's exact derivatives, 12 Newton iterations a step, the +; same backward Euler grid): the largest |difference| of any node at any step is 15.2 mV at 2 mV (bound 15.4) and +; 19.0 mV at 0.5 mV (bound 19.2). This is the largest error here: the 62.5 mV cells, the knee and the slopes +; that are constant over a cell. In time it is -0.14 %, -0.20 % (decision) and +0.47 %, +0.53 % (resolution). +; model the analytic square law against ngspice LEVEL=1 with the same W, L, m, KP, VTO and capacitors: 1.44e-9 V at the +; eight recorded steps and eight nodes at 2 mV (1.38e-9 V at 0.5 mV), the times agreeing to the 6 printed digits. +; Over every node at every step 1..140 the largest difference is 1.47e-9 V at 2 mV and 1.42e-9 V at 0.5 mV (a +; float64 reader of the raw file, not in the repo). The deck is the spec's circuit. +; The discretization (backward Euler at 0.1 ps) is shared with ngspice by construction: both are the same recurrence on the +; same grid, so it is not an error of this comparison. It is an error of the circuit: the decision times are the grid's, not the +; continuous circuit's. Nothing here claims otherwise. +; - +; THE LOSER. The output that loses rises with the winner and falls back. ngspice's loser peaks 3 mV under Vdd / 2 (0.7473 V at 2 mV, +; 0.8221 V at 0.5 mV); the table model's peaks at 0.7530 V (+5.7 mV), a few mV OVER Vdd / 2. Whether the loser "crosses" Vdd / 2 at +; 2 mV is therefore below this model's resolution (15 mV) and is not used. The winner's crossing and the sign are. +; - +; THE METASTABILITY, reported. Each quartering of the input (2 mV, 0.5, 0.125, 0.031) adds 7.67 steps (0.77 ps) to the resolution +; time: one e-fold of the input costs 5.54 steps = 0.554 ps, the regeneration time constant. At 0 V of input the circuit does not +; decide: in real arithmetic it stays exactly symmetric (dp = dn at every step); in fixed point one unit of rounding at step 60 +; (|dp - dn| = 1) grows to 2027 at step 100 and 2659643 (0.158 V) at step 140, short of Vdd / 2 and set by rounding, not by an +; input; left to run, that noise would decide it later, by chance. The decision time grows without bound as the input goes to +; 0 V, and the spec says so rather than hiding it. +; - +; THE NEGATIVE CONTROLS, each a test that passes by seeing the failure: +; K = 1 the residual reaches 29361 units (1.75 mV) against the bound of 2; K = 2 reaches 55 (3.3 uV) and changes no time; +; K = 4 reaches the same 2 and the same times as K = 3, which is therefore enough. +; inputs swapped, or the input pair's m swapped when they differ (38, 34 against 34, 38): the decision follows the mismatch, not +; the input; the sign check fails for all four stimuli in one orientation and the opposite one in the other. A +; 5.6 % mismatch (35, 37) already outweighs 2 mV. +; Vth = 0.30 a table built with the wrong threshold is 525686 units off the law at a node (bound 0.5), decides the right way +; but 16.7 % early (decision) and 13.7 % early (resolution) against ngspice (bound 5 %), and 0.39 V off at the +; recorded steps (bound 15 mV): the sign alone would not catch it, the comparison with ngspice does. +; And two symmetries: in real arithmetic the circuit with VDP and VDN exchanged gives the node-mirrored run to the last bit (a +; device on the wrong side, a gate on the wrong net would break it; an unequal pair does, by more than 1e6 units), and the +; Jacobian is the derivative of the residual (a central difference, 1e-4 of a pure number; the transposed matrix is off by more +; than 0.5). +; - +; THE DECKS. specs/analog/decks/sar_cmp_p2.cir, _m2, _p05, _m05: one per stimulus, plain netlists. The recipe of bench.t27 for a +; fixed-step backward-Euler reference: `.options method=gear maxord=1 trtol=1e7` (gear with maxord 1 is backward Euler and the +; truncation test cannot shrink a step), a clock as a B source with no breakpoint, `nint(time/dt)` picking the grid index whatever +; the offset of ngspice's opening ramp (here 4 fs before a multiple of dt: step j is the row with time (j + 10) 0.1 ps - 4 fs), +; reltol, vntol, abstol and gmin tightened so that Newton ends far below the 1e-6 V compared here. The circuit rests until the +; clock moves, from .ic. Probe sources mx = max(dp, dn) and diff = |dp - dn| give the two times. Run `ngspice -b `: it +; prints tck, tfirst, tres, dpmax, dnmax, dpend, dnend and the .print table (7 digits). The decision time is tfirst - tck and +; the resolution time tres - tck, in steps: (7.29229e-12 - 1.79600e-12) / 1e-13 = 54.9629 at 2 mV, 69.3083; 55.6936, 76.9385 at +; 0.5 mV. The recorded checkpoints SC_NG_P2 and SC_NG_P05 are the nodes ap an net1 net2 net3 dp dn nclk at steps 20 40 52 60 68 +; 76 90 130, from `ngspice -b -r deck.raw ` (a float64 raw file; ngspice evaluates no .meas with -r), 10 digits. This +; file does not run ngspice: it holds one run's numbers (ngspice-47 from Homebrew, SPARSE 1.3, 2026-10-19). +; - +; THE HARDWARE. Everything from sc_clk_at to sc_sign takes and returns bool and unsigned words and has no float: it is what a +; board would run. `fold` and `fabs64` come from M1 (`use analog::linear`), the table reader, `mag_diff` and the u64 blend from +; M3a (`use analog::mos_table`); the f64 side (the table builder, the real twin of the iteration, the analytic law, the +; measurements on real runs, the ngspice numbers and every comparison) is reached only by the tests. The device table, the node +; capacitances and the multipliers are consts (a ROM on a board), the multipliers also a slice argument so that the tests can +; mismatch the pair. No board top and no clock are in this file; the cost of the 8 x 9 solve (the divisions, the multipliers, +; the table memory) is measured when there is one. sc_board_step is the step with a constant trip count. +; - +; WHAT A PASS SAYS. A comparator of eight nodes and fifteen square-law devices, each read from M3a's bilinear table, is solved by +; backward Euler with three fixed Newton iterations of an 8 x 8 pivoted elimination in unsigned sign-magnitude words, last +; residual at most 2 units (0.12 uV) at all 840 steps of six runs; it decides +-2 mV and +-0.5 mV the right way, with the first +; output at 5.489 and 5.558 ps and the resolution at 6.963 and 7.734 ps, within 0.25 % and 0.6 % of ngspice MOSFET LEVEL=1; its +; arithmetic is within 68 and 251 units of a real run of the same iteration, its table within 15 and 19 mV of the analytic law, +; and the analytic law is ngspice to 1.4e-9 V. Each of K = 1, an inverted pair, a mismatched pair and a wrong threshold fails. +; WHAT IT DOES NOT SAY. Nothing about a transistor: the devices are the square law, not PSP103 (M3b), and nothing is a PDK +; corner; nothing about the continuous circuit (the times are the 0.1 ps backward Euler grid's); nothing about the tracking +; switches, the DAC, the NAND "Valid" or the inverters x5, x6 (the outputs are dp and dn, and Valid is taken at Vdd / 2); nothing +; about noise, kickback, offset from mismatch beyond the controls above, or a clock other than one 1.6 ps rise. Nothing about a +; board: no bitstream, no LUT, DSP or BRAM, no division cost, no UART line. The loser's peak is wrong by 5.7 mV and the soft knee +; moves the rest state of net1 and the early steps; both are the table's, and both are measured above. Nothing about power or +; about speed against SPICE. +; ASCII only, English only. +; phi^2 + 1/phi^2 = 3 | TRINITY + +module AnalogSarComparator; + +use analog::linear::{fold, fabs64}; +use analog::mos_table::{id_of, gm_cell, gds_cell, mag_diff, idf, gdsf, cellf, nodef, fmax64, monotone, law_amps, law_gm, BETA, CELL, CELL_BITS, CELL_F, NODES, VDD, VDD_F}; + +pub const KIND : str = "analog-cell"; +pub const ID : str = "analog/sar_comparator"; +pub const REFERENCE : str = "uberWAVE wp2-parser/golden/top.spice, subckt dynamic_comparator (GPL-2.0: read at run time, never vendored); specs/analog/mos_table.t27"; + +; --- the devices: the square law, one table for both types --------------------------------------- +pub const SC_KP_N : f64 = 0.0003; +pub const SC_KP_P : f64 = 0.0001; +pub const SC_W_N : f64 = 1.0; +pub const SC_W_P : f64 = 3.0; +pub const SC_L : f64 = 0.45; +/// beta of a device with m = 1, A/V^2: KP W / L = 300 uA/V^2 * 1 / 0.45 = 100 uA/V^2 * 3 / 0.45, the same for both types. +pub const SC_BETA : f64 = 0.000666666666666667; +pub const SC_VTH_V : f64 = 0.35; +pub const SC_VTH_WRONG_V : f64 = 0.30; + +; --- the units ----------------------------------------------------------------------------------- +/// The step, in femtoseconds: dt = 0.1 ps. +pub const SC_DT_FS : u32 = 100; +/// The charge unit, in attofarads: 0.125 fF. +pub const SC_CUNIT_AF : u32 = 125; +/// dt / C_unit = 800 V/A, so a table entry is the voltage a 0.125 fF node moves in one step. +pub const SC_STEP_V_PER_A : f64 = 800.0; +/// Voltage units per volt: 2^24. The M3a grid of 24 cells then spans 1.5 V, a cell is 62.5 mV. +pub const SC_UPV : f64 = 16777216.0; +/// Vdd / 2 = 0.75 V, the common mode and the decision level. +pub const SC_VCM : u32 = 12582912; +/// round(0.35 * 2^24): the threshold in units (the table puts it inside cell 5). +pub const SC_VTH_U : u32 = 5872026; +/// Half the differential input, in units: the inputs are VCM + h and VCM - h. 16777 * 2 = 33554 units = 2 mV less 0.4 uV. +pub const SC_H2 : u32 = 16777; +/// 0.5 mV. +pub const SC_H05 : u32 = 4194; +/// 0.125 mV. +pub const SC_H0125 : u32 = 1049; +/// 31 uV. +pub const SC_H003 : u32 = 262; +pub const SC_TBL_MAX : u32 = 5916765; + +; --- the circuit: eight internal nodes, fifteen devices -------------------------------------------- +pub const SC_NN : u32 = 8; +pub const SC_AP : u32 = 0; +pub const SC_AN : u32 = 1; +pub const SC_NET1 : u32 = 2; +pub const SC_NET2 : u32 = 3; +pub const SC_NET3 : u32 = 4; +pub const SC_DP : u32 = 5; +pub const SC_DN : u32 = 6; +pub const SC_NCLK : u32 = 7; +pub const SC_T_VDD : u32 = 8; +pub const SC_T_GND : u32 = 9; +pub const SC_T_CLK : u32 = 10; +pub const SC_T_VDP : u32 = 11; +pub const SC_T_VDN : u32 = 12; +/// Node capacitance in 0.125 fF: 4, 4, 6, 2, 2, 4, 4, 2 fF. +pub const SC_CAP : [8]u32 = [32, 32, 48, 16, 16, 32, 32, 16]; +pub const SC_NDEV : u32 = 15; +/// 0 = NMOS, 1 = PMOS. XM1 XM2 XM3 XM4 XM5 XM6 XM7 XM8 XM9 XM10 XM11 XM12 XM13, then the inverter's PMOS and NMOS. +pub const SC_DV_KIND : [15]u32 = [1, 1, 0, 0, 0, 1, 1, 1, 1, 0, 0, 0, 0, 1, 0]; +pub const SC_DV_D : [15]u32 = [0, 1, 0, 1, 2, 3, 4, 5, 6, 6, 5, 6, 5, 7, 7]; +pub const SC_DV_G : [15]u32 = [10, 10, 12, 11, 10, 5, 6, 1, 0, 5, 6, 7, 7, 10, 10]; +pub const SC_DV_S : [15]u32 = [8, 8, 2, 2, 9, 8, 8, 4, 3, 9, 9, 9, 9, 8, 9]; +pub const SC_DV_M : [15]u32 = [9, 9, 36, 36, 18, 30, 30, 30, 30, 12, 12, 6, 6, 1, 1]; + +; --- the run ------------------------------------------------------------------------------------- +/// Steps of the clock's rise. +pub const SC_R : u32 = 16; +/// Vdd / 16 = 1572864 exactly. +pub const SC_CLK_STEP : u32 = 1572864; +pub const SC_N : u32 = 140; +/// Words of a run's states: 141 states of 8 nodes. +pub const SC_LEN : usize = 1128; +pub const SC_K : u32 = 3; +pub const SC_RES_BOUND : u32 = 2; +/// Fractional bits of a time in steps. +pub const SC_QB : u32 = 12; +pub const SC_NONE : u32 = 0xFFFFFFFF; +/// Newton iterations of the analytic reference. +pub const SC_K_ANA : u32 = 12; + +; --- the registers of one stamp: positive and negative terms of J and G --------------------------- +pub const SC_RJP : u32 = 0; +pub const SC_RJN : u32 = 64; +pub const SC_RGP : u32 = 128; +pub const SC_RGN : u32 = 136; +pub const SC_REG_LEN : usize = 144; +/// The real stamp's output: J row-major in words 0..63, G in words 64..71. +pub const SC_FG : u32 = 64; +/// The step, in units, of the central differences that check the Jacobian: 12 uV. +pub const SC_FD_DELTA : f64 = 200.0; +/// Sign bit of a sign-magnitude u64. +pub const SC_SGN : u64 = 0x8000000000000000; +pub const SC_MAGM : u64 = 0x7FFFFFFFFFFFFFFF; + +; --- the hardware side: unsigned words only --------------------------------------------------------- + +/// The clock at step j: a ramp of SC_R steps from 0 to Vdd, then Vdd. +pub fn sc_clk_at(j: u32) -> u32 { + if (j >= SC_R) { return VDD; } + return j * SC_CLK_STEP; +} + +/// The voltage on terminal tm: an internal node from v, or the supply, ground, the clock or an input. +pub fn sc_term_v(tm: u32, v: []const u32, clk: u32, vdp: u32, vdn: u32) -> u32 { + if (tm < SC_NN) { return v[tm as usize]; } + if (tm == SC_T_VDD) { return VDD; } + if (tm == SC_T_GND) { return 0; } + if (tm == SC_T_CLK) { return clk; } + if (tm == SC_T_VDP) { return vdp; } + return vdn; +} + +/// The gate-drive argument of the table: an NMOS reads vg - vlo, a PMOS vhi - vg (its source is the higher channel end); 0 below. +pub fn sc_drive(kind: u32, vg: u32, vhi: u32, vlo: u32) -> u32 { + if (kind == 0) { + if (vg > vlo) { return vg - vlo; } + return 0; + } + if (vhi > vg) { return vhi - vg; } + return 0; +} + +/// Add the signed term (mag, neg) to the Jacobian entry (r, c): the negative register when neg. Rails and inputs have no row or column. +pub fn sc_add_j(reg: []u64, r: u32, c: u32, mag: u64, neg: bool) { + if (r >= SC_NN or c >= SC_NN) { return; } + if (neg) { + const k : usize = (SC_RJN + r * SC_NN + c) as usize; + reg[k] = reg[k] + mag; + } else { + const k : usize = (SC_RJP + r * SC_NN + c) as usize; + reg[k] = reg[k] + mag; + } +} + +/// The partial p of the channel current with respect to column c: +p into the row of the high end, -p into the row of the low end. +pub fn sc_add_pair(reg: []u64, hi: u32, lo: u32, c: u32, mag: u64, neg: bool) { + sc_add_j(reg, hi, c, mag, neg); + sc_add_j(reg, lo, c, mag, !neg); +} + +/// Stamp device e: its channel current m * T(x, y) leaves the high end and enters the low end, and its partials go into J. +pub fn sc_stamp_dev(t: []const u32, e: usize, v: []const u32, clk: u32, vdp: u32, vdn: u32, m: u32, reg: []u64) { + const kind : u32 = SC_DV_KIND[e]; + const d : u32 = SC_DV_D[e]; + const g : u32 = SC_DV_G[e]; + const s : u32 = SC_DV_S[e]; + const vd : u32 = sc_term_v(d, v, clk, vdp, vdn); + const vg : u32 = sc_term_v(g, v, clk, vdp, vdn); + const vs : u32 = sc_term_v(s, v, clk, vdp, vdn); + var hi : u32 = d; + var lo : u32 = s; + var vhi : u32 = vd; + var vlo : u32 = vs; + if (vd < vs) { + hi = s; + lo = d; + vhi = vs; + vlo = vd; + } + const x : u32 = sc_drive(kind, vg, vhi, vlo); + const y : u32 = vhi - vlo; + const cur : u64 = (m as u64) * (id_of(t, x, y) as u64); + const a : u64 = (m as u64) * (gm_cell(t, x, y) as u64); + const b : u64 = (m as u64) * (gds_cell(t, x, y) as u64); + if (hi < SC_NN) { + const k : usize = (SC_RGP + hi) as usize; + reg[k] = reg[k] + cur; + } + if (lo < SC_NN) { + const k : usize = (SC_RGN + lo) as usize; + reg[k] = reg[k] + cur; + } + if (kind == 0) { + sc_add_pair(reg, hi, lo, hi, b, false); + sc_add_pair(reg, hi, lo, g, a, false); + sc_add_pair(reg, hi, lo, lo, a + b, true); + } else { + sc_add_pair(reg, hi, lo, hi, a + b, false); + sc_add_pair(reg, hi, lo, g, a, true); + sc_add_pair(reg, hi, lo, lo, b, true); + } +} + +/// All registers of the Newton system at state v, for the step from v0: G = C (v - v0) + (currents out - currents in), J = dG/dv. +pub fn sc_stamp(t: []const u32, v: []const u32, v0: []const u32, clk: u32, vdp: u32, vdn: u32, mult: []const u32, reg: []u64) { + var i : usize = 0; + while (i < SC_REG_LEN) { + reg[i] = 0; + i = i + 1; + } + var n : u32 = 0; + while (n < SC_NN) { + const c : u64 = (SC_CAP[n as usize]) as u64; + reg[(SC_RGP + n) as usize] = c * (v[n as usize] as u64); + reg[(SC_RGN + n) as usize] = c * (v0[n as usize] as u64); + reg[(SC_RJP + n * SC_NN + n) as usize] = c * (CELL as u64); + n = n + 1; + } + var e : usize = 0; + while (e < (SC_NDEV as usize)) { + sc_stamp_dev(t, e, v, clk, vdp, vdn, mult[e], reg); + e = e + 1; + } +} + +/// The residual of node n in voltage units: round(|G| 2^20 / J_nn), the Newton move that is still owed there. +pub fn sc_resid_at(reg: []const u64, n: u32) -> u32 { + const gp : u64 = reg[(SC_RGP + n) as usize]; + const gn : u64 = reg[(SC_RGN + n) as usize]; + var g : u64 = 0; + if (gp >= gn) { g = gp - gn; } else { g = gn - gp; } + const jp : u64 = reg[(SC_RJP + n * SC_NN + n) as usize]; + const jn : u64 = reg[(SC_RJN + n * SC_NN + n) as usize]; + var jd : u64 = 1; + if (jp > jn) { + if (jp - jn > 1) { jd = jp - jn; } + } + return (((g << CELL_BITS) + (jd >> 1)) / jd) as u32; +} + +/// The largest residual over the eight nodes. +pub fn sc_resid(reg: []const u64) -> u32 { + var worst : u32 = 0; + var n : u32 = 0; + while (n < SC_NN) { + const r : u32 = sc_resid_at(reg, n); + if (r > worst) { worst = r; } + n = n + 1; + } + return worst; +} + +; --- sign-magnitude words: bit 63 is the direction, bits 0..62 the magnitude --------------------------- + +pub fn sc_pack(m: u64, neg: bool) -> u64 { + if (m == 0) { return 0; } + if (neg) { return m | SC_SGN; } + return m; +} + +pub fn sc_neg(w: u64) -> bool { + return (w & SC_SGN) != 0; +} + +pub fn sc_mag(w: u64) -> u64 { + return w & SC_MAGM; +} + +/// a + b. +pub fn sc_sadd(a: u64, b: u64) -> u64 { + const ma : u64 = sc_mag(a); + const mb : u64 = sc_mag(b); + if (sc_neg(a) == sc_neg(b)) { return sc_pack(ma + mb, sc_neg(a)); } + if (ma >= mb) { return sc_pack(ma - mb, sc_neg(a)); } + return sc_pack(mb - ma, sc_neg(b)); +} + +/// -w. +pub fn sc_flip(w: u64) -> u64 { + return sc_pack(sc_mag(w), !sc_neg(w)); +} + +/// p - n for two unsigned registers, as a word. +pub fn sc_sub_reg(p: u64, n: u64) -> u64 { + if (p >= n) { return sc_pack(p - n, false); } + return sc_pack(n - p, true); +} + +/// The 8 x 9 augmented matrix [J | -G] from the registers, row-major, as words. +pub fn sc_load(reg: []const u64, a: []u64) { + var i : u32 = 0; + while (i < SC_NN) { + var j : u32 = 0; + while (j < SC_NN) { + a[(i * 9 + j) as usize] = sc_sub_reg(reg[(SC_RJP + i * SC_NN + j) as usize], reg[(SC_RJN + i * SC_NN + j) as usize]); + j = j + 1; + } + a[(i * 9 + 8) as usize] = sc_sub_reg(reg[(SC_RGN + i) as usize], reg[(SC_RGP + i) as usize]); + i = i + 1; + } +} + +/// Forward elimination with partial pivoting, in place. The multiplier is round(|a_rc| 2^20 / |pivot|) <= 2^20. +pub fn sc_eliminate(a: []u64) { + var c : u32 = 0; + while (c < SC_NN) { + var p : u32 = c; + var r : u32 = c + 1; + while (r < SC_NN) { + if (sc_mag(a[(r * 9 + c) as usize]) > sc_mag(a[(p * 9 + c) as usize])) { p = r; } + r = r + 1; + } + if (p != c) { + var k : u32 = 0; + while (k < 9) { + const tmp : u64 = a[(c * 9 + k) as usize]; + a[(c * 9 + k) as usize] = a[(p * 9 + k) as usize]; + a[(p * 9 + k) as usize] = tmp; + k = k + 1; + } + } + const piv : u64 = a[(c * 9 + c) as usize]; + var pm : u64 = sc_mag(piv); + if (pm == 0) { pm = 1; } + r = c + 1; + while (r < SC_NN) { + const arc : u64 = a[(r * 9 + c) as usize]; + const am : u64 = sc_mag(arc); + if (am != 0) { + const f : u64 = ((am << CELL_BITS) + (pm >> 1)) / pm; + const fneg : bool = sc_neg(arc) != sc_neg(piv); + var k : u32 = c; + while (k < 9) { + const ck : u64 = a[(c * 9 + k) as usize]; + const pr : u64 = sc_pack((f * sc_mag(ck) + ((CELL >> 1) as u64)) >> CELL_BITS, fneg != sc_neg(ck)); + a[(r * 9 + k) as usize] = sc_sadd(a[(r * 9 + k) as usize], sc_flip(pr)); + k = k + 1; + } + } + r = r + 1; + } + c = c + 1; + } +} + +/// Back substitution on the eliminated matrix: x_i = round((b_i 2^20 - sum_{j>i} A_ij x_j) / A_ii), as words. +pub fn sc_backsub(a: []const u64, x: []u64) { + var i : u32 = SC_NN; + while (i > 0) { + i = i - 1; + var pos : u64 = 0; + var neg : u64 = 0; + const b : u64 = a[(i * 9 + 8) as usize]; + if (sc_neg(b)) { neg = neg + (sc_mag(b) << CELL_BITS); } else { pos = pos + (sc_mag(b) << CELL_BITS); } + var j : u32 = i + 1; + while (j < SC_NN) { + const aij : u64 = a[(i * 9 + j) as usize]; + const xj : u64 = x[j as usize]; + const pr : u64 = sc_mag(aij) * sc_mag(xj); + if (sc_neg(aij) != sc_neg(xj)) { pos = pos + pr; } else { neg = neg + pr; } + j = j + 1; + } + const num : u64 = sc_sub_reg(pos, neg); + const aii : u64 = a[(i * 9 + i) as usize]; + var den : u64 = sc_mag(aii); + if (den == 0) { den = 1; } + const q : u64 = (sc_mag(num) + (den >> 1)) / den; + x[i as usize] = sc_pack(q, sc_neg(num) != sc_neg(aii)); + } +} + +/// v + x for the eight nodes, held to [0, Vdd]. +pub fn sc_apply(v: []u32, x: []const u64) { + var i : u32 = 0; + while (i < SC_NN) { + const m : u64 = sc_mag(x[i as usize]); + const cur : u64 = (v[i as usize]) as u64; + if (sc_neg(x[i as usize])) { + if (m >= cur) { v[i as usize] = 0; } else { v[i as usize] = (cur - m) as u32; } + } else { + if (m >= (VDD as u64) - cur) { v[i as usize] = VDD; } else { v[i as usize] = (cur + m) as u32; } + } + i = i + 1; + } +} + +/// One Newton iteration on the eight nodes at once: stamp, load, eliminate, back-substitute, apply. +pub fn sc_iter(t: []const u32, v: []u32, v0: []const u32, clk: u32, vdp: u32, vdn: u32, mult: []const u32) { + var reg : [144]u64 = [0] ** 144; + var a : [72]u64 = [0] ** 72; + var x : [8]u64 = [0] ** 8; + sc_stamp(t, v, v0, clk, vdp, vdn, mult, ®); + sc_load(®, &a); + sc_eliminate(&a); + sc_backsub(&a, &x); + sc_apply(v, &x); +} + +/// One backward Euler step to step j from the state in v (updated in place): k Newton iterations, then the residual. +pub fn sc_step(t: []const u32, v: []u32, j: u32, vdp: u32, vdn: u32, k: u32, mult: []const u32) -> u32 { + var v0 : [8]u32 = [0] ** 8; + var n : u32 = 0; + while (n < SC_NN) { + v0[n as usize] = v[n as usize]; + n = n + 1; + } + const clk : u32 = sc_clk_at(j); + var i : u32 = 0; + while (i < k) { + sc_iter(t, v, &v0, clk, vdp, vdn, mult); + i = i + 1; + } + var reg : [144]u64 = [0] ** 144; + sc_stamp(t, v, &v0, clk, vdp, vdn, mult, ®); + return sc_resid(®); +} + +/// The step the board runs: SC_K Newton iterations. The trip count is a constant, so the loop unrolls (`t27b verilog` takes a +/// `while` only with a static bound); sc_step takes k as an argument so that the tests can run K = 1 and K = 2. +pub fn sc_board_step(t: []const u32, v: []u32, j: u32, vdp: u32, vdn: u32, mult: []const u32) -> u32 { + var v0 : [8]u32 = [0] ** 8; + var n : u32 = 0; + while (n < SC_NN) { + v0[n as usize] = v[n as usize]; + n = n + 1; + } + const clk : u32 = sc_clk_at(j); + var i : u32 = 0; + while (i < SC_K) { + sc_iter(t, v, &v0, clk, vdp, vdn, mult); + i = i + 1; + } + var reg : [144]u64 = [0] ** 144; + sc_stamp(t, v, &v0, clk, vdp, vdn, mult, ®); + return sc_resid(®); +} + +/// The state at rest with the clock low: Vdd on ap, an, net2, net3 and n_clkc, 0 on dp and dn, net1 one threshold below the higher input. +pub fn sc_rest(vdp: u32, vdn: u32, v: []u32) { + var top : u32 = vdp; + if (vdn > vdp) { top = vdn; } + v[0] = VDD; + v[1] = VDD; + v[2] = top - SC_VTH_U; + v[3] = VDD; + v[4] = VDD; + v[5] = 0; + v[6] = 0; + v[7] = VDD; +} + +/// The run: SC_N steps from rest, inputs VCM + h and VCM - h (exchanged when swap). out gets the 8 nodes of every state, step 0 first; +/// res the residual of every step. Returns the order-sensitive fold of every sample. +pub fn sc_run(t: []const u32, h: u32, swap: bool, k: u32, mult: []const u32, out: []u32, res: []u32) -> u64 { + var vdp : u32 = SC_VCM + h; + var vdn : u32 = SC_VCM - h; + if (swap) { + vdp = SC_VCM - h; + vdn = SC_VCM + h; + } + var v : [8]u32 = [0] ** 8; + sc_rest(vdp, vdn, &v); + var s : u64 = 0; + var n : u32 = 0; + while (n < SC_NN) { + out[n as usize] = v[n as usize]; + n = n + 1; + } + var j : u32 = 1; + while (j <= SC_N) { + res[(j - 1) as usize] = sc_step(t, &v, j, vdp, vdn, k, mult); + n = 0; + while (n < SC_NN) { + out[(j * SC_NN + n) as usize] = v[n as usize]; + s = fold(s, v[n as usize]); + n = n + 1; + } + j = j + 1; + } + return s; +} + +/// The series of one node over the run (sel < 8), or of |dp - dn| (sel = 8): SC_N + 1 samples. +pub fn sc_pick(out: []const u32, sel: u32, s: []u32) { + var j : u32 = 0; + while (j <= SC_N) { + if (sel < SC_NN) { + s[j as usize] = out[(j * SC_NN + sel) as usize]; + } else { + s[j as usize] = mag_diff(out[(j * SC_NN + SC_DP) as usize], out[(j * SC_NN + SC_DN) as usize]); + } + j = j + 1; + } +} + +/// The time, in steps with SC_QB fractional bits, of the first rising crossing of lvl by the series, linear between samples; +/// SC_NONE when there is none. +pub fn sc_cross(s: []const u32) -> u32 { + const lvl : u32 = VDD >> 1; + var j : u32 = 1; + while (j <= SC_N) { + const y0 : u32 = s[(j - 1) as usize]; + const y1 : u32 = s[j as usize]; + if (y0 < lvl and y1 >= lvl) { + const frac : u64 = (((lvl - y0) as u64) << SC_QB) / ((y1 - y0) as u64); + return ((j - 1) << SC_QB) + (frac as u32); + } + j = j + 1; + } + return SC_NONE; +} + +/// The decision time of a run: the first output through Vdd / 2, in steps (SC_QB fractional bits) from step 0. +pub fn sc_t_first(out: []const u32) -> u32 { + var s : [141]u32 = [0] ** 141; + sc_pick(out, SC_DP, &s); + const tp : u32 = sc_cross(&s); + sc_pick(out, SC_DN, &s); + const tn : u32 = sc_cross(&s); + if (tp < tn) { return tp; } + return tn; +} + +/// The resolution time of a run: |dp - dn| through Vdd / 2. +pub fn sc_t_res(out: []const u32) -> u32 { + var s : [141]u32 = [0] ** 141; + sc_pick(out, 8, &s); + return sc_cross(&s); +} + +/// A time measured from step 0, moved to the clock's 50 % point (step SC_R / 2), still in steps with SC_QB fractional bits. +pub fn sc_after_clock(tq: u32) -> u32 { + if (tq == SC_NONE) { return SC_NONE; } + return tq - ((SC_R >> 1) << SC_QB); +} + +/// Is dp above dn in the last state: the decision "VDP > VDN". +pub fn sc_sign(out: []const u32) -> bool { + return out[(SC_N * SC_NN + SC_DP) as usize] > out[(SC_N * SC_NN + SC_DN) as usize]; +} + +; --- the f64 side: references that only the tests reach ------------------------------------------------------ + +/// Fill the table from the law with threshold vth: entry (ix, iy) = round(Id(ix * 62.5 mV, iy * 62.5 mV) * dt / C_unit) in units. +pub fn sc_build_table(t: []u32, vth: f64) { + var ix : u32 = 0; + while (ix < NODES) { + var iy : u32 = 0; + while (iy < NODES) { + const i : f64 = law_amps((ix as f64) / 16.0, (iy as f64) / 16.0, vth) * (SC_BETA / BETA); + t[(ix * NODES + iy) as usize] = (i * SC_STEP_V_PER_A * SC_UPV + 0.5) as u32; + iy = iy + 1; + } + ix = ix + 1; + } +} + +/// The unit device's law in table units: T(x, y) for x, y in voltage units. +pub fn sc_law_t(x: f64, y: f64) -> f64 { + return law_amps(x / SC_UPV, y / SC_UPV, SC_VTH_V) * (SC_BETA / BETA) * SC_STEP_V_PER_A * SC_UPV; +} + +/// dT/dx of the law, a pure number. +pub fn sc_law_gm_t(x: f64, y: f64) -> f64 { + return law_gm(x / SC_UPV, y / SC_UPV, SC_VTH_V) * (SC_BETA / BETA) * SC_STEP_V_PER_A; +} + +/// dT/dy of the law, a pure number: beta (ov - y) in the triode region, y = 0 included (the limit from above, which M3a's law_gds +/// leaves at 0), and 0 in saturation or below threshold. +pub fn sc_law_gds_t(x: f64, y: f64) -> f64 { + const ov : f64 = x / SC_UPV - SC_VTH_V; + const yv : f64 = y / SC_UPV; + if (ov <= 0.0 or yv >= ov) { return 0.0; } + return SC_BETA * (ov - yv) * SC_STEP_V_PER_A; +} + +/// gm_cell in real arithmetic: how much the interpolant moves across one cell of the gate axis. +pub fn sc_gmf(t: []const u32, vg: f64, vd: f64) -> f64 { + const cx : u32 = cellf(vg); + const cy : u32 = cellf(vd); + const fy : f64 = vd - (cy as f64) * CELL_F; + const d0 : f64 = fmax64(nodef(t, cx + 1, cy) - nodef(t, cx, cy), 0.0); + const d1 : f64 = fmax64(nodef(t, cx + 1, cy + 1) - nodef(t, cx, cy + 1), 0.0); + return (d0 * (CELL_F - fy) + d1 * fy) / CELL_F; +} + +pub fn sc_term_vf(tm: u32, v: []const f64, clk: f64, vdp: f64, vdn: f64) -> f64 { + if (tm < SC_NN) { return v[tm as usize]; } + if (tm == SC_T_VDD) { return VDD_F; } + if (tm == SC_T_GND) { return 0.0; } + if (tm == SC_T_CLK) { return clk; } + if (tm == SC_T_VDP) { return vdp; } + return vdn; +} + +pub fn sc_addf(jg: []f64, r: u32, c: u32, val: f64) { + if (r >= SC_NN or c >= SC_NN) { return; } + jg[(r * SC_NN + c) as usize] = jg[(r * SC_NN + c) as usize] + val; +} + +/// +p into the row of the high end, -p into the row of the low end. +pub fn sc_addf_pair(jg: []f64, hi: u32, lo: u32, c: u32, val: f64) { + sc_addf(jg, hi, c, val); + sc_addf(jg, lo, c, 0.0 - val); +} + +/// The Newton system in real arithmetic. mode 0 reads the table (the real-valued twin of sc_stamp), mode 1 the analytic law. +/// drv is the clock, VDP and VDN in units, in that order. out gets J (row-major, a pure number per entry) in words 0..63 and G +/// in words 64..71. +pub fn sc_stamp_f(mode: u32, t: []const u32, v: []const f64, v0: []const f64, drv: []const f64, mult: []const u32, jg: []f64) { + var i : u32 = 0; + while (i < SC_NN) { + const c : f64 = (SC_CAP[i as usize]) as f64; + jg[(SC_FG + i) as usize] = c * (v[i as usize] - v0[i as usize]); + var j : u32 = 0; + while (j < SC_NN) { + jg[(i * SC_NN + j) as usize] = 0.0; + j = j + 1; + } + jg[(i * SC_NN + i) as usize] = c; + i = i + 1; + } + var e : usize = 0; + while (e < (SC_NDEV as usize)) { + const kind : u32 = SC_DV_KIND[e]; + const d : u32 = SC_DV_D[e]; + const g : u32 = SC_DV_G[e]; + const s : u32 = SC_DV_S[e]; + const m : f64 = (mult[e]) as f64; + const vd : f64 = sc_term_vf(d, v, drv[0], drv[1], drv[2]); + const vg : f64 = sc_term_vf(g, v, drv[0], drv[1], drv[2]); + const vs : f64 = sc_term_vf(s, v, drv[0], drv[1], drv[2]); + var hi : u32 = d; + var lo : u32 = s; + var vhi : f64 = vd; + var vlo : f64 = vs; + if (vd < vs) { + hi = s; + lo = d; + vhi = vs; + vlo = vd; + } + var x : f64 = 0.0; + if (kind == 0) { + if (vg > vlo) { x = vg - vlo; } + } else { + if (vhi > vg) { x = vhi - vg; } + } + if (x > VDD_F) { x = VDD_F; } + var y : f64 = vhi - vlo; + if (y > VDD_F) { y = VDD_F; } + var cur : f64 = 0.0; + var gmv : f64 = 0.0; + var gdv : f64 = 0.0; + if (mode == 0) { + cur = m * idf(t, x, y); + gmv = m * sc_gmf(t, x, y) / CELL_F; + gdv = m * gdsf(t, x, y) / CELL_F; + } else { + cur = m * sc_law_t(x, y); + gmv = m * sc_law_gm_t(x, y); + gdv = m * sc_law_gds_t(x, y); + } + if (hi < SC_NN) { jg[(SC_FG + hi) as usize] = jg[(SC_FG + hi) as usize] + cur; } + if (lo < SC_NN) { jg[(SC_FG + lo) as usize] = jg[(SC_FG + lo) as usize] - cur; } + if (kind == 0) { + sc_addf_pair(jg, hi, lo, hi, gdv); + sc_addf_pair(jg, hi, lo, g, gmv); + sc_addf_pair(jg, hi, lo, lo, 0.0 - (gmv + gdv)); + } else { + sc_addf_pair(jg, hi, lo, hi, gmv + gdv); + sc_addf_pair(jg, hi, lo, g, 0.0 - gmv); + sc_addf_pair(jg, hi, lo, lo, 0.0 - gdv); + } + e = e + 1; + } +} + +/// Solve the 8 x 8 system in the 8 x 9 augmented matrix a (row-major) by elimination with partial pivoting; the solution goes to x. +pub fn sc_solve_f(a: []f64, x: []f64) { + var c : u32 = 0; + while (c < SC_NN) { + var p : u32 = c; + var r : u32 = c + 1; + while (r < SC_NN) { + if (fabs64(a[(r * 9 + c) as usize]) > fabs64(a[(p * 9 + c) as usize])) { p = r; } + r = r + 1; + } + if (p != c) { + var k : u32 = 0; + while (k < 9) { + const tmp : f64 = a[(c * 9 + k) as usize]; + a[(c * 9 + k) as usize] = a[(p * 9 + k) as usize]; + a[(p * 9 + k) as usize] = tmp; + k = k + 1; + } + } + r = c + 1; + while (r < SC_NN) { + const f : f64 = a[(r * 9 + c) as usize] / a[(c * 9 + c) as usize]; + var k : u32 = c; + while (k < 9) { + a[(r * 9 + k) as usize] = a[(r * 9 + k) as usize] - f * a[(c * 9 + k) as usize]; + k = k + 1; + } + r = r + 1; + } + c = c + 1; + } + var i : u32 = SC_NN; + while (i > 0) { + i = i - 1; + var acc : f64 = a[(i * 9 + 8) as usize]; + var j : u32 = i + 1; + while (j < SC_NN) { + acc = acc - a[(i * 9 + j) as usize] * x[j as usize]; + j = j + 1; + } + x[i as usize] = acc / a[(i * 9 + i) as usize]; + } +} + +/// One real Newton iteration on v: the same system, the same clamp to [0, Vdd], no rounding. +pub fn sc_iter_f(mode: u32, t: []const u32, v: []f64, v0: []const f64, drv: []const f64, mult: []const u32) { + var jg : [72]f64 = [0.0] ** 72; + var a : [72]f64 = [0.0] ** 72; + var x : [8]f64 = [0.0] ** 8; + sc_stamp_f(mode, t, v, v0, drv, mult, &jg); + var i : u32 = 0; + while (i < SC_NN) { + var j : u32 = 0; + while (j < SC_NN) { + a[(i * 9 + j) as usize] = jg[(i * SC_NN + j) as usize]; + j = j + 1; + } + a[(i * 9 + 8) as usize] = 0.0 - jg[(SC_FG + i) as usize]; + i = i + 1; + } + sc_solve_f(&a, &x); + i = 0; + while (i < SC_NN) { + var z : f64 = v[i as usize] + x[i as usize]; + if (z < 0.0) { z = 0.0; } + if (z > VDD_F) { z = VDD_F; } + v[i as usize] = z; + i = i + 1; + } +} + +/// One real backward Euler step to step j: k iterations; returns the residual in voltage units, max_n |G_n| / J_nn. +pub fn sc_step_f(mode: u32, t: []const u32, v: []f64, j: u32, vdp: f64, vdn: f64, k: u32, mult: []const u32) -> f64 { + var v0 : [8]f64 = [0.0] ** 8; + var n : u32 = 0; + while (n < SC_NN) { + v0[n as usize] = v[n as usize]; + n = n + 1; + } + var drv : [3]f64 = [0.0] ** 3; + drv[0] = (sc_clk_at(j)) as f64; + drv[1] = vdp; + drv[2] = vdn; + var i : u32 = 0; + while (i < k) { + sc_iter_f(mode, t, v, &v0, &drv, mult); + i = i + 1; + } + var jg : [72]f64 = [0.0] ** 72; + sc_stamp_f(mode, t, v, &v0, &drv, mult, &jg); + var worst : f64 = 0.0; + n = 0; + while (n < SC_NN) { + worst = fmax64(worst, fabs64(jg[(SC_FG + n) as usize]) / jg[(n * SC_NN + n) as usize]); + n = n + 1; + } + return worst; +} + +/// The real run, mode 0 (table) or 1 (law), k iterations a step: out gets the 8 nodes of every state in units, res the residual. +pub fn sc_run_f(mode: u32, t: []const u32, h: u32, swap: bool, k: u32, mult: []const u32, out: []f64, res: []f64) { + var vdp : u32 = SC_VCM + h; + var vdn : u32 = SC_VCM - h; + if (swap) { + vdp = SC_VCM - h; + vdn = SC_VCM + h; + } + var vi : [8]u32 = [0] ** 8; + sc_rest(vdp, vdn, &vi); + var v : [8]f64 = [0.0] ** 8; + var n : u32 = 0; + while (n < SC_NN) { + v[n as usize] = (vi[n as usize]) as f64; + out[n as usize] = v[n as usize]; + n = n + 1; + } + var j : u32 = 1; + while (j <= SC_N) { + res[(j - 1) as usize] = sc_step_f(mode, t, &v, j, (vdp) as f64, (vdn) as f64, k, mult); + n = 0; + while (n < SC_NN) { + out[(j * SC_NN + n) as usize] = v[n as usize]; + n = n + 1; + } + j = j + 1; + } +} + +/// The real-valued series of one node over a real run (sel < 8), or of |dp - dn| (sel = 8). +pub fn sc_pick_f(out: []const f64, sel: u32, s: []f64) { + var j : u32 = 0; + while (j <= SC_N) { + if (sel < SC_NN) { + s[j as usize] = out[(j * SC_NN + sel) as usize]; + } else { + s[j as usize] = fabs64(out[(j * SC_NN + SC_DP) as usize] - out[(j * SC_NN + SC_DN) as usize]); + } + j = j + 1; + } +} + +/// The first rising crossing of Vdd / 2 by a real series, in steps, linear between samples; -1 when there is none. +pub fn sc_cross_f(s: []const f64) -> f64 { + const lvl : f64 = VDD_F / 2.0; + var j : u32 = 1; + while (j <= SC_N) { + const y0 : f64 = s[(j - 1) as usize]; + const y1 : f64 = s[j as usize]; + if (y0 < lvl and y1 >= lvl) { return ((j - 1) as f64) + (lvl - y0) / (y1 - y0); } + j = j + 1; + } + return 0.0 - 1.0; +} + +/// Decision time of a real run: the first output through Vdd / 2, in steps after the clock's 50 % point. +pub fn sc_t_first_f(out: []const f64) -> f64 { + var s : [141]f64 = [0.0] ** 141; + sc_pick_f(out, SC_DP, &s); + const tp : f64 = sc_cross_f(&s); + sc_pick_f(out, SC_DN, &s); + const tn : f64 = sc_cross_f(&s); + var r : f64 = tp; + if (tp < 0.0 or (tn >= 0.0 and tn < tp)) { r = tn; } + return r - ((SC_R >> 1) as f64); +} + +/// Resolution time of a real run: |dp - dn| through Vdd / 2, in steps after the clock's 50 % point. +pub fn sc_t_res_f(out: []const f64) -> f64 { + var s : [141]f64 = [0.0] ** 141; + sc_pick_f(out, 8, &s); + return sc_cross_f(&s) - ((SC_R >> 1) as f64); +} + +; --- ngspice, MOSFET LEVEL=1, the same square law: what the decks in specs/analog/decks print ------------------- +/// The steps at which the decks' .print tables are recorded below. +pub const SC_CHK : [8]u32 = [20, 40, 52, 60, 68, 76, 90, 130]; +/// ap an net1 net2 net3 dp dn nclk in volts at the eight steps, deck sar_cmp_p2.cir (VDP - VDN = +2 mV), 7 significant digits. +pub const SC_NG_P2 : [64]f64 = [ + 1.373046736e+00, 1.370477235e+00, 1.441352199e-01, 1.500000000e+00, 1.500000000e+00, 6.736590897e-11, 6.736590897e-11, 1.310805003e+00, + 9.507006870e-01, 9.417399819e-01, 1.322801241e-01, 1.486894963e+00, 1.485730242e+00, 3.084394696e-02, 2.777715565e-02, 8.760007962e-01, + 6.945813029e-01, 6.817653523e-01, 1.322691549e-01, 1.423827061e+00, 1.420885309e+00, 2.741662624e-01, 2.571980653e-01, 6.447477966e-01, + 5.238323059e-01, 5.084461711e-01, 1.322689964e-01, 1.317786861e+00, 1.317543276e+00, 6.212632187e-01, 5.776557689e-01, 5.135600229e-01, + 3.536222315e-01, 3.368194066e-01, 1.311135087e-01, 1.066957943e+00, 1.166696564e+00, 9.197645140e-01, 7.471906692e-01, 4.027723340e-01, + 2.110072685e-01, 2.011838989e-01, 1.073065865e-01, 7.966138827e-01, 1.289248736e+00, 1.169177576e+00, 5.296883274e-01, 3.118478457e-01, + 7.774419066e-02, 7.630726171e-02, 4.863912853e-02, 5.704152208e-01, 1.484195818e+00, 1.471793117e+00, 8.285936318e-02, 1.946202492e-01, + 4.038136293e-03, 4.027760198e-03, 2.722788440e-03, 4.070910078e-01, 1.499997711e+00, 1.499995923e+00, 3.815979111e-03, 4.615473271e-02 +]; +/// The same, deck sar_cmp_p05.cir (+0.5 mV). +pub const SC_NG_P05 : [64]f64 = [ + 1.371813195e+00, 1.371169184e+00, 1.440960714e-01, 1.500000000e+00, 1.500000000e+00, 6.736590897e-11, 6.736590897e-11, 1.310805003e+00, + 9.470585626e-01, 9.448168128e-01, 1.322791900e-01, 1.486424242e+00, 1.486132487e+00, 2.978133758e-02, 2.901243643e-02, 8.760007962e-01, + 6.894950650e-01, 6.862895581e-01, 1.322682567e-01, 1.422619525e+00, 1.421883678e+00, 2.681461432e-01, 2.638977182e-01, 6.447477966e-01, + 5.177833270e-01, 5.139353102e-01, 1.322680986e-01, 1.317469010e+00, 1.317412421e+00, 6.054538344e-01, 5.945322536e-01, 5.135600229e-01, + 3.469996517e-01, 3.428000900e-01, 1.311497970e-01, 1.097060536e+00, 1.123555727e+00, 8.580290551e-01, 8.142634308e-01, 4.027723340e-01, + 2.070600595e-01, 2.046084482e-01, 1.072664318e-01, 9.148061367e-01, 1.085241154e+00, 9.495221142e-01, 7.704152909e-01, 3.118478457e-01, + 7.710041741e-02, 7.674220608e-02, 4.858147968e-02, 5.932671856e-01, 1.424987297e+00, 1.374568049e+00, 2.466656408e-01, 1.946202492e-01, + 4.028570114e-03, 4.025980268e-03, 2.718979635e-03, 4.077966063e-01, 1.499988215e+00, 1.499979013e+00, 4.011192655e-03, 4.615473271e-02 +]; +/// Decision time and resolution time of the four decks, in steps (0.1 ps) after the clock through 0.75 V: (tfirst - tck) and +/// (tres - tck) of the .meas lines. The -2 mV and -0.5 mV decks print the same times as +2 mV and +0.5 mV, to the 6 digits. +pub const SC_NG_FIRST_2 : f64 = 54.9629; +pub const SC_NG_RES_2 : f64 = 69.3083; +pub const SC_NG_FIRST_05 : f64 = 55.6936; +pub const SC_NG_RES_05 : f64 = 76.9385; +/// The loser's peak (dnmax of the +2 mV and +0.5 mV decks, dpmax of the others), in volts. +pub const SC_NG_PEAK_2 : f64 = 0.747264; +pub const SC_NG_PEAK_05 : f64 = 0.822056; + +/// The largest |model - ngspice| over the eight nodes at the eight recorded steps, in volts; the model's states are real, in units. +pub fn sc_dev_ng_f(out: []const f64, ng: []const f64) -> f64 { + var worst : f64 = 0.0; + var c : u32 = 0; + while (c < 8) { + var n : u32 = 0; + while (n < SC_NN) { + const model : f64 = out[(SC_CHK[c as usize] * SC_NN + n) as usize] / SC_UPV; + worst = fmax64(worst, fabs64(model - ng[(c * SC_NN + n) as usize])); + n = n + 1; + } + c = c + 1; + } + return worst; +} + +/// The same for the fixed-point run's unsigned states. +pub fn sc_dev_ng_u(out: []const u32, ng: []const f64) -> f64 { + var worst : f64 = 0.0; + var c : u32 = 0; + while (c < 8) { + var n : u32 = 0; + while (n < SC_NN) { + const model : f64 = ((out[(SC_CHK[c as usize] * SC_NN + n) as usize]) as f64) / SC_UPV; + worst = fmax64(worst, fabs64(model - ng[(c * SC_NN + n) as usize])); + n = n + 1; + } + c = c + 1; + } + return worst; +} + +/// The largest |a - b| between a fixed-point run and a real run over every node and step, in units. +pub fn sc_dev_fix_f(a: []const u32, b: []const f64) -> f64 { + var worst : f64 = 0.0; + var i : usize = 0; + while (i < SC_LEN) { + worst = fmax64(worst, fabs64((a[i] as f64) - b[i])); + i = i + 1; + } + return worst; +} + +/// The largest |a - b| between two real runs, in units. +pub fn sc_dev_f_f(a: []const f64, b: []const f64) -> f64 { + var worst : f64 = 0.0; + var i : usize = 0; + while (i < SC_LEN) { + worst = fmax64(worst, fabs64(a[i] - b[i])); + i = i + 1; + } + return worst; +} + +/// The largest |entry - law| over the 625 nodes, in table units, against the law with the true threshold. +pub fn sc_node_err(t: []const u32) -> f64 { + var worst : f64 = 0.0; + var ix : u32 = 0; + while (ix < NODES) { + var iy : u32 = 0; + while (iy < NODES) { + const want : f64 = sc_law_t(((ix as f64) / 16.0) * SC_UPV, ((iy as f64) / 16.0) * SC_UPV); + worst = fmax64(worst, fabs64((t[(ix * NODES + iy) as usize] as f64) - want)); + iy = iy + 1; + } + ix = ix + 1; + } + return worst; +} + +/// A copy of the device multipliers with the input pair's m set to (a, b): XM3 gets a, XM4 gets b. +pub fn sc_mult_pair(a: u32, b: u32, mult: []u32) { + var e : usize = 0; + while (e < (SC_NDEV as usize)) { + mult[e] = SC_DV_M[e]; + e = e + 1; + } + mult[2] = a; + mult[3] = b; +} + +/// The node that takes the place of node n when the two halves of the comparator are exchanged. +pub fn sc_mirror_node(n: u32) -> u32 { + if (n == SC_AP) { return SC_AN; } + if (n == SC_AN) { return SC_AP; } + if (n == SC_NET2) { return SC_NET3; } + if (n == SC_NET3) { return SC_NET2; } + if (n == SC_DP) { return SC_DN; } + if (n == SC_DN) { return SC_DP; } + return n; +} + +/// The largest |a(step, n) - b(step, mirror(n))| over a pair of fixed-point runs, in units: 0 for a circuit that is symmetric. +pub fn sc_mirror_dev_u(a: []const u32, b: []const u32) -> u32 { + var worst : u32 = 0; + var j : u32 = 0; + while (j <= SC_N) { + var n : u32 = 0; + while (n < SC_NN) { + const d : u32 = mag_diff(a[(j * SC_NN + n) as usize], b[(j * SC_NN + sc_mirror_node(n)) as usize]); + if (d > worst) { worst = d; } + n = n + 1; + } + j = j + 1; + } + return worst; +} + +/// The same for two real runs. +pub fn sc_mirror_dev_f(a: []const f64, b: []const f64) -> f64 { + var worst : f64 = 0.0; + var j : u32 = 0; + while (j <= SC_N) { + var n : u32 = 0; + while (n < SC_NN) { + worst = fmax64(worst, fabs64(a[(j * SC_NN + n) as usize] - b[(j * SC_NN + sc_mirror_node(n)) as usize])); + n = n + 1; + } + j = j + 1; + } + return worst; +} + +/// The state at step j of a real run as eight reals. +pub fn sc_state_f(out: []const f64, j: u32, v: []f64) { + var n : u32 = 0; + while (n < SC_NN) { + v[n as usize] = out[(j * SC_NN + n) as usize]; + n = n + 1; + } +} + +/// The state at step j of a fixed-point run as eight words. +pub fn sc_state_u(out: []const u32, j: u32, v: []u32) { + const base : u32 = j * SC_NN; + var n : u32 = 0; + while (n < SC_NN) { + v[n as usize] = out[(base + n) as usize]; + n = n + 1; + } +} + +/// A sign-magnitude word as a real number. +pub fn sc_word_f(w: u64) -> f64 { + const m : f64 = (sc_mag(w)) as f64; + if (sc_neg(w)) { return 0.0 - m; } + return m; +} + +/// The largest |dG_i / dv_j (central difference, +-SC_FD_DELTA units) - J_ij| of the real system at v, for the step from v0 with +/// the clock, VDP and VDN in drv. With transposed set it compares against J_ji instead, which is wrong on purpose. +pub fn sc_fd_err(mode: u32, t: []const u32, v: []const f64, v0: []const f64, drv: []const f64, mult: []const u32, transposed: bool) -> f64 { + var j0 : [72]f64 = [0.0] ** 72; + var jp : [72]f64 = [0.0] ** 72; + var jq : [72]f64 = [0.0] ** 72; + var vp : [8]f64 = [0.0] ** 8; + var vq : [8]f64 = [0.0] ** 8; + sc_stamp_f(mode, t, v, v0, drv, mult, &j0); + var worst : f64 = 0.0; + var j : u32 = 0; + while (j < SC_NN) { + var n : u32 = 0; + while (n < SC_NN) { + vp[n as usize] = v[n as usize]; + vq[n as usize] = v[n as usize]; + n = n + 1; + } + vp[j as usize] = v[j as usize] + SC_FD_DELTA; + vq[j as usize] = v[j as usize] - SC_FD_DELTA; + sc_stamp_f(mode, t, &vp, v0, drv, mult, &jp); + sc_stamp_f(mode, t, &vq, v0, drv, mult, &jq); + var i : u32 = 0; + while (i < SC_NN) { + const fd : f64 = (jp[(SC_FG + i) as usize] - jq[(SC_FG + i) as usize]) / (2.0 * SC_FD_DELTA); + var an : f64 = j0[(i * SC_NN + j) as usize]; + if (transposed) { an = j0[(j * SC_NN + i) as usize]; } + worst = fmax64(worst, fabs64(fd - an)); + i = i + 1; + } + j = j + 1; + } + return worst; +} + +/// The largest |G_fixed - G_real| over the eight nodes at the fixed-point state v (a table system, both stamps on the same words), units. +pub fn sc_gap_g(t: []const u32, v: []const u32, v0: []const u32, clk: u32, vdp: u32, vdn: u32, mult: []const u32) -> f64 { + var reg : [144]u64 = [0] ** 144; + sc_stamp(t, v, v0, clk, vdp, vdn, mult, ®); + var vf : [8]f64 = [0.0] ** 8; + var v0f : [8]f64 = [0.0] ** 8; + var n : u32 = 0; + while (n < SC_NN) { + vf[n as usize] = (v[n as usize]) as f64; + v0f[n as usize] = (v0[n as usize]) as f64; + n = n + 1; + } + var jg : [72]f64 = [0.0] ** 72; + var drv : [3]f64 = [0.0] ** 3; + drv[0] = (clk) as f64; + drv[1] = (vdp) as f64; + drv[2] = (vdn) as f64; + sc_stamp_f(0, t, &vf, &v0f, &drv, mult, &jg); + var worst : f64 = 0.0; + n = 0; + while (n < SC_NN) { + const gx : f64 = (reg[(SC_RGP + n) as usize] as f64) - (reg[(SC_RGN + n) as usize] as f64); + worst = fmax64(worst, fabs64(gx - jg[(SC_FG + n) as usize])); + n = n + 1; + } + return worst; +} + +/// The largest |J_fixed / 2^20 - J_real| over the 64 entries at the same state. +pub fn sc_gap_j(t: []const u32, v: []const u32, v0: []const u32, clk: u32, vdp: u32, vdn: u32, mult: []const u32) -> f64 { + var reg : [144]u64 = [0] ** 144; + sc_stamp(t, v, v0, clk, vdp, vdn, mult, ®); + var vf : [8]f64 = [0.0] ** 8; + var v0f : [8]f64 = [0.0] ** 8; + var n : u32 = 0; + while (n < SC_NN) { + vf[n as usize] = (v[n as usize]) as f64; + v0f[n as usize] = (v0[n as usize]) as f64; + n = n + 1; + } + var jg : [72]f64 = [0.0] ** 72; + var drv : [3]f64 = [0.0] ** 3; + drv[0] = (clk) as f64; + drv[1] = (vdp) as f64; + drv[2] = (vdn) as f64; + sc_stamp_f(0, t, &vf, &v0f, &drv, mult, &jg); + var worst : f64 = 0.0; + var k : usize = 0; + while (k < 64) { + const jx : f64 = ((reg[(SC_RJP as usize) + k] as f64) - (reg[(SC_RJN as usize) + k] as f64)) / CELL_F; + worst = fmax64(worst, fabs64(jx - jg[k])); + k = k + 1; + } + return worst; +} + +/// The largest |fixed - real| over the first steps 0..last, in units. +pub fn sc_dev_prefix(a: []const u32, b: []const f64, last: u32) -> f64 { + var worst : f64 = 0.0; + var i : usize = 0; + while (i < (((last + 1) * SC_NN) as usize)) { + worst = fmax64(worst, fabs64((a[i] as f64) - b[i])); + i = i + 1; + } + return worst; +} + +/// The highest value node sel takes over a real run, in units. +pub fn sc_peak_f(out: []const f64, sel: u32) -> f64 { + var m : f64 = 0.0; + var j : u32 = 0; + while (j <= SC_N) { + m = fmax64(m, out[(j * SC_NN + sel) as usize]); + j = j + 1; + } + return m; +} + +/// The largest magnitude among the first n sign-magnitude words. +pub fn sc_max_mag(w: []const u64, n: u32) -> u64 { + var m : u64 = 0; + var i : u32 = 0; + while (i < n) { + if (sc_mag(w[i as usize]) > m) { m = sc_mag(w[i as usize]); } + i = i + 1; + } + return m; +} + +/// The widest word the fixed-point solve forms in a whole run, as a bound: the largest register, the largest matrix word +/// before and after elimination, and the largest sum back substitution can form, 8 |a|max |x|max + |b|max 2^20. A u64 holds +/// the result when it is below 2^63; the tests ask for 2^62. +pub fn sc_headroom(t: []const u32, h: u32, swap: bool, k: u32, mult: []const u32) -> u64 { + var vdp : u32 = SC_VCM + h; + var vdn : u32 = SC_VCM - h; + if (swap) { + vdp = SC_VCM - h; + vdn = SC_VCM + h; + } + var v : [8]u32 = [0] ** 8; + var v0 : [8]u32 = [0] ** 8; + sc_rest(vdp, vdn, &v); + var wide : u64 = 0; + var j : u32 = 1; + while (j <= SC_N) { + var n : u32 = 0; + while (n < SC_NN) { + v0[n as usize] = v[n as usize]; + n = n + 1; + } + const clk : u32 = sc_clk_at(j); + var i : u32 = 0; + while (i < k) { + var reg : [144]u64 = [0] ** 144; + var a : [72]u64 = [0] ** 72; + var x : [8]u64 = [0] ** 8; + sc_stamp(t, &v, &v0, clk, vdp, vdn, mult, ®); + var r : u64 = 0; + var q : usize = 0; + while (q < SC_REG_LEN) { + if (reg[q] > r) { r = reg[q]; } + q = q + 1; + } + sc_load(®, &a); + const before : u64 = sc_max_mag(&a, 72); + sc_eliminate(&a); + sc_backsub(&a, &x); + const amax : u64 = sc_max_mag(&a, 72); + const xmax : u64 = sc_max_mag(&x, 8); + var bmax : u64 = 0; + n = 0; + while (n < SC_NN) { + if (sc_mag(a[(n * 9 + 8) as usize]) > bmax) { bmax = sc_mag(a[(n * 9 + 8) as usize]); } + n = n + 1; + } + const sums : u64 = 8 * amax * xmax + (bmax << CELL_BITS); + if (r > wide) { wide = r; } + if (before > wide) { wide = before; } + if (amax > wide) { wide = amax; } + if (sums > wide) { wide = sums; } + sc_apply(&v, &x); + i = i + 1; + } + j = j + 1; + } + return wide; +} + +/// How many of the first n words equal x. +pub fn sc_count_eq(r: []const u32, n: u32, x: u32) -> u32 { + var c : u32 = 0; + var i : u32 = 0; + while (i < n) { + if (r[i as usize] == x) { c = c + 1; } + i = i + 1; + } + return c; +} + +/// The highest value node sel takes over a fixed-point run. +pub fn sc_peak(out: []const u32, sel: u32) -> u32 { + var m : u32 = 0; + var j : u32 = 0; + while (j <= SC_N) { + if (out[(j * SC_NN + sel) as usize] > m) { m = out[(j * SC_NN + sel) as usize]; } + j = j + 1; + } + return m; +} + +/// Does node sel never rise (never exceed its previous sample) from step 1 on in a fixed-point run? +pub fn sc_never_rises(out: []const u32, sel: u32) -> bool { + var j : u32 = 2; + while (j <= SC_N) { + if (out[(j * SC_NN + sel) as usize] > out[((j - 1) * SC_NN + sel) as usize]) { return false; } + j = j + 1; + } + return true; +} + +/// The first crossing of Vdd / 2 by one node, as sc_t_first reports it (SC_NONE when the node never gets there). +pub fn sc_t_node(out: []const u32, sel: u32) -> u32 { + var s : [141]u32 = [0] ** 141; + sc_pick(out, sel, &s); + return sc_cross(&s); +} + +/// The largest of the first n words. +pub fn sc_max_u(r: []const u32, n: u32) -> u32 { + var m : u32 = 0; + var i : u32 = 0; + while (i < n) { + if (r[i as usize] > m) { m = r[i as usize]; } + i = i + 1; + } + return m; +} + +/// The largest of the first n reals. +pub fn sc_max_f(r: []const f64, n: u32) -> f64 { + var m : f64 = 0.0; + var i : u32 = 0; + while (i < n) { + m = fmax64(m, r[i as usize]); + i = i + 1; + } + return m; +} + +; --- the tests ---------------------------------------------------------------------------------- + + +test the_units_the_devices_and_the_circuit_are_as_stated { + // Variables, not constants: a constant assert is folded by the compiler and counts as no assert at all. + var cell : u32 = CELL; + var vdd : u32 = VDD; + var vcm : u32 = SC_VCM; + var r : u32 = SC_R; + var cstep : u32 = SC_CLK_STEP; + var upv : f64 = SC_UPV; + var vddf : f64 = VDD_F; + var vthu : u32 = SC_VTH_U; + var tmax : u32 = SC_TBL_MAX; + var h2 : u32 = SC_H2; + var h05 : u32 = SC_H05; + var h0125 : u32 = SC_H0125; + var h003 : u32 = SC_H003; + // The M3a grid, read as 24 cells of 62.5 mV: Vdd = 1.5 V = 1.5 * 2^24 units, the common mode is half of it. + assert_eq(vdd, 24 * cell); + assert_eq(vddf, 1.5 * upv); + assert_eq(vcm * 2, vdd); + // The clock is a ramp of 16 steps of exactly Vdd / 16, and is at Vdd / 2 at step 8. + assert_eq(cstep * r, vdd); + assert_eq(sc_clk_at(0), 0); + assert_eq(sc_clk_at(r >> 1), vcm); + assert_eq(sc_clk_at(r), vdd); + assert_eq(sc_clk_at(r + 50), vdd); + // dt = 0.1 ps over the 0.125 fF charge unit is 800 V/A; 4, 4, 6, 2, 2, 4, 4, 2 fF. + const dtc : f64 = (SC_DT_FS as f64) * 0.000000000000001 / ((SC_CUNIT_AF as f64) * 0.000000000000000001) - SC_STEP_V_PER_A; + assert dtc * dtc < 0.000001 * 0.000001; + var capsum : u32 = 0; + var i : u32 = 0; + while (i < SC_NN) { + capsum = capsum + SC_CAP[i as usize]; + i = i + 1; + } + assert_eq(capsum, 224); + assert_eq(SC_CAP[5], 32); + assert_eq(SC_CAP[6], 32); + // beta = KP W / L is 666.7 uA/V^2 per unit m for both types, which is why ONE table serves both transistors. + const bn : f64 = SC_KP_N * SC_W_N / SC_L - SC_BETA; + const bp : f64 = SC_KP_P * SC_W_P / SC_L - SC_BETA; + assert bn * bn < 0.00000000000000001 * 0.00000000000000001; + assert bp * bp < 0.00000000000000001 * 0.00000000000000001; + // The threshold 0.35 V is 5.6 cells: inside cell 5, so the table has a soft knee (the next test measures it). + const vthd : f64 = (vthu as f64) - SC_VTH_V * upv; + assert vthd * vthd < 0.5 * 0.5; + assert vthu > 5 * cell; + assert vthu < 6 * cell; + // The inputs are Vcm +- h; 2 h is the stated differential to within half a unit each side (59.6 nV). + const d2 : f64 = 2.0 * (h2 as f64) / upv - 0.002; + const d05 : f64 = 2.0 * (h05 as f64) / upv - 0.0005; + const d0125 : f64 = 2.0 * (h0125 as f64) / upv - 0.000125; + const d003 : f64 = 2.0 * (h003 as f64) / upv - 0.00003125; + assert d2 * d2 < 0.0000001 * 0.0000001; + assert d05 * d05 < 0.0000001 * 0.0000001; + assert d0125 * d0125 < 0.0000001 * 0.0000001; + assert d003 * d003 < 0.0000001 * 0.0000001; + // The device table: fifteen devices, eight NMOS (kind 0) and seven PMOS, every terminal one of the 13 names. + var ndev : u32 = SC_NDEV; + var pmos : u32 = 0; + var msum : u32 = 0; + var ok : bool = true; + var e : usize = 0; + while (e < (ndev as usize)) { + pmos = pmos + SC_DV_KIND[e]; + msum = msum + SC_DV_M[e]; + if (SC_DV_D[e] >= 13 or SC_DV_G[e] >= 13 or SC_DV_S[e] >= 13) { ok = false; } + if (SC_DV_D[e] == SC_DV_S[e] or SC_DV_G[e] == SC_DV_D[e] or SC_DV_G[e] == SC_DV_S[e]) { ok = false; } + if (SC_DV_D[e] >= SC_NN and SC_DV_S[e] >= SC_NN) { ok = false; } + e = e + 1; + } + assert(ok); + assert_eq(pmos, 7); + assert_eq(msum, 266); + // The input pair: XM3 (gate VDN) and XM4 (gate VDP), both m = 36, on ap and an, sharing net1 with the tail XM5 (m = 18). + assert_eq(SC_DV_G[2], SC_T_VDN); + assert_eq(SC_DV_G[3], SC_T_VDP); + assert_eq(SC_DV_M[2], SC_DV_M[3]); + assert_eq(SC_DV_S[2], SC_DV_S[3]); + assert_eq(SC_DV_S[2], SC_DV_D[4]); + // The cross-coupled NMOS pair XM10, XM11 (m = 12) and PMOS pair XM8, XM9 (m = 30): each output gates the other's pull-down. + assert_eq(SC_DV_G[9], SC_DV_D[10]); + assert_eq(SC_DV_G[10], SC_DV_D[9]); + assert_eq(SC_DV_D[7], SC_DV_D[10]); + assert_eq(SC_DV_D[8], SC_DV_D[9]); + // The widest product of the interpolation: TBL_MAX * 2^40 and the rounding half fit in a u64 (the same rule as M3a). + assert tmax < 8388608; + assert (tmax as u64) * 1099511627776 + 549755813888 < 9223372036854775808; +} + +test the_table_is_the_law_at_the_nodes_and_the_knee_sits_inside_a_cell { + var t : [625]u32 = [0] ** 625; + sc_build_table(&t, SC_VTH_V); + // At the 625 nodes the entry is the law to within half a unit (59.6 nV of node voltage, in the charge unit). + assert sc_node_err(&t) <= 0.5 + 0.0000001; + assert sc_node_err(&t) > 0.3; + assert(monotone(&t)); + assert_eq(t[624], SC_TBL_MAX); + assert_eq(t[0], 0); + // Vgs on node 5 (0.3125 V) and below: no current. Vds = 0: no current. + var ix : u32 = 0; + var iy : u32 = 0; + while (ix <= 5) { + iy = 0; + while (iy < NODES) { + assert_eq(t[(ix * NODES + iy) as usize], 0); + iy = iy + 1; + } + ix = ix + 1; + } + ix = 0; + while (ix < NODES) { + assert_eq(t[(ix * NODES) as usize], 0); + ix = ix + 1; + } + // Node 6 (0.375 V) has 25 mV of overdrive: T = beta ov^2 / 2 * 800 * 2^24 = 2796 units, flat in Vds from there on. + assert_eq(t[(6 * NODES + 24) as usize], 2796); + assert_eq(t[(6 * NODES + 2) as usize], t[(6 * NODES + 24) as usize]); + // Saturation is flat in Vds: Vgs = 0.75 V has ov = 0.4 V, so Vds >= 0.4375 V (node 7) gives one value, node 6 less. + assert_eq(t[(12 * NODES + 7) as usize], t[(12 * NODES + 24) as usize]); + assert t[(12 * NODES + 6) as usize] < t[(12 * NODES + 7) as usize]; + // THE KNEE. The law is exactly 0 at Vgs = 0.33 V, a chord from node 5 (0) to node 6 (2796) is not: the table reads + // 783 units there, about 28 % of the next node's value. The table's threshold is soft inside cell 5; it is the price of + // reading M3a's 24-cell grid as 62.5 mV cells, and PSP103 tables (M3b) will have a subthreshold tail of their own. + const leak : u32 = id_of(&t, 5536481, SC_VCM); + assert_eq(leak, 783); + assert leak > 0; + assert sc_law_t(5536481.0, 12582912.0) == 0.0; + // N and P read this one table: the unit device's law at the top node is beta (1.5 - 0.35)^2 / 2 = 0.44083 mA. + const top : f64 = sc_law_t(VDD_F, VDD_F) / (SC_STEP_V_PER_A * SC_UPV) - 0.000440833333333; + assert top * top < 0.000000001 * 0.000000001; +} + +test a_device_conducts_from_its_high_end_to_its_low_end { + var t : [625]u32 = [0] ** 625; + sc_build_table(&t, SC_VTH_V); + var v : [8]u32 = [0] ** 8; + var reg : [144]u64 = [0] ** 144; + // XM3 (index 2): NMOS, drain ap, gate VDN = 0.75 V, source net1, m = 36. Forward: ap at 1.0 V, net1 at 0.2 V. + v[0] = 16777216; + v[2] = 3355443; + sc_stamp_dev(&t, 2, &v, 0, SC_VCM, SC_VCM, 36, ®); + const x : u32 = SC_VCM - 3355443; + const y : u32 = 16777216 - 3355443; + const cur : u64 = 36 * (id_of(&t, x, y) as u64); + assert cur > 0; + assert_eq(reg[(SC_RGP + 0) as usize], cur); + assert_eq(reg[(SC_RGN + 2) as usize], cur); + assert_eq(reg[(SC_RGP + 2) as usize], 0); + assert_eq(reg[(SC_RGN + 0) as usize], 0); + // The table's current against the law: m T(x, y), here 0.55 V of drive and 0.8 V across the channel. + const law : f64 = 36.0 * sc_law_t((x as f64), (y as f64)); + const rel : f64 = ((cur as f64) - law) / law; + assert rel * rel < 0.02 * 0.02; + assert rel * rel > 0.01 * 0.01; + // Reverse: net1 at 1.0 V, ap at 0.2 V. The channel is symmetric, so the same current now leaves net1 and enters ap. + var w : [8]u32 = [0] ** 8; + var reg2 : [144]u64 = [0] ** 144; + w[0] = 3355443; + w[2] = 16777216; + sc_stamp_dev(&t, 2, &w, 0, SC_VCM, SC_VCM, 36, ®2); + assert_eq(reg2[(SC_RGP + 2) as usize], cur); + assert_eq(reg2[(SC_RGN + 0) as usize], cur); + assert_eq(reg2[(SC_RGP + 0) as usize], 0); + assert_eq(reg2[(SC_RGN + 2) as usize], 0); + // XM8 (index 7): PMOS, drain dp, gate an, source net3, m = 30; all three are internal nodes. Forward: net3 at 1.4 V. + var p : [8]u32 = [0] ** 8; + var reg3 : [144]u64 = [0] ** 144; + p[4] = 23488102; + p[5] = 5033165; + p[1] = 6710886; + sc_stamp_dev(&t, 7, &p, 0, SC_VCM, SC_VCM, 30, ®3); + const px : u32 = 23488102 - 6710886; + const py : u32 = 23488102 - 5033165; + const pcur : u64 = 30 * (id_of(&t, px, py) as u64); + assert pcur > 0; + assert_eq(reg3[(SC_RGP + 4) as usize], pcur); + assert_eq(reg3[(SC_RGN + 5) as usize], pcur); + // Reverse: dp at 1.4 V and net3 at 0.3 V, the gate unchanged: the same drive (1.4 - 0.4 V) and span, the roles exchanged. + var q : [8]u32 = [0] ** 8; + var reg4 : [144]u64 = [0] ** 144; + q[5] = 23488102; + q[4] = 5033165; + q[1] = 6710886; + sc_stamp_dev(&t, 7, &q, 0, SC_VCM, SC_VCM, 30, ®4); + assert_eq(reg4[(SC_RGP + 5) as usize], pcur); + assert_eq(reg4[(SC_RGN + 4) as usize], pcur); + // A common shift of all three terminals changes no current, so every row of the Jacobian sums to zero (all three are + // internal here): row sums of (J+ - J-) over the eight columns. + var row : u32 = 0; + var worst : f64 = 0.0; + while (row < SC_NN) { + var sum : f64 = 0.0; + var col : u32 = 0; + while (col < SC_NN) { + sum = sum + (reg3[(SC_RJP + row * SC_NN + col) as usize] as f64) - (reg3[(SC_RJN + row * SC_NN + col) as usize] as f64); + sum = sum + (reg4[(SC_RJP + row * SC_NN + col) as usize] as f64) - (reg4[(SC_RJN + row * SC_NN + col) as usize] as f64); + col = col + 1; + } + worst = fmax64(worst, fabs64(sum)); + row = row + 1; + } + assert worst < 0.5; + // And the Jacobian of the PMOS is not empty: the gm entry (gate column an) of the high end's row is negative for a PMOS. + assert reg3[(SC_RJN + 4 * SC_NN + 1) as usize] > 0; + assert_eq(reg3[(SC_RJP + 4 * SC_NN + 1) as usize], 0); +} + +test the_jacobian_is_the_derivative_of_the_residual { + var t : [625]u32 = [0] ** 625; + sc_build_table(&t, SC_VTH_V); + var outa : [1128]f64 = [0.0] ** 1128; + var resa : [140]f64 = [0.0] ** 140; + var outt : [1128]f64 = [0.0] ** 1128; + var rest : [140]f64 = [0.0] ** 140; + sc_run_f(1, &t, SC_H2, false, SC_K_ANA, &SC_DV_M, &outa, &resa); + sc_run_f(0, &t, SC_H2, false, SC_K, &SC_DV_M, &outt, &rest); + var drv : [3]f64 = [0.0] ** 3; + drv[1] = ((SC_VCM + SC_H2) as f64); + drv[2] = ((SC_VCM - SC_H2) as f64); + var steps : [7]u32 = [1, 12, 24, 40, 62, 70, 110]; + var v : [8]f64 = [0.0] ** 8; + var v0 : [8]f64 = [0.0] ** 8; + var worst_law : f64 = 0.0; + var worst_tab : f64 = 0.0; + var wrong : f64 = 0.0; + var k : usize = 0; + while (k < 7) { + const j : u32 = steps[k]; + drv[0] = (sc_clk_at(j)) as f64; + sc_state_f(&outa, j, &v); + sc_state_f(&outa, j - 1, &v0); + worst_law = fmax64(worst_law, sc_fd_err(1, &t, &v, &v0, &drv, &SC_DV_M, false)); + wrong = fmax64(wrong, sc_fd_err(1, &t, &v, &v0, &drv, &SC_DV_M, true)); + sc_state_f(&outt, j, &v); + sc_state_f(&outt, j - 1, &v0); + worst_tab = fmax64(worst_tab, sc_fd_err(0, &t, &v, &v0, &drv, &SC_DV_M, false)); + k = k + 1; + } + // The analytic law is C1, so a central difference of the residual (+-200 units, 12 uV) is its Jacobian to 1e-4 of a pure + // number, entries being O(1) to O(50). What is left is the kink in the second derivative where a device sits exactly at + // Vds = 0, as several do at rest (steps 1, 12, 24); from step 40 on the difference is below 1e-6. + assert worst_law < 0.0002; + assert worst_law > 0.00005; + // The table's slopes ARE the derivative of its bilinear interpolant inside a cell; the same kinks (and the cell edges) leave 2e-5. + assert worst_tab < 0.0001; + assert worst_tab > 0.00001; + // The check has teeth: the transposed Jacobian (a row and a column exchanged, the commonest stamping slip) is far off. + assert wrong > 0.5; +} + +test the_fixed_point_system_is_the_real_system_to_rounding { + var t : [625]u32 = [0] ** 625; + sc_build_table(&t, SC_VTH_V); + var out : [1128]u32 = [0] ** 1128; + var res : [140]u32 = [0] ** 140; + const s : u64 = sc_run(&t, SC_H2, false, SC_K, &SC_DV_M, &out, &res); + var steps : [5]u32 = [12, 40, 62, 70, 110]; + var v : [8]u32 = [0] ** 8; + var v0 : [8]u32 = [0] ** 8; + var wg : f64 = 0.0; + var wj : f64 = 0.0; + var k : usize = 0; + while (k < 5) { + const j : u32 = steps[k]; + sc_state_u(&out, j, &v); + sc_state_u(&out, j - 1, &v0); + wg = fmax64(wg, sc_gap_g(&t, &v, &v0, sc_clk_at(j), SC_VCM + SC_H2, SC_VCM - SC_H2, &SC_DV_M)); + wj = fmax64(wj, sc_gap_j(&t, &v, &v0, sc_clk_at(j), SC_VCM + SC_H2, SC_VCM - SC_H2, &SC_DV_M)); + k = k + 1; + } + // G: each device's current is read, rounded to a unit, and multiplied by m (at most 36): half of m per device, 266 / 2 in all. + assert wg < 133.0; + assert wg > 1.0; + // J: each slope is read and rounded to 2^-20, times m (at most 36): well under 1e-4 of a pure number. + assert wj < 0.0001; + assert wj > 0.000001; + assert_eq(s, 858149221550211914); +} + +test the_gauss_solver_pivots_and_agrees_with_the_real_solver { + // (a) A permutation system with signs: A has 1.0 (2^20) at (i, 7 - i), so every diagonal entry is zero and each column + // needs a row swap. The solution is exact: x = (100, -200, 300, -400, 500, -600, 700, -800) units. + var a : [72]u64 = [0] ** 72; + var x : [8]u64 = [0] ** 8; + var xt : [8]u64 = [0] ** 8; + var i : u32 = 0; + while (i < 8) { + xt[i as usize] = sc_pack(((i + 1) * 100) as u64, (i % 2) == 1); + i = i + 1; + } + i = 0; + while (i < 8) { + const cneg : bool = (i % 3) == 1; + const j : u32 = 7 - i; + a[(i * 9 + j) as usize] = sc_pack(CELL as u64, cneg); + a[(i * 9 + 8) as usize] = sc_pack(sc_mag(xt[j as usize]), cneg != sc_neg(xt[j as usize])); + i = i + 1; + } + sc_eliminate(&a); + sc_backsub(&a, &x); + i = 0; + while (i < 8) { + assert_eq(x[i as usize], xt[i as usize]); + i = i + 1; + } + // (b) The real matrix of the run at step 60, from the previous state as the first guess: the full first Newton move. + var t : [625]u32 = [0] ** 625; + sc_build_table(&t, SC_VTH_V); + var out : [1128]u32 = [0] ** 1128; + var res : [140]u32 = [0] ** 140; + const s : u64 = sc_run(&t, SC_H2, false, SC_K, &SC_DV_M, &out, &res); + var v0 : [8]u32 = [0] ** 8; + sc_state_u(&out, 59, &v0); + var reg : [144]u64 = [0] ** 144; + sc_stamp(&t, &v0, &v0, sc_clk_at(60), SC_VCM + SC_H2, SC_VCM - SC_H2, &SC_DV_M, ®); + var b : [72]u64 = [0] ** 72; + var y : [8]u64 = [0] ** 8; + sc_load(®, &b); + var af : [72]f64 = [0.0] ** 72; + var yf : [8]f64 = [0.0] ** 8; + i = 0; + while (i < 72) { + af[i as usize] = sc_word_f(b[i as usize]); + if ((i % 9) < 8) { af[i as usize] = af[i as usize] / CELL_F; } + i = i + 1; + } + sc_eliminate(&b); + sc_backsub(&b, &y); + sc_solve_f(&af, &yf); + var worst : f64 = 0.0; + var biggest : f64 = 0.0; + i = 0; + while (i < 8) { + worst = fmax64(worst, fabs64(sc_word_f(y[i as usize]) - yf[i as usize])); + biggest = fmax64(biggest, fabs64(yf[i as usize])); + i = i + 1; + } + // The first move at step 60 is thousands of units; the fixed-point solve agrees with the real one to a few. + assert biggest > 1000.0; + assert worst < 4.0; + assert worst > 0.1; + assert_eq(s, 858149221550211914); +} + +test the_circuit_rests_until_the_clock_moves { + var t : [625]u32 = [0] ** 625; + sc_build_table(&t, SC_VTH_V); + var v : [8]u32 = [0] ** 8; + var v0 : [8]u32 = [0] ** 8; + const vdp : u32 = SC_VCM + SC_H2; + const vdn : u32 = SC_VCM - SC_H2; + sc_rest(vdp, vdn, &v); + sc_rest(vdp, vdn, &v0); + assert_eq(v[0], VDD); + assert_eq(v[1], VDD); + assert_eq(v[2], vdp - SC_VTH_U); + assert_eq(v[3], VDD); + assert_eq(v[4], VDD); + assert_eq(v[5], 0); + assert_eq(v[6], 0); + assert_eq(v[7], VDD); + // At the rest state with the clock low five nodes are in balance, exactly: net2, net3, dp, dn and n_clkc. + var reg : [144]u64 = [0] ** 144; + sc_stamp(&t, &v, &v0, 0, vdp, vdn, &SC_DV_M, ®); + assert_eq(reg[(SC_RGP + SC_NET2) as usize], reg[(SC_RGN + SC_NET2) as usize]); + assert_eq(reg[(SC_RGP + SC_NET3) as usize], reg[(SC_RGN + SC_NET3) as usize]); + assert_eq(reg[(SC_RGP + SC_DP) as usize], reg[(SC_RGN + SC_DP) as usize]); + assert_eq(reg[(SC_RGP + SC_DN) as usize], reg[(SC_RGN + SC_DN) as usize]); + assert_eq(reg[(SC_RGP + SC_NCLK) as usize], reg[(SC_RGN + SC_NCLK) as usize]); + // Three are not: both input devices sit at the threshold 0.35 V, which is inside the table's cell 5, so their chords + // conduct. The current leaves ap and an (G > 0: they sag a little under their precharge PMOS) and charges net1 (G < 0) + // until net1 is a whole 0.3125 V under the higher input. The real circuit (the decks) rests at once; this is the knee. + assert reg[(SC_RGP + SC_AP) as usize] > reg[(SC_RGN + SC_AP) as usize]; + assert reg[(SC_RGP + SC_AN) as usize] > reg[(SC_RGN + SC_AN) as usize]; + assert reg[(SC_RGN + SC_NET1) as usize] > reg[(SC_RGP + SC_NET1) as usize]; + var w : [8]u32 = [0] ** 8; + var n : u32 = 0; + while (n < SC_NN) { + w[n as usize] = v[n as usize]; + n = n + 1; + } + const r : u32 = sc_step(&t, &w, 0, vdp, vdn, SC_K, &SC_DV_M); + assert r <= SC_RES_BOUND; + // One step at clock 0 (0.1 ps): the five balanced nodes do not move; ap and an sag by 91 and 96 uV, net1 rises by 145 uV. + assert_eq(w[SC_NET2 as usize], v[SC_NET2 as usize]); + assert_eq(w[SC_NET3 as usize], v[SC_NET3 as usize]); + assert_eq(w[SC_DP as usize], v[SC_DP as usize]); + assert_eq(w[SC_DN as usize], v[SC_DN as usize]); + assert_eq(w[SC_NCLK as usize], v[SC_NCLK as usize]); + assert_eq(v[0] - w[0], 1524); + assert_eq(v[1] - w[1], 1610); + assert_eq(w[2] - v[2], 2440); +} + +test the_plus_2_mv_run_is_pinned { + var t : [625]u32 = [0] ** 625; + sc_build_table(&t, SC_VTH_V); + var out : [1128]u32 = [0] ** 1128; + var res : [140]u32 = [0] ** 140; + const s : u64 = sc_run(&t, SC_H2, false, SC_K, &SC_DV_M, &out, &res); + assert_eq(s, 858149221550211914); + // VDP > VDN: dp wins. It is the first output through Vdd / 2, at 62.8876 steps (the clock is at 50 % at step 8, so + // 54.8876 steps = 5.489 ps after the clock edge). The loser dn rises with it, crosses Vdd / 2 3.7 steps later and falls + // back: ngspice's loser peaks 3 mV UNDER Vdd / 2 (0.7473 V), the table model's peak is a few mV over, so whether the loser + // "crosses" is below the model's resolution and is not used; the winner's crossing and the final sign are. + assert(sc_sign(&out)); + assert_eq(sc_t_first(&out), 257587); + assert_eq(sc_t_node(&out, SC_DP), 257587); + assert_eq(sc_t_node(&out, SC_DN), 272643); + // |dp - dn| reaches Vdd / 2 at 77.6312 steps: 69.6312 steps = 6.963 ps after the clock edge. + assert_eq(sc_t_res(&out), 317977); + const tf : f64 = ((sc_after_clock(sc_t_first(&out))) as f64) / 4096.0; + const tr : f64 = ((sc_after_clock(sc_t_res(&out))) as f64) / 4096.0; + assert tf > 54.887; + assert tf < 54.888; + assert tr > 69.631; + assert tr < 69.632; + // At step 140 the latch has settled: dp within 10 units of Vdd, dn at 47375 units = 2.8 mV. + assert_eq(out[(140 * 8 + 5) as usize], 25165814); + assert_eq(out[(140 * 8 + 6) as usize], 47375); + // The loser's peak, in units: ngspice 0.747264 V (the decks' dnmax). + assert_eq(sc_peak(&out, SC_DN), 12632942); + // The shape of the evaluation: the precharged nodes ap and an and the inverted clock only fall. + assert(sc_never_rises(&out, SC_AP)); + assert(sc_never_rises(&out, SC_AN)); + assert(sc_never_rises(&out, SC_NCLK)); +} + +test all_four_stimuli_decide_the_right_way_and_the_halves_mirror { + var t : [625]u32 = [0] ** 625; + sc_build_table(&t, SC_VTH_V); + var a : [1128]u32 = [0] ** 1128; + var b : [1128]u32 = [0] ** 1128; + var c : [1128]u32 = [0] ** 1128; + var d : [1128]u32 = [0] ** 1128; + var res : [140]u32 = [0] ** 140; + const sa : u64 = sc_run(&t, SC_H2, false, SC_K, &SC_DV_M, &a, &res); + const sb : u64 = sc_run(&t, SC_H2, true, SC_K, &SC_DV_M, &b, &res); + const sc : u64 = sc_run(&t, SC_H05, false, SC_K, &SC_DV_M, &c, &res); + const sd : u64 = sc_run(&t, SC_H05, true, SC_K, &SC_DV_M, &d, &res); + // +2 mV and +0.5 mV decide "VDP > VDN" (dp high); -2 mV and -0.5 mV decide it the other way. + assert(sc_sign(&a)); + assert(!sc_sign(&b)); + assert(sc_sign(&c)); + assert(!sc_sign(&d)); + // The first output to cross is the winner: dp for the positive inputs, dn for the negative ones, at the same time. + assert_eq(sc_t_node(&a, SC_DP), sc_t_node(&b, SC_DN)); + assert_eq(sc_t_node(&c, SC_DP), sc_t_node(&d, SC_DN)); + assert_eq(sc_t_first(&a), sc_t_first(&b)); + assert_eq(sc_t_first(&c), 260428); + assert_eq(sc_t_first(&d), sc_t_first(&c)); + assert_eq(sc_t_res(&a), sc_t_res(&b)); + assert_eq(sc_t_res(&c), 349568); + assert_eq(sc_t_res(&d), sc_t_res(&c)); + // At 0.5 mV the loser (ngspice: peak 0.822 V) also gets through Vdd / 2, after the winner. + assert sc_t_node(&c, SC_DN) != SC_NONE; + assert sc_t_node(&c, SC_DN) > sc_t_node(&c, SC_DP); + // The slower decision at the smaller input: first output (commonly driven) 0.69 steps later, resolution 7.71 steps later. + assert sc_t_first(&c) > sc_t_first(&a); + assert sc_t_res(&c) - sc_t_res(&a) > 31000; + assert sc_t_res(&c) - sc_t_res(&a) < 32000; + assert_eq(sa, 858149221550211914); + assert_eq(sb, 7832980476838438025); + assert_eq(sc, 7639752914144809206); + assert_eq(sd, 42403085586441860); + // The comparator is its own mirror. In real arithmetic the exchanged run is the mirrored run (see the f64 test); in + // fixed point the two halves round differently (the pivot order differs), by at most 2 units at 2 mV and 61 at 0.5 mV, + // which regeneration has amplified from one. + assert_eq(sc_mirror_dev_u(&a, &b), 2); + assert_eq(sc_mirror_dev_u(&c, &d), 61); +} + +test the_residual_is_bounded_at_every_step { + var t : [625]u32 = [0] ** 625; + sc_build_table(&t, SC_VTH_V); + var out : [1128]u32 = [0] ** 1128; + var res : [140]u32 = [0] ** 140; + var hs : [6]u32 = [SC_H2, SC_H05, SC_H0125, SC_H003, SC_H2, SC_H05]; + var n : u32 = 0; + var worst : u32 = 0; + var twos : u32 = 0; + while (n < 6) { + const s : u64 = sc_run(&t, hs[n as usize], n >= 4, SC_K, &SC_DV_M, &out, &res); + var j : u32 = 0; + while (j < SC_N) { + assert res[j as usize] <= SC_RES_BOUND; + j = j + 1; + } + if (sc_max_u(&res, 140) > worst) { worst = sc_max_u(&res, 140); } + twos = twos + sc_count_eq(&res, 140, 2); + n = n + 1; + } + // From both sides, over six runs of 140 steps (2 mV, 0.5 mV, 0.125 mV, 31 uV, and the two mirrored): the worst residual is + // exactly 2 units, so the bound is tight; it is the rounding of the table reads times m, over the node's own slope, and it + // is reached at 31 of the 840 steps and exceeded at none. + assert_eq(worst, 2); + assert_eq(twos, 31); + // The real twin converges further: after K = 3 real Newton iterations the residual is 1.5e-4 of a unit at 2 mV and 4.6e-5 + // at 0.5 mV (bound 2e-4, from both sides). + var outf : [1128]f64 = [0.0] ** 1128; + var resf : [140]f64 = [0.0] ** 140; + sc_run_f(0, &t, SC_H2, false, SC_K, &SC_DV_M, &outf, &resf); + const r2 : f64 = sc_max_f(&resf, 140); + sc_run_f(0, &t, SC_H05, false, SC_K, &SC_DV_M, &outf, &resf); + const r05 : f64 = sc_max_f(&resf, 140); + assert r2 < 0.0002; + assert r2 > 0.0001; + assert r05 < 0.00005; + assert r05 > 0.00004; +} + +test k1_and_k2_are_the_negative_controls_they_break_the_residual_bound { + var t : [625]u32 = [0] ** 625; + sc_build_table(&t, SC_VTH_V); + var out : [1128]u32 = [0] ** 1128; + var res : [140]u32 = [0] ** 140; + var out1 : [1128]u32 = [0] ** 1128; + var res1 : [140]u32 = [0] ** 140; + var out2 : [1128]u32 = [0] ** 1128; + var res2 : [140]u32 = [0] ** 140; + const s3 : u64 = sc_run(&t, SC_H2, false, SC_K, &SC_DV_M, &out, &res); + const s1 : u64 = sc_run(&t, SC_H2, false, 1, &SC_DV_M, &out1, &res1); + const s2 : u64 = sc_run(&t, SC_H2, false, 2, &SC_DV_M, &out2, &res2); + // One Newton iteration leaves 29361 units (1.75 mV) at the worst step, 14680 times the bound. + assert sc_max_u(&res1, 140) > SC_RES_BOUND; + assert_eq(sc_max_u(&res1, 140), 29361); + // Two leave 55 units (3.3 uV), 27 times the bound; three leave 2. + assert sc_max_u(&res2, 140) > SC_RES_BOUND; + assert_eq(sc_max_u(&res2, 140), 55); + assert_eq(sc_max_u(&res, 140), 2); + assert s1 != s3; + assert s2 != s3; + // K = 1 also moves the decision: the first output 88 counts of 1/4096 step late (0.021 steps), the resolution 870 counts (0.21 steps). + // K = 2 does not move it at this resolution: the bound is about the solve, not about the answer. + assert_eq(sc_t_first(&out1), 257675); + assert_eq(sc_t_res(&out1), 318847); + assert_eq(sc_t_first(&out2), sc_t_first(&out)); + assert_eq(sc_t_res(&out2), sc_t_res(&out)); + // And K = 3 is enough: a fourth iteration buys nothing. Its worst residual is also 2 units and its times are the same. + var out4 : [1128]u32 = [0] ** 1128; + var res4 : [140]u32 = [0] ** 140; + const s4 : u64 = sc_run(&t, SC_H2, false, 4, &SC_DV_M, &out4, &res4); + assert_eq(sc_max_u(&res4, 140), 2); + assert_eq(sc_t_first(&out4), sc_t_first(&out)); + assert_eq(sc_t_res(&out4), sc_t_res(&out)); + assert s4 != 0; +} + +test the_arithmetic_error_against_the_real_twin { + var t : [625]u32 = [0] ** 625; + sc_build_table(&t, SC_VTH_V); + var out : [1128]u32 = [0] ** 1128; + var res : [140]u32 = [0] ** 140; + var outf : [1128]f64 = [0.0] ** 1128; + var resf : [140]f64 = [0.0] ** 140; + const s2 : u64 = sc_run(&t, SC_H2, false, SC_K, &SC_DV_M, &out, &res); + sc_run_f(0, &t, SC_H2, false, SC_K, &SC_DV_M, &outf, &resf); + // The fixed-point run against the real run of the SAME iteration (the same table words, K = 3, the same clamps), the largest + // |difference| over the eight nodes and 141 states, in units of 59.6 nV: 68.05 at 2 mV (bound 70, from both sides). + const a2 : f64 = sc_dev_fix_f(&out, &outf); + assert a2 < 70.0; + assert a2 > 66.0; + // Until step 50 it is under 3 units (0.2 uV): the arithmetic itself is tight. The regeneration then amplifies it, e-fold per + // 5.3 steps, to the 68 above. (2.915 at 2 mV, bound 3.) + const p2 : f64 = sc_dev_prefix(&out, &outf, 50); + assert p2 < 3.0; + assert p2 > 2.8; + const f2 : f64 = (sc_after_clock(sc_t_first(&out)) as f64) / 4096.0 - sc_t_first_f(&outf); + const r2 : f64 = (sc_after_clock(sc_t_res(&out)) as f64) / 4096.0 - sc_t_res_f(&outf); + // The decision times differ by under 2.5e-4 steps, which is the 1/4096 step the fixed-point time is truncated to. + assert f2 * f2 < 0.00025 * 0.00025; + assert r2 * r2 < 0.00025 * 0.00025; + const s05 : u64 = sc_run(&t, SC_H05, false, SC_K, &SC_DV_M, &out, &res); + sc_run_f(0, &t, SC_H05, false, SC_K, &SC_DV_M, &outf, &resf); + // At 0.5 mV the decision takes longer, so the amplification has longer to act: 250.8 units (15 uV; bound 255, from both sides). + const a05 : f64 = sc_dev_fix_f(&out, &outf); + assert a05 < 255.0; + assert a05 > 245.0; + const p05 : f64 = sc_dev_prefix(&out, &outf, 50); + assert p05 < 3.4; + assert p05 > 3.2; + const f05 : f64 = (sc_after_clock(sc_t_first(&out)) as f64) / 4096.0 - sc_t_first_f(&outf); + const r05 : f64 = (sc_after_clock(sc_t_res(&out)) as f64) / 4096.0 - sc_t_res_f(&outf); + assert f05 * f05 < 0.00025 * 0.00025; + assert r05 * r05 < 0.00035 * 0.00035; + assert r05 * r05 > 0.0003 * 0.0003; + assert_eq(s2, 858149221550211914); + assert_eq(s05, 7639752914144809206); +} + +test the_table_error_against_the_analytic_law { + var t : [625]u32 = [0] ** 625; + sc_build_table(&t, SC_VTH_V); + var tab : [1128]f64 = [0.0] ** 1128; + var law : [1128]f64 = [0.0] ** 1128; + var res : [140]f64 = [0.0] ** 140; + sc_run_f(0, &t, SC_H2, false, SC_K, &SC_DV_M, &tab, &res); + sc_run_f(1, &t, SC_H2, false, SC_K_ANA, &SC_DV_M, &law, &res); + // The real table run against the analytic square law solved to convergence (12 Newton iterations a step) on the same + // backward Euler grid: the largest |difference| of any node at any step, in volts. This is the table's error, the 62.5 mV + // cells and the knee, and it is the largest error here: 15.2 mV at 2 mV (bound 15.4 mV, from both sides). + const e2 : f64 = sc_dev_f_f(&tab, &law) / SC_UPV; + assert e2 < 0.0154; + assert e2 > 0.0150; + const f2 : f64 = sc_t_first_f(&tab) - sc_t_first_f(&law); + const r2 : f64 = sc_t_res_f(&tab) - sc_t_res_f(&law); + sc_run_f(0, &t, SC_H05, false, SC_K, &SC_DV_M, &tab, &res); + sc_run_f(1, &t, SC_H05, false, SC_K_ANA, &SC_DV_M, &law, &res); + const e05 : f64 = sc_dev_f_f(&tab, &law) / SC_UPV; + // 19.0 mV at 0.5 mV (bound 19.2 mV). + assert e05 < 0.0192; + assert e05 > 0.0188; + const f05 : f64 = sc_t_first_f(&tab) - sc_t_first_f(&law); + const r05 : f64 = sc_t_res_f(&tab) - sc_t_res_f(&law); + // In time the table is 0.075 and 0.112 steps early to the first output, 0.323 and 0.406 steps late to the resolution + // (0.14 %, 0.20 %, 0.47 %, 0.53 %): the first output is the common-mode rise, the resolution is the regeneration, which + // the table's slopes (piecewise constant over a cell) slow a little. + assert f2 < 0.0 - 0.07; + assert f2 > 0.0 - 0.08; + assert f05 < 0.0 - 0.11; + assert f05 > 0.0 - 0.12; + assert r2 > 0.32; + assert r2 < 0.33; + assert r05 > 0.40; + assert r05 < 0.41; +} + +test the_analytic_run_is_the_ngspice_deck { + var t : [625]u32 = [0] ** 625; + sc_build_table(&t, SC_VTH_V); + var law : [1128]f64 = [0.0] ** 1128; + var law20 : [1128]f64 = [0.0] ** 1128; + var res : [140]f64 = [0.0] ** 140; + sc_run_f(1, &t, SC_H2, false, SC_K_ANA, &SC_DV_M, &law, &res); + sc_run_f(1, &t, SC_H2, false, 20, &SC_DV_M, &law20, &res); + // Newton is converged: 12 iterations and 20 give the same run to a millionth of a unit. + assert sc_dev_f_f(&law, &law20) < 0.000001; + // The analytic square law, backward Euler on the grid, K = 12, against ngspice MOSFET LEVEL=1 on the same grid: the largest + // |difference| at the eight nodes and eight recorded steps, in volts. It is 1.44e-9 V at 2 mV (bound 1.6e-9, from both + // sides): the deck IS the spec's circuit, device for device and farad for farad. Wrong by a device, a node cap or a + // wire, it would be millivolts. + const e2 : f64 = sc_dev_ng_f(&law, &SC_NG_P2); + assert e2 < 0.0000000016; + assert e2 > 0.0000000013; + const d1 : f64 = sc_t_first_f(&law) - SC_NG_FIRST_2; + const d2 : f64 = sc_t_res_f(&law) - SC_NG_RES_2; + assert d1 * d1 < 0.00005 * 0.00005; + assert d2 * d2 < 0.00005 * 0.00005; + // The loser's peak: ngspice 0.747264 V (printed to 6 digits). + const pk : f64 = sc_peak_f(&law, SC_DN) / SC_UPV - SC_NG_PEAK_2; + assert pk * pk < 0.000001 * 0.000001; + sc_run_f(1, &t, SC_H05, false, SC_K_ANA, &SC_DV_M, &law, &res); + const e05 : f64 = sc_dev_ng_f(&law, &SC_NG_P05); + assert e05 < 0.0000000015; + assert e05 > 0.0000000012; + const d3 : f64 = sc_t_first_f(&law) - SC_NG_FIRST_05; + const d4 : f64 = sc_t_res_f(&law) - SC_NG_RES_05; + assert d3 * d3 < 0.00005 * 0.00005; + assert d4 * d4 < 0.00005 * 0.00005; + // The four time differences are the 6 printed digits of the .meas lines (a step is 0.1 ps, 7 digits of ps). + assert d1 * d1 > 0.00001 * 0.00001; + const pk05 : f64 = sc_peak_f(&law, SC_DN) / SC_UPV - SC_NG_PEAK_05; + assert pk05 * pk05 < 0.000001 * 0.000001; +} + +test the_fixed_point_run_against_ngspice { + var t : [625]u32 = [0] ** 625; + sc_build_table(&t, SC_VTH_V); + var a : [1128]u32 = [0] ** 1128; + var b : [1128]u32 = [0] ** 1128; + var c : [1128]u32 = [0] ** 1128; + var d : [1128]u32 = [0] ** 1128; + var res : [140]u32 = [0] ** 140; + const sa : u64 = sc_run(&t, SC_H2, false, SC_K, &SC_DV_M, &a, &res); + const sb : u64 = sc_run(&t, SC_H2, true, SC_K, &SC_DV_M, &b, &res); + const sc : u64 = sc_run(&t, SC_H05, false, SC_K, &SC_DV_M, &c, &res); + const sd : u64 = sc_run(&t, SC_H05, true, SC_K, &SC_DV_M, &d, &res); + // THE DECISION. ngspice LEVEL=1 resolves +2 mV and +0.5 mV with dp high and -2 mV and -0.5 mV with dn high (the decks' + // dpend and dnend); so does the fixed-point run. + assert(sc_sign(&a)); + assert(!sc_sign(&b)); + assert(sc_sign(&c)); + assert(!sc_sign(&d)); + // THE DECISION TIME: the first output through Vdd / 2, in steps after the clock's 50 % point, against the decks' tfirst - tck. + const fa : f64 = (sc_after_clock(sc_t_first(&a)) as f64) / 4096.0; + const fb : f64 = (sc_after_clock(sc_t_first(&b)) as f64) / 4096.0; + const fc : f64 = (sc_after_clock(sc_t_first(&c)) as f64) / 4096.0; + const fd : f64 = (sc_after_clock(sc_t_first(&d)) as f64) / 4096.0; + const ra : f64 = (fa - SC_NG_FIRST_2) / SC_NG_FIRST_2; + const rb : f64 = (fb - SC_NG_FIRST_2) / SC_NG_FIRST_2; + const rc : f64 = (fc - SC_NG_FIRST_05) / SC_NG_FIRST_05; + const rd : f64 = (fd - SC_NG_FIRST_05) / SC_NG_FIRST_05; + // 5.489 ps against 5.496 ps at 2 mV (-0.14 %), 5.558 ps against 5.569 ps at 0.5 mV (-0.20 %): within 0.25 %, from both sides. + assert ra * ra < 0.0025 * 0.0025; + assert rb * rb < 0.0025 * 0.0025; + assert rc * rc < 0.0025 * 0.0025; + assert rd * rd < 0.0025 * 0.0025; + assert rc * rc > 0.0018 * 0.0018; + assert ra * ra > 0.0012 * 0.0012; + // THE RESOLUTION TIME: |dp - dn| through Vdd / 2 against the decks' tres - tck. 6.963 ps against 6.931 ps (+0.47 %) and + // 7.734 ps against 7.694 ps (+0.53 %): within 0.6 %, from both sides. + const ga : f64 = (sc_after_clock(sc_t_res(&a)) as f64) / 4096.0; + const gb : f64 = (sc_after_clock(sc_t_res(&b)) as f64) / 4096.0; + const gc : f64 = (sc_after_clock(sc_t_res(&c)) as f64) / 4096.0; + const gd : f64 = (sc_after_clock(sc_t_res(&d)) as f64) / 4096.0; + const qa : f64 = (ga - SC_NG_RES_2) / SC_NG_RES_2; + const qb : f64 = (gb - SC_NG_RES_2) / SC_NG_RES_2; + const qc : f64 = (gc - SC_NG_RES_05) / SC_NG_RES_05; + const qd : f64 = (gd - SC_NG_RES_05) / SC_NG_RES_05; + assert qa * qa < 0.006 * 0.006; + assert qb * qb < 0.006 * 0.006; + assert qc * qc < 0.006 * 0.006; + assert qd * qd < 0.006 * 0.006; + assert qc * qc > 0.005 * 0.005; + assert qa * qa > 0.004 * 0.004; + // The issue's 5 % is met by a factor 9 at the worst of the four. + assert qc * qc < 0.05 * 0.05 / 81.0; + // THE NODES: the fixed-point states against the 8 x 8 recorded ngspice values, the largest |difference| in volts. This is + // the table's error (15 mV) seen at the checkpoints: 14.66 mV at 2 mV (bound 14.8 mV), 15.40 mV at 0.5 mV (bound 15.6 mV). + const n2 : f64 = sc_dev_ng_u(&a, &SC_NG_P2); + const n05 : f64 = sc_dev_ng_u(&c, &SC_NG_P05); + assert n2 < 0.0148; + assert n2 > 0.0145; + assert n05 < 0.0156; + assert n05 > 0.0152; + // THE REAL TWIN (f64, the same iteration, no rounding) against ngspice: the model's own error with the arithmetic taken out. + // 14.67 mV and 15.42 mV at the recorded steps; -0.137 % and +0.466 % in the two times at 2 mV, -0.202 % and +0.527 % at 0.5 mV. + var tw : [1128]f64 = [0.0] ** 1128; + var rs : [140]f64 = [0.0] ** 140; + sc_run_f(0, &t, SC_H2, false, SC_K, &SC_DV_M, &tw, &rs); + const w2 : f64 = sc_dev_ng_f(&tw, &SC_NG_P2); + const wf2 : f64 = (sc_t_first_f(&tw) - SC_NG_FIRST_2) / SC_NG_FIRST_2; + const wr2 : f64 = (sc_t_res_f(&tw) - SC_NG_RES_2) / SC_NG_RES_2; + assert w2 < 0.0148; + assert w2 > 0.0145; + assert wf2 < 0.0 - 0.0013; + assert wf2 > 0.0 - 0.0014; + assert wr2 > 0.0046; + assert wr2 < 0.0047; + sc_run_f(0, &t, SC_H05, false, SC_K, &SC_DV_M, &tw, &rs); + const w05 : f64 = sc_dev_ng_f(&tw, &SC_NG_P05); + const wf05 : f64 = (sc_t_first_f(&tw) - SC_NG_FIRST_05) / SC_NG_FIRST_05; + const wr05 : f64 = (sc_t_res_f(&tw) - SC_NG_RES_05) / SC_NG_RES_05; + assert w05 < 0.0156; + assert w05 > 0.0152; + assert wf05 < 0.0 - 0.0020; + assert wf05 > 0.0 - 0.0021; + assert wr05 > 0.0052; + assert wr05 < 0.0053; + // THE LOSER'S PEAK: 0.7530 V against ngspice's 0.7473 V at 2 mV (+5.7 mV), and ngspice's loser at 0.5 mV peaks at 0.822 V. + const k2 : f64 = (sc_peak(&a, SC_DN) as f64) / SC_UPV - SC_NG_PEAK_2; + assert k2 > 0.005; + assert k2 < 0.006; + const k05 : f64 = (sc_peak(&c, SC_DN) as f64) / SC_UPV - SC_NG_PEAK_05; + assert k05 * k05 < 0.01 * 0.01; + assert_eq(sa, 858149221550211914); + assert_eq(sb, 7832980476838438025); + assert_eq(sc, 7639752914144809206); + assert_eq(sd, 42403085586441860); +} + +test swapped_inputs_and_a_mismatched_pair_are_the_negative_controls { + var t : [625]u32 = [0] ** 625; + sc_build_table(&t, SC_VTH_V); + var out : [1128]u32 = [0] ** 1128; + var res : [140]u32 = [0] ** 140; + var mm : [15]u32 = [0] ** 15; + // The right answer for a positive input (VDP > VDN) is dp high. The nominal circuit gives it for 2 mV and 0.5 mV. + const s0 : u64 = sc_run(&t, SC_H2, false, SC_K, &SC_DV_M, &out, &res); + assert(sc_sign(&out)); + const s1 : u64 = sc_run(&t, SC_H05, false, SC_K, &SC_DV_M, &out, &res); + assert(sc_sign(&out)); + // CONTROL 1, the inputs wired the wrong way round (the stimulus +2 mV with VDP and VDN exchanged): dn goes high, the + // check "dp high for VDP > VDN" fails, at both sizes. + const s2 : u64 = sc_run(&t, SC_H2, true, SC_K, &SC_DV_M, &out, &res); + assert(!sc_sign(&out)); + const s3 : u64 = sc_run(&t, SC_H05, true, SC_K, &SC_DV_M, &out, &res); + assert(!sc_sign(&out)); + // The pair of equal m = 36 gives the nominal run, bit for bit, through the same table of multipliers. + sc_mult_pair(36, 36, &mm); + const s4 : u64 = sc_run(&t, SC_H2, false, SC_K, &mm, &out, &res); + assert_eq(s4, s0); + // CONTROL 2, the pair's m swapped when they differ: m3 = 38, m4 = 34 against m3 = 34, m4 = 38. The decision follows the + // mismatch, not the input: (38, 34) gives dn high for every stimulus, (34, 38) dp high for every stimulus, and the two are + // opposite for all four stimuli (+-2 mV, +-0.5 mV). A 5.6 % mismatch (35, 37) already outweighs 2 mV. + var hs : [4]u32 = [SC_H2, SC_H2, SC_H05, SC_H05]; + var k : u32 = 0; + while (k < 4) { + sc_mult_pair(38, 34, &mm); + const sx : u64 = sc_run(&t, hs[k as usize], (k % 2) == 1, SC_K, &mm, &out, &res); + const hi_lo : bool = sc_sign(&out); + sc_mult_pair(34, 38, &mm); + const sy : u64 = sc_run(&t, hs[k as usize], (k % 2) == 1, SC_K, &mm, &out, &res); + const lo_hi : bool = sc_sign(&out); + assert(!hi_lo); + assert(lo_hi); + k = k + 1; + } + sc_mult_pair(35, 37, &mm); + const s5 : u64 = sc_run(&t, SC_H2, true, SC_K, &mm, &out, &res); + assert(sc_sign(&out)); + assert s5 != s0; + assert s2 != s0; + assert s3 != s1; +} + +test a_wrong_threshold_table_is_caught { + var good : [625]u32 = [0] ** 625; + var bad : [625]u32 = [0] ** 625; + sc_build_table(&good, SC_VTH_V); + sc_build_table(&bad, SC_VTH_WRONG_V); + // The law check: the right table is within half a unit of the law, the wrong one (Vth = 0.30 V) is 525686 units off at a node. + assert sc_node_err(&good) <= 0.5 + 0.0000001; + assert sc_node_err(&bad) > 525000.0; + assert sc_node_err(&bad) < 526000.0; + // In the circuit: the run on the wrong table still decides the right way (the sign alone would not catch it) ... + var out : [1128]u32 = [0] ** 1128; + var res : [140]u32 = [0] ** 140; + const s : u64 = sc_run(&bad, SC_H2, false, SC_K, &SC_DV_M, &out, &res); + assert(sc_sign(&out)); + // ... but the decision time is 16.7 % early against ngspice (bound 5 %), the resolution 13.7 % early, and the states are + // 0.39 V off at the recorded steps (bound 15 mV). + const f : f64 = (sc_after_clock(sc_t_first(&out)) as f64) / 4096.0; + const r : f64 = (sc_after_clock(sc_t_res(&out)) as f64) / 4096.0; + const rf : f64 = (f - SC_NG_FIRST_2) / SC_NG_FIRST_2; + const rr : f64 = (r - SC_NG_RES_2) / SC_NG_RES_2; + assert rf * rf > 0.05 * 0.05; + assert rr * rr > 0.05 * 0.05; + assert rf < 0.0 - 0.16; + assert rf > 0.0 - 0.17; + assert rr < 0.0 - 0.13; + assert rr > 0.0 - 0.14; + const n : f64 = sc_dev_ng_u(&out, &SC_NG_P2); + assert n > 0.39; + assert n < 0.40; + assert n > 20.0 * 0.0148; + assert s != 858149221550211914; +} + +test metastability_is_reported_not_hidden { + var t : [625]u32 = [0] ** 625; + sc_build_table(&t, SC_VTH_V); + var out : [1128]u32 = [0] ** 1128; + var res : [140]u32 = [0] ** 140; + // The input quartered three times, 2 mV down to 31 uV: every one is decided the right way within 140 steps, and the + // resolution time grows by the same 7.67 steps (0.77 ps) at every quartering, the logarithm of the input. + var hs : [4]u32 = [SC_H2, SC_H05, SC_H0125, SC_H003]; + var prev : u32 = 0; + var k : u32 = 0; + while (k < 4) { + const s : u64 = sc_run(&t, hs[k as usize], false, SC_K, &SC_DV_M, &out, &res); + assert(sc_sign(&out)); + assert sc_max_u(&res, 140) <= SC_RES_BOUND; + const tr : u32 = sc_t_res(&out); + assert tr != SC_NONE; + if (k > 0) { + assert tr - prev > 31300; + assert tr - prev < 31700; + } + prev = tr; + k = k + 1; + } + assert_eq(prev, 412373); + // So one e-fold of the input costs 7.67 / ln 4 = 5.53 steps = 0.553 ps: the regeneration time constant. + const tau : f64 = (((prev - 317977) as f64) / 4096.0) / 3.0 / 1.3862943611198906; + assert tau > 5.5; + assert tau < 5.6; + // At 0 V of input the SAME circuit does not decide. In real arithmetic it stays exactly symmetric: dp = dn at every step + // to a millionth of a unit, and |dp - dn| never reaches Vdd / 2. + var outf : [1128]f64 = [0.0] ** 1128; + var resf : [140]f64 = [0.0] ** 140; + sc_run_f(0, &t, 0, false, SC_K, &SC_DV_M, &outf, &resf); + var j : u32 = 0; + while (j <= SC_N) { + const d : f64 = fabs64(outf[(j * SC_NN + SC_DP) as usize] - outf[(j * SC_NN + SC_DN) as usize]); + assert d < 0.000001; + j = j + 1; + } + // In fixed point the same run is NOT symmetric to the last unit: one unit of rounding at step 60 (|dp - dn| = 1) grows, + // 2027 units at step 100 and 2659643 at step 140 (0.158 V), still short of Vdd / 2: unresolved, the side it leans to + // set by rounding, not by an input. This is the metastable state, reported as it is. + const s0 : u64 = sc_run(&t, 0, false, SC_K, &SC_DV_M, &out, &res); + assert_eq(sc_t_res(&out), SC_NONE); + assert_eq(mag_diff(out[(60 * 8 + 5) as usize], out[(60 * 8 + 6) as usize]), 1); + assert_eq(mag_diff(out[(100 * 8 + 5) as usize], out[(100 * 8 + 6) as usize]), 2027); + assert_eq(mag_diff(out[(140 * 8 + 5) as usize], out[(140 * 8 + 6) as usize]), 2659643); + assert mag_diff(out[(140 * 8 + 5) as usize], out[(140 * 8 + 6) as usize]) < (VDD >> 1); +} + +test the_halves_mirror_exactly_in_real_arithmetic { + var t : [625]u32 = [0] ** 625; + sc_build_table(&t, SC_VTH_V); + var a : [1128]f64 = [0.0] ** 1128; + var b : [1128]f64 = [0.0] ** 1128; + var res : [140]f64 = [0.0] ** 140; + // Exchanging VDP and VDN gives the run with every node replaced by its mirror (ap <-> an, net2 <-> net3, dp <-> dn): to the + // last bit in real arithmetic, for the table and for the law. A device on the wrong side, or a gate on the wrong net, breaks it. + sc_run_f(0, &t, SC_H2, false, SC_K, &SC_DV_M, &a, &res); + sc_run_f(0, &t, SC_H2, true, SC_K, &SC_DV_M, &b, &res); + assert sc_mirror_dev_f(&a, &b) < 0.000000001; + sc_run_f(1, &t, SC_H05, false, SC_K_ANA, &SC_DV_M, &a, &res); + sc_run_f(1, &t, SC_H05, true, SC_K_ANA, &SC_DV_M, &b, &res); + assert sc_mirror_dev_f(&a, &b) < 0.000000001; + // And the check can fail: a run with the input pair's m unequal is not its own mirror. + var mm : [15]u32 = [0] ** 15; + sc_mult_pair(38, 34, &mm); + sc_run_f(0, &t, SC_H2, false, SC_K, &mm, &a, &res); + sc_run_f(0, &t, SC_H2, true, SC_K, &mm, &b, &res); + assert sc_mirror_dev_f(&a, &b) > 1000000.0; +} + +test the_board_step_is_the_k3_step_and_the_words_have_headroom { + var t : [625]u32 = [0] ** 625; + sc_build_table(&t, SC_VTH_V); + // The step a board would run (a static trip count) is sc_step at K = 3, state for state and residual for residual. + var v : [8]u32 = [0] ** 8; + var w : [8]u32 = [0] ** 8; + sc_rest(SC_VCM + SC_H2, SC_VCM - SC_H2, &v); + sc_rest(SC_VCM + SC_H2, SC_VCM - SC_H2, &w); + var j : u32 = 1; + while (j <= 80) { + const ra : u32 = sc_step(&t, &v, j, SC_VCM + SC_H2, SC_VCM - SC_H2, SC_K, &SC_DV_M); + const rb : u32 = sc_board_step(&t, &w, j, SC_VCM + SC_H2, SC_VCM - SC_H2, &SC_DV_M); + assert_eq(ra, rb); + var n : u32 = 0; + while (n < SC_NN) { + assert_eq(v[n as usize], w[n as usize]); + n = n + 1; + } + j = j + 1; + } + // HEADROOM. Every register, matrix word and move of the solve, and the widest sum back substitution can form (8 |a| |x| + // + |b| 2^20), over the whole run at 2 mV, 0.5 mV and 0 V: the widest is 4.95e14, under 2^49. A u64 holds 2^63: a factor + // of more than 10^4 to spare, and the test asks for 2^62. + const w2 : u64 = sc_headroom(&t, SC_H2, false, SC_K, &SC_DV_M); + const w05 : u64 = sc_headroom(&t, SC_H05, false, SC_K, &SC_DV_M); + const w0 : u64 = sc_headroom(&t, 0, false, SC_K, &SC_DV_M); + var wide : u64 = w2; + if (w05 > wide) { wide = w05; } + if (w0 > wide) { wide = w0; } + assert wide < 4611686018427387904; + assert wide < 562949953421312; + assert wide > 281474976710656; + assert_eq(wide, 495357036992680); +}