admin/LithosAnanake
1
1
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

feat(v4.0.0): Q.* from four UM* cell products

Browse Source 
Q.* is floor(a*b / 2^16), signed (D-8), cut to two cells. Cells 0..2 of
the product come from the unsigned cell products a0*b0, a0*b1, a1*b0 and
the low cell of a1*b1; reading a1 and b1 as signed takes
(a1<0 ? b0 : 0) + (b1<0 ? a0 : 0) off cell 2. The result is cells 0..2
shifted right 16 with Q.TO-INT twice. Clobbers A.

Checked at 32- and 64-bit cells, optimised and ASan+UBSan, on every pair
of 17 edge Q values plus 20000 pseudo-random pairs:
- against an independent reference (16-bit limbs, magnitudes multiplied
  schoolbook, negated, shifted), and
- against v3's q48_mul for non-negative operands.

The first draft applied only a1's sign correction; the reference caught
the missing b1 term (9860 failures at 32-bit cells).

v3 bug found, not fixed: q48_mul's portable branch (used where there is
no __int128) returns (result_hi << 16) | (result_lo >> 16); the high
part must move up 48, so it is wrong whenever the product passes bit
64, e.g. 0.5 * Q max gives 0000ffffffffffff instead of 3fffffffffffffff.
The parity reference here is that branch with the shift corrected,
which matches v3's __int128 branch.

