From 8f4f92db8faf9ded279d9cbaa10bfccdad022145 Mon Sep 17 00:00:00 2001 From: Simon Branch Date: Wed, 5 Aug 2026 20:43:50 -0700 Subject: [PATCH] Fix runtime quaternion division The runtime implementations of quaternion division are incorrect when calculating q/r with nonzero jmag(r) and either nonzero imag(q) or kmag(q). Notably, the inverse 1/r is still calculated correctly because imag(1) = kmag(1) = 0. The mistake can be verified by calculating (a/b)*b which will be very different from both a and a * (1/b) * b. The compile-time constant folding is correct, see exact_value.cpp in function exact_binary_operator_value -> ExactValue_Quaternion -> Token_Quo. quo256 :: proc(q, r: quaternion256) -> quaternion256 { return q/r } a: quaternion256 : 3 + 5i + 7j + 11k b: quaternion256 : 2 + 7i + 3j + 5k c: quaternion256 : a/b fmt.println("comp", abs((c*b) - a)) fmt.println("run ", abs(quo256(a,b)*b - a)) // both values should be very close to zero; // without fix, only the first is --- base/runtime/internal.odin | 12 ++++++------ 1 file changed, 6 insertions(+), 6 deletions(-) diff --git a/base/runtime/internal.odin b/base/runtime/internal.odin index e08d0e01d..a682f7e76 100644 --- a/base/runtime/internal.odin +++ b/base/runtime/internal.odin @@ -1032,9 +1032,9 @@ quo_quaternion64 :: proc "contextless" (q, r: quaternion64) -> quaternion64 { invmag2 := 1.0 / (r0*r0 + r1*r1 + r2*r2 + r3*r3) t0 := (r0*q0 + r1*q1 + r2*q2 + r3*q3) * invmag2 - t1 := (r0*q1 - r1*q0 - r2*q3 - r3*q2) * invmag2 + t1 := (r0*q1 - r1*q0 + r2*q3 - r3*q2) * invmag2 t2 := (r0*q2 - r1*q3 - r2*q0 + r3*q1) * invmag2 - t3 := (r0*q3 + r1*q2 + r2*q1 - r3*q0) * invmag2 + t3 := (r0*q3 + r1*q2 - r2*q1 - r3*q0) * invmag2 return quaternion(w=f16(t0), x=f16(t1), y=f16(t2), z=f16(t3)) } @@ -1046,9 +1046,9 @@ quo_quaternion128 :: proc "contextless" (q, r: quaternion128) -> quaternion128 { invmag2 := 1.0 / (r0*r0 + r1*r1 + r2*r2 + r3*r3) t0 := (r0*q0 + r1*q1 + r2*q2 + r3*q3) * invmag2 - t1 := (r0*q1 - r1*q0 - r2*q3 - r3*q2) * invmag2 + t1 := (r0*q1 - r1*q0 + r2*q3 - r3*q2) * invmag2 t2 := (r0*q2 - r1*q3 - r2*q0 + r3*q1) * invmag2 - t3 := (r0*q3 + r1*q2 + r2*q1 - r3*q0) * invmag2 + t3 := (r0*q3 + r1*q2 - r2*q1 - r3*q0) * invmag2 return quaternion(w=t0, x=t1, y=t2, z=t3) } @@ -1060,9 +1060,9 @@ quo_quaternion256 :: proc "contextless" (q, r: quaternion256) -> quaternion256 { invmag2 := 1.0 / (r0*r0 + r1*r1 + r2*r2 + r3*r3) t0 := (r0*q0 + r1*q1 + r2*q2 + r3*q3) * invmag2 - t1 := (r0*q1 - r1*q0 - r2*q3 - r3*q2) * invmag2 + t1 := (r0*q1 - r1*q0 + r2*q3 - r3*q2) * invmag2 t2 := (r0*q2 - r1*q3 - r2*q0 + r3*q1) * invmag2 - t3 := (r0*q3 + r1*q2 + r2*q1 - r3*q0) * invmag2 + t3 := (r0*q3 + r1*q2 - r2*q1 - r3*q0) * invmag2 return quaternion(w=t0, x=t1, y=t2, z=t3) }