Geometric Algebra 4D ❡

Addition STA Format e₂₃ e₃₁ e₁₂: It turns out in STA Heidelberg teaches to go e₂₃ e₃₁ e₁₂ BECAUSE you can then say "go around X,Y,Z and for 4D to X to W, Y to W, Z to W". Thats sane and i thus adapted to that. But i derived using e₁₂, e₂₃, e₃₁, e₃₄, e₁₄, e₂₄ and swapped the code which is the same.

e₂₃     is orthogonal to    X           e₁₄(!)
e₃₁     is orthogonal to    Y           e₂₄
e₁₂     is orthogonal to    Z           e₃₄

e₂₃ e₁₄ = e₂₃₁₄ = -e₂₁₃₄ = e₁₂₃₄ e₃₁ e₂₄ = e₃₁₂₄ = -e₁₃₂₄ = e₁₂₃₄ e₁₂ e₃₄ = e₁₂₃₄

all have ring 1-2-3-4-1...

e₁ ∧ e₂₃₄ = +e₁₂₃₄      0 swaps
e₂ ∧ e₁₃₄ = -e₁₂₃₄      1 swap
e₃ ∧ e₁₂₄ = +e₁₂₃₄      2 swaps
e₄ ∧ e₁₂₃ = -e₁₂₃₄      3 swaps

all good hypervolumes. These need to alternate similar to how a 2D face needs to be of alternating orientation at the face edges in 3D.

@TBH I THINK this has to do with all 8 cubes needing to be summed 0 volume. To be checked!!

Shorthand for R̃vR ❡

Doing vR then R̃(vR) sounds tedious? It is!

v = v₁e₁ v₂e₂ v₃e₃ v₄e₄

B = R[1]e₂₃ + R[2]e₃₁ + R[3]e₁₂ + R[4]e₁₄ + R[5]e₂₄ + R[6]e₃₄

I = R[7]e₁₂₃₄

s = R[0]

R̃vR = (s - B + I) v (s + B + I)

Now you can group s + I to term lets say J

R̃vR = (J - B) v (J + B)

You can do that because in Addition you can shuffle around, you only need to consider commutativity in multiplication

Next move was new to me

R̃vR = (Jv - Bv) (J + B)

That's distributivity and to the left and its ok because we put the v AFTER those terms

Now lets expand that, then tackle one by one

R̃vR = JvJ + JvB - BvJ - BvB

JvJ

Jv = (s + I)v = sv + Iv

JvJ is thus (sv + Iv) * (s + I) and thats svs + svI + Ivs + IvI

as I anticommutes with vectors putting -vII makes that -vI²

I stumbled here but you can see any odd and grade 4 will anticommute, even will commute if you do it on paper as an exercise. e₁ e₁₂₃₄ or e₂₃ e₁₂₃ or e₁ e₁₂₃ or e₁₂ e₂₃

so

svs + svI + Ivs - vI²

s is a scalar so we can do svI + sIv then svI - svI

leaves s²v - vI² which can be distributed to (s² - I²)v

JvB - BvJ

A mathmatician will see this can be done together as 3 terms are the same so lets do them in one go

JvB - BvJ

    = (s + I)vB - Bv(s + I)

expand

    = svB + IvB - Bvs - BvI

we now want s and I grouped out. Remember s commutes with anything and I and B are both even graded so they commute with each other and s ofc. But v is grade 1 so moving it will flip signs of B or I

So you reorder first and in case you did not know, over addition the terms always commute

    = svB - Bvs + IvB - BvI

    // remember v is grade one and B and I will commute so -BvI = BIv = IBv

    = s(vB - Bv) + I(vB + Bv)

GA Identities ❡

I will admit i know nothing much about all of these tricks but still do all by hand. So i got hinted now these include the standard GA identities:

    vB - Bv = 2(v · B)
    vB + Bv = 2(v ∧ B)

So far I can tell from Doran Chapter 2 you can do two vectors(! grade 1 !) a and b like

ba = b · a + b ∧ a = a · b - a ∧ b

so now when you do ab + ba you get twice the dot product but due to anticommutivity of the wedge the wedge will cancel out

ab + ba = 2 a · b !

and for ab - ba now your dot product comes to 0 while the wedge will double (as its - - so +)

ab - ba = 2 a ∧ b

Wedge ∧ ❡

But we dont deal with grade 1 vectors here where e₁ and e₂ etc would anticommute. B and v are grade 2 and 1 though so when you have v and B now

vB = v · B + v ∧ B

when you reorder these to Bv the sign of the wedge product will not change. e₃ ∧ e₁₂ and e₁₂ ∧ e₃ are the same!

e₁ e₂₃ = e₁₂₃

e₂₃ e₁ = e₂₃₁ = -e₂₁₃ = e₁₂₃

Dot · ❡

However ... the dot product will now sign flip and cancel out :D Why? See that a plane and a vector will not dot like 2 vectors

v · B = e₁ · (e₁ ∧ e₂)

The reason is not explained in Hollmeier and noted its the "hodge adjoint"

but it turns out you can derive the hack of "bringing the same basis close at the paranthesis" with known geometric product rules:

Because the dot product is defined to LOWER any geometric product by the amount these 2 differ in their grade

Aᵣ · Bₛ = ⟨AᵣBₛ⟩(|r-s|)

Imagine B · v

B · v = (e₁ e₂) · e₁

when you dot two vectors you want to measure how much they overlap and 2 orthogonal bases eₘ eₙ will always be 0 and the same eₙ² is 1 in our 4D Cl(4,0). But that's just telling you how much eₙ overlaps with eₘ. When you want to know how much a line eₙ overlaps with an area/plane eₘₙ you will get back the other line that makes that area. But direction matters now. Because (e₁ e₂) · e₁ leaves you -e₂ and e₁ · (e₁ e₂) leaves you e₂!

When you know or think about or test some

Aᵣ · Bₛ = ⟨AᵣBₛ⟩(|r-s|)

you will understand that and how this works.

e3 · e₁ ∧ e₂

thats just e₃ (e₁e₂) and we want all grade 1 parts. what are they?

e3 · e1 · e2 thats a scalar and its 0. none overlap!

e3 ∧ e₁ ∧ e₂ = e₁ ∧ e₂ ∧ e₃

problem is ... these are only grade 0 and 3 - so the result is zero.

e₂₁ · (e₁₂₃)

that is

e₂₁₁₂₃ so thats just e₃ left you can spare the rest

leaves

e₂₁ · (e₁₂₃) = e₃

try the other way

e₁₂ · (e₁₂₃) = e₁₂₁₂₃ = -e₃

another very simple example

e₁ · (e₂ ∧ e₁) = ⟨e₁₂₁ = -e₂₁₁ = -e₂⟩(1) = -e₂

back to R̃vR = JvJ + JvB - BvJ - BvB part JvB - BvJ ❡

vB - Bv = 2(v · B)
vB + Bv = 2(v ∧ B)

so to prove this is correct lets see for vB + Bv

vB = v · B + v ∧ B
Bv = B · v + B ∧ v

we know the wedges commute as these are known grade 2 bivectors and v is grade 1

f.e. e₁ ∧ e₂₃ and e₂₃ ∧ e₁ seem e₁₂₃ and e₂₃₁ thats 2 swaps. They commute

and the dot products give

v · B that will be some grade 1 vector left

B · v that will be the same grade 1 but negated

so it makes sense we get that twice on vB - Bv and none on vB + Bv

and as B ∧ v == v ∧ B it makes sense they go twice on vB + Bv and 0 on vB - Bv

JvB - BvJ == s(vB - Bv) + I(vB + Bv) ❡

so we are at

(s² - I²)v + (s(vB - Bv) + I(vB + Bv)) - BvB

-BvB ❡

Ok so as departure is getting closer ill short this one a bit

we know

vB - Bv = 2(v · B)

so

vB = 2(v · B) + Bv

• 2(v · B)

  vB - 2(v · B) = Bv

flip

Bv = vB - 2(v · B)

Now add B to the right

BvB = vB² - 2(v · B)B

Minus

-BvB = - vB² + 2(v · B) B

-BvB = 2(v · B)B - vB²

(s² - I²)v + s(vB - Bv) + I(vB + Bv) + 2(v · B)B - vB² ❡

(s² - I²)v  + s(vB - Bv) + I(vB + Bv) + 2(v · B)B  - vB²

puh

Remember

vB - Bv = 2(v · B)
vB + Bv = 2(v ∧ B)

So

(s² - I²)v  + 2s(v · B) + 2I(v ∧ B) + 2(v · B)B  - vB²

2(v · B)B - vB² === -BvB ❡

• 2(v · B)B will give 1-vectors (v · B) but with ∧ B also when eₙ ∧ eₘₐ is lin. independent a trivector

  • Grade 1: (v · B) · B
  • Grade 3: (v · B) ∧ B

Now this trivector will cancel out with whatever comes from vB²

- B² has 6 quadvectors in its results. e₁₂ * e₂₃ = e₁₃ but e₁₂ * e₃₄ = e₁₂₃₄

- (B₁₂² + B₂₃² + B₃₁² + B₃₄² + B₁₄² + B₂₄²) + 2(B₁₂B₃₄ + B₂₃B₁₄ + B₃₁B₂₄)I

when you multiply this with a vector this will return trivectors

Trivectors from 2(v · B)B

(v · B)

(v₁e₁ + v₂e₂ + v₃e₃ + v₄e₄) · (B₂₃e₂₃ + B₃₁e₃₁ + B₁₂e₁₂ + B₁₄e₁₄ + B₂₄e₂₄ + B₃₄e₃₄)

=

    v₁(0 - B₃₁e₃ + B₁₂e₂ + B₁₄e₄ + 0 + 0) +
    v₂(B₂₃e₃ + 0 - B₁₂e₁ + 0 + B₂₄e₄ + 0) +
    v₃(-B₂₃e₂ + B₃₁e₁ + 0 + 0 + 0 + B₃₄e₄) +
    v₄(0 + 0 + 0 - B₁₄e₁ - B₂₄e₂ - B₃₄e₃) +

--------------

=   v₁(-B₃₁e₃ + B₁₂e₂ + B₁₄e₄)
    v₂(+B₂₃e₃ - B₁₂e₁ + B₂₄e₄)
    v₃(-B₂₃e₂ + B₃₁e₁ + B₃₄e₄)
    v₄(-B₁₄e₁ - B₂₄e₂ - B₃₄e₃)

(v · B)B