Stack use is tight: Q.* leaves its caller 2 data cells and 1 return
entry, because UM* parks three values on the return stack. Recorded in
5.26; to be improved next.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
This commit is contained in:
rajames
2026-10-02 23:27:57 -04:00
co-authored by Claude Opus 5.5
parent d0accdbb88
commit ab3b708540
2 changed files with 160 additions and 4 deletions
Download Patch File Download Diff File Expand all files Collapse all files
docs/v4.0.0/DECOMPOSITION.md
+2 -2
View File
@@ -147,7 +147,7 @@ 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 SEND RECV`. Words that clobber `B`:
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`.
**Return-stack words** (`>R R> R@ 2>R 2R> 2R@ I J UNLOOP` and the loop runtimes) are always IN.
@@ -647,7 +647,7 @@ word. Per D-10 it occupies two cells on a 64-bit node too. Per D-8, v4 Q values
| Word | Fate | Notes |
| --- | --- | --- |
| `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 | 64×64 product from four `UM*` partial products, shifted right 16. |
| `Q.*` | CAP | `2OVER drop over 0< over and push over UM* pop - push push push over pop UM* drop pop pop push SWAP pop + push push over 0< over and pop pop push SWAP pop SWAP - push push over over UM* pop pop D+ push push push drop pop UM* 0 pop pop D+ push dup push Q.TO-INT pop pop Q.TO-INT` — `floor(a*b / 2^16)`, signed (D-8), cut to two cells. Cells 0–2 of the product come from the four unsigned cell products `a0*b0`, `a0*b1`, `a1*b0`, `a1*b1` (low cell only); reading `a1` and `b1` as signed takes `(a1<0 ? b0 : 0) + (b1<0 ? a0 : 0)` off cell 2; 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. **Stack use is tight:** it leaves its caller 2 data cells under its arguments and 1 return entry, because `UM*` itself parks three values on the return stack; to be improved. Clobbers `A`. |
| `Q./` | CAP | Shifted long division. **Division by zero sets an error** (v3 returned 0 silently). |
| `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. |
v4/tests/test_foundation.c
+158 -2
View File
@@ -22,7 +22,8 @@
* against v3's q48_from_u64 and q48_to_u64 where v4 keeps v3's meaning.
* <, =, 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. Every word is
* 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
* also probed for how much of the 10- and 9-deep circular stacks (D-2) it
* leaves to its caller.
*
@@ -53,7 +54,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_dmax, w_dmin, w_qgt_doc, w_qgt, w_dltkeep, w_qstar;
#define O(name) v4_asm_op(&as, V4_OP_##name)
#define LIT(v) v4_asm_lit(&as, (v4_cell)(v))
@@ -516,6 +517,42 @@ static void build(void)
w_qgt = v4_asm_label(&as);
CALL(w_2swap); CALL(w_dless); O(SEMI);
/* : Q.* ( a b -- c ) 5.26
* a = a0 a1, b = b0 b1. Cells 0..2 of the signed product a*b, then
* shifted right 16. The cell products are unsigned; reading a1 and b1
* as signed takes (a1<0 ? b0 : 0) + (b1<0 ? a0 : 0) off cell 2, the
* only cell above the unsigned product's that the result reaches.
* 2OVER drop over 0< over and push R: b1<0 ? a0 : 0
* over UM* pop - a0 a1 b0 b1 l01 h01 a0*b1, corrected
* push push R: h01 l01
* push over pop UM* drop a0 a1 b0 t3 a1*b1 low
* pop pop push SWAP pop + a0 a1 b0 x0 x1
* push push over 0< over and a0 a1 b0 c R: x1 x0
* pop pop push SWAP pop SWAP - a0 a1 b0 x0 x1-c
* push push over over UM* a0 a1 b0 l10 h10 a1*b0
* pop pop D+ a0 a1 b0 X0 X1
* push push push drop pop a0 b0 R: X1 X0
* UM* 0 pop pop D+ c0 c1 c2
* push dup push Q.TO-INT pop pop Q.TO-INT ;
* SWAP in line. Clobbers A. */
#define SWAP_INLINE() do { O(OVER); O(PUSH); O(PUSH); O(DROP); O(RPOP); O(RPOP); } while (0)
w_qstar = v4_asm_label(&as);
CALL(w_2over); O(DROP);
O(OVER); CALL(w_zless); O(OVER); O(AND); O(PUSH);
O(OVER); CALL(w_umstar); O(RPOP); CALL(w_minus);
O(PUSH); O(PUSH);
O(PUSH); O(OVER); O(RPOP); CALL(w_umstar); O(DROP);
O(RPOP); O(RPOP); O(PUSH); SWAP_INLINE(); O(RPOP); O(ADD);
O(PUSH); O(PUSH); O(OVER); CALL(w_zless); O(OVER); O(AND);
O(RPOP); O(RPOP); O(PUSH); SWAP_INLINE(); O(RPOP); SWAP_INLINE(); CALL(w_minus);
O(PUSH); O(PUSH); O(OVER); O(OVER); CALL(w_umstar);
O(RPOP); O(RPOP); CALL(w_dplus);
O(PUSH); O(PUSH); O(PUSH); O(DROP); O(RPOP);
CALL(w_umstar); LIT(0); O(RPOP); O(RPOP); CALL(w_dplus);
O(PUSH); O(DUP); O(PUSH); CALL(w_qtoint); O(RPOP); O(RPOP); CALL(w_qtoint);
O(SEMI);
#undef SWAP_INLINE
CHECK(v4_asm_ok(&as), "foundation words assemble");
}
@@ -727,6 +764,82 @@ static int64_t qs(uint64_t q)
return (q >> 63) ? -(int64_t)(~q) - 1 : (int64_t)q;
}
/* Reference Q.* for any cell width, by a different method from the FORTH:
* the Q values (two cells each) as signed integers in 16-bit limbs, the
* magnitudes multiplied schoolbook, the product negated if the signs differ,
* then shifted right 16 (one limb) and cut to two cells. Shifting the
* two's-complement product rounds toward minus infinity. */
#define QLIMBS (2 * V4_CELL_BITS / 16)
static void q_to_limbs(v4_ucell lo, v4_ucell hi, uint32_t *l, int *neg)
{
unsigned i;
*neg = (hi & V4_MSB) != 0;
if (*neg) dneg(lo, hi, &lo, &hi);
for (i = 0; i < QLIMBS / 2; i++) {
l[i] = (uint32_t)((lo >> (16 * i)) & 0xFFFFu);
l[i + QLIMBS / 2] = (uint32_t)((hi >> (16 * i)) & 0xFFFFu);
}
}
static void qmul_ref(v4_ucell al, v4_ucell ah, v4_ucell bl, v4_ucell bh,
v4_ucell *rl, v4_ucell *rh)
{
uint32_t a[QLIMBS], b[QLIMBS], p[2 * QLIMBS];
int an, bn;
unsigned i, j;
q_to_limbs(al, ah, a, &an);
q_to_limbs(bl, bh, b, &bn);
for (i = 0; i < 2 * QLIMBS; i++) p[i] = 0;
for (i = 0; i < QLIMBS; i++) {
uint32_t carry = 0;
for (j = 0; j < QLIMBS; j++) {
uint32_t t = a[i] * b[j] + p[i + j] + carry; /* < 2^32 */
p[i + j] = t & 0xFFFFu;
carry = t >> 16;
}
p[i + QLIMBS] = carry;
}
if (an != bn) { /* negate, 2's complement */
uint32_t carry = 1;
for (i = 0; i < 2 * QLIMBS; i++) {
uint32_t t = (~p[i] & 0xFFFFu) + carry;
p[i] = t & 0xFFFFu;
carry = t >> 16;
}
}
*rl = 0; *rh = 0; /* limbs 1 .. QLIMBS */
for (i = 0; i < QLIMBS / 2; i++) {
*rl |= (v4_ucell)p[1 + i] << (16 * i);
*rh |= (v4_ucell)p[1 + QLIMBS / 2 + i] << (16 * i);
}
}
/* v3's q48_mul: unsigned (a * b) >> 16, cut to 64 bits, as its __int128
* branch computes it on 64-bit gcc/clang hosts. This is the portable
* branch of v3/src/word_source/q48_16_words.c with its last line corrected:
* v3 returns (result_hi << 16) | (result_lo >> 16), which is wrong whenever
* the product reaches past bit 64 (the high cell must move up 48 bits, not
* 16). v3 uses that branch only on hosts without __int128. */
static uint64_t v3_q48_mul(uint64_t a, uint64_t b)
{
uint64_t a_hi = a >> 32;
uint64_t a_lo = a & 0xFFFFFFFFULL;
uint64_t b_hi = b >> 32;
uint64_t b_lo = b & 0xFFFFFFFFULL;
uint64_t p_ll = a_lo * b_lo;
uint64_t p_lh = a_lo * b_hi;
uint64_t p_hl = a_hi * b_lo;
uint64_t p_hh = a_hi * b_hi;
uint64_t carry = 0;
uint64_t mid = p_lh + p_hl;
uint64_t result_hi, result_lo;
if (mid < p_lh) carry++;
result_hi = p_hh + (mid >> 32) + (carry << 32);
result_lo = p_ll + ((mid & 0xFFFFFFFFULL) << 32);
if (result_lo < p_ll) result_hi++;
return (result_hi << 48) | (result_lo >> 16);
}
/* 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
@@ -1174,6 +1287,43 @@ int main(void)
printf(" Q.> as written (SWAP D<): %u wrong of %u\n", doc_wrong, doc_runs);
}
/* Q.* against the limb reference (signed, D-8), on edge Q values and
* pseudo-random ones, and against v3's q48_mul for non-negative
* operands, where v3's unsigned product agrees. */
{
static const uint64_t qm[] = {
0, 1, 0x8000u, 0x10000u, 0x18000u, 0x20000u, 0xFFFFu,
(uint64_t)0 - 0x10000u, (uint64_t)0 - 1u, (uint64_t)0 - 0x18000u,
0xFFFFFFFFu, 0x100000000u, 0x7FFFFFFFu, (uint64_t)0 - 0x100000000u,
((uint64_t)1 << 63) - 1u, (uint64_t)1 << 63, (uint64_t)1 << 47
};
const unsigned nm = sizeof qm / sizeof qm[0];
uint64_t x = 0xD1B54A32D192ED03u;
for (unsigned i = 0; i < nm * nm + 20000u; i++) {
uint64_t a, b;
v4_cell al, ah, bl, bh;
v4_ucell rl, rh;
if (i < nm * nm) { a = qm[i / nm]; b = qm[i % nm]; }
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);
}
qcells(a, &al, &ah);
qcells(b, &bl, &bh);
qmul_ref((v4_ucell)al, (v4_ucell)ah, (v4_ucell)bl, (v4_ucell)bh, &rl, &rh);
CHECK(call4(w_qstar, al, ah, bl, bh) && left2((v4_cell)rl, (v4_cell)rh),
"Q.* [%u]", i);
if (!(a >> 63) && !(b >> 63)) {
v4_cell vl, vh;
qcells(v3_q48_mul(a, b), &vl, &vh);
CHECK(call4(w_qstar, al, ah, bl, bh) && n.ds.s == vl
&& (V4_CELL_BITS == 64 || n.ds.t == vh), "Q.* v3 [%u]", i);
}
}
}
/* Stack headroom of the two longest definitions. */
{
static const v4_cell umstar_args[][4] = {
@@ -1238,6 +1388,12 @@ int main(void)
headroom("DMIN", w_dmin, 4, cmp_args, 9, 2, &dh, &rh);
CHECK(dh >= 2 && rh >= 2, "DMIN leaves room");
headroom("Q.> (2SWAP D<)", w_qgt, 4, cmp_args, 9, 1, &dh, &rh);
static const v4_cell qstar_args[][4] = {
{ 0, 1, 0, 1 }, { -1, -1, -1, -1 }, { 0, -1, 0, 1 }, { 5, 7, 9, 11 },
{ (v4_cell)V4_MSB, -1, 3, (v4_cell)V4_MSB }
};
headroom("Q.*", w_qstar, 4, qstar_args, 5, 2, &dh, &rh);
CHECK(dh >= 0 && rh >= 0, "Q.* runs at all");
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);