mirrors/Odin
1
0
Fork 0
mirror of https://github.com/odin-lang/Odin.git synced 2026-08-14 01:34:35 +00:00
Code Issues Packages Projects Releases Wiki Activity

Fix runtime quaternion division

Browse Source 
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
This commit is contained in:
Simon Branch
2026-08-05 20:43:50 -07:00
parent 9caff6370f
commit 8f4f92db8f
1 changed files with 6 additions and 6 deletions
Download Patch File Download Diff File Expand all files Collapse all files

12
base/runtime/internal.odin
View File

@@ -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)
}