(all the v · B) (B₂₃e₂₃ + B₃₁e₃₁ + B₁₂e₁₂ + B₁₄e₁₄ + B₂₄e₂₄ + B₃₄e₃₄)

Dot grade 1:

=   v₁( B₃₁B₂₃e₂ - B₃₁²e₁ + 0 + 0 + 0 - B₃₁B₃₄e₄ +
        B₁₂B₂₃e₃ + 0 - B₁₂²e₁ + 0 + B₁₂B₂₄e₄ + 0 +
        0 + 0 + 0 - B₁₄²e₁ - B₁₄B₂₄e₂ - B₁₄B₃₄e₃
    ) +
    v₂(-B₂₃²e₂ + B₂₃B₃₁e₁ + 0 + 0 + 0 + B₂₃B₃₄e₄ +
        0 + B₁₂B₃₁e₃ - B₁₂²e₂ - B₁₂B₁₄e₄ + 0 + 0 +
        0 + 0 + 0 - B₂₄B₁₄e₁ - B₂₄²e₂ - B₂₄B₃₄e₃
    ) +
    v₃(-B₂₃²e₃ + 0 + B₂₃B₁₂e₁ + 0 - B₂₃B₂₄e₄ + 0 +
        0 - B₃₁²e₃ + B₃₁B₁₂e₂ + B₃₁B₁₄e₄ + 0 + 0 +
        0 + 0 + 0 - B₃₄B₁₄e₁ - B₃₄B₂₄e₂ - B₃₄²e₃
    ) +
    v₄( 0 + B₁₄B₃₁e₃ - B₁₄B₁₂e₂ - B₁₄²e₄ + 0 + 0
       -B₂₄B₂₃e₃ + 0 + B₁₂B₂₄e₁ + 0 - B₂₄²e₄ + 0
        B₃₄B₂₃e₂ - B₃₄B₃₁e₁ + 0 + 0 + 0 - B₃₄²e₄
    )

=   v₁(+ B₃₁B₂₃e₂ - B₃₁²e₁ - B₃₁B₃₄e₄
       + B₁₂B₂₃e₃ - B₁₂²e₁ + B₁₂B₂₄e₄
       - B₁₄²e₁ - B₁₄B₂₄e₂ - B₁₄B₃₄e₃
    ) +
    v₂(- B₂₃²e₂ + B₂₃B₃₁e₁ + B₂₃B₃₄e₄
       + B₁₂B₃₁e₃ - B₁₂²e₂ - B₁₂B₁₄e₄
       - B₂₄B₁₄e₁ - B₂₄²e₂ - B₂₄B₃₄e₃
    ) +
    v₃(- B₂₃²e₃ + B₂₃B₁₂e₁ - B₂₃B₂₄e₄
       - B₃₁²e₃ + B₃₁B₁₂e₂ + B₃₁B₁₄e₄
       - B₃₄B₁₄e₁ - B₃₄B₂₄e₂ - B₃₄²e₃
    ) +
    v₄(+ B₁₄B₃₁e₃ - B₁₄B₁₂e₂ - B₁₄²e₄
       - B₂₄B₂₃e₃ + B₁₂B₂₄e₁ - B₂₄²e₄
       + B₃₄B₂₃e₂ - B₃₄B₃₁e₁ - B₃₄²e₄
    )

Wedge grade 3 - there are only grade 3 here as sharing an axis in a wedge is an immediate 0

=   v₁( 0 + 0 - B₃₁B₁₂e₁₂₃ + B₃₁B₁₄e₁₃₄ + B₃₁B₂₄e₂₃₄ + 0 +
        0 + B₁₂B₃₁e₁₂₃ + 0 - B₁₂B₁₄e₁₂₄ + 0 + B₁₂B₃₄e₂₃₄ +
        B₁₄B₂₃e₂₃₄ - B₁₄B₃₁e₁₃₄ + B₁₄B₁₂e₁₂₄ + 0 + 0 + 0
    ) +
    v₂( 0 + 0 + B₂₃B₁₂e₁₂₃ - B₂₃B₁₄e₁₃₄ - B₂₃B₂₄e₂₃₄ + 0 +
       -B₁₂B₂₃e₁₂₃ + 0 + 0 + 0 - B₁₂B₂₄e₁₂₄ - B₁₂B₃₄e₁₃₄ +
        B₂₄B₂₃e₂₃₄ - B₂₄B₃₁e₁₃₄ + B₂₄B₁₂e₁₂₄ + 0 + 0 + 0
    ) + 
    v₃( 0 - B₂₃B₃₁e₁₂₃ + 0 + B₂₃B₁₄e₁₂₄ + 0 - B₂₃B₃₄e₂₃₄ +
        B₃₁B₂₃e₁₂₃ + 0 + 0 + 0 + B₃₁B₂₄e₁₂₄ + B₃₁B₃₄e₁₃₄ +
        B₃₄B₂₃e₂₃₄ - B₃₄B₃₁e₁₃₄ + B₃₄B₁₂e₁₂₄ + 0 + 0 + 0
    ) +
    v₄(-B₁₄B₂₃e₁₂₃ + 0 + 0 + 0 - B₁₄B₂₄e₁₂₄ - B₁₄B₃₄e₁₃₄ +
        0 - B₂₄B₃₁e₁₂₃ + 0 + B₂₄B₁₄e₁₂₄ + 0 - B₂₄B₃₄e₂₃₄ +
        0 + 0 - B₃₄B₁₂e₁₂₃ + B₃₄B₁₄e₁₃₄ + B₃₄B₂₄e₂₃₄ + 0
    )

=   v₁(+B₃₁B₂₄e₂₃₄ + B₁₂B₃₄e₂₃₄ + B₁₄B₂₃e₂₃₄) +
    v₂(-B₂₃B₁₄e₁₃₄ - B₁₂B₃₄e₁₃₄ - B₂₄B₃₁e₁₃₄) + 
    v₃(+B₂₃B₁₄e₁₂₄ + B₃₁B₂₄e₁₂₄ + B₃₄B₁₂e₁₂₄) +
    v₄(-B₁₄B₂₃e₁₂₃ - B₂₄B₃₁e₁₂₃ - B₃₄B₁₂e₁₂₃)

So now as you can see these trivectors all share the same 3 factors B₁₂B₃₄, B₂₃B₁₄, B₃₁B₂₄ so we can just call these like F

v₁(+F e₂₃₄)
v₂(-F e₁₃₄)
v₃(+F e₁₂₄)
v₄(-F e₁₂₃)

and as we need to know 2(v · B) ∧ B multiply by 2:

2 v₁(+F e₂₃₄)
2 v₂(-F e₁₃₄)
2 v₃(+F e₁₂₄)
2 v₄(-F e₁₂₃)

Trivectors from vB²

(v₁e₁ + v₂e₂ + v₃e₃ + v₄e₄) * ((B₁₂² + B₂₃² + B₃₁² + B₃₄² + B₁₄² + B₂₄²) + 2(B₁₂B₃₄ + B₂₃B₁₄ + B₃₁B₂₄)e₁₂₃₄)

we can ignore all the grade 0 components like B₁₂² for now they are just v₁B₁₂²e₁ etc. all grade 1

(v₁e₁ + v₂e₂ + v₃e₃ + v₄e₄) * (2 B₁₂B₃₄ e₁₂₃₄ + 2 B₂₃B₁₄ e₁₂₃₄ + 2 B₃₁B₂₄ e₁₂₃₄)

now that's rather easy as you know a wedge is immediately 0 when an axis is shared so these are all wedged to 0

Dot is left to get grade 3 trivectors

(v₁e₁ + v₂e₂ + v₃e₃ + v₄e₄) · (2 B₁₂B₃₄ e₁₂₃₄ + 2 B₂₃B₁₄ e₁₂₃₄ + 2 B₃₁B₂₄ e₁₂₃₄)

= 
    v₁(2 B₁₂B₃₄ e₂₃₄ + 2 B₁₂B₁₄ e₂₃₄ + 2 B₃₁B₂₄ e₂₃₄) +
    v₂(-2 B₁₂B₃₄ e₁₃₄ - 2 B₂₃B₁₄ e₁₃₄ - 2 B₃₁B₂₄ e₁₃₄) +
    v₃(2 B₁₂B₃₄ e₁₂₄ + 2 B₂₃B₁₄ e₁₂₄ + 2 B₃₁B₂₄ e₁₂₄) +
    v₄(-2 B₁₂B₃₄ e₁₂₃ - 2 B₂₃B₁₄ e₁₂₃ - 2 B₃₁B₂₄ e₁₂₃)

now that wasnt hard, its just 0,1,2,3 swaps.

cleanup with the same F = B₁₂B₃₄, B₂₃B₁₄, B₃₁B₂₄

= 
    2 v₁(+F e₂₃₄) +
    2 v₂(-F e₁₃₄) +
    2 v₃(+F e₁₂₄) +
    2 v₄(-F e₁₂₃)

Tadaaaa................ they all cancel out!

So the quadvector from B² will cancel out by going to trivectors with v while the trivectors that cancel these come from wedging grade 1 with grade 2 so trivectors again in the first part 2(v · B)B

Means -BvB will always be a vector / grade 1 so we could do either

⟨-BvB⟩₁ = ⟨2(v · B)B⟩₁ - v⟨B²⟩₀

This means take only grade 1 (from 1 + 3) then minus vector * grade 0 (vector again grade 1) from grade 1 + 3

=> leaves only grade 1 vector part

also correct is to use · not × as this will only return ther grade 1 dot product (see above)

⟨-BvB⟩₁ = (2(v · B) · B) - v⟨B²⟩₀

Bonus Level for -BvB ❡

Now ofc you may not like the pain of tracking all trivectors. So why do they disappear?

We know that a rotor reversed will flip all planes. So when

R = s + B + I

R̃ = s - B + I

Similar each multivector in general will have rules in which grade a sign flips:

Grade 0 reversed is Grade 0
Grade 1 reversed is Grade 1
Grade 2 reversed is -Grade 2 because you need to swap once

    e₁e₂ = -e₂e₁

Grade 3 reversed is -Grade 3 because it takes 3 swaps

    e₁e₂e₃ = -e₂e₁e₃ = e₂e₃e₁ = -e₃e₂e₁ !

Grade 4 reversed is Grade 4 because it takes 6 swaps

            1        2       3        4       5        6
    e₁₂₃₄ = -e₁₂₄₃ = e₁₄₂₃ = -e₄₁₂₃ = e₄₁₃₂ = -e₄₃₁₂ = e₄₃₂₁

General Rule: (-1) ^ (k (k-1)) / 2 (this works as oddeven is always even)

For B this means B̃ = -B because they are all grade 2!

Now to check whether this sandwich would change under reversion

BvB = (BvB)~ ?

(BvB)~ = B̃ ṽ B̃

(BvB)~ = (-1)B (+1)v (-1)B

Thus you get -1² and 1 which are scalars and commute:

(BvB)~ = (-1)(+1)(-1) BvB = BvB

Tadaa

Because this means we can reverse the whole sandwich, parting it into its grades will hold reversion too:

(BvB)~ = (⟨BvB⟩₁ + ⟨BvB⟩₃)~

Remember grade 1 is 1 reversed, but grade 3 is -1 reversed

(BvB)~ = ⟨BvB⟩₁~ + ⟨BvB⟩₃~

so

⟨BvB⟩₁~ = ⟨BvB⟩₁
⟨BvB⟩₃~ = - ⟨BvB⟩₃

Now we already showed BvB == (BvB)~ !!

thus

BvB = (BvB)~

can be replaced with ⟨BvB⟩₁ + ⟨BvB⟩₃ and ⟨BvB⟩₁~ + ⟨BvB⟩₃~

and we know

⟨BvB⟩₁~ + ⟨BvB⟩₃~ = ⟨BvB⟩₁ - ⟨BvB⟩₃

So hold on to your seat because that means

⟨BvB⟩₁ + ⟨BvB⟩₃ = ⟨BvB⟩₁ - ⟨BvB⟩₃

you can already see it but for completeness

subtract ⟨BvB⟩₁

⟨BvB⟩₃ = - ⟨BvB⟩₃

add ⟨BvB⟩₃

2 ⟨BvB⟩₃ = 0

The only solution is to make ⟨BvB⟩₃ == 0

- vB² ❡

Caution here. when you square B that will always introduce a minus as all planes square to -1 and we only want the scalar parts then.

(Just means instead of (B₂₃² - B₃₁² - B₁₂² - B₁₄² - B₂₄² - B₃₄²) we go -(... + ... + ...))

So

- v⟨B²⟩₀ == -v (- |B|²) == |B|²v

you cannot do that in STA though as B² there is a mix of plus and minus planes. Which is why i will stay with ⟨B²⟩₀

As we later want to know the planes we can

(s² - I²)v - vB² == v(s² - I² - ⟨B²⟩₀)

B² is thus -(B₂₃² + B₃₁² + B₁₂² + B₁₄² + B₂₄² + B₃₄²) + 2(B₂₃B₁₄ + B₃₁B₂₄ + B₁₂B₃₄)I (grade 4) + all other mixed terms cancel out

- Pairs sharing an axis anticommute and cancel (e₁₂e₂₃ + e₂₃e₁₂ = 0)
- orthogonal pairs commute and form the quadvector (2(B₂₃B₁₄ + B₃₁B₂₄ + B₁₂B₃₄)I)

For more fun on this part see the GA STA Doc in 705-GA-STA

Some do

From the formula to the algo ❡

Ok so have this now as a usable vector multiplication function Rotor4.rotateVector4(v)

Easy parts

s = R[0]

p = R[7]

Bivector parts

B₁₂ = R[1]
B₂₃ = R[2]
B₃₁ = R[3]
B₃₄ = R[4]
B₁₄ = R[5]
B₂₄ = R[6]

Squares makes it easier and faster

s_2  = s²
xy_2 = B₁₂²
yz_2 = B₂₃²
zx_2 = B₃₁²
zw_2 = B₃₄²
xw_2 = B₁₄²
yw_2 = B₂₄²
p_2  = I²

v = e₁ ❡

We set v = e₁, start on v(s² - I² - ⟨B²⟩₀)

v'  = (s² - p² + B₂₃² + B₃₁² + B₁₂² + B₁₄² + B₂₄² + B₃₄²) e₁ 
    + ...

dot all bivectors with e₁ for 2s(v · B). gives

e₁ · B₁₂e₁₂ = +B₁₂e₂
e₁ · B₃₁e₃₁ = -B₃₁e₃
e₁ · B₁₄e₁₄ = +B₁₄e₄

The other 3 bivectors/planes have no e₁ in them, thus will be 0

v'  = (s² - p² + B₂₃² + B₃₁² + B₁₂² + B₁₄² + B₂₄² + B₃₄²) e₁ 
    + 2 * s * (B₁₂e₂ - B₃₁e₃ + B₁₄e₄)
    + ...

wedge all bivectors with e₁ for 2I(v ∧ B) AND multiplying by I from the left !! gives

e₁₂₃₄ * (e₁ ∧ B₂₃e₂₃) = e₁₂₃₄ * B₂₃e₁₂₃ = B₂₃ e₁₂₃₄₁₂₃ = -e₁₂₃₁₄₂₃ = e₁₂₁₃₄₂₃ = -e₂₃₄₂₃ = e₂₃₄₃₂ = -e₂₄₂ =          B₂₃e₄
e₁₂₃₄ * (e₁ ∧ B₂₄e₂₄) = e₁₂₃₄ * B₂₄e₁₂₄ = B₂₄ e₁₂₃₄₁₂₄ = -e₁₂₃₁₄₂₄ = e₁₂₃₁₂ = -e₁₂₃₂₁ = e₁₃₁ =                     -B₂₄e₃
e₁₂₃₄ * (e₁ ∧ B₃₄e₃₄) = e₁₂₃₄ * B₃₄e₁₃₄ = B₃₄ e₁₂₃₄₁₃₄ = -e₁₂₃₁₄₃₄ = e₁₂₃₁₃ = -e₁₂₁ =                               B₃₄e₂

We still need to multiply by the pseudoscalars scalar component which ill just use p so you know its p times I

The other 3 planes/bivectors have e₁ in them, thus the wedge will be 0.

v'  = (s² - p² + B₂₃² + B₃₁² + B₁₂² + B₁₄² + B₂₄² + B₃₄²) e₁ 
    + 2 * s * (B₁₂e₂ - B₃₁e₃ + B₁₄e₄)
    + 2 * p * (B₂₃e₄ - B₂₄e₃ + B₃₄e₂)
    + ...

Last term is 2(v · B) · B

We already know v · B is B₁₂e₂ - B₃₁e₃ + B₁₄e₄

Now this is a vector and we know its components e₂ e₃ e₄ will have dot 0 with all planes that are not that basis vector

another dot of this vector with B gives

(B₁₂e₂ - B₃₁e₃ + B₁₄e₄) · (B₂₃e₂₃ + B₃₁e₃₁ + B₁₂e₁₂ + B₁₄e₁₄ + B₂₄e₂₄ + B₃₄e₃₄)

= -B₁₂²e₁ + B₁₂B₂₃e₃ + 0 + 0 + 0 + B₁₂B₂₄e₄
+ 0 + B₃₁B₂₃e₂ - B₃₁²e₁ -B₃₁B₃₄e₄ - 0 - 0
+ 0 + 0 + 0 - B₁₄B₃₄e₃ - B₁₄²e₁ - B₁₄B₂₄e₂

= -B₁₂²e₁ + B₁₂B₂₃e₃ + B₁₂B₂₄e₄ + B₃₁B₂₃e₂ - B₃₁²e₁ -B₃₁B₃₄e₄ - B₁₄B₃₄e₃ - B₁₄²e₁ - B₁₄B₂₄e₂

Now don't forget it's 2(v · B) · B so 2 *

    2 * (-B₁₂²e₁ + B₁₂B₂₃e₃ + B₁₂B₂₄e₄ + B₃₁B₂₃e₂ - B₃₁²e₁ -B₃₁B₃₄e₄ - B₁₄B₃₄e₃ - B₁₄²e₁ - B₁₄B₂₄e₂)

Almost there ❡

v'  = (s² - p² + B₂₃² + B₃₁² + B₁₂² + B₁₄² + B₂₄² + B₃₄²) e₁ 
    + 2 * s * (B₁₂e₂ - B₃₁e₃ + B₁₄e₄)
    + 2 * p * (B₂₃e₄ - B₂₄e₃ + B₃₄e₂)
    + 2 * (-B₁₂²e₁ + B₁₂B₂₃e₃ + B₁₂B₂₄e₄ + B₃₁B₂₃e₂ - B₃₁²e₁ -B₃₁B₃₄e₄ - B₁₄B₃₄e₃ - B₁₄²e₁ - B₁₄B₂₄e₂)

This would already work ... but we can simplify this ofc

e₁ Parts

    (s² - p² + B₂₃² + B₃₁² + B₁₂² + B₁₄² + B₂₄² + B₃₄²) e₁
    2 * (-B₁₂²e₁ - B₃₁²e₁ - B₁₄²e₁)

    = (s² - p² - B₁₂² + B₂₃² - B₃₁² + B₃₄² - B₁₄² + B₂₄²) e₁

e₂ Parts

    2 * s * B₁₂e₂
    2 * p * B₃₄e₂
    2 * (B₃₁B₂₃e₂ - B₁₄B₂₄e₂)

    = 2 * (s * B₁₂ + p * B₃₄ + B₃₁B₂₃ - B₁₄B₂₄) e₂

e₃ Parts

    2 * s * -B₃₁e₃
    2 * p * -B₂₄e₃
    2 * (B₁₂B₂₃e₃ - B₁₄B₃₄e₃)

    = 2 * (s * -B₃₁ + p * -B₂₄ + B₁₂B₂₃ - B₁₄B₃₄) e₃

e₄ Parts

    2 * s * B₁₄e₄
    2 * p * B₂₃e₄
    2 * (B₁₂B₂₄e₄ - B₃₁B₃₄e₄)

    = 2 (s * B₁₄ + p * B₂₃ + B₁₂B₂₄ - B₃₁B₃₄) e₄

not bad!

Full v = e₁ ❡

= (s² - p² + B₂₃² + B₃₁² + B₁₂² + B₁₄² + B₂₄² + B₃₄²) e₁
+ 2 * (s * B₁₂ + p * B₃₄ + B₃₁B₂₃ - B₁₄B₂₄) e₂
+ 2 * (s * -B₃₁ + p * -B₂₄ + B₁₂B₂₃ - B₁₄B₃₄) e₃
+ 2 (s * B₁₄ + p * B₂₃ + B₁₂B₂₄ - B₃₁B₃₄) e₄

Ok so now i have to admit thats actually only the COLUMNS so this is not part v₂, v₃, etc. but actually go there so that's one part of each!

v'₁ = (s² - p² + B₂₃² + B₃₁² + B₁₂² + B₁₄² + B₂₄² + B₃₄²) * v₁ + ... v₂ + ... v₃ + ... v₄
v'₂ = 2 * (s * B₁₂ + p * B₃₄ + B₃₁B₂₃ - B₁₄B₂₄) * v₁ + ... v₂ + ... v₃ + ... v₄
v'₃ = 2 * (s * -B₃₁ + p * -B₂₄ + B₁₂B₂₃ - B₁₄B₃₄) * v₁ + ... v₂ + ... v₃ + ... v₄
v'₄ = 2 (s * B₁₄ + p * B₂₃ + B₁₂B₂₄ - B₃₁B₃₄) * v₁ + ... v₂ + ... v₃ + ... v₄

v = e₂ ❡

v(s² - I² + B²) ❡

v'  = (s² - p² + B₁₂² + B₂₃² + B₃₁² + B₃₄² + B₁₄² + B₂₄²) e₂
    + ...

2s(v · B) ❡

e₂ · B₁₂e₁₂ = -B₁₂e₁
e₂ · B₂₃e₂₃ = B₂₃e₃
e₂ · B₂₄e₂₄ = B₂₄e₄

=> 2 * s * (-B₁₂e₁ + B₂₃e₃ + B₂₄e₄)

2I(v ∧ B) ❡

e₁₂₃₄ * (e₂ ∧ B₃₁e₃₁) = ... e₁₂₃₄₂₃₁ = -e₁₂₃₂₄₃₁ = e₁₃₄₃₁ = -e₁₄₁ =         B₃₁e₄
e₁₂₃₄ * (e₂ ∧ B₃₄e₃₄) = ... e₁₂₃₄₂₃₄ = -e₁₂₃₂₄₃₄ = e₁₂₃₂₃ =                 -B₃₄e₁ 
e₁₂₃₄ * (e₂ ∧ B₁₄e₁₄) = ... e₁₂₃₄₂₁₄ = -e₁₂₃₂₄₁₄ = e₁₂₃₂₁ = -e₁₃₁ =         B₁₄e₃

=> 2 * p * (B₃₁e₄ - B₃₄e₁ + B₁₄e₃)

2(v · B) · B ❡

(-B₁₂e₁ + B₂₃e₃ + B₂₄e₄) · (B₂₃e₂₃ + B₃₁e₃₁ + B₁₂e₁₂ + B₁₄e₁₄ + B₂₄e₂₄ + B₃₄e₃₄)

= -B₁₂²e₂ + 0 + B₁₂B₃₁e₃ + 0 - B₁₂B₁₄e₄ + 0
+ 0 - B₂₃²e₂ + B₂₃B₃₁e₁ + B₂₃B₃₄e₄ + 0 + 0
+ 0 + 0 + 0 - B₂₄B₃₄e₃ - B₂₄B₁₄e₁ - B₂₄²e₂

=> 2 * (-B₁₂²e₂ + B₁₂B₃₁e₃ - B₁₂B₁₄e₄ - B₂₃²e₂ + B₂₃B₃₁e₁ + B₂₃B₃₄e₄ - B₂₄B₃₄e₃ - B₂₄B₁₄e₁ - B₂₄²e₂)

Simplify ❡

v'  = (s² - p² + B₁₂² + B₂₃² + B₃₁² + B₃₄² + B₁₄² + B₂₄²) e₂
    + 2 * s * (-B₁₂e₁ + B₂₃e₃ + B₂₄e₄)
    + 2 * p * (B₃₁e₄ - B₃₄e₁ + B₁₄e₃)
    + 2 * (-B₁₂²e₂ + B₁₂B₃₁e₃ - B₁₂B₁₄e₄ - B₂₃²e₂ + B₂₃B₃₁e₁ + B₂₃B₃₄e₄ - B₂₄B₃₄e₃ - B₂₄B₁₄e₁ - B₂₄²e₂)

e₂ Parts

(s² - p² - B₁₂² - B₂₃² + B₃₁² + B₃₄² + B₁₄² - B₂₄²) e₂

e₁ Parts

2 * (s * -B₁₂ + p * -B₃₄ + B₂₃B₃₁ - B₂₄B₁₄) e₁

e₃ Parts

2 * (s * B₂₃ + p * B₁₄ + B₁₂B₃₁ - B₂₄B₃₄) e₃

e₄ Parts

2 * (s * B₂₄ + p * B₃₁ - B₁₂B₁₄ + B₂₃B₃₄) e₄

Combine with v' e₁ ❡

v'₁ = (s² - p² + B₂₃² + B₃₁² + B₁₂² + B₁₄² + B₂₄² + B₃₄²) * v₁
    + 2 * (s * -B₁₂ + p * -B₃₄ + B₂₃B₃₁ - B₂₄B₁₄) * v₂
    +
    +

v'₂ = 2 * (s * B₁₂ + p * B₃₄ + B₃₁B₂₃ - B₁₄B₂₄) * v₁
    + (s² - p² - B₁₂² - B₂₃² + B₃₁² + B₃₄² + B₁₄² - B₂₄²) * v₂
    + 
    +

v'₃ = 2 * (s * -B₃₁ + p * -B₂₄ + B₁₂B₂₃ - B₁₄B₃₄) * v₁
    + 2 * (s * B₂₃ + p * B₁₄ + B₁₂B₃₁ - B₂₄B₃₄) * v₂
    +
    +

v'₄ = 2 (s * B₁₄ + p * B₂₃ + B₁₂B₂₄ - B₃₁B₃₄) * v₁
    + 2 * (s * B₂₄ + p * B₃₁ - B₁₂B₁₄ + B₂₃B₃₄) * v₂
    +
    +

v = e₃ ❡

v(s² - I² + B²) ❡

v'  = (s² - p² + B₂₃² + B₃₁² + B₁₂² + B₁₄² + B₂₄² + B₃₄²) e₃
    + ...

2s(v · B) ❡

e₃ · B₂₃e₂₃ = -B₂₃e₂
e₃ · B₃₁e₃₁ = B₃₁e₁
e₃ · B₃₄e₃₄ = B₃₄e₄

=> 2 * s * (-B₂₃e₂ + B₃₁e₁ + B₃₄e₄)

2I(v ∧ B) ❡

e₁₂₃₄ * (e₃ ∧ B₁₂e₁₂) = e₁₂₃₄₃₁₂ = -e₁₂₄₁₂ = e₁₂₁₄₂ = -e₂₄₂ =               B₁₂e₄
e₁₂₃₄ * (e₃ ∧ B₁₄e₁₄) = e₁₂₃₄₃₁₄ = -e₁₂₄₁₄ = e₁₂₁ =                         -B₁₄e₂
e₁₂₃₄ * (e₃ ∧ B₂₄e₂₄) = e₁₂₃₄₃₂₄ = -e₁₂₄₂₄ =                                B₂₄e₁

=> 2 * p * (B₁₂e₄ - B₁₄e₂ + B₂₄e₁)

2(v · B) · B ❡

(-B₂₃e₂ + B₃₁e₁ + B₃₄e₄) · (B₂₃e₂₃ + B₃₁e₃₁ + B₁₂e₁₂ + B₁₄e₁₄ + B₂₄e₂₄ + B₃₄e₃₄)

= B₂₃B₁₂e₁ - B₂₃²e₃ + 0 + 0 + 0 - B₂₃B₂₄e₄
+ B₃₁B₁₂e₂ + 0 - B₃₁²e₃ + 0 + B₃₁B₁₄e₄ + 0
+ 0 + 0 + 0 - B₃₄²e₃ - B₃₄B₁₄e₁ - B₃₄B₂₄e₂

=> 2 * (B₂₃B₁₂e₁ - B₂₃²e₃ - B₂₃B₂₄e₄ + B₃₁B₁₂e₂ - B₃₁²e₃ + B₃₁B₁₄e₄ - B₃₄²e₃ - B₃₄B₁₄e₁ - B₃₄B₂₄e₂)

Simplify ❡

v'  = (s² - p² + B₂₃² + B₃₁² + B₁₂² + B₁₄² + B₂₄² + B₃₄²) e₃
    + 2 * s * (-B₂₃e₂ + B₃₁e₁ + B₃₄e₄)
    + 2 * p * (B₁₂e₄ - B₁₄e₂ + B₂₄e₁)
    + 2 * (B₂₃B₁₂e₁ - B₂₃²e₃ - B₂₃B₂₄e₄ + B₃₁B₁₂e₂ - B₃₁²e₃ + B₃₁B₁₄e₄ - B₃₄²e₃ - B₃₄B₁₄e₁ - B₃₄B₂₄e₂)

e₃ Parts

    (s² - p² + B₁₂² - B₂₃² - B₃₁² - B₃₄² + B₁₄² + B₂₄²) e₃

e₁ Parts

    2 * (s * B₃₁ + p * B₂₄ + B₂₃B₁₂ - B₃₄B₁₄) e₁

e₂ Parts

    2 * (s * -B₂₃ + p * -B₁₄ + B₃₁B₁₂ - B₃₄B₂₄) e₂

e₄ Parts

    2 * (s * B₃₄ + p * B₁₂ - B₂₃B₂₄ + B₃₁B₁₄) e₄

Combine with v' e₁ e₂ ❡

v'₁ = (s² - p² + B₂₃² + B₃₁² + B₁₂² + B₁₄² + B₂₄² + B₃₄²) * v₁
    + 2 * (s * -B₁₂ + p * -B₃₄ + B₂₃B₃₁ - B₂₄B₁₄) * v₂
    + 2 * (s * B₃₁ + p * B₂₄ + B₂₃B₁₂ - B₃₄B₁₄) * v₃
    +

v'₂ = 2 * (s * B₁₂ + p * B₃₄ + B₃₁B₂₃ - B₁₄B₂₄) * v₁
    + (s² - p² - B₁₂² - B₂₃² + B₃₁² + B₃₄² + B₁₄² - B₂₄²) * v₂
    + 2 * (s * -B₂₃ + p * -B₁₄ + B₃₁B₁₂ - B₃₄B₂₄) * v₃
    +

v'₃ = 2 * (s * -B₃₁ + p * -B₂₄ + B₁₂B₂₃ - B₁₄B₃₄) * v₁
    + 2 * (s * B₂₃ + p * B₁₄ + B₁₂B₃₁ - B₂₄B₃₄) * v₂
    + (s² - p² + B₁₂² - B₂₃² - B₃₁² - B₃₄² + B₁₄² + B₂₄²) * v₃
    +

v'₄ = 2 (s * B₁₄ + p * B₂₃ + B₁₂B₂₄ - B₃₁B₃₄) * v₁
    + 2 * (s * B₂₄ + p * B₃₁ - B₁₂B₁₄ + B₂₃B₃₄) * v₂
    + 2 * (s * B₃₄ + p * B₁₂ - B₂₃B₂₄ + B₃₁B₁₄) * v₃
    +

v = e₄ ❡

v(s² - I² + B²) ❡

v'  = (s² - p² + B₂₃² + B₃₁² + B₁₂² + B₁₄² + B₂₄² + B₃₄²) e₄
    + ...

2s(v · B) ❡

e₄ · B₃₄e₃₄ = -B₃₄e₃
e₄ · B₁₄e₁₄ = -B₁₄e₁
e₄ · B₂₄e₂₄ = -B₂₄e₂

=> 2 * s * (-B₃₄e₃ -B₁₄e₁ -B₂₄e₂)

2I(v ∧ B) ❡

e₁₂₃₄ * (e₄ ∧ B₁₂e₁₂) = e₁₂₃₁₂ = -e₁₂₁₃₂ = e₂₃₂ =                           -B₁₂e₃
e₁₂₃₄ * (e₄ ∧ B₂₃e₂₃) = e₁₂₃₂₃ =                                            -B₂₃e₁
e₁₂₃₄ * (e₄ ∧ B₃₁e₃₁) = e₁₂₁ =                                              -B₃₁e₂

=> 2 * p * (-B₁₂e₃ - B₂₃e₁ - B₃₁e₂)

2(v · B) · B ❡

(-B₃₄e₃ -B₁₄e₁ -B₂₄e₂) · (B₂₃e₂₃ + B₃₁e₃₁ + B₁₂e₁₂ + B₁₄e₁₄ + B₂₄e₂₄ + B₃₄e₃₄)

= 0 + B₃₄B₂₃e₂ - B₃₄B₃₁e₁ - B₃₄²e₄ + 0 + 0
+ -B₁₄B₁₂e₂ + 0 + B₁₄B₃₁e₃ + 0 - B₁₄²e₄ + 0
+ B₂₄B₁₂e₁ - B₂₄B₂₃e₃ + 0 + 0 + 0 - B₂₄²e₄

=> 2 * (B₃₄B₂₃e₂ - B₃₄B₃₁e₁ - B₃₄²e₄ - B₁₄B₁₂e₂ + B₁₄B₃₁e₃ - B₁₄²e₄ + B₂₄B₁₂e₁ - B₂₄B₂₃e₃ - B₂₄²e₄)

Simplify ❡

v'  = (s² - p² + B₂₃² + B₃₁² + B₁₂² + B₁₄² + B₂₄² + B₃₄²) e₄
    + 2 * s * (-B₃₄e₃ -B₁₄e₁ -B₂₄e₂)
    + 2 * p * (-B₁₂e₃ - B₂₃e₁ - B₃₁e₂)
    + 2 * (B₃₄B₂₃e₂ - B₃₄B₃₁e₁ - B₃₄²e₄ - B₁₄B₁₂e₂ + B₁₄B₃₁e₃ - B₁₄²e₄ + B₂₄B₁₂e₁ - B₂₄B₂₃e₃ - B₂₄²e₄)

e₄ Parts

(s² - p² + B₁₂² + B₂₃² + B₃₁² - B₃₄² - B₁₄² - B₂₄²) e₄

e₁ Parts

2 * (s * -B₁₄ + p * -B₂₃ -B₃₄B₃₁ + B₂₄B₁₂) e₁

e₂ Parts

2 * (s * -B₂₄ + p * -B₃₁ + B₃₄B₂₃ - B₁₄B₁₂) e₂

e₃ Parts

2 * (s * -B₃₄ + p * -B₁₂ + B₁₄B₃₁ - B₂₄B₂₃) e₃

Combine with v' e₁ e₂ e₃ ❡

v'₁ = (s² - p² - B₁₂² + B₂₃² - B₃₁² + B₃₄² - B₁₄² + B₂₄²) * v₁
    + 2 * (s * -B₁₂ + p * -B₃₄ + B₂₃B₃₁ - B₂₄B₁₄) * v₂
    + 2 * (s * B₃₁ + p * B₂₄ + B₂₃B₁₂ - B₃₄B₁₄) * v₃
    + 2 * (s * -B₁₄ + p * -B₂₃ -B₃₄B₃₁ + B₂₄B₁₂) * v₄

v'₂ = 2 * (s * B₁₂ + p * B₃₄ + B₃₁B₂₃ - B₁₄B₂₄) * v₁
    + (s² - p² - B₁₂² - B₂₃² + B₃₁² + B₃₄² + B₁₄² - B₂₄²) * v₂
    + 2 * (s * -B₂₃ + p * -B₁₄ + B₃₁B₁₂ - B₃₄B₂₄) * v₃
    + 2 * (s * -B₂₄ + p * -B₃₁ + B₃₄B₂₃ - B₁₄B₁₂) * v₄

v'₃ = 2 * (s * -B₃₁ + p * -B₂₄ + B₁₂B₂₃ - B₁₄B₃₄) * v₁
    + 2 * (s * B₂₃ + p * B₁₄ + B₁₂B₃₁ - B₂₄B₃₄) * v₂
    + (s² - p² + B₁₂² - B₂₃² - B₃₁² - B₃₄² + B₁₄² + B₂₄²) * v₃
    + 2 * (s * -B₃₄ + p * -B₁₂ + B₁₄B₃₁ - B₂₄B₂₃) * v₄

v'₄ = 2 * (s * B₁₄ + p * B₂₃ + B₁₂B₂₄ - B₃₁B₃₄) * v₁
    + 2 * (s * B₂₄ + p * B₃₁ - B₁₂B₁₄ + B₂₃B₃₄) * v₂
    + 2 * (s * B₃₄ + p * B₁₂ - B₂₃B₂₄ + B₃₁B₁₄) * v₃
    + (s² - p² + B₁₂² + B₂₃² + B₃₁² - B₃₄² - B₁₄² - B₂₄²) * v₄

Code (old e₁₂!) ❡

static rotateVector4HB(r: Rot4, v: Vec4): Vec4 {
    const s = r[0], xy = r[1], yz = r[2], zx = r[3],
        zw = r[4], xw = r[5], yw = r[6], p = r[7];

    // Helper squares
    const s_2 = s * s;
    const xy_2 = xy * xy, yz_2 = yz * yz, zx_2 = zx * zx;
    const zw_2 = zw * zw, xw_2 = xw * xw, yw_2 = yw * yw;
    const p_2 = p * p;

    const out = new Float32Array(4);

    // (s² - I²)v  + 2s(v · B) + 2I(v ∧ B) + 2(v · B)B  - vB²

    // ((())) looks stupid but makes it WGSL-safe (no PEMDAS)

    out[0] = ((s_2 - p_2 - xy_2 + yz_2 - zx_2 + zw_2 - xw_2 + yw_2) * v[0])
        + (2 * ((s * -xy) + (p * -zw) + (yz * zx) - (yw * xw)) * v[1])
        + (2 * ((s * zx) + (p * yw) + (yz * xy) - (zw * xw)) * v[2])
        + (2 * ((s * -xw) + (p * -yz) - (zw * zx) + (yw * xy)) * v[3]);

    out[1] = (2 * ((s * xy) + (p * zw) + (zx * yz) - (xw * yw)) * v[0])
        + ((s_2 - p_2 - xy_2 - yz_2 + zx_2 + zw_2 + xw_2 - yw_2) * v[1])
        + (2 * ((s * -yz) + (p * -xw) + (zx * xy) - (zw * yw)) * v[2])
        + (2 * ((s * -yw) + (p * -zx) + (zw * yz) - (xw * xy)) * v[3]);

    out[2] = (2 * ((s * -zx) + (p * -yw) + (xy * yz) - (xw * zw)) * v[0])
        + (2 * ((s * yz) + (p * xw) + (xy * zx) - (yw * zw)) * v[1])
        + ((s_2 - p_2 + xy_2 - yz_2 - zx_2 - zw_2 + xw_2 + yw_2) * v[2])
        + (2 * ((s * -zw) + (p * -xy) + (xw * zx) - (yw * yz)) * v[3]);

    out[3] = (2 * ((s * xw) + (p * yz) + (xy * yw) - (zx * zw)) * v[0])
        + (2 * ((s * yw) + (p * zx) - (xy * xw) + (yz * zw)) * v[1])
        + (2 * ((s * zw) + (p * xy) - (yz * yw) + (zx * xw)) * v[2])
        + ((s_2 - p_2 + xy_2 + yz_2 + zx_2 - zw_2 - xw_2 - yw_2) * v[3]);

    return out;
}

To export this as a matrix thats really easy, just remember you want COLUMN mayor.

So 0-3 are from v₁ to v₄ onto v'₁

static toMat4HB(r: Rot4): Mat4 {
    const s = r[0], xy = r[1], yz = r[2], zx = r[3],
        zw = r[4], xw = r[5], yw = r[6], p = r[7];

    // Helper squares
    const s_2 = s * s;
    const xy_2 = xy * xy, yz_2 = yz * yz, zx_2 = zx * zx;
    const zw_2 = zw * zw, xw_2 = xw * xw, yw_2 = yw * yw;
    const p_2 = p * p;

    const out = MAT4_IDENTITY.slice();

    // all vₙ parts IN TO OUT v'₁ (0-3) to v'₄ (12-15) - COLUMN MAYOR!!

    // COLUMN 0 (Basis X: v=[1, 0, 0, 0])  
    out[0] = (s_2 - p_2 - xy_2 + yz_2 - zx_2 + zw_2 - xw_2 + yw_2);
    out[1] = 2 * ((s * xy) + (p * zw) + (zx * yz) - (xw * yw));
    out[2] = 2 * ((s * -zx) + (p * -yw) + (xy * yz) - (xw * zw));
    out[3] = 2 * ((s * xw) + (p * yz) + (xy * yw) - (zx * zw));

    // COLUMN 1 (Basis Y: v=[0, 1, 0, 0])
    out[4] = 2 * ((s * -xy) + (p * -zw) + (yz * zx) - (yw * xw));
    out[5] = (s_2 - p_2 - xy_2 - yz_2 + zx_2 + zw_2 + xw_2 - yw_2);
    out[6] = 2 * ((s * yz) + (p * xw) + (xy * zx) - (yw * zw));
    out[7] = 2 * ((s * yw) + (p * zx) - (xy * xw) + (yz * zw));

    // COLUMN 2 (Basis Z: v=[0, 0, 1, 0])
    out[8] = 2 * ((s * zx) + (p * yw) + (yz * xy) - (zw * xw));
    out[9] = 2 * ((s * -yz) + (p * -xw) + (zx * xy) - (zw * yw));
    out[10] = (s_2 - p_2 + xy_2 - yz_2 - zx_2 - zw_2 + xw_2 + yw_2);
    out[11] = 2 * ((s * zw) + (p * xy) - (yz * yw) + (zx * xw));

    // COLUMN 3 (Basis W: v=[0, 0, 0, 1])
    out[12] = 2 * ((s * -xw) + (p * -yz) - (zw * zx) + (yw * xy));
    out[13] = 2 * ((s * -yw) + (p * -zx) + (zw * yz) - (xw * xy));
    out[14] = 2 * ((s * -zw) + (p * -xy) + (xw * zx) - (yw * yz));
    out[15] = (s_2 - p_2 + xy_2 + yz_2 + zx_2 - zw_2 - xw_2 - yw_2);

    return out;
}

R×R ❡

R1 × R2 - the SANE Heidelberg way (i think)

Im sure you can short this somehow but after all its expanding to 64 terms and will cancel 8 i think .. dono. ill just do it manually

R1*R2 = (R1[0] + R1[1] e₁₂ + R1[2] e₂₃ + R1[3] e₃₁ + R1[4] e₃₄ + R1[5] e₁₄ + R1[6] e₂₄ + R1[7] e₁₂₃₄)
      * (R2[0] + R2[1] e₁₂ + R2[2] e₂₃ + R2[3] e₃₁ + R2[4] e₃₄ + R2[5] e₁₄ + R2[6] e₂₄ + R2[7] e₁₂₃₄)

      = R1[0]     R2[0] + R1[0]     R2[1] e₁₂ + R1[0]     R2[2] e₂₃ + R1[0]     R2[3] e₃₁ + R1[0]     R2[4] e₃₄ + R1[0]     R2[5] e₁₄ + R1[0]     R2[6] e₂₄ + R1[0]     R2[7] e₁₂₃₄

      + R1[1] e₁₂ R2[0] + R1[1] e₁₂ R2[1] e₁₂ + R1[1] e₁₂ R2[2] e₂₃ + R1[1] e₁₂ R2[3] e₃₁ + R1[1] e₁₂ R2[4] e₃₄ + R1[1] e₁₂ R2[5] e₁₄ + R1[1] e₁₂ R2[6] e₂₄ + R1[1] e₁₂ R2[7] e₁₂₃₄
      + R1[2] e₂₃ R2[0] + R1[2] e₂₃ R2[1] e₁₂ + R1[2] e₂₃ R2[2] e₂₃ + R1[2] e₂₃ R2[3] e₃₁ + R1[2] e₂₃ R2[4] e₃₄ + R1[2] e₂₃ R2[5] e₁₄ + R1[2] e₂₃ R2[6] e₂₄ + R1[2] e₂₃ R2[7] e₁₂₃₄
      + R1[3] e₃₁ R2[0] + R1[3] e₃₁ R2[1] e₁₂ + R1[3] e₃₁ R2[2] e₂₃ + R1[3] e₃₁ R2[3] e₃₁ + R1[3] e₃₁ R2[4] e₃₄ + R1[3] e₃₁ R2[5] e₁₄ + R1[3] e₃₁ R2[6] e₂₄ + R1[3] e₃₁ R2[7] e₁₂₃₄
      + R1[4] e₃₄ R2[0] + R1[4] e₃₄ R2[1] e₁₂ + R1[4] e₃₄ R2[2] e₂₃ + R1[4] e₃₄ R2[3] e₃₁ + R1[4] e₃₄ R2[4] e₃₄ + R1[4] e₃₄ R2[5] e₁₄ + R1[4] e₃₄ R2[6] e₂₄ + R1[4] e₃₄ R2[7] e₁₂₃₄
      + R1[5] e₁₄ R2[0] + R1[5] e₁₄ R2[1] e₁₂ + R1[5] e₁₄ R2[2] e₂₃ + R1[5] e₁₄ R2[3] e₃₁ + R1[5] e₁₄ R2[4] e₃₄ + R1[5] e₁₄ R2[5] e₁₄ + R1[5] e₁₄ R2[6] e₂₄ + R1[5] e₁₄ R2[7] e₁₂₃₄
      + R1[6] e₂₄ R2[0] + R1[6] e₂₄ R2[1] e₁₂ + R1[6] e₂₄ R2[2] e₂₃ + R1[6] e₂₄ R2[3] e₃₁ + R1[6] e₂₄ R2[4] e₃₄ + R1[6] e₂₄ R2[5] e₁₄ + R1[6] e₂₄ R2[6] e₂₄ + R1[6] e₂₄ R2[7] e₁₂₃₄

      + R1[7] e₁₂₃₄ R2[0] + R1[7] e₁₂₃₄ R2[1] e₁₂ + R1[7] e₁₂₃₄ R2[2] e₂₃ + R1[7] e₁₂₃₄ R2[3] e₃₁ + R1[7] e₁₂₃₄ R2[4] e₃₄ + R1[7] e₁₂₃₄ R2[5] e₁₄ + R1[7] e₁₂₃₄ R2[6] e₂₄ + R1[7] e₁₂₃₄ R2[7] e₁₂₃₄

First lets solve all same planes to -1 / clean all eₘₙ² to -1

R1*R2 = (R1[0] + R1[1] e₁₂ + R1[2] e₂₃ + R1[3] e₃₁ + R1[4] e₃₄ + R1[5] e₁₄ + R1[6] e₂₄ + R1[7] e₁₂₃₄)
      * (R2[0] + R2[1] e₁₂ + R2[2] e₂₃ + R2[3] e₃₁ + R2[4] e₃₄ + R2[5] e₁₄ + R2[6] e₂₄ + R2[7] e₁₂₃₄)

      = R1[0]     R2[0]     + R1[0]     R2[1] e₁₂   + R1[0]     R2[2] e₂₃   + R1[0]     R2[3] e₃₁   + R1[0]     R2[4] e₃₄   + R1[0]     R2[5] e₁₄   + R1[0]     R2[6] e₂₄   + R1[0]     R2[7] e₁₂₃₄

      + R1[1] e₁₂ R2[0]     - R1[1]     R2[1]       + R1[1] e₁₂ R2[2] e₂₃   + R1[1] e₁₂ R2[3] e₃₁   + R1[1] e₁₂ R2[4] e₃₄   + R1[1] e₁₂ R2[5] e₁₄   + R1[1] e₁₂ R2[6] e₂₄   + R1[1] e₁₂ R2[7] e₁₂₃₄
      + R1[2] e₂₃ R2[0]     + R1[2] e₂₃ R2[1] e₁₂   - R1[2]     R2[2]       + R1[2] e₂₃ R2[3] e₃₁   + R1[2] e₂₃ R2[4] e₃₄   + R1[2] e₂₃ R2[5] e₁₄   + R1[2] e₂₃ R2[6] e₂₄   + R1[2] e₂₃ R2[7] e₁₂₃₄
      + R1[3] e₃₁ R2[0]     + R1[3] e₃₁ R2[1] e₁₂   + R1[3] e₃₁ R2[2] e₂₃   - R1[3]     R2[3]       + R1[3] e₃₁ R2[4] e₃₄   + R1[3] e₃₁ R2[5] e₁₄   + R1[3] e₃₁ R2[6] e₂₄   + R1[3] e₃₁ R2[7] e₁₂₃₄
      + R1[4] e₃₄ R2[0]     + R1[4] e₃₄ R2[1] e₁₂   + R1[4] e₃₄ R2[2] e₂₃   + R1[4] e₃₄ R2[3] e₃₁   - R1[4]     R2[4]       + R1[4] e₃₄ R2[5] e₁₄   + R1[4] e₃₄ R2[6] e₂₄   + R1[4] e₃₄ R2[7] e₁₂₃₄
      + R1[5] e₁₄ R2[0]     + R1[5] e₁₄ R2[1] e₁₂   + R1[5] e₁₄ R2[2] e₂₃   + R1[5] e₁₄ R2[3] e₃₁   + R1[5] e₁₄ R2[4] e₃₄   - R1[5]     R2[5]       + R1[5] e₁₄ R2[6] e₂₄   + R1[5] e₁₄ R2[7] e₁₂₃₄
      + R1[6] e₂₄ R2[0]     + R1[6] e₂₄ R2[1] e₁₂   + R1[6] e₂₄ R2[2] e₂₃   + R1[6] e₂₄ R2[3] e₃₁   + R1[6] e₂₄ R2[4] e₃₄   + R1[6] e₂₄ R2[5] e₁₄   - R1[6]     R2[6]       + R1[6] e₂₄ R2[7] e₁₂₃₄

      + R1[7] e₁₂₃₄ R2[0]   + R1[7] e₁₂₃₄ R2[1] e₁₂ + R1[7] e₁₂₃₄ R2[2] e₂₃ + R1[7] e₁₂₃₄ R2[3] e₃₁ + R1[7] e₁₂₃₄ R2[4] e₃₄ + R1[7] e₁₂₃₄ R2[5] e₁₄ + R1[7] e₁₂₃₄ R2[6] e₂₄ + R1[7]     R2[7]

Thats 6 so 3 + 2 + 1 swaps for I² then how the axis square

e₁₂₃₄² = e₁₂₃₄₁₂₃₄ = -e₁₂₃₁₄₂₃₄ = e₁₂₁₃₄₂₃₄ = -e₂₃₄₂₃₄ = e₂₃₂₄₃₄ = -e₃₄₃₄ = 1

Lets do all plane multipications sharing an axis

R1*R2 = (R1[0] + R1[1] e₁₂ + R1[2] e₂₃ + R1[3] e₃₁ + R1[4] e₃₄ + R1[5] e₁₄ + R1[6] e₂₄ + R1[7] e₁₂₃₄)
      * (R2[0] + R2[1] e₁₂ + R2[2] e₂₃ + R2[3] e₃₁ + R2[4] e₃₄ + R2[5] e₁₄ + R2[6] e₂₄ + R2[7] e₁₂₃₄)

      = R1[0]     R2[0]     + R1[0]     R2[1] e₁₂   + R1[0]     R2[2] e₂₃   + R1[0]     R2[3] e₃₁   + R1[0]     R2[4] e₃₄   + R1[0]     R2[5] e₁₄   + R1[0]     R2[6] e₂₄   + R1[0]     R2[7] e₁₂₃₄
      //                                                                                                                                                                      e₁₂₁₂₃₄ = -e₃₄
      + R1[1] e₁₂ R2[0]     - R1[1]     R2[1]       - R1[1]     R2[2] e₃₁   + R1[1]     R2[3] e₂₃   + R1[1]     R2[4] e₁₂₃₄ - R1[1]     R2[5] e₂₄   + R1[1]     R2[6] e₁₄   - R1[1]     R2[7] e₃₄
      //                                                                                                                      e₂₃₁₄ = -e₂₁₃₄ = e₁₂₃₄                          e₂₃₁₂₃₄ = -e₃₂₁₂₃₄ = e₃₁₃₄ = -e₁₄
      + R1[2] e₂₃ R2[0]     + R1[2]     R2[1] e₃₁   - R1[2]     R2[2]       - R1[2]     R2[3] e₁₂   + R1[2]     R2[4] e₂₄   + R1[2]     R2[5] e₁₂₃₄ - R1[2]     R2[6] e₃₄   - R1[2]     R2[7] e₁₄
      //                                                                                                                                              e₃₁₂₄ = -e₁₃₂₄ = e₁₂₃₄  e₃₂₃₄ = -e₂₄
      + R1[3] e₃₁ R2[0]     - R1[3]     R2[1] e₂₃   + R1[3]     R2[2] e₁₂   - R1[3]     R2[3]       - R1[3]     R2[4] e₁₄   + R1[3]     R2[5] e₃₄   + R1[3]     R2[6] e₁₂₃₄ - R1[3]     R2[7] e₂₄
      //                      e₃₄₁₂ = -e₃₁₄₂ = e₁₃₄₂ = -e₁₃₂₄ = e₁₂₃₄                                                                                                         e₃₄₁₂₃₄ = -e₃₁₄₂₃₄ = e₁₃₄₂₃₄ = -e₁₃₂₄₃₄ = e₁₂₃₄₃₄ = -e₁₂
      + R1[4] e₃₄ R2[0]     + R1[4]     R2[1] e₁₂₃₄ - R1[4]     R2[2] e₂₄   + R1[4]     R2[3] e₁₄   - R1[4]     R2[4]       - R1[4]     R2[5] e₃₁   + R1[4]     R2[6] e₂₃   - R1[4]     R2[7] e₁₂
      //                                              e₁₄₂₃ = -e₁₂₄₃ = e₁₂₃₄                                                                                                  e₁₄₁₂₃₄ = -e₄₂₃₄ = e₂₄₃₄ = -e₂₃
      + R1[5] e₁₄ R2[0]     + R1[5]     R2[1] e₂₄   + R1[5]     R2[2] e₁₂₃₄ - R1[5]     R2[3] e₃₄   + R1[5]     R2[4] e₃₁   - R1[5]     R2[5]       - R1[5]     R2[6] e₁₂   - R1[5]     R2[7] e₂₃
      //                                                                      e₂₄₃₁ = -e₂₃₄₁ = e₂₃₁₄ = -e₂₁₃₄ = e₁₂₃₄                                                         e₂₄₁₂₃₄ = -e₂₁₄₂₃₄ = e₂₁₂₄₃₄ = -e₁₄₃₄ = e₁₃ = -e₃₁
      + R1[6] e₂₄ R2[0]     - R1[6]     R2[1] e₁₄   + R1[6]     R2[2] e₃₄   + R1[6]     R2[3] e₁₂₃₄ - R1[6]     R2[4] e₂₃   + R1[6]     R2[5] e₁₂   - R1[6]     R2[6]       - R1[6]     R2[7] e₃₁
      //
      //                      e₁₂₃₄₁₂ = -e₁₂₃₁₄₂ = e₁₂₃₁₂₄ = -e₁₂₁₃₂₄ = e₁₂₁₂₃₄ = -e₃₄
      //                                              e₁₂₃₄₂₃ = -e₁₂₃₂₄₃ = e₁₂₃₂₃₄ = -e₁₄
      //                                                                      e₁₂₃₄₃₁ = -e₁₂₄₁ = e₂₁₄₁ = -e₂₄
      //                                                                                            e₁₂₃₄₃₄ = -e₁₂            e₁₂₃₄₁₄ = -e₁₂₃₁ = e₁₂₁₃ = -e₂₃
      //                                                                                                                                              e₁₂₃₄₂₄ = -e₁₂₃₂ = e₁₃ = -e₃₁
      + R1[7] e₁₂₃₄ R2[0]   - R1[7]     R2[1] e₃₄   - R1[7]     R2[2] e₁₄   - R1[7]     R2[3] e₂₄   - R1[7]     R2[4] e₁₂   - R1[7]     R2[5] e₂₃   - R1[7]     R2[6] e₃₁   + R1[7]     R2[7]

clean up and no notes (btw i made 4 mistakes here before having it checked for correctness)

R1*R2 = (R1[0] + R1[1] e₁₂ + R1[2] e₂₃ + R1[3] e₃₁ + R1[4] e₃₄ + R1[5] e₁₄ + R1[6] e₂₄ + R1[7] e₁₂₃₄)
      * (R2[0] + R2[1] e₁₂ + R2[2] e₂₃ + R2[3] e₃₁ + R2[4] e₃₄ + R2[5] e₁₄ + R2[6] e₂₄ + R2[7] e₁₂₃₄)

      = R1[0] R2[0]     + R1[0] R2[1] e₁₂   + R1[0] R2[2] e₂₃   + R1[0] R2[3] e₃₁   + R1[0] R2[4] e₃₄   + R1[0] R2[5] e₁₄   + R1[0] R2[6] e₂₄   + R1[0] R2[7] e₁₂₃₄
      + R1[1] R2[0] e₁₂ - R1[1] R2[1]       - R1[1] R2[2] e₃₁   + R1[1] R2[3] e₂₃   + R1[1] R2[4] e₁₂₃₄ - R1[1] R2[5] e₂₄   + R1[1] R2[6] e₁₄   - R1[1] R2[7] e₃₄
      + R1[2] R2[0] e₂₃ + R1[2] R2[1] e₃₁   - R1[2] R2[2]       - R1[2] R2[3] e₁₂   + R1[2] R2[4] e₂₄   + R1[2] R2[5] e₁₂₃₄ - R1[2] R2[6] e₃₄   - R1[2] R2[7] e₁₄
      + R1[3] R2[0] e₃₁ - R1[3] R2[1] e₂₃   + R1[3] R2[2] e₁₂   - R1[3] R2[3]       - R1[3] R2[4] e₁₄   + R1[3] R2[5] e₃₄   + R1[3] R2[6] e₁₂₃₄ - R1[3] R2[7] e₂₄
      + R1[4] R2[0] e₃₄ + R1[4] R2[1] e₁₂₃₄ - R1[4] R2[2] e₂₄   + R1[4] R2[3] e₁₄   - R1[4] R2[4]       - R1[4] R2[5] e₃₁   + R1[4] R2[6] e₂₃   - R1[4] R2[7] e₁₂
      + R1[5] R2[0] e₁₄ + R1[5] R2[1] e₂₄   + R1[5] R2[2] e₁₂₃₄ - R1[5] R2[3] e₃₄   + R1[5] R2[4] e₃₁   - R1[5] R2[5]       - R1[5] R2[6] e₁₂   - R1[5] R2[7] e₂₃
      + R1[6] R2[0] e₂₄ - R1[6] R2[1] e₁₄   + R1[6] R2[2] e₃₄   + R1[6] R2[3] e₁₂₃₄ - R1[6] R2[4] e₂₃   + R1[6] R2[5] e₁₂   - R1[6] R2[6]       - R1[6] R2[7] e₃₁

      + R1[7] R2[0] e₁₂₃₄ - R1[7] R2[1] e₃₄   - R1[7] R2[2] e₁₄   - R1[7] R2[3] e₂₄   - R1[7] R2[4] e₁₂   - R1[7] R2[5] e₂₃   - R1[7] R2[6] e₃₁   + R1[7] R2[7]

Collect

R'[0] s     = R1[0] R2[0] - R1[1] R2[1] - R1[2] R2[2] - R1[3] R2[3] - R1[4] R2[4] - R1[5] R2[5] - R1[6] R2[6] + R1[7] R2[7]
R'[1] e₁₂   = R1[0] R2[1] + R1[1] R2[0] - R1[2] R2[3] + R1[3] R2[2] - R1[4] R2[7] - R1[5] R2[6] + R1[6] R2[5] - R1[7] R2[4]
R'[2] e₂₃   = R1[0] R2[2] + R1[1] R2[3] + R1[2] R2[0] - R1[3] R2[1] + R1[4] R2[6] - R1[5] R2[7] - R1[6] R2[4] - R1[7] R2[5]
R'[3] e₃₁   = R1[0] R2[3] - R1[1] R2[2] + R1[2] R2[1] + R1[3] R2[0] - R1[4] R2[5] + R1[5] R2[4] - R1[6] R2[7] - R1[7] R2[6]
R'[4] e₃₄   = R1[0] R2[4] - R1[1] R2[7] - R1[2] R2[6] + R1[3] R2[5] + R1[4] R2[0] - R1[5] R2[3] + R1[6] R2[2] - R1[7] R2[1]
R'[5] e₁₄   = R1[0] R2[5] + R1[1] R2[6] - R1[2] R2[7] - R1[3] R2[4] + R1[4] R2[3] + R1[5] R2[0] - R1[6] R2[1] - R1[7] R2[2]
R'[6] e₂₄   = R1[0] R2[6] - R1[1] R2[5] + R1[2] R2[4] - R1[3] R2[7] - R1[4] R2[2] + R1[5] R2[1] + R1[6] R2[0] - R1[7] R2[3]
R'[7] e₁₂₃₄ = R1[0] R2[7] + R1[1] R2[4] + R1[2] R2[5] + R1[3] R2[6] + R1[4] R2[1] + R1[5] R2[2] + R1[6] R2[3] + R1[7] R2[0]

Code (old e₁₂!) ❡

static multiplyHB(a: Rot4, b: Rot4): Rot4 {

    const s1 = a[0], xy1 = a[1], yz1 = a[2], zx1 = a[3], zw1 = a[4], xw1 = a[5], yw1 = a[6], p1 = a[7];
    const s2 = b[0], xy2 = b[1], yz2 = b[2], zx2 = b[3], zw2 = b[4], xw2 = b[5], yw2 = b[6], p2 = b[7];

    const s = (s1 * s2) - (xy1 * xy2) - (yz1 * yz2) - (zx1 * zx2) - (zw1 * zw2) - (xw1 * xw2) - (yw1 * yw2) + (p1 * p2);
    const exy = (s1 * xy2) + (xy1 * s2) - (yz1 * zx2) + (zx1 * yz2) - (zw1 * p2) - (xw1 * yw2) + (yw1 * xw2) - (p1 * zw2);
    const eyz = (s1 * yz2) + (xy1 * zx2) + (yz1 * s2) - (zx1 * xy2) + (zw1 * yw2) - (xw1 * p2) - (yw1 * zw2) - (p1 * xw2);
    const ezx = (s1 * zx2) - (xy1 * yz2) + (yz1 * xy2) + (zx1 * s2) - (zw1 * xw2) + (xw1 * zw2) - (yw1 * p2) - (p1 * yw2);
    const ezw = (s1 * zw2) - (xy1 * p2) - (yz1 * yw2) + (zx1 * xw2) + (zw1 * s2) - (xw1 * zx2) + (yw1 * yz2) - (p1 * xy2);
    const exw = (s1 * xw2) + (xy1 * yw2) - (yz1 * p2) - (zx1 * zw2) + (zw1 * zx2) + (xw1 * s2) - (yw1 * xy2) - (p1 * yz2);
    const eyw = (s1 * yw2) - (xy1 * xw2) + (yz1 * zw2) - (zx1 * p2) - (zw1 * yz2) + (xw1 * xy2) + (yw1 * s2) - (p1 * zx2);
    const exyzw = (s1 * p2) + (xy1 * zw2) + (yz1 * xw2) + (zx1 * yw2) + (zw1 * xy2) + (xw1 * yz2) + (yw1 * zx2) + (p1 * s2);

    const out = new Float32Array([s, exy, eyz, ezx, ezw, exw, eyw, exyzw]);

    return out;

}

New Code for e₂₃ e₃₁ e₁₂ e₁₄ e₂₄ e₃₄ ❡

// MARK: Rotate

/**
 * Rotates a vector 'v' by rotor 'r' using the Sandwich Product R̃vR
 * Derived manually
 * 
 */
static rotateVector4(r: Rot4, v: Vec4): Vec4 {
    const s = r[0], yz = r[1], zx = r[2], xy = r[3],
        xw = r[4], yw = r[5], zw = r[6], p = r[7];

    // Helper squares
    const s_2 = s * s;
    const yz_2 = yz * yz, zx_2 = zx * zx, xy_2 = xy * xy;
    const xw_2 = xw * xw, yw_2 = yw * yw, zw_2 = zw * zw;
    const p_2 = p * p;

    const out = new Float32Array(4);

    // (s² - I²)v  + 2s(v · B) + 2I(v ∧ B) + 2(v · B)B  - vB²

    // ((())) looks stupid but makes it wgsl-safe

    out[0] = ((s_2 - p_2 - xy_2 + yz_2 - zx_2 + zw_2 - xw_2 + yw_2) * v[0])
        + (2 * ((s * -xy) + (p * -zw) + (yz * zx) - (yw * xw)) * v[1])
        + (2 * ((s * zx) + (p * yw) + (yz * xy) - (zw * xw)) * v[2])
        + (2 * ((s * -xw) + (p * -yz) - (zw * zx) + (yw * xy)) * v[3]);

    out[1] = (2 * ((s * xy) + (p * zw) + (zx * yz) - (xw * yw)) * v[0])
        + ((s_2 - p_2 - xy_2 - yz_2 + zx_2 + zw_2 + xw_2 - yw_2) * v[1])
        + (2 * ((s * -yz) + (p * -xw) + (zx * xy) - (zw * yw)) * v[2])
        + (2 * ((s * -yw) + (p * -zx) + (zw * yz) - (xw * xy)) * v[3]);

    out[2] = (2 * ((s * -zx) + (p * -yw) + (xy * yz) - (xw * zw)) * v[0])
        + (2 * ((s * yz) + (p * xw) + (xy * zx) - (yw * zw)) * v[1])
        + ((s_2 - p_2 + xy_2 - yz_2 - zx_2 - zw_2 + xw_2 + yw_2) * v[2])
        + (2 * ((s * -zw) + (p * -xy) + (xw * zx) - (yw * yz)) * v[3]);

    out[3] = (2 * ((s * xw) + (p * yz) + (xy * yw) - (zx * zw)) * v[0])
        + (2 * ((s * yw) + (p * zx) - (xy * xw) + (yz * zw)) * v[1])
        + (2 * ((s * zw) + (p * xy) - (yz * yw) + (zx * xw)) * v[2])
        + ((s_2 - p_2 + xy_2 + yz_2 + zx_2 - zw_2 - xw_2 - yw_2) * v[3]);

    return out;
}

// MARK: Export

static toMat4(r: Rot4): Mat4 {
    const s = r[0], yz = r[1], zx = r[2], xy = r[3],
        xw = r[4], yw = r[5], zw = r[6], p = r[7];

    // Helper squares
    const s_2 = s * s;
    const yz_2 = yz * yz, zx_2 = zx * zx, xy_2 = xy * xy;
    const xw_2 = xw * xw, yw_2 = yw * yw, zw_2 = zw * zw;
    const p_2 = p * p;

    const out = MAT4_IDENTITY.slice();

    // all vₙ parts IN TO OUT v'₁ to v'₄ - COLUMN MAYOR!!

    // COLUMN 0 (Basis X: v=[1, 0, 0, 0])  
    out[0] = (s_2 - p_2 - xy_2 + yz_2 - zx_2 + zw_2 - xw_2 + yw_2);
    out[1] = 2 * ((s * xy) + (p * zw) + (zx * yz) - (xw * yw));
    out[2] = 2 * ((s * -zx) + (p * -yw) + (xy * yz) - (xw * zw));
    out[3] = 2 * ((s * xw) + (p * yz) + (xy * yw) - (zx * zw));

    // COLUMN 1 (Basis Y: v=[0, 1, 0, 0])
    out[4] = 2 * ((s * -xy) + (p * -zw) + (yz * zx) - (yw * xw));
    out[5] = (s_2 - p_2 - xy_2 - yz_2 + zx_2 + zw_2 + xw_2 - yw_2);
    out[6] = 2 * ((s * yz) + (p * xw) + (xy * zx) - (yw * zw));
    out[7] = 2 * ((s * yw) + (p * zx) - (xy * xw) + (yz * zw));

    // COLUMN 2 (Basis Z: v=[0, 0, 1, 0])
    out[8] = 2 * ((s * zx) + (p * yw) + (yz * xy) - (zw * xw));
    out[9] = 2 * ((s * -yz) + (p * -xw) + (zx * xy) - (zw * yw));
    out[10] = (s_2 - p_2 + xy_2 - yz_2 - zx_2 - zw_2 + xw_2 + yw_2);
    out[11] = 2 * ((s * zw) + (p * xy) - (yz * yw) + (zx * xw));

    // COLUMN 3 (Basis W: v=[0, 0, 0, 1])
    out[12] = 2 * ((s * -xw) + (p * -yz) - (zw * zx) + (yw * xy));
    out[13] = 2 * ((s * -yw) + (p * -zx) + (zw * yz) - (xw * xy));
    out[14] = 2 * ((s * -zw) + (p * -xy) + (xw * zx) - (yw * yz));
    out[15] = (s_2 - p_2 + xy_2 + yz_2 + zx_2 - zw_2 - xw_2 - yw_2);

    return out;
}

// MARK: R×R

// Rotor4-Rotor4 Product
static multiply(a: Rot4, b: Rot4): Rot4 {

    const s1 = a[0], yz1 = a[1], zx1 = a[2], xy1 = a[3], xw1 = a[4], yw1 = a[5], zw1 = a[6], p1 = a[7];
    const s2 = b[0], yz2 = b[1], zx2 = b[2], xy2 = b[3], xw2 = b[4], yw2 = b[5], zw2 = b[6], p2 = b[7];

    const s = (s1 * s2) - (xy1 * xy2) - (yz1 * yz2) - (zx1 * zx2) - (zw1 * zw2) - (xw1 * xw2) - (yw1 * yw2) + (p1 * p2);
    const eyz = (s1 * yz2) + (xy1 * zx2) + (yz1 * s2) - (zx1 * xy2) + (zw1 * yw2) - (xw1 * p2) - (yw1 * zw2) - (p1 * xw2);
    const exy = (s1 * xy2) + (xy1 * s2) - (yz1 * zx2) + (zx1 * yz2) - (zw1 * p2) - (xw1 * yw2) + (yw1 * xw2) - (p1 * zw2);
    const ezx = (s1 * zx2) - (xy1 * yz2) + (yz1 * xy2) + (zx1 * s2) - (zw1 * xw2) + (xw1 * zw2) - (yw1 * p2) - (p1 * yw2);
    const exw = (s1 * xw2) + (xy1 * yw2) - (yz1 * p2) - (zx1 * zw2) + (zw1 * zx2) + (xw1 * s2) - (yw1 * xy2) - (p1 * yz2);
    const eyw = (s1 * yw2) - (xy1 * xw2) + (yz1 * zw2) - (zx1 * p2) - (zw1 * yz2) + (xw1 * xy2) + (yw1 * s2) - (p1 * zx2);
    const ezw = (s1 * zw2) - (xy1 * p2) - (yz1 * yw2) + (zx1 * xw2) + (zw1 * s2) - (xw1 * zx2) + (yw1 * yz2) - (p1 * xy2);
    const exyzw = (s1 * p2) + (xy1 * zw2) + (yz1 * xw2) + (zx1 * yw2) + (zw1 * xy2) + (xw1 * yz2) + (yw1 * zx2) + (p1 * s2);

    const out = new Float32Array([s, eyz, ezx, exy, exw, eyw, ezw, exyzw]);

    return out;

}