diff --git a/FABRIC.md b/FABRIC.md index dd2925e..7ac4539 100644 --- a/FABRIC.md +++ b/FABRIC.md @@ -3728,11 +3728,35 @@ document and committing that amendment as its own item.* > > Not wired into `init.4th`'s boot chain, per plan — Console isn't a fleet VM yet. -- [ ] **4.3.3a — Q48.16 trigonometry.** `Q.SIN`/`Q.COS` (radian input) added to `q48_16.c`, +- [x] **4.3.3a — Q48.16 trigonometry.** `Q.SIN`/`Q.COS` (radian input) added to `q48_16.c`, Taylor series after range-reducing into `[-π, π]` — same pattern as this file's existing `Q.LOG`/`Q.EXP`/`Q.SQRT`, not a new precedent. Raised 2026-08-07 while scoping 4.3.3: needed by 4.3.3b, does not exist anywhere in this codebase today (checked). *Refs:* §27.3. + > **Done, 2026-08-07.** `q48_reduce_angle()` (single integer division on the raw Q48.16 + > representations to strip full `2·PI_Q48` turns, then a bounded fix-up loop) plus + > `q48_sin_approx`/`q48_cos_approx` (Taylor series, terms `n=3,5,7,9,11` / `n=2,4,6,8,10`, + > early exit below 10). `PI_Q48 = 205887`; `TWO_PI_Q48` is *derived* as `2·PI_Q48` rather + > than independently rounded, so the ±π reduction boundary has no seam. + > + > **A real duplication, not previously flagged:** this codebase has *two* independent + > Q48.16 implementations — `src/word_source/q48_16_words.c` (vendored/hosted) and + > `src/starkernel/math/q48_16.c` (kernel-only; the kernel build does not compile the + > former at all). Found the hard way — the hosted build linked fine, the kernel build + > failed with `undefined reference to q48_sin_approx` until the same two functions were + > added to both files (plus both `q48_16.h` headers — `include/q48_16.h` and + > `include/starkernel/q48_16.h`, which also declare the same functions independently). + > Not fixed at the root (de-duplicating the two implementations is a much larger change + > than this item), just navigated correctly — `Q.LOG`/`Q.EXP`/`Q.SQRT` already had this + > same four-file duplication, unremarked until now. + > + > **Verified live on amd64** via serial injection: `sin(0)=0`, `cos(0)=65536` (exact), + > `sin(π/2)=65536`, `cos(π/2)=0` (exact), `sin(-π/2)=-65536` (exact, confirms the odd- + > function sign handling), `sin(π)≈-27` (residual from `PI_Q48` rounding, ~0.04%), + > `cos(π)≈-65656` (Taylor truncation near the interval edge, ~0.18%), and `sin(3π)` reduces + > to the same `-27` as `sin(π)`, confirming range reduction across multiple turns. Both + > hosted and kernel (amd64) builds clean. + - [ ] **4.3.3b — Geometry drawing primitive wordset.** `LINE`, `CIRCLE`, `ARC`, `ELLIPSE` in `capsules/fabric.4th`, built on 4.3.3's `PLOT`/`CART-PLOT` and 4.3.3a's `Q.SIN`/`Q.COS`. Raised 2026-08-07. Q48.16 throughout; resolution-agnostic (48 integer bits comfortably diff --git a/capsules/BLOCK_MAP.md b/capsules/BLOCK_MAP.md index 7b73e60..b17dbed 100644 --- a/capsules/BLOCK_MAP.md +++ b/capsules/BLOCK_MAP.md @@ -1,5 +1,5 @@ # Capsule Block Manifest — Auto-generated - + diff --git a/include/q48_16.h b/include/q48_16.h index d57f1ee..cf33dce 100644 --- a/include/q48_16.h +++ b/include/q48_16.h @@ -246,6 +246,28 @@ q48_16_t q48_exp_approx(q48_16_t q); */ q48_16_t q48_sqrt_approx(q48_16_t q); +/* + * @brief Approximate sin(q) in Q48.16 (integer-only, Taylor series) + * + * Purpose: Console drawing fabric's CIRCLE/ARC/ELLIPSE (FABRIC.md item 4.3.3b) + * Method: Range-reduce into [-pi, pi], then Taylor series + * + * @param q Angle in Q48.16 format (radians, any magnitude) + * @return sin(q) in Q48.16 format + */ +q48_16_t q48_sin_approx(q48_16_t q); + +/* + * @brief Approximate cos(q) in Q48.16 (integer-only, Taylor series) + * + * Purpose: Console drawing fabric's CIRCLE/ARC/ELLIPSE (FABRIC.md item 4.3.3b) + * Method: Range-reduce into [-pi, pi], then Taylor series + * + * @param q Angle in Q48.16 format (radians, any magnitude) + * @return cos(q) in Q48.16 format + */ +q48_16_t q48_cos_approx(q48_16_t q); + /* ============================================================================ * Diagnostic / Testing Utilities * ============================================================================ diff --git a/include/starkernel/q48_16.h b/include/starkernel/q48_16.h index 3c2c4e9..237ea4c 100644 --- a/include/starkernel/q48_16.h +++ b/include/starkernel/q48_16.h @@ -145,5 +145,19 @@ q48_16_t q48_exp_approx(q48_16_t q); */ q48_16_t q48_sqrt_approx(q48_16_t q); +/* + * Approximate sin(q) in Q48.16 (integer-only, Taylor series, radians) + * Input: q in Q48.16 format (any magnitude) + * Output: sin(q) in Q48.16 format + */ +q48_16_t q48_sin_approx(q48_16_t q); + +/* + * Approximate cos(q) in Q48.16 (integer-only, Taylor series, radians) + * Input: q in Q48.16 format (any magnitude) + * Output: cos(q) in Q48.16 format + */ +q48_16_t q48_cos_approx(q48_16_t q); + #endif /* !Q48_16_H */ #endif /* STARKERNEL_Q48_16_H */ diff --git a/lfs/amd64/starforth b/lfs/amd64/starforth index ac22db8..5079f34 100755 Binary files a/lfs/amd64/starforth and b/lfs/amd64/starforth differ diff --git a/src/starkernel/math/q48_16.c b/src/starkernel/math/q48_16.c index 91bd76f..026703e 100644 --- a/src/starkernel/math/q48_16.c +++ b/src/starkernel/math/q48_16.c @@ -258,3 +258,85 @@ q48_16_t q48_sqrt_approx(q48_16_t q) return x; } + +/* ============================================================================ + * Approximation: Sine / Cosine (Taylor Series, Integer-Only) + * ============================================================================ + * + * FABRIC.md item 4.3.3a -- needed by the Console drawing fabric's + * CIRCLE/ARC/ELLIPSE (item 4.3.3b). Radian input. + * + * PI_Q48 = 205887 (pi * 65536, rounded). TWO_PI_Q48 is derived as + * 2 * PI_Q48 rather than independently rounded, so the range-reduction + * boundary at +-PI_Q48 is self-consistent (no seam). + */ + +/* Reduce a signed Q48.16 angle into [-PI_Q48, PI_Q48]. A single integer + * division on the raw Q48.16 representations gives the correct (scale- + * independent) quotient of how many full 2*pi turns to remove, then a + * bounded fix-up loop (at most one or two iterations) handles the + * remainder landing just outside the target interval. */ +static int64_t q48_reduce_angle(int64_t x) +{ + const int64_t PI_Q48 = 205887; /* pi * 65536, rounded */ + const int64_t TWO_PI_Q48 = 2 * PI_Q48; /* derived, not independently rounded */ + + int64_t k = x / TWO_PI_Q48; + x -= k * TWO_PI_Q48; + + while (x > PI_Q48) x -= TWO_PI_Q48; + while (x < -PI_Q48) x += TWO_PI_Q48; + + return x; +} + +q48_16_t q48_sin_approx(q48_16_t q) +{ + int64_t x = q48_reduce_angle((int64_t)q); + + int negative = (x < 0); + if (negative) x = -x; + + q48_16_t xu = (q48_16_t)x; + q48_16_t x2 = q48_mul(xu, xu); + q48_16_t term = xu; + q48_16_t result = xu; + int subtract = 1; + + /* sin(x) = x - x^3/3! + x^5/5! - x^7/7! + x^9/9! - x^11/11! ... */ + for (int n = 3; n <= 11; n += 2) { + term = q48_mul(term, x2); + term = q48_div(term, q48_from_u64((uint64_t)(n * (n - 1)))); + + result = subtract ? q48_sub(result, term) : q48_add(result, term); + subtract = !subtract; + + if (term < 10) break; + } + + return negative ? (q48_16_t)(0ULL - result) : result; +} + +q48_16_t q48_cos_approx(q48_16_t q) +{ + int64_t x = q48_reduce_angle((int64_t)q); + + q48_16_t xu = (q48_16_t)(x < 0 ? -x : x); + q48_16_t x2 = q48_mul(xu, xu); + q48_16_t term = Q48_ONE; + q48_16_t result = term; + int subtract = 1; + + /* cos(x) = 1 - x^2/2! + x^4/4! - x^6/6! + x^8/8! - x^10/10! ... */ + for (int n = 2; n <= 10; n += 2) { + term = q48_mul(term, x2); + term = q48_div(term, q48_from_u64((uint64_t)(n * (n - 1)))); + + result = subtract ? q48_sub(result, term) : q48_add(result, term); + subtract = !subtract; + + if (term < 10) break; + } + + return result; +} diff --git a/src/word_source/q48_16_words.c b/src/word_source/q48_16_words.c index 0e5994a..892520c 100644 --- a/src/word_source/q48_16_words.c +++ b/src/word_source/q48_16_words.c @@ -414,6 +414,117 @@ q48_16_t q48_sqrt_approx(q48_16_t q) return x; } +/* ============================================================================ + * Approximation: Sine / Cosine (Integer-Only, Taylor Series) + * ============================================================================ + * + * FABRIC.md item 4.3.3a -- needed by the Console drawing fabric's + * CIRCLE/ARC/ELLIPSE (item 4.3.3b). Radian input, same Taylor-series + * approach as q48_exp_approx above. + * + * PI_Q48 = 205887 (pi * 65536, rounded). TWO_PI_Q48 is derived as + * 2 * PI_Q48 rather than independently rounded, so the range-reduction + * boundary at +-PI_Q48 is self-consistent (no seam). + */ + +/** + * @brief Reduce a signed Q48.16 angle into [-PI_Q48, PI_Q48]. + * + * A single integer division on the raw Q48.16 representations gives the + * correct (scale-independent) quotient of how many full 2*pi turns to + * remove, then a bounded fix-up loop (at most one or two iterations) + * handles the remainder landing just outside the target interval. + * + * @param x Signed angle in Q48.16 format (any magnitude) + * @return Equivalent angle in Q48.16 format, within [-PI_Q48, PI_Q48] + */ +static int64_t q48_reduce_angle(int64_t x) +{ + const int64_t PI_Q48 = 205887; /* pi * 65536, rounded */ + const int64_t TWO_PI_Q48 = 2 * PI_Q48; /* derived, not independently rounded */ + + int64_t k = x / TWO_PI_Q48; + x -= k * TWO_PI_Q48; + + while (x > PI_Q48) x -= TWO_PI_Q48; + while (x < -PI_Q48) x += TWO_PI_Q48; + + return x; +} + +/** + * @brief Compute an approximation of @c sin(q) in Q48.16 format (radians). + * + * Range-reduces into [-pi, pi] via @c q48_reduce_angle(), then sums the + * Taylor series @c x - x^3/3! + x^5/5! - ... for terms n = 3,5,7,9,11, + * with an early exit when a term falls below 10 (sub-LSB in Q48.16). Sine + * is odd (@c sin(-x) = -sin(x)), so the sign is factored out before the + * series and reapplied to the result. + * + * @param q Angle in Q48.16 format (radians, any magnitude) + * @return @c sin(q) in Q48.16 format + */ +q48_16_t q48_sin_approx(q48_16_t q) +{ + int64_t x = q48_reduce_angle((int64_t)q); + + int negative = (x < 0); + if (negative) x = -x; + + q48_16_t xu = (q48_16_t)x; + q48_16_t x2 = q48_mul(xu, xu); + q48_16_t term = xu; + q48_16_t result = xu; + int subtract = 1; + + for (int n = 3; n <= 11; n += 2) { + term = q48_mul(term, x2); + term = q48_div(term, q48_from_u64((uint64_t)(n * (n - 1)))); + + result = subtract ? q48_sub(result, term) : q48_add(result, term); + subtract = !subtract; + + if (term < 10) break; + } + + return negative ? (q48_16_t)(0ULL - result) : result; +} + +/** + * @brief Compute an approximation of @c cos(q) in Q48.16 format (radians). + * + * Range-reduces into [-pi, pi] via @c q48_reduce_angle(), then sums the + * Taylor series @c 1 - x^2/2! + x^4/4! - ... for terms n = 2,4,6,8,10, + * with an early exit when a term falls below 10. Cosine is even + * (@c cos(-x) = cos(x)), so no sign correction is needed -- @c x^2 already + * absorbs the sign of @c x. + * + * @param q Angle in Q48.16 format (radians, any magnitude) + * @return @c cos(q) in Q48.16 format + */ +q48_16_t q48_cos_approx(q48_16_t q) +{ + int64_t x = q48_reduce_angle((int64_t)q); + + q48_16_t xu = (q48_16_t)(x < 0 ? -x : x); + q48_16_t x2 = q48_mul(xu, xu); + q48_16_t term = 65536; /* 1.0 in Q48.16 */ + q48_16_t result = term; + int subtract = 1; + + for (int n = 2; n <= 10; n += 2) { + term = q48_mul(term, x2); + term = q48_div(term, q48_from_u64((uint64_t)(n * (n - 1)))); + + result = subtract ? q48_sub(result, term) : q48_add(result, term); + subtract = !subtract; + + if (term < 10) break; + } + + return result; +} + /* ============================================================================ * Diagnostic / Testing Utilities * ============================================================================ diff --git a/src/word_source/q48_words.c b/src/word_source/q48_words.c index 3afb8a6..084e73b 100644 --- a/src/word_source/q48_words.c +++ b/src/word_source/q48_words.c @@ -24,6 +24,8 @@ * Q.LOG ( q -- ln_q ) natural log, Newton-Raphson * Q.EXP ( q -- e^q ) Taylor series * Q.SQRT ( q -- sqrt_q ) Newton-Raphson + * Q.SIN ( q -- sin_q ) Taylor series, radians, any magnitude + * Q.COS ( q -- cos_q ) Taylor series, radians, any magnitude * Q.FROM-INT ( n -- q ) integer → Q48.16 (n << 16) * Q.TO-INT ( q -- n ) Q48.16 → integer (q >> 16, truncate) * Q.1 ( -- 65536 ) 1.0 in Q48.16 @@ -115,6 +117,16 @@ static void q48_word_sqrt(VM *vm) q48_push(vm, q48_sqrt_approx(q48_pop(vm))); } +static void q48_word_sin(VM *vm) +{ + q48_push(vm, q48_sin_approx(q48_pop(vm))); +} + +static void q48_word_cos(VM *vm) +{ + q48_push(vm, q48_cos_approx(q48_pop(vm))); +} + /* ── conversions ─────────────────────────────────────────────────────── */ static void q48_word_from_int(VM *vm) @@ -208,6 +220,8 @@ void register_q48_words(VM *vm) register_word(vm, "Q.LOG", q48_word_log); register_word(vm, "Q.EXP", q48_word_exp); register_word(vm, "Q.SQRT", q48_word_sqrt); + register_word(vm, "Q.SIN", q48_word_sin); + register_word(vm, "Q.COS", q48_word_cos); register_word(vm, "Q.FROM-INT", q48_word_from_int); register_word(vm, "Q.TO-INT", q48_word_to_int); register_word(vm, "Q.1", q48_word_one);