From 9e03bdab192a7cc634fbf28689c3aced74199a8e Mon Sep 17 00:00:00 2001 From: rajames Date: Sat, 3 Oct 2026 00:17:11 -0400 Subject: [PATCH] feat(v4.0.0): Q./ under D-11, with (UQ/), D2*C and the (Q/) variable D-11 (ruled 2026-10-02): Q./ rounds toward zero, saturates to Q max / Q min on overflow, and on division by zero returns Q max / Q min by the dividend's sign (0 for 0/0) and sets the new NODE-ERROR register (7), which VM-ERROR? reads. - (UQ/): unsigned floor(a * 2^16 / b) for b != 0 and no overflow. Restoring division, 2N+16 steps, on one shifting register (quotient above remainder). Quotient register and divisor live in a 5-cell variable (Q/), like BASE and the hold buffer, so the stacks carry only the remainder, the quotient bit and the loop count. A first version kept the divisor on the return stack and overflowed it when called from Q./ (it never returned). - D2*C: double shift left with carry in and out. - Q./: zero-divisor handling, signs (kept in (Q/)), an overflow test that needs no division (|a| * 2^16 >= |b| * 2^(2N-1), possible only for |b| < 2^17), then (UQ/) and the sign. Checked at 32- and 64-bit cells, optimised and ASan+UBSan, against an independent limb-by-limb long division (including NODE-ERROR): every pair of 18 edge Q values, 13 width-native edge values (true Q max/min at either width), the overflow boundary, divisors in [2^(2N-2), 2^(2N-1)), and pseudo-random cases; and against v3's q48_div for non-negative a < 2^48. Q.* now also runs on the width-native edge values. Mutations of the compare paths, the take path, the error flag and the overflow mask are all caught (the mask matters only for Q min / -1.0). Headroom (data under args / return): Q./ 3/2. Co-Authored-By: Claude Opus 5.5 --- docs/v4.0.0/DECOMPOSITION.md | 11 +- v4/tests/test_foundation.c | 469 ++++++++++++++++++++++++++++++++++- 2 files changed, 475 insertions(+), 5 deletions(-) diff --git a/docs/v4.0.0/DECOMPOSITION.md b/docs/v4.0.0/DECOMPOSITION.md index 1a160c5f..a858a7a5 100644 --- a/docs/v4.0.0/DECOMPOSITION.md +++ b/docs/v4.0.0/DECOMPOSITION.md @@ -147,8 +147,8 @@ IF body1 ELSE body2 THEN → if L1 drop body1 jump L2 times, as on the F18. `FOR ... UNEXT` is the same but the body must fit in one instruction word. **Register conventions.** `A` and `B` are caller-saved. A word that uses them says so. Words in this -document that clobber `A`: `@ ! +! -! 2@ 2! C@ C! UM* * UM/MOD Q.FROM-INT Q.TO-INT Q.* SEND RECV`. Words that clobber `B`: -`SEND RECV`. +document that clobber `A`: `@ ! +! -! 2@ 2! C@ C! UM* * UM/MOD Q.FROM-INT Q.TO-INT Q.* Q./ SEND RECV`. Words that clobber `B`: +`Q./ SEND RECV`. **Return-stack words** (`>R R> R@ 2>R 2R> 2R@ I J UNLOOP` and the loop runtimes) are always IN. @@ -172,6 +172,7 @@ definition below depends on one, it says so. | **D-8** | Q48.16 signedness. | **Signed.** v3's unsigned comparisons and `Q.FROM-INT` clamping are retired. | | **D-9** | Instruction word width on a 64-bit-cell host (ruled 2026-10-02). | **32 bits at every cell width.** Six 5-bit slots plus 2 spare bits, as in §1.2. On a 64-bit host the instruction word is the low 32 bits of the cell and the high half is ignored, so compiled code is identical at both widths. Only data and `@p` literals are a full cell wide. | | **D-10** | Q48.16 width on a 64-bit-cell node (ruled 2026-10-02). | **Two cells at every cell width.** A Q value is a signed double on 32- and 64-bit nodes alike, so every Q word is the same double word at both widths. At 32-bit cells this is bit-for-bit v3's 64-bit Q. At 64-bit cells the low cell is v3's value, and where v3 wraps on Q48.16 overflow (a sum past Q max, `ABS` or `NEG` of Q min) the high cell carries the true result instead. | +| **D-11** | `Q./` semantics: division by zero, rounding, overflow (ruled 2026-10-02). | **Saturate and flag; round toward zero; saturate on overflow.** The quotient of `a * 2^16 / b` is rounded toward zero, like `SM/REM`. When it does not fit a Q value it is clamped to Q max or Q min by its sign. Division by zero returns Q max or Q min by the sign of the dividend (`0 / 0` gives 0) and sets the node's `NODE-ERROR` register (§7), which `VM-ERROR?` reads. v3 returned 0 on division by zero and saturated whenever the dividend was 2^48 or more, even when the quotient would have fit. | **Consequences of D-2 that every definition must respect.** The data stack holds 10 items and the return stack 9, and every `call`, `FOR`, `DO` loop frame and `push` uses return-stack slots. Nesting @@ -657,7 +658,10 @@ word. Per D-10 it occupies two cells on a 64-bit node too. Per D-8, v4 Q values | --- | --- | --- | | `Q.+` `Q.-` | CAP | `D+`, `D-` — executed on the golden model (2026-10-02) against v3's `q48_add`/`q48_sub` (`uint64_t` wrapping): bit-for-bit at 32-bit cells; at 64-bit cells the low cell is v3's, per D-10. | | `Q.*` | CAP | `SWAP push over over UM* drop push push over pop -if L1 SWAP jump L2 L1: SWAP drop 0 L2: pop SWAP - push push over pop UM* pop + ROT pop dup push over -if L3 drop dup jump L4 L3: drop 0 L4: push UM* pop - D+ ROT pop UM* SWAP push 0 D+ pop push over pop SWAP Q.TO-INT push Q.TO-INT pop SWAP` (`SWAP` and `ROT` in line) — `floor(a*b / 2^16)`, signed (D-8), cut to two cells. Cells 0–2 of the product come from the unsigned cell products `a1*b1` (low cell), `a0*b1`, `a1*b0`, `a0*b0`, in that order, each input dropped after its last use and `b0` waiting on the return stack. Reading `a1` and `b1` as signed takes `(b1<0 ? a0 : 0)` and `(a1<0 ? b0 : 0)` off cell 2; each is folded in as soon as its operands are adjacent. The result is cells 0–2 shifted right 16 with `Q.TO-INT` twice. Executed on the golden model (2026-10-02) against an independent limb-by-limb reference, and against v3's `q48_mul` for non-negative operands. Leaves its caller 3 data cells under its arguments and 3 return entries. Clobbers `A`. | -| `Q./` | CAP | Shifted long division. **Division by zero sets an error** (v3 returned 0 silently). | +| `Q./` | CAP | `over over OR if ZERO drop push over pop SWAP over xor (Q/) 4 + b! !b -if L1 DNEGATE L1: push push -if L2 DNEGATE L2: pop pop dup if CHK drop jump DIV CHK: drop over -131072 and if CHK2 drop jump DIV CHK2: drop push push dup pop dup push N-18 FOR 2* UNEXT over over xor -if SAME drop drop -if NOOV jump OV SAME: drop push inv pop + inv -if OV NOOV: drop pop pop jump DIV OV: drop drop drop pop pop drop drop (Q/) 4 + b! @b -if OVP drop 0 MSB ; OVP: drop -1 MAXHI ; DIV: (UQ/) (Q/) 4 + b! @b -if L3 drop DNEGATE ; L3: drop ; ZERO: drop drop drop NODE-ERROR b! -1 !b over over OR if Z0 drop -if ZP drop drop 0 MSB ; ZP: drop drop -1 MAXHI ; Z0: drop ;` — D-11: `a * 2^16 / b`, rounded toward zero; saturated to Q max / Q min on overflow; on division by zero, Q max / Q min by the dividend's sign (0 for `0 / 0`) and `NODE-ERROR` set. The sign of the result waits in `(Q/)`. Overflow needs no division: the quotient reaches `2^(2N-1)` exactly when `|a| * 2^16 >= |b| * 2^(2N-1)`, possible only for `|b| < 2^17`, where it is `a1 >= b0 << (N-17)` unsigned. `MSB` and `MAXHI` are the high cells of Q min and Q max. Executed on the golden model (2026-10-02) against an independent limb-by-limb reference (including `NODE-ERROR`), on width-native edge values, and against v3's `q48_div` for non-negative `a < 2^48`. Leaves its caller 3 data cells under its arguments and 2 return entries. Clobbers `A` and `B`. | +| `(UQ/)` | CAP | Internal to `Q./`: unsigned `floor(a * 2^16 / b)` for `b != 0` and no overflow. Restoring division, `2N+16` steps, on one shifting register — quotient (starting as `a`) above remainder — with the quotient and divisor in `(Q/)` so the stacks carry only the remainder, the quotient bit and the loop count. Each step's quotient bit enters at the next step's shift; with no overflow the bits leaving the quotient in the last 16 steps are zeros. Trial subtraction by in-line sign tests (`x - y` with `y` on top is `push inv pop + inv`). Executed on the golden model (2026-10-02). Clobbers `A` and `B`. | +| `D2*C` | CAP | `( lo hi cin -- lo' hi' cout )`: the double shifted left one bit, `cin` entering at the bottom and `cout` the bit leaving the top (both 0 or 1). `ROT dup -if L1 drop 1 jump L2 L1: drop 0 L2: push 2* + SWAP dup -if L3 drop 1 jump L4 L3: drop 0 L4: push 2* pop pop ROT + SWAP`, `ROT` and `SWAP` in line. Executed on the golden model (2026-10-02). | +| `(Q/)` | CAP | Variable, 5 cells: `Q./`'s quotient register, divisor and result sign. Like `BASE` and the hold buffer, it lives in node memory. | | `Q.ABS` `Q.NEG` | CAP | `DABS`, `DNEGATE` — executed on the golden model (2026-10-02) against v3's `q48_abs` and `0 - q`: bit-for-bit at 32-bit cells, including v3's wrap of Q min to itself; at 64-bit cells the low cell is v3's and the high cell is the true sign, per D-10. | | `Q.LOG` `Q.EXP` `Q.SQRT` `Q.SIN` `Q.COS` | CAP | Algorithms ported from `q48_words.c`; the hosted C versions are the golden model. | | `Q.FROM-INT` | CAP | `push 0 a! 0 pop 15 FOR +* UNEXT push drop a pop` — `n * 2^16` as a signed double. `+*` with `S = 0` never adds, so each step is an exact arithmetic right shift of `T:A`; starting from `T:A = n:0` (that is, `n * 2^N`) and shifting `N-16` bits leaves `n * 2^16`. The count is `N-17`: 15 at 32-bit cells, 47 at 64. Negative values are no longer clamped to 0. Executed on the golden model (2026-10-02): matches v3's `q48_from_u64` for `n >= 0` (at 64-bit cells in the low cell, where v3 wraps for `n >= 2^47`, per D-10). Clobbers `A`. | @@ -765,3 +769,4 @@ Addresses are assigned in the node memory map (D-4). Names only here. | `GOV-*` | R/W | Governor parameters and Jacquard selector state (`L8-*`, decay rate). | | `REC-ENABLE` | R/W | DoE recorder on/off (`HB-ON` / `HB-OFF`). | | `PORT-STATUS` | R | Per-port ready flags, for non-blocking polls. | +| `NODE-ERROR` | R/W | Arithmetic error flag: set to −1 by `Q./` on division by zero (D-11), cleared by writing 0. `VM-ERROR?` reads it. | diff --git a/v4/tests/test_foundation.c b/v4/tests/test_foundation.c index 1bcabfa7..5d9f94bb 100644 --- a/v4/tests/test_foundation.c +++ b/v4/tests/test_foundation.c @@ -23,7 +23,10 @@ * <, =, D<, 2OVER, the call-free (D<), DMAX, DMIN and 2SWAP, and the Q * comparisons on them, are checked against C, signed (D-8), and against v3 * where v3's unsigned Q comparisons agree. Q.* is checked against an independent limb-by-limb - * reference, signed, and against v3's q48_mul for non-negative operands. Every word is + * reference, signed, and against v3's q48_mul for non-negative operands. + * Q./ (D-11) and its core (UQ/) are checked against an independent + * limb-by-limb long division, including NODE-ERROR, on width-native edge + * values, and against v3's q48_div where v4 keeps v3's meaning. Every word is * also probed for how much of the 10- and 9-deep circular stacks (D-2) it * leaves to its caller. * @@ -41,6 +44,11 @@ static int failures = 0, checks = 0; #define CHECK(c,...) do{checks++; if(!(c)){failures++; printf("FAIL %s:%d: ",__FILE__,__LINE__); printf(__VA_ARGS__); printf("\n");}}while(0) #define CANARY ((v4_cell)0x0C0FFEE5) +/* NODE-ERROR (section 7); D-4 leaves its address open, so the test puts it + * in the word next to the halt address. */ +#define NODE_ERROR ((v4_cell)(V4_NODE_WORDS - 2u)) +/* (Q/), Q./'s 5-cell scratch variable: q0 q1 b0 b1 sign. */ +#define QS ((v4_cell)(V4_NODE_WORDS - 8u)) #define MAXU ((v4_ucell)~(v4_ucell)0) #define FLAG(c) ((c) ? V4_ALL_ONES : (v4_cell)0) @@ -54,7 +62,7 @@ static v4_cell w_nip, w_swap, w_or, w_negate, w_rot, w_zless, w_zequal, w_ugreater, w_abs, w_s2d, w_dplus, w_dnegate, w_dabs, w_smrem, w_slashmod, w_mplus, w_dminus, w_d0equal, w_dequal, w_qfromint, w_qtoint, w_less, w_equal, w_dless, w_2swap, w_2over, - w_dmax, w_dmin, w_qgt_doc, w_qgt, w_dltkeep, w_qstar; + w_dmax, w_dmin, w_qgt_doc, w_qgt, w_dltkeep, w_qstar, w_d2starc, w_uqdiv, w_qslash; #define O(name) v4_asm_op(&as, V4_OP_##name) #define LIT(v) v4_asm_lit(&as, (v4_cell)(v)) @@ -593,6 +601,225 @@ static void build(void) O(SEMI); } #undef ROT_INLINE +#undef SWAP_INLINE + + /* : D2*C ( lo hi cin -- lo' hi' cout ) 5.26, for Q./ + * The double shifted left one bit, cin (0 or 1) entering at the bottom, + * cout (0 or 1) the bit that left the top. + * ROT dup -if L1 drop 1 jump L2 L1: drop 0 L2: push hi cin lo R: m + * 2* + SWAP lo' hi + * dup -if L3 drop 1 jump L4 L3: drop 0 L4: push R: m cout + * 2* pop pop ROT + SWAP ; lo' hi' cout + * ROT and SWAP in line. */ +#define SWAP_INLINE() do { O(OVER); O(PUSH); O(PUSH); O(DROP); O(RPOP); O(RPOP); } while (0) +#define ROT_INLINE() do { O(PUSH); SWAP_INLINE(); O(RPOP); SWAP_INLINE(); } while (0) + { + v4_asm_ref l1, l2, l3, l4; + w_d2starc = v4_asm_label(&as); + ROT_INLINE(); O(DUP); + l1 = v4_asm_branch_fwd(&as, V4_OP_MINUS_IF); + O(DROP); LIT(1); + l2 = v4_asm_branch_fwd(&as, V4_OP_JUMP); + v4_asm_resolve(&as, l1, v4_asm_label(&as)); + O(DROP); LIT(0); + v4_asm_resolve(&as, l2, v4_asm_label(&as)); + O(PUSH); O(TWO_STAR); O(ADD); SWAP_INLINE(); + O(DUP); + l3 = v4_asm_branch_fwd(&as, V4_OP_MINUS_IF); + O(DROP); LIT(1); + l4 = v4_asm_branch_fwd(&as, V4_OP_JUMP); + v4_asm_resolve(&as, l3, v4_asm_label(&as)); + O(DROP); LIT(0); + v4_asm_resolve(&as, l4, v4_asm_label(&as)); + O(PUSH); O(TWO_STAR); O(RPOP); O(RPOP); ROT_INLINE(); O(ADD); SWAP_INLINE(); + O(SEMI); + } +#undef ROT_INLINE +#undef SWAP_INLINE + + /* : (UQ/) ( a0 a1 b0 b1 -- q0 q1 ) 5.26, core of Q./ + * Unsigned floor(a * 2^16 / b), for b != 0 and a quotient below + * 2^(2N-1) (Q./ checks both first). Restoring division, 2N+16 steps, + * on one shifting register: quotient q (starts as a) above remainder r. + * The bit leaving q's top enters r; each step's quotient bit enters q at + * the next step's shift, and a final shift takes in the last one. With + * no overflow the bits leaving q in the last 16 steps are zeros. + * q and b live in the scratch variable (Q/) = QS: q0 q1 b0 b1, so the + * stacks hold only r, the quotient bit and the loop count. + * QS 2 + a! SWAP !+ ! QS a! SWAP !+ ! 0 0 0 r0 r1 qb + * 2N+15 FOR + * QS a! @+ @ ROT D2*C r0 r1 q0 q1 c + * push SWAP QS a! !+ ! pop D2*C r0 r1 t + * if CMP drop jump TAKE + * CMP: drop QS 3 + b! @b r0 r1 b1 + * over over xor if EQH -if SAMEH + * drop drop -if NT1 jump TAKE NT1: jump NOTAKE + * SAMEH: drop over SWAP push inv pop + inv r0 r1 r1-b1 + * -if T2 drop jump NOTAKE T2: drop jump TAKE + * EQH: drop drop over QS 2 + b! @b r0 r1 r0 b0 + * over over xor -if SAMEL + * drop drop -if NT3 drop jump TAKEL NT3: drop jump NOTAKEL + * SAMEL: drop push inv pop + inv r0 r1 r0-b0 + * -if T4 drop jump NOTAKEL T4: drop jump TAKEL + * TAKEL: drop QS 2 + b! @b push inv pop + inv 0 1 jump END + * NOTAKEL: 0 jump END + * TAKE: QS 2 + a! @+ @ DNEGATE D+ 1 jump END + * NOTAKE: 0 + * END: + * NEXT + * QS a! @+ @ ROT D2*C drop push push drop drop pop pop ; + * x - y with y on top is push inv pop + inv. SWAP, ROT in line. + * Clobbers A and B. */ +#define SWAP_INLINE() do { O(OVER); O(PUSH); O(PUSH); O(DROP); O(RPOP); O(RPOP); } while (0) +#define ROT_INLINE() do { O(PUSH); SWAP_INLINE(); O(RPOP); SWAP_INLINE(); } while (0) +#define SUB_INLINE() do { O(PUSH); O(INV); O(RPOP); O(ADD); O(INV); } while (0) +#define FWD(op) v4_asm_branch_fwd(&as, V4_OP_##op) +#define HERE_(r) v4_asm_resolve(&as, (r), v4_asm_label(&as)) + { + v4_cell loop, l_take, l_notake, l_takel, l_notakel, l_end; + v4_asm_ref cmp, eqh, sameh, nt1, t2, samel, nt3, t4; + v4_asm_ref j_take[4], j_notake[3], j_takel[3], j_notakel[3], j_end[3]; + unsigned nt = 0, nn = 0, ntl = 0, nnl = 0, ne = 0, k; + + w_uqdiv = v4_asm_label(&as); + LIT(QS + 2); O(BANG_A); SWAP_INLINE(); O(STORE_INC); O(STORE_A); + LIT(QS); O(BANG_A); SWAP_INLINE(); O(STORE_INC); O(STORE_A); + LIT(0); LIT(0); LIT(0); + LIT(2 * V4_CELL_BITS + 15); O(PUSH); + loop = v4_asm_label(&as); + LIT(QS); O(BANG_A); O(FETCH_INC); O(FETCH_A); ROT_INLINE(); CALL(w_d2starc); + O(PUSH); SWAP_INLINE(); LIT(QS); O(BANG_A); O(STORE_INC); O(STORE_A); O(RPOP); + CALL(w_d2starc); + cmp = FWD(IF); + O(DROP); j_take[nt++] = FWD(JUMP); + HERE_(cmp); /* CMP */ + O(DROP); LIT(QS + 3); O(BANG_B); O(FETCH_B); + O(OVER); O(OVER); O(XOR); eqh = FWD(IF); + sameh = FWD(MINUS_IF); + O(DROP); O(DROP); nt1 = FWD(MINUS_IF); + j_take[nt++] = FWD(JUMP); + HERE_(nt1); j_notake[nn++] = FWD(JUMP); + HERE_(sameh); /* SAMEH */ + O(DROP); O(OVER); SWAP_INLINE(); SUB_INLINE(); t2 = FWD(MINUS_IF); + O(DROP); j_notake[nn++] = FWD(JUMP); + HERE_(t2); O(DROP); j_take[nt++] = FWD(JUMP); + HERE_(eqh); /* EQH */ + O(DROP); O(DROP); O(OVER); LIT(QS + 2); O(BANG_B); O(FETCH_B); + O(OVER); O(OVER); O(XOR); samel = FWD(MINUS_IF); + O(DROP); O(DROP); nt3 = FWD(MINUS_IF); + O(DROP); j_takel[ntl++] = FWD(JUMP); + HERE_(nt3); O(DROP); j_notakel[nnl++] = FWD(JUMP); + HERE_(samel); /* SAMEL */ + O(DROP); SUB_INLINE(); t4 = FWD(MINUS_IF); + O(DROP); j_notakel[nnl++] = FWD(JUMP); + HERE_(t4); O(DROP); j_takel[ntl++] = FWD(JUMP); + l_takel = v4_asm_label(&as); /* TAKEL */ + O(DROP); LIT(QS + 2); O(BANG_B); O(FETCH_B); SUB_INLINE(); LIT(0); LIT(1); + j_end[ne++] = FWD(JUMP); + l_notakel = v4_asm_label(&as); /* NOTAKEL */ + LIT(0); j_end[ne++] = FWD(JUMP); + l_take = v4_asm_label(&as); /* TAKE */ + LIT(QS + 2); O(BANG_A); O(FETCH_INC); O(FETCH_A); + CALL(w_dnegate); CALL(w_dplus); LIT(1); j_end[ne++] = FWD(JUMP); + l_notake = v4_asm_label(&as); /* NOTAKE */ + LIT(0); + l_end = v4_asm_label(&as); /* END */ + v4_asm_branch(&as, V4_OP_NEXT, loop); + LIT(QS); O(BANG_A); O(FETCH_INC); O(FETCH_A); ROT_INLINE(); CALL(w_d2starc); O(DROP); + O(PUSH); O(PUSH); O(DROP); O(DROP); O(RPOP); O(RPOP); O(SEMI); + + for (k = 0; k < nt; k++) v4_asm_resolve(&as, j_take[k], l_take); + for (k = 0; k < nn; k++) v4_asm_resolve(&as, j_notake[k], l_notake); + for (k = 0; k < ntl; k++) v4_asm_resolve(&as, j_takel[k], l_takel); + for (k = 0; k < nnl; k++) v4_asm_resolve(&as, j_notakel[k], l_notakel); + for (k = 0; k < ne; k++) v4_asm_resolve(&as, j_end[k], l_end); + } +#undef HERE_ +#undef FWD +#undef SUB_INLINE +#undef ROT_INLINE +#undef SWAP_INLINE + + /* : Q./ ( a b -- q ) 5.26, D-11 + * over over OR if ZERO drop + * push over pop SWAP over xor QS 4 + b! !b a0 a1 b0 b1 (Q/)+4: sign + * -if L1 DNEGATE L1: push push -if L2 DNEGATE L2: pop pop |a| |b| + * dup if CHK drop jump DIV b1 != 0: no overflow + * CHK: drop over -131072 and if CHK2 drop jump DIV b0 >= 2^17: none + * CHK2: drop push push dup pop dup push N-18 FOR 2* UNEXT + * over over xor -if SAME drop drop -if NOOV jump OV + * SAME: drop push inv pop + inv -if OV + * NOOV: drop pop pop jump DIV a0 a1 b0 0 + * OV: drop drop drop pop pop drop drop + * QS 4 + b! @b -if OVP drop 0 MSB ; OVP: drop -1 MAXHI ; + * DIV: (UQ/) QS 4 + b! @b -if L3 drop DNEGATE ; L3: drop ; + * ZERO: drop drop drop NODE-ERROR b! -1 !b + * over over OR if Z0 drop -if ZP drop drop 0 MSB ; + * ZP: drop drop -1 MAXHI ; Z0: drop ; + * MSB and MAXHI are the high cells of Q min and Q max. Clobbers A and B. */ +#define SWAP_INLINE() do { O(OVER); O(PUSH); O(PUSH); O(DROP); O(RPOP); O(RPOP); } while (0) +#define SUB_INLINE() do { O(PUSH); O(INV); O(RPOP); O(ADD); O(INV); } while (0) +#define FWD(op) v4_asm_branch_fwd(&as, V4_OP_##op) +#define HERE_(r) v4_asm_resolve(&as, (r), v4_asm_label(&as)) + { + v4_asm_ref zero, l1, l2, chk, chk2, same, noov1, ov1, ov2, ovp, l3, z0, zp; + v4_asm_ref d1, d2, d3; + v4_cell l_div, l_noov, l_ov; + + w_qslash = v4_asm_label(&as); + O(OVER); O(OVER); CALL(w_or); zero = FWD(IF); + O(DROP); + O(PUSH); O(OVER); O(RPOP); SWAP_INLINE(); O(OVER); O(XOR); + LIT(QS + 4); O(BANG_B); O(STORE_B); + l1 = FWD(MINUS_IF); CALL(w_dnegate); HERE_(l1); + O(PUSH); O(PUSH); + l2 = FWD(MINUS_IF); CALL(w_dnegate); HERE_(l2); + O(RPOP); O(RPOP); + O(DUP); chk = FWD(IF); + O(DROP); d1 = FWD(JUMP); + HERE_(chk); + O(DROP); O(OVER); LIT(-131072); O(AND); chk2 = FWD(IF); + O(DROP); d2 = FWD(JUMP); + HERE_(chk2); + O(DROP); O(PUSH); O(PUSH); O(DUP); O(RPOP); O(DUP); O(PUSH); + LIT(V4_CELL_BITS - 18); O(PUSH); + (void)v4_asm_label(&as); + O(TWO_STAR); O(UNEXT); + O(OVER); O(OVER); O(XOR); same = FWD(MINUS_IF); + O(DROP); O(DROP); noov1 = FWD(MINUS_IF); + ov1 = FWD(JUMP); + HERE_(same); + O(DROP); SUB_INLINE(); ov2 = FWD(MINUS_IF); + l_noov = v4_asm_label(&as); + O(DROP); O(RPOP); O(RPOP); d3 = FWD(JUMP); + l_ov = v4_asm_label(&as); + O(DROP); O(DROP); O(DROP); O(RPOP); O(RPOP); O(DROP); O(DROP); + LIT(QS + 4); O(BANG_B); O(FETCH_B); ovp = FWD(MINUS_IF); + O(DROP); LIT(0); LIT((v4_cell)V4_MSB); O(SEMI); + HERE_(ovp); O(DROP); LIT(-1); LIT((v4_cell)(MAXU >> 1)); O(SEMI); + l_div = v4_asm_label(&as); + CALL(w_uqdiv); LIT(QS + 4); O(BANG_B); O(FETCH_B); l3 = FWD(MINUS_IF); + O(DROP); CALL(w_dnegate); O(SEMI); + HERE_(l3); O(DROP); O(SEMI); + HERE_(zero); + O(DROP); O(DROP); O(DROP); + LIT(NODE_ERROR); O(BANG_B); LIT(-1); O(STORE_B); + O(OVER); O(OVER); CALL(w_or); z0 = FWD(IF); + O(DROP); zp = FWD(MINUS_IF); + O(DROP); O(DROP); LIT(0); LIT((v4_cell)V4_MSB); O(SEMI); + HERE_(zp); O(DROP); O(DROP); LIT(-1); LIT((v4_cell)(MAXU >> 1)); O(SEMI); + HERE_(z0); O(DROP); O(SEMI); + + v4_asm_resolve(&as, d1, l_div); + v4_asm_resolve(&as, d2, l_div); + v4_asm_resolve(&as, d3, l_div); + v4_asm_resolve(&as, noov1, l_noov); + v4_asm_resolve(&as, ov1, l_ov); + v4_asm_resolve(&as, ov2, l_ov); + } +#undef HERE_ +#undef FWD +#undef SUB_INLINE #undef SWAP_INLINE CHECK(v4_asm_ok(&as), "foundation words assemble"); @@ -882,6 +1109,69 @@ static uint64_t v3_q48_mul(uint64_t a, uint64_t b) return (result_hi << 48) | (result_lo >> 16); } +/* Reference Q./ (D-11), by a different method from the FORTH: magnitudes + * in 16-bit limbs, |a| * 2^16 divided by |b| one bit at a time, then the + * sign applied (rounding toward zero) and the result saturated to Q max or + * Q min. Division by zero gives Q max / Q min by the dividend's sign, 0 for + * 0 / 0, and *err = 1. */ +static void qdiv_ref(v4_ucell al, v4_ucell ah, v4_ucell bl, v4_ucell bh, + v4_ucell *rl, v4_ucell *rh, int *err) +{ + uint32_t a[QLIMBS], b[QLIMBS], num[QLIMBS + 1], quo[QLIMBS + 1], rem[QLIMBS + 1]; + int an, bn, neg, big = 0; + unsigned i, k; + + *err = 0; + if (bl == 0 && bh == 0) { + *err = 1; + if (al == 0 && ah == 0) { *rl = 0; *rh = 0; } + else if (ah & V4_MSB) { *rl = 0; *rh = V4_MSB; } + else { *rl = MAXU; *rh = MAXU >> 1; } + return; + } + q_to_limbs(al, ah, a, &an); + q_to_limbs(bl, bh, b, &bn); + neg = an != bn; + num[0] = 0; /* |a| * 2^16 */ + for (i = 0; i < QLIMBS; i++) num[i + 1] = a[i]; + for (i = 0; i <= QLIMBS; i++) { quo[i] = 0; rem[i] = 0; } + for (k = 16 * (QLIMBS + 1); k-- > 0; ) { + uint32_t bit = (num[k / 16] >> (k % 16)) & 1u, c = bit; + int ge = 1; + for (i = 0; i <= QLIMBS; i++) { /* rem = 2 rem + bit */ + uint32_t t = (rem[i] << 1) | c; + rem[i] = t & 0xFFFFu; + c = t >> 16; + } + for (i = QLIMBS + 1; i-- > 0; ) { /* rem >= |b| ? */ + uint32_t bv = i < QLIMBS ? b[i] : 0; + if (rem[i] != bv) { ge = rem[i] > bv; break; } + } + if (ge) { + uint32_t br = 0; + for (i = 0; i <= QLIMBS; i++) { + uint32_t bv = i < QLIMBS ? b[i] : 0; + uint32_t t = rem[i] - bv - br; + br = (t >> 16) & 1u; + rem[i] = t & 0xFFFFu; + } + quo[k / 16] |= 1u << (k % 16); + } + } + if (quo[QLIMBS] || (quo[QLIMBS - 1] & 0x8000u)) big = 1; /* >= 2^(2N-1) */ + if (big) { + if (neg) { *rl = 0; *rh = V4_MSB; } + else { *rl = MAXU; *rh = MAXU >> 1; } + return; + } + *rl = 0; *rh = 0; + for (i = 0; i < QLIMBS / 2; i++) { + *rl |= (v4_ucell)quo[i] << (16 * i); + *rh |= (v4_ucell)quo[QLIMBS / 2 + i] << (16 * i); + } + if (neg) dneg(*rl, *rh, rl, rh); +} + /* Stack headroom. Run `word` with `dfill` marked cells under the canary and * `rfill` marked cells under its return address, and report whether the * results, the canary and every marked cell come back intact. The F18 stacks @@ -1366,6 +1656,174 @@ int main(void) } } + /* (UQ/) against the Q./ reference, for non-negative a and b where the + * quotient fits (|a| * 2^16 < |b| * 2^(2N-1), i.e. not (b < 2^17 and + * a1 >= b0 << (N-17))). */ + { + uint64_t x = 0x94D049BB133111EBu; + unsigned ran = 0; + for (unsigned i = 0; i < 5000; i++) { + v4_ucell al, ah, bl, bh, rl, rh; + int err; + x ^= x << 13; x ^= x >> 7; x ^= x << 17; al = (v4_ucell)x; + x ^= x << 13; x ^= x >> 7; x ^= x << 17; ah = (v4_ucell)x >> 1; + x ^= x << 13; x ^= x >> 7; x ^= x << 17; bl = (v4_ucell)x; + x ^= x << 13; x ^= x >> 7; x ^= x << 17; bh = (v4_ucell)x >> 1; + switch (i & 7u) { + case 1: bh = 0; break; + case 2: bh = 0; bl >>= (unsigned)(x >> 59); break; + case 3: ah >>= (unsigned)(x >> 59) % V4_CELL_BITS; break; + case 4: bh >>= (unsigned)(x >> 58) % V4_CELL_BITS; break; + case 5: ah = 0; bh = 0; break; + case 6: al = bl; ah = bh; break; /* a = b */ + default: break; + } + if (bl == 0 && bh == 0) bl = 1; + if (bh == 0 && (bl >> 17) == 0 && ah >= (bl << (V4_CELL_BITS - 17))) continue; + qdiv_ref(al, ah, bl, bh, &rl, &rh, &err); + ran++; + CHECK(call4(w_uqdiv, (v4_cell)al, (v4_cell)ah, (v4_cell)bl, (v4_cell)bh) + && left2((v4_cell)rl, (v4_cell)rh), "(UQ/) [%u]", i); + } + /* Divisors in [2^(2N-2), 2^(2N-1)), where the remainder's and the + * divisor's top bits can differ. (Q./ passes magnitudes, so b never + * exceeds 2^(2N-1); that one value is covered through Q./.) */ + for (unsigned i = 0; i < 4000; i++) { + v4_ucell al, ah, bl, bh, rl, rh; + int err; + x ^= x << 13; x ^= x >> 7; x ^= x << 17; al = (v4_ucell)x; + x ^= x << 13; x ^= x >> 7; x ^= x << 17; ah = (v4_ucell)x >> 1; + x ^= x << 13; x ^= x >> 7; x ^= x << 17; bl = (v4_ucell)x; + x ^= x << 13; x ^= x >> 7; x ^= x << 17; bh = (v4_ucell)x; + bh = (bh >> 2) | (V4_MSB >> 1); /* b in [2^(2N-2), 2^(2N-1)) */ + if (i & 1u) bh |= (V4_MSB >> 1) - 1u; /* near the top of it */ + qdiv_ref(al, ah, bl, bh, &rl, &rh, &err); + ran++; + CHECK(call4(w_uqdiv, (v4_cell)al, (v4_cell)ah, (v4_cell)bl, (v4_cell)bh) + && left2((v4_cell)rl, (v4_cell)rh), "(UQ/) large b [%u]", i); + } + printf(" (UQ/) cases run: %u\n", ran); + } + /* D2*C against C. */ + for (unsigned i = 0; i < NVEC; i++) + for (unsigned j = 0; j < NVEC; j++) + for (unsigned c = 0; c < 2; c++) { + v4_ucell lo = (v4_ucell)vec[i], hi = (v4_ucell)vec[j]; + v4_cell w[3]; + w[0] = (v4_cell)((lo << 1) | c); + w[1] = (v4_cell)((hi << 1) | (lo >> (V4_CELL_BITS - 1))); + w[2] = (v4_cell)(hi >> (V4_CELL_BITS - 1)); + CHECK(call(w_d2starc, 3, vec[i], vec[j], (v4_cell)c) && leftk(w, 3), + "D2*C [%u,%u,%u]", i, j, c); + } + + /* Q./ against the reference (D-11), and NODE-ERROR set exactly on + * division by zero; against v3's q48_div for non-negative a < 2^48 and + * b != 0 where v3's quotient is below 2^63. */ + { + static const uint64_t qd[] = { + 0, 1, 3, 0x8000u, 0x10000u, 0x18000u, 0x30000u, 0xFFFFu, + (uint64_t)0 - 0x10000u, (uint64_t)0 - 1u, (uint64_t)0 - 0x30000u, + 0xFFFFFFFFu, 0x100000000u, (uint64_t)0 - 0x100000000u, + ((uint64_t)1 << 47) - 1u, (uint64_t)1 << 47, + ((uint64_t)1 << 63) - 1u, (uint64_t)1 << 63 + }; + const unsigned nd = sizeof qd / sizeof qd[0]; + uint64_t x = 0xBF58476D1CE4E5B9u; + for (unsigned i = 0; i < nd * nd + 5000u; i++) { + uint64_t a, b; + v4_cell al, ah, bl, bh; + v4_ucell rl, rh; + int err; + if (i < nd * nd) { a = qd[i / nd]; b = qd[i % nd]; } + else { + x ^= x << 13; x ^= x >> 7; x ^= x << 17; a = x; + x ^= x << 13; x ^= x >> 7; x ^= x << 17; b = x; + if (i & 1u) a = asr64(a, (unsigned)(x >> 58)); + if (i & 2u) b = asr64(b, (unsigned)(x >> 52) & 63u); + if ((i & 12u) == 12u) b = 0; + } + qcells(a, &al, &ah); + qcells(b, &bl, &bh); + qdiv_ref((v4_ucell)al, (v4_ucell)ah, (v4_ucell)bl, (v4_ucell)bh, &rl, &rh, &err); + v4_node_store(&n, NODE_ERROR, 0); + CHECK(call4(w_qslash, al, ah, bl, bh) && left2((v4_cell)rl, (v4_cell)rh), + "Q./ [%u]", i); + CHECK(v4_node_load(&n, NODE_ERROR) == (err ? -1 : 0), "Q./ NODE-ERROR [%u]", i); + if (!(a >> 63) && !(b >> 63) && b != 0 && a < ((uint64_t)1 << 48)) { + uint64_t v3 = (a << 16) / b; + if (!(v3 >> 63)) { + v4_cell vl, vh; + qcells(v3, &vl, &vh); + CHECK(call4(w_qslash, al, ah, bl, bh) && n.ds.s == vl + && (V4_CELL_BITS == 64 || n.ds.t == vh), "Q./ v3 [%u]", i); + } + } + } + } + + /* Q.* and Q./ on edge values built from cells, so that they are the true + * extremes at either cell width: 0, ulp, -ulp, 0.5, 1.0, -1.0, 2.0, Q max, + * Q max - ulp, Q min, Q min + ulp, and the largest Q below 1.0 * 2^(N-16). */ + { + v4_ucell ce[13][2]; + unsigned k = 0, ne; + ce[k][0] = 0; ce[k][1] = 0; k++; + ce[k][0] = 1; ce[k][1] = 0; k++; + ce[k][0] = MAXU; ce[k][1] = MAXU; k++; + ce[k][0] = 0x8000u; ce[k][1] = 0; k++; + ce[k][0] = 0x10000u; ce[k][1] = 0; k++; + ce[k][0] = (v4_ucell)0 - 0x10000u; ce[k][1] = MAXU; k++; + ce[k][0] = 0x20000u; ce[k][1] = 0; k++; + ce[k][0] = MAXU; ce[k][1] = MAXU >> 1; k++; + ce[k][0] = MAXU - 1u; ce[k][1] = MAXU >> 1; k++; + ce[k][0] = 0; ce[k][1] = V4_MSB; k++; + ce[k][0] = 1; ce[k][1] = V4_MSB; k++; + ce[k][0] = 0; ce[k][1] = 1; k++; + ce[k][0] = MAXU; ce[k][1] = 0; k++; + ne = k; + for (unsigned i = 0; i < ne; i++) + for (unsigned j = 0; j < ne; j++) { + v4_ucell rl, rh; + int err; + qmul_ref(ce[i][0], ce[i][1], ce[j][0], ce[j][1], &rl, &rh); + CHECK(call4(w_qstar, (v4_cell)ce[i][0], (v4_cell)ce[i][1], + (v4_cell)ce[j][0], (v4_cell)ce[j][1]) + && left2((v4_cell)rl, (v4_cell)rh), "Q.* cell edge [%u,%u]", i, j); + qdiv_ref(ce[i][0], ce[i][1], ce[j][0], ce[j][1], &rl, &rh, &err); + v4_node_store(&n, NODE_ERROR, 0); + CHECK(call4(w_qslash, (v4_cell)ce[i][0], (v4_cell)ce[i][1], + (v4_cell)ce[j][0], (v4_cell)ce[j][1]) + && left2((v4_cell)rl, (v4_cell)rh) + && v4_node_load(&n, NODE_ERROR) == (err ? -1 : 0), + "Q./ cell edge [%u,%u]", i, j); + } + } + + /* Q./ at the overflow boundary: |a| * 2^16 against |b| * 2^(2N-1), for + * divisors around 2^16 and 2^17, every sign. */ + { + static const v4_ucell bv[] = { 0x10000u, 0x10001u, 0x18000u, 0x1FFFFu, 0x20000u, 0x20001u, 0xFFFFu }; + for (unsigned i = 0; i < sizeof bv / sizeof bv[0]; i++) + for (unsigned j = 0; j < 4; j++) + for (unsigned sg = 0; sg < 4; sg++) { + v4_ucell thr = bv[i] << (V4_CELL_BITS - 17), al, ah, bl = bv[i], bh = 0, rl, rh; + int err; + switch (j) { + case 0: al = 0; ah = thr; break; /* at the threshold */ + case 1: al = MAXU; ah = thr - 1u; break; /* just below */ + case 2: al = 1; ah = thr; break; /* just above */ + default: al = MAXU; ah = thr; break; + } + if (ah & V4_MSB) continue; /* |a| must be a Q value */ + if (sg & 1u) dneg(al, ah, &al, &ah); + if (sg & 2u) dneg(bl, bh, &bl, &bh); + qdiv_ref(al, ah, bl, bh, &rl, &rh, &err); + CHECK(call4(w_qslash, (v4_cell)al, (v4_cell)ah, (v4_cell)bl, (v4_cell)bh) + && left2((v4_cell)rl, (v4_cell)rh), "Q./ boundary [%u,%u,%u]", i, j, sg); + } + } + /* Stack headroom of the two longest definitions. */ { static const v4_cell umstar_args[][4] = { @@ -1436,6 +1894,13 @@ int main(void) }; headroom("Q.*", w_qstar, 4, qstar_args, 5, 2, &dh, &rh); CHECK(dh >= 0 && rh >= 0, "Q.* runs at all"); + static const v4_cell qdiv_args[][4] = { + { 0, 1, 0, 1 }, { 0, 3, 0, 1 }, { 0, -1, 0, 1 }, { 5, 7, 9, 11 }, + { 0, 1, 0, 0 }, { 0, 0x7FFF, 1, 0 }, { -1, 0x7FFF, 0, -3 } + }; + headroom("(UQ/)", w_uqdiv, 4, qdiv_args, 2, 2, &dh, &rh); + headroom("Q./", w_qslash, 4, qdiv_args, 7, 2, &dh, &rh); + CHECK(dh >= 2 && rh >= 2, "Q./ leaves room"); headroom("SM/REM", w_smrem, 3, smrem_args, 5, 2, &dh, &rh); CHECK(rh >= 1, "SM/REM leaves room for /MOD's return address"); headroom("/MOD", w_slashmod, 2, slashmod_args, 4, 2, &dh, &rh);