3D PGA Reference Implementation in Wolfram Notebook

Usage examples:

  1. Line from point A (xa, ya, za) to point B (xb, yb, zb)
    MatrixForm[Join[Dual[e0 + xa*e1 + ya*e2 + za*e3], Dual[e0 + xb*e1 + yb*e2 + zb*e3]]]

  2. Distance from plane 2x + 3y + 4z + 5 = 0 to point (10, 11, 12)
    MatrixForm[Join[Normalized[2e1 + 3e2 + 4e3 + 5e0], Dual[e0 + 10e1 + 11e2 + 12e3]]]

  3. Plane parallel to plane P and containing normalized point X
    MatrixForm[(plane[p].pointN[x])**Reverse[pointN[x]]]

  4. Apply normalized motor to ideal point
    mn=motorN[m]
    ip=pointI[p]
    MatrixForm[SimplifyReals[mn**ip**Reverse[mn]]]

Wolfram code:

{"1","e0","e1","e2","e3","e01","e02","e03","e12","e31","e23","e021","e013","e032","e123","e0123"}

e0=tPGA3[0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0];
e1=tPGA3[0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0];
e2=tPGA3[0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0];
e3=tPGA3[0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0];
e01=tPGA3[0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0];
e02=tPGA3[0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0];
e03=tPGA3[0,0,0,0,0,0,0,1,0,0,0,0,0,0,0,0];
e12=tPGA3[0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0];
e31=tPGA3[0,0,0,0,0,0,0,0,0,1,0,0,0,0,0,0];
e23=tPGA3[0,0,0,0,0,0,0,0,0,0,1,0,0,0,0,0];
e021=tPGA3[0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,0];
e013=tPGA3[0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0];
e032=tPGA3[0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0];
e123=tPGA3[0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0];
e0123=tPGA3[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1];

(*Returns dualized multivector.*)
tPGA3/:Dual[a_tPGA3]:=tPGA3[a[[16]],-a[[15]],-a[[14]],-a[[13]],-a[[12]],a[[11]],a[[10]],a[[9]],a[[8]],a[[7]],a[[6]],a[[5]],a[[4]],a[[3]],a[[2]],a[[1]]]
tPieceValue/:Dual[a_tPieceValue]:=MapAt[Dual,a,{All,1}]
Dual/:Dual[scalar:Except[_tPGA3|_tPieceValue|_Dual]]:=Dual[ScalarPGA3[scalar]]

(*Returns undualized multivector.*)
tPGA3/:UnDual[a_tPGA3]:=tPGA3[a[[16]],a[[15]],a[[14]],a[[13]],a[[12]],a[[11]],a[[10]],a[[9]],a[[8]],a[[7]],a[[6]],-a[[5]],-a[[4]],-a[[3]],-a[[2]],a[[1]]]
tPieceValue/:UnDual[a_tPieceValue]:=MapAt[UnDual,a,{All,1}]
UnDual/:UnDual[scalar:Except[_tPGA3|_tPieceValue|_UnDual]]:=UnDual[ScalarPGA3[scalar]]

(*Returns reversed multivector.*)
tPGA3/:Reverse[a_tPGA3]:=tPGA3[a[[1]],a[[2]],a[[3]],a[[4]],a[[5]],-a[[6]],-a[[7]],-a[[8]],-a[[9]],-a[[10]],-a[[11]],-a[[12]],-a[[13]],-a[[14]],-a[[15]],a[[16]]]
tPieceValue/:Reverse[a_tPieceValue]:=MapAt[Reverse,a,{All,1}]

(*Returns the involute of multivector.*)
tPGA3/:Involute[a_tPGA3]:=tPGA3[a[[1]],-a[[2]],-a[[3]],-a[[4]],-a[[5]],a[[6]],a[[7]],a[[8]],a[[9]],a[[10]],a[[11]],-a[[12]],-a[[13]],-a[[14]],-a[[15]],a[[16]]]
tPieceValue/:Involute[a_tPieceValue]:=MapAt[Involute,a,{All,1}]
Involute/:Involute[scalar:Except[_tPGA3|_tPieceValue|_Involute]]:=Involute[ScalarPGA3[scalar]]

(*Returns the conjugate of multivector.*)
tPGA3/:Conjugate[a_tPGA3]:=tPGA3[a[[1]],-a[[2]],-a[[3]],-a[[4]],-a[[5]],-a[[6]],-a[[7]],-a[[8]],-a[[9]],-a[[10]],-a[[11]],a[[12]],a[[13]],a[[14]],a[[15]],a[[16]]]
tPieceValue/:Conjugate[a_tPieceValue]:=MapAt[Conjugate,a,{All,1}]

(*Returns the inverse of multivector.*)
tPGA3/:Inverse[a_tPGA3]:=Simplify[
COND[(a[[5]]*a[[9]]+a[[4]]*a[[10]]+a[[3]]*a[[11]]-a[[1]]*a[[15]]!=0)||(a[[1]]^2-a[[3]]^2-a[[4]]^2-a[[5]]^2+a[[9]]^2+a[[10]]^2+a[[11]]^2-a[[15]]^2!=0),
With[{
p=2*(a[[5]]*a[[9]]+a[[4]]*a[[10]]+a[[3]]*a[[11]]-a[[1]]*a[[15]]),
s=a[[1]]^2-a[[3]]^2-a[[4]]^2-a[[5]]^2+a[[9]]^2+a[[10]]^2+a[[11]]^2-a[[15]]^2,
k1=2*(a[[1]]*a[[14]]+a[[2]]*a[[11]]+a[[3]]*a[[16]]-a[[4]]*a[[8]]+a[[5]]*a[[7]]-a[[6]]*a[[15]]+a[[9]]*a[[13]]-a[[10]]*a[[12]]),
k2=2*(a[[1]]*a[[13]]+a[[2]]*a[[10]]+a[[3]]*a[[8]]+a[[4]]*a[[16]]-a[[5]]*a[[6]]-a[[7]]*a[[15]]-a[[9]]*a[[14]]+a[[11]]*a[[12]]),
k3=2*(a[[1]]*a[[12]]+a[[2]]*a[[9]]-a[[3]]*a[[7]]+a[[4]]*a[[6]]+a[[5]]*a[[16]]-a[[8]]*a[[15]]+a[[10]]*a[[14]]-a[[11]]*a[[13]]),
k4=2*(a[[1]]*a[[16]]+a[[2]]*a[[15]]+a[[3]]*a[[14]]+a[[4]]*a[[13]]+a[[5]]*a[[12]]-a[[6]]*a[[11]]-a[[7]]*a[[10]]-a[[8]]*a[[9]])},
With[{dr=1/(p^2+s^2)},tPGA3[
(a[[1]]*s-a[[15]]*p)*dr,
(-a[[2]]*s-a[[16]]*p+a[[11]]*k1+a[[10]]*k2+a[[9]]*k3-a[[15]]*k4)*dr,
(-a[[3]]*s+a[[11]]*p)*dr,
(-a[[4]]*s+a[[10]]*p)*dr,
(-a[[5]]*s+a[[9]]*p)*dr,
(-a[[6]]*s+a[[14]]*p+a[[15]]*k1+a[[5]]*k2-a[[4]]*k3-a[[11]]*k4)*dr,
(-a[[7]]*s+a[[13]]*p-a[[5]]*k1+a[[15]]*k2+a[[3]]*k3-a[[10]]*k4)*dr,
(-a[[8]]*s+a[[12]]*p+a[[4]]*k1-a[[3]]*k2+a[[15]]*k3-a[[9]]*k4)*dr,
(-a[[9]]*s-a[[5]]*p)*dr,
(-a[[10]]*s-a[[4]]*p)*dr,
(-a[[11]]*s-a[[3]]*p)*dr,
(a[[12]]*s+a[[8]]*p+a[[10]]*k1-a[[11]]*k2-a[[1]]*k3+a[[5]]*k4)*dr,
(a[[13]]*s+a[[7]]*p-a[[9]]*k1-a[[1]]*k2+a[[11]]*k3+a[[4]]*k4)*dr,
(a[[14]]*s+a[[6]]*p-a[[1]]*k1+a[[9]]*k2-a[[10]]*k3+a[[3]]*k4)*dr,
(a[[15]]*s+a[[1]]*p)*dr,
(a[[16]]*s-a[[2]]*p+a[[3]]*k1+a[[4]]*k2+a[[5]]*k3-a[[1]]*k4)*dr]]]]]
tPieceValue/:Inverse[a_tPieceValue]:=Simplify[MapAt[Inverse,a,{All,1}]]

(*Returns the sum of 2 multivectors.*)
tPGA3/:a_tPGA3+b_tPGA3:=Simplify[tPGA3[
a[[1]]+b[[1]],
a[[2]]+b[[2]],
a[[3]]+b[[3]],
a[[4]]+b[[4]],
a[[5]]+b[[5]],
a[[6]]+b[[6]],
a[[7]]+b[[7]],
a[[8]]+b[[8]],
a[[9]]+b[[9]],
a[[10]]+b[[10]],
a[[11]]+b[[11]],
a[[12]]+b[[12]],
a[[13]]+b[[13]],
a[[14]]+b[[14]],
a[[15]]+b[[15]],
a[[16]]+b[[16]]]]
tPieceValue/:Plus[a_tPieceValue,b_]:=Simplify[MapAt[Function[{x},Plus[x,b]],a,{All,1}]]
tPGA3/:Plus[a_tPGA3,scalar:Except[_tPGA3|_tPieceValue|_Plus]]:=Plus[a,ScalarPGA3[scalar]]
tPGA3/:Times[_tPGA3,_tPGA3]:=Echo[Undefined,"ERROR: Times[_tPGA3,_tPGA3]"]
tPieceValue/:Times[a_tPieceValue,b_]:=Simplify[MapAt[Function[{x},Times[x,b]],a,{All,1}]]
tPGA3/:Times[a_tPGA3,scalar:Except[_tPGA3|_tPieceValue|_Times]]:=NonCommutativeMultiply[a,ScalarPGA3[scalar]]

(*Returns the geometric product of 2 multivectors.*)
tPGA3/:a_tPGA3**b_tPGA3:=Simplify[tPGA3[
b[[1]]*a[[1]]+b[[3]]*a[[3]]+b[[4]]*a[[4]]+b[[5]]*a[[5]]-b[[9]]*a[[9]]-b[[10]]*a[[10]]-b[[11]]*a[[11]]-b[[15]]*a[[15]],
b[[2]]*a[[1]]+b[[1]]*a[[2]]-b[[6]]*a[[3]]-b[[7]]*a[[4]]-b[[8]]*a[[5]]+b[[3]]*a[[6]]+b[[4]]*a[[7]]+b[[5]]*a[[8]]+b[[12]]*a[[9]]+b[[13]]*a[[10]]+b[[14]]*a[[11]]+b[[9]]*a[[12]]+b[[10]]*a[[13]]+b[[11]]*a[[14]]+b[[16]]*a[[15]]-b[[15]]*a[[16]],
b[[3]]*a[[1]]+b[[1]]*a[[3]]-b[[9]]*a[[4]]+b[[10]]*a[[5]]+b[[4]]*a[[9]]-b[[5]]*a[[10]]-b[[15]]*a[[11]]-b[[11]]*a[[15]],
b[[4]]*a[[1]]+b[[9]]*a[[3]]+b[[1]]*a[[4]]-b[[11]]*a[[5]]-b[[3]]*a[[9]]-b[[15]]*a[[10]]+b[[5]]*a[[11]]-b[[10]]*a[[15]],b[[5]]*a[[1]]-b[[10]]*a[[3]]+b[[11]]*a[[4]]+b[[1]]*a[[5]]-b[[15]]*a[[9]]+b[[3]]*a[[10]]-b[[4]]*a[[11]]-b[[9]]*a[[15]],
b[[6]]*a[[1]]+b[[3]]*a[[2]]-b[[2]]*a[[3]]-b[[12]]*a[[4]]+b[[13]]*a[[5]]+b[[1]]*a[[6]]-b[[9]]*a[[7]]+b[[10]]*a[[8]]+b[[7]]*a[[9]]-b[[8]]*a[[10]]-b[[16]]*a[[11]]-b[[4]]*a[[12]]+b[[5]]*a[[13]]+b[[15]]*a[[14]]-b[[14]]*a[[15]]-b[[11]]*a[[16]],
b[[7]]*a[[1]]+b[[4]]*a[[2]]+b[[12]]*a[[3]]-b[[2]]*a[[4]]-b[[14]]*a[[5]]+b[[9]]*a[[6]]+b[[1]]*a[[7]]-b[[11]]*a[[8]]-b[[6]]*a[[9]]-b[[16]]*a[[10]]+b[[8]]*a[[11]]+b[[3]]*a[[12]]+b[[15]]*a[[13]]-b[[5]]*a[[14]]-b[[13]]*a[[15]]-b[[10]]*a[[16]],
b[[8]]*a[[1]]+b[[5]]*a[[2]]-b[[13]]*a[[3]]+b[[14]]*a[[4]]-b[[2]]*a[[5]]-b[[10]]*a[[6]]+b[[11]]*a[[7]]+b[[1]]*a[[8]]-b[[16]]*a[[9]]+b[[6]]*a[[10]]-b[[7]]*a[[11]]+b[[15]]*a[[12]]-b[[3]]*a[[13]]+b[[4]]*a[[14]]-b[[12]]*a[[15]]-b[[9]]*a[[16]],
b[[9]]*a[[1]]+b[[4]]*a[[3]]-b[[3]]*a[[4]]+b[[15]]*a[[5]]+b[[1]]*a[[9]]+b[[11]]*a[[10]]-b[[10]]*a[[11]]+b[[5]]*a[[15]],
b[[10]]*a[[1]]-b[[5]]*a[[3]]+b[[15]]*a[[4]]+b[[3]]*a[[5]]-b[[11]]*a[[9]]+b[[1]]*a[[10]]+b[[9]]*a[[11]]+b[[4]]*a[[15]],b[[11]]*a[[1]]+b[[15]]*a[[3]]+b[[5]]*a[[4]]-b[[4]]*a[[5]]+b[[10]]*a[[9]]-b[[9]]*a[[10]]+b[[1]]*a[[11]]+b[[3]]*a[[15]],
b[[12]]*a[[1]]-b[[9]]*a[[2]]+b[[7]]*a[[3]]-b[[6]]*a[[4]]+b[[16]]*a[[5]]-b[[4]]*a[[6]]+b[[3]]*a[[7]]-b[[15]]*a[[8]]-b[[2]]*a[[9]]+b[[14]]*a[[10]]-b[[13]]*a[[11]]+b[[1]]*a[[12]]+b[[11]]*a[[13]]-b[[10]]*a[[14]]+b[[8]]*a[[15]]-b[[5]]*a[[16]],
b[[13]]*a[[1]]-b[[10]]*a[[2]]-b[[8]]*a[[3]]+b[[16]]*a[[4]]+b[[6]]*a[[5]]+b[[5]]*a[[6]]-b[[15]]*a[[7]]-b[[3]]*a[[8]]-b[[14]]*a[[9]]-b[[2]]*a[[10]]+b[[12]]*a[[11]]-b[[11]]*a[[12]]+b[[1]]*a[[13]]+b[[9]]*a[[14]]+b[[7]]*a[[15]]-b[[4]]*a[[16]],
b[[14]]*a[[1]]-b[[11]]*a[[2]]+b[[16]]*a[[3]]+b[[8]]*a[[4]]-b[[7]]*a[[5]]-b[[15]]*a[[6]]-b[[5]]*a[[7]]+b[[4]]*a[[8]]+b[[13]]*a[[9]]-b[[12]]*a[[10]]-b[[2]]*a[[11]]+b[[10]]*a[[12]]-b[[9]]*a[[13]]+b[[1]]*a[[14]]+b[[6]]*a[[15]]-b[[3]]*a[[16]],
b[[15]]*a[[1]]+b[[11]]*a[[3]]+b[[10]]*a[[4]]+b[[9]]*a[[5]]+b[[5]]*a[[9]]+b[[4]]*a[[10]]+b[[3]]*a[[11]]+b[[1]]*a[[15]],
b[[16]]*a[[1]]+b[[15]]*a[[2]]+b[[14]]*a[[3]]+b[[13]]*a[[4]]+b[[12]]*a[[5]]+b[[11]]*a[[6]]+b[[10]]*a[[7]]+b[[9]]*a[[8]]+b[[8]]*a[[9]]+b[[7]]*a[[10]]+b[[6]]*a[[11]]-b[[5]]*a[[12]]-b[[4]]*a[[13]]-b[[3]]*a[[14]]-b[[2]]*a[[15]]+b[[1]]*a[[16]]]]
tPieceValue/:NonCommutativeMultiply[a_tPieceValue,b_]:=Simplify[MapAt[Function[{x},NonCommutativeMultiply[x,b]],a,{All,1}]]
tPieceValue/:NonCommutativeMultiply[a_,b_tPieceValue]:=Simplify[MapAt[Function[{x},NonCommutativeMultiply[a,x]],b,{All,1}]]
tPGA3/:NonCommutativeMultiply[a_tPGA3,scalar:Except[_tPGA3|_tPieceValue|_NonCommutativeMultiply]]:=NonCommutativeMultiply[a,ScalarPGA3[scalar]]
tPGA3/:NonCommutativeMultiply[scalar:Except[_tPGA3|_tPieceValue|_NonCommutativeMultiply],b_tPGA3]:=NonCommutativeMultiply[ScalarPGA3[scalar],b]

(*Returns the outer product of 2 multivectors.*)
tPGA3/:Meet[a_tPGA3,b_tPGA3]:=Simplify[tPGA3[
b[[1]]*a[[1]],
b[[2]]*a[[1]]+b[[1]]*a[[2]],
b[[3]]*a[[1]]+b[[1]]*a[[3]],
b[[4]]*a[[1]]+b[[1]]*a[[4]],
b[[5]]*a[[1]]+b[[1]]*a[[5]],
b[[6]]*a[[1]]+b[[3]]*a[[2]]-b[[2]]*a[[3]]+b[[1]]*a[[6]],
b[[7]]*a[[1]]+b[[4]]*a[[2]]-b[[2]]*a[[4]]+b[[1]]*a[[7]],
b[[8]]*a[[1]]+b[[5]]*a[[2]]-b[[2]]*a[[5]]+b[[1]]*a[[8]],
b[[9]]*a[[1]]+b[[4]]*a[[3]]-b[[3]]*a[[4]]+b[[1]]*a[[9]],
b[[10]]*a[[1]]-b[[5]]*a[[3]]+b[[3]]*a[[5]]+b[[1]]*a[[10]],
b[[11]]*a[[1]]+b[[5]]*a[[4]]-b[[4]]*a[[5]]+b[[1]]*a[[11]],
b[[12]]*a[[1]]-b[[9]]*a[[2]]+b[[7]]*a[[3]]-b[[6]]*a[[4]]-b[[4]]*a[[6]]+b[[3]]*a[[7]]-b[[2]]*a[[9]]+b[[1]]*a[[12]],
b[[13]]*a[[1]]-b[[10]]*a[[2]]-b[[8]]*a[[3]]+b[[6]]*a[[5]]+b[[5]]*a[[6]]-b[[3]]*a[[8]]-b[[2]]*a[[10]]+b[[1]]*a[[13]],
b[[14]]*a[[1]]-b[[11]]*a[[2]]+b[[8]]*a[[4]]-b[[7]]*a[[5]]-b[[5]]*a[[7]]+b[[4]]*a[[8]]-b[[2]]*a[[11]]+b[[1]]*a[[14]],
b[[15]]*a[[1]]+b[[11]]*a[[3]]+b[[10]]*a[[4]]+b[[9]]*a[[5]]+b[[5]]*a[[9]]+b[[4]]*a[[10]]+b[[3]]*a[[11]]+b[[1]]*a[[15]],
b[[16]]*a[[1]]+b[[15]]*a[[2]]+b[[14]]*a[[3]]+b[[13]]*a[[4]]+b[[12]]*a[[5]]+b[[11]]*a[[6]]+b[[10]]*a[[7]]+b[[9]]*a[[8]]+b[[8]]*a[[9]]+b[[7]]*a[[10]]+b[[6]]*a[[11]]-b[[5]]*a[[12]]-b[[4]]*a[[13]]-b[[3]]*a[[14]]-b[[2]]*a[[15]]+b[[1]]*a[[16]]]]
SetAttributes[Meet,Flat]
tPieceValue/:Meet[a_tPieceValue,b_]:=Simplify[MapAt[Function[{x},Meet[x,b]],a,{All,1}]]
tPieceValue/:Meet[a_,b_tPieceValue]:=Simplify[MapAt[Function[{x},Meet[a,x]],b,{All,1}]]
Meet/:Meet[a_tPGA3,scalar:Except[_tPGA3|_tPieceValue|_Meet]]:=Meet[a,ScalarPGA3[scalar]]
Meet/:Meet[scalar:Except[_tPGA3|_tPieceValue|_Meet],b_tPGA3]:=Meet[ScalarPGA3[scalar],b]

(*Returns the regressive product of 2 multivectors.*)
tPGA3/:Join[a_tPGA3,b_tPGA3]:=Dual[Meet[UnDual[a],UnDual[b]]]
tPieceValue/:Join[a_tPieceValue,b_]:=Simplify[MapAt[Function[{x},Join[x,b]],a,{All,1}]]
tPieceValue/:Join[a_,b_tPieceValue]:=Simplify[MapAt[Function[{x},Join[a,x]],b,{All,1}]]
tPGA3/:Join[a_tPGA3,scalar:Except[_tPGA3|_tPieceValue|_Join]]:=Join[a,ScalarPGA3[scalar]]
tPGA3/:Join[scalar:Except[_tPGA3|_tPieceValue|_Join],b_tPGA3]:=Join[ScalarPGA3[scalar],b]

(*Returns the inner product of 2 multivectors.*)
tPGA3/:a_tPGA3.b_tPGA3:=Simplify[tPGA3[
b[[1]]*a[[1]]+b[[3]]*a[[3]]+b[[4]]*a[[4]]+b[[5]]*a[[5]]-b[[9]]*a[[9]]-b[[10]]*a[[10]]-b[[11]]*a[[11]]-b[[15]]*a[[15]],
b[[2]]*a[[1]]+b[[1]]*a[[2]]-b[[6]]*a[[3]]-b[[7]]*a[[4]]-b[[8]]*a[[5]]+b[[3]]*a[[6]]+b[[4]]*a[[7]]+b[[5]]*a[[8]]+b[[12]]*a[[9]]+b[[13]]*a[[10]]+b[[14]]*a[[11]]+b[[9]]*a[[12]]+b[[10]]*a[[13]]+b[[11]]*a[[14]]+b[[16]]*a[[15]]-b[[15]]*a[[16]],
b[[3]]*a[[1]]+b[[1]]*a[[3]]-b[[9]]*a[[4]]+b[[10]]*a[[5]]+b[[4]]*a[[9]]-b[[5]]*a[[10]]-b[[15]]*a[[11]]-b[[11]]*a[[15]],
b[[4]]*a[[1]]+b[[9]]*a[[3]]+b[[1]]*a[[4]]-b[[11]]*a[[5]]-b[[3]]*a[[9]]-b[[15]]*a[[10]]+b[[5]]*a[[11]]-b[[10]]*a[[15]],
b[[5]]*a[[1]]-b[[10]]*a[[3]]+b[[11]]*a[[4]]+b[[1]]*a[[5]]-b[[15]]*a[[9]]+b[[3]]*a[[10]]-b[[4]]*a[[11]]-b[[9]]*a[[15]],
b[[6]]*a[[1]]-b[[12]]*a[[4]]+b[[13]]*a[[5]]+b[[1]]*a[[6]]-b[[16]]*a[[11]]-b[[4]]*a[[12]]+b[[5]]*a[[13]]-b[[11]]*a[[16]],
b[[7]]*a[[1]]+b[[12]]*a[[3]]-b[[14]]*a[[5]]+b[[1]]*a[[7]]-b[[16]]*a[[10]]+b[[3]]*a[[12]]-b[[5]]*a[[14]]-b[[10]]*a[[16]],
b[[8]]*a[[1]]-b[[13]]*a[[3]]+b[[14]]*a[[4]]+b[[1]]*a[[8]]-b[[16]]*a[[9]]-b[[3]]*a[[13]]+b[[4]]*a[[14]]-b[[9]]*a[[16]],
b[[9]]*a[[1]]+b[[15]]*a[[5]]+b[[1]]*a[[9]]+b[[5]]*a[[15]],
b[[10]]*a[[1]]+b[[15]]*a[[4]]+b[[1]]*a[[10]]+b[[4]]*a[[15]],
b[[11]]*a[[1]]+b[[15]]*a[[3]]+b[[1]]*a[[11]]+b[[3]]*a[[15]],
b[[12]]*a[[1]]+b[[16]]*a[[5]]+b[[1]]*a[[12]]-b[[5]]*a[[16]],
b[[13]]*a[[1]]+b[[16]]*a[[4]]+b[[1]]*a[[13]]-b[[4]]*a[[16]],
b[[14]]*a[[1]]+b[[16]]*a[[3]]+b[[1]]*a[[14]]-b[[3]]*a[[16]],
b[[15]]*a[[1]]+b[[1]]*a[[15]],
b[[16]]*a[[1]]+b[[1]]*a[[16]]]]
ClearAttributes[Dot,Flat]
tPieceValue/:Dot[a_tPieceValue,b_]:=Simplify[MapAt[Function[{x},Dot[x,b]],a,{All,1}]]
tPieceValue/:Dot[a_,b_tPieceValue]:=Simplify[MapAt[Function[{x},Dot[a,x]],b,{All,1}]]
tPGA3/:Dot[a_tPGA3,scalar:Except[_tPGA3|_tPieceValue|_Dot]]:=Dot[a,ScalarPGA3[scalar]]
tPGA3/:Dot[scalar:Except[_tPGA3|_tPieceValue|_Dot],b_tPGA3]:=Dot[ScalarPGA3[scalar],b]

(*Returns the commutator product of 2 multivectors.*)
tPGA3/:Commutator[a_tPGA3,b_tPGA3]:=(a**b-b**a)/2
tPieceValue/:Commutator[a_tPieceValue,b_]:=Simplify[MapAt[Function[{x},Commutator[x,b]],a,{All,1}]]
tPieceValue/:Commutator[a_,b_tPieceValue]:=Simplify[MapAt[Function[{x},Commutator[a,x]],b,{All,1}]]
Commutator/:Commutator[a_tPGA3,scalar:Except[_tPGA3|_tPieceValue|_Commutator]]:=Commutator[a,ScalarPGA3[scalar]]
Commutator/:Commutator[scalar:Except[_tPGA3|_tPieceValue|_Commutator],b_tPGA3]:=Commutator[ScalarPGA3[scalar],b]

(*Returns the square root of normalized motor.*)
tPGA3/:Sqrt[a_tPGA3]:=Simplify[
COND[IsMotor[a]&&IsNormalized[a]&&(a[[1]]!=-1),
Normalized[1+a]]]
tPieceValue/:Sqrt[a_tPieceValue]:=Simplify[MapAt[Sqrt,a,{All,1}]]

(*Returns the exponentiation of bivector.*)
tPGA3/:Exp[bi_tPGA3]:=Simplify[
COND[IsBivector[bi],
With[{k=bi[[9]]*bi[[9]]+bi[[10]]*bi[[10]]+bi[[11]]*bi[[11]]},
IF[k==0,
tPGA3[1,0,0,0,0,bi[[6]],bi[[7]],bi[[8]],0,0,0,0,0,0,0,0],
With[{m=bi[[6]]*bi[[11]]+bi[[7]]*bi[[10]]+bi[[8]]*bi[[9]],a=Sqrt[k]},
With[{c=Cos[a],s=Sin[a]/a},
With[{t=(m/k)*(c-s)},
tPGA3[c,0,0,0,0,s*bi[[6]]+t*bi[[11]],s*bi[[7]]+t*bi[[10]],s*bi[[8]]+t*bi[[9]],s*bi[[9]],s*bi[[10]],s*bi[[11]],0,0,0,0,m*s]]]]]]]]
tPieceValue/:Exp[a_tPieceValue]:=Simplify[MapAt[Exp,a,{All,1}]]

(*Returns the logarithm of normalized motor.*)
tPGA3/:Log[r_tPGA3]:=Simplify[
COND[IsMotor[r]&&IsNormalized[r]&&(r[[1]]!=-1),
IF[r[[1]]==1,
tPGA3[0,0,0,0,0,r[[6]],r[[7]],r[[8]],0,0,0,0,0,0,0,0],
With[{a=1/(1-r[[1]]*r[[1]])},
With[{b=ArcCos[r[[1]]]*Sqrt[a]},
With[{c=a*r[[16]]*(1-r[[1]]*b)},
tPGA3[0,0,0,0,0,b*r[[6]]+c*r[[11]],b*r[[7]]+c*r[[10]],b*r[[8]]+c*r[[9]],b*r[[9]],b*r[[10]],b*r[[11]],0,0,0,0,0]]]]]]]
tPieceValue/:Log[a_tPieceValue]:=Simplify[MapAt[Log,a,{All,1}]]

(*Returns the inertia mapping of bivector in eigenbasis body frame with inertia {i1,i2,i3,m}.*)
tPGA3/:DualInertia[inertia_List,v_tPGA3]:=Simplify[
COND[IsBivector[v],
With[{diagonal=Join[inertia,{inertia[[4]],inertia[[4]]}],in=Apply[List,Dual[v]][[6;;11]]},
With[{out=in*diagonal},
Apply[tPGA3,Join[{0,0,0,0,0},out,{0,0,0,0,0}]]]]]]
tPieceValue/:DualInertia[a_tPieceValue,b_]:=Simplify[MapAt[Function[{x},DualInertia[x,b]],a,{All,1}]]
tPieceValue/:DualInertia[a_,b_tPieceValue]:=Simplify[MapAt[Function[{x},DualInertia[a,x]],b,{All,1}]]
DualInertia/:DualInertia[inertia_List,scalar:Except[_tPGA3|_tPieceValue|_DualInertia]]:=DualInertia[inertia,ScalarPGA3[scalar]]

(*Returns the inverse inertia mapping of bivector in eigenbasis body frame with inertia {i1,i2,i3,m}.*)
tPGA3/:UnDualInertia[inertia_List,v_tPGA3]:=Simplify[
COND[IsBivector[v],
With[{diagonal=Join[inertia,{inertia[[4]],inertia[[4]]}],in=Apply[List,v][[6;;11]]},
With[{out=in/diagonal},
UnDual[Apply[tPGA3,Join[{0,0,0,0,0},out,{0,0,0,0,0}]]]]]]]
tPieceValue/:UnDualInertia[a_tPieceValue,b_]:=Simplify[MapAt[Function[{x},UnDualInertia[x,b]],a,{All,1}]]
tPieceValue/:UnDualInertia[a_,b_tPieceValue]:=Simplify[MapAt[Function[{x},UnDualInertia[a,x]],b,{All,1}]]
UnDualInertia/:UnDualInertia[inertia_List,scalar:Except[_tPGA3|_tPieceValue|_UnDualInertia]]:=UnDualInertia[inertia,ScalarPGA3[scalar]]

(*Returns the grade selection of multivector.*)
tPGA3/:Grade[grade_Integer,v_tPGA3]:=
If[TrueQ[grade==0],
tPGA3[v[[1]],0,0,0,0,0,0,0,0,0,0,0,0,0,0,0],
If[TrueQ[grade==1],
tPGA3[0,v[[2]],v[[3]],v[[4]],v[[5]],0,0,0,0,0,0,0,0,0,0,0],
If[TrueQ[grade==2],
tPGA3[0,0,0,0,0,v[[6]],v[[7]],v[[8]],v[[9]],v[[10]],v[[11]],0,0,0,0,0],
If[TrueQ[grade==3],
tPGA3[0,0,0,0,0,0,0,0,0,0,0,v[[12]],v[[13]],v[[14]],v[[15]],0],
If[TrueQ[grade==4],
tPGA3[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,v[[16]]],
Echo[Undefined,"ERROR: grade_Integer"]]]]]]
tPieceValue/:Grade[grade_Integer,v_tPieceValue]:=MapAt[Function[{x},Grade[grade,x]],v,{All,1}]
Grade/:Grade[grade_Integer,scalar:Except[_tPGA3|_tPieceValue|_Grade]]:=Grade[grade,ScalarPGA3[scalar]]

(*Returns the ideal part selection of multivector.*)
tPGA3/:IdealPart[a_tPGA3]:=tPGA3[0,a[[2]],0,0,0,a[[6]],a[[7]],a[[8]],0,0,0,a[[12]],a[[13]],a[[14]],0,a[[16]]]
tPieceValue/:IdealPart[a_tPieceValue]:=MapAt[IdealPart,a,{All,1}]
IdealPart/:IdealPart[scalar:Except[_tPGA3|_tPieceValue|_IdealPart]]:=IdealPart[ScalarPGA3[scalar]]

(*Returns the origin part selection of multivector.*)
tPGA3/:OriginPart[a_tPGA3]:=tPGA3[a[[1]],0,a[[3]],a[[4]],a[[5]],0,0,0,a[[9]],a[[10]],a[[11]],0,0,0,a[[15]],0]
tPieceValue/:OriginPart[a_tPieceValue]:=MapAt[OriginPart,a,{All,1}]
OriginPart/:OriginPart[scalar:Except[_tPGA3|_tPieceValue|_OriginPart]]:=OriginPart[ScalarPGA3[scalar]]

(*Returns the motor part selection of multivector.*)
tPGA3/:MotorPart[a_tPGA3]:=tPGA3[a[[1]],0,0,0,0,a[[6]],a[[7]],a[[8]],a[[9]],a[[10]],a[[11]],0,0,0,0,a[[16]]]
tPieceValue/:MotorPart[a_tPieceValue]:=MapAt[MotorPart,a,{All,1}]
MotorPart/:MotorPart[scalar:Except[_tPGA3|_tPieceValue|_MotorPart]]:=MotorPart[ScalarPGA3[scalar]]

(*Returns the study number part selection of multivector.*)
tPGA3/:StudyPart[a_tPGA3]:=tPGA3[a[[1]],0,0,0,0,0,0,0,0,0,0,0,0,0,0,a[[16]]]
tPieceValue/:StudyPart[a_tPieceValue]:=MapAt[StudyPart,a,{All,1}]
StudyPart/:StudyPart[scalar:Except[_tPGA3|_tPieceValue|_StudyPart]]:=StudyPart[ScalarPGA3[scalar]]

(*Returns the square root of study number.*)
tPGA3/:StudySqrt[a_tPGA3]:=Simplify[With[{first=First[a],last=Last[a]},
COND[IsStudy[a],
If[TrueQ[last==0],
COND[0<=first,
tPGA3[Sqrt[first],0,0,0,0,0,0,0,0,0,0,0,0,0,0,0]],
IF[0<first,
tPGA3[Sqrt[first],0,0,0,0,0,0,0,0,0,0,0,0,0,0,last/(2*Sqrt[first])],
COND[(0<=first)&&(first==0)&&(last==0),
tPGA3[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0]]]]]]]
tPieceValue/:StudySqrt[a_tPieceValue]:=Simplify[MapAt[StudySqrt,a,{All,1}]]
StudySqrt/:StudySqrt[scalar:Except[_tPGA3|_tPieceValue|_StudySqrt]]:=StudySqrt[ScalarPGA3[scalar]]

(*Returns the inverse of study number.*)
tPGA3/:StudyInverse[a_tPGA3]:=Simplify[With[{first=First[a],last=Last[a]},
COND[IsStudy[a]&&(first!=0),
tPGA3[1/first,0,0,0,0,0,0,0,0,0,0,0,0,0,0,-last/(first^2)]]]]
tPieceValue/:StudyInverse[a_tPieceValue]:=Simplify[MapAt[StudyInverse,a,{All,1}]]
StudyInverse/:StudyInverse[scalar:Except[_tPGA3|_tPieceValue|_StudyInverse]]:=StudyInverse[ScalarPGA3[scalar]]

(*Returns the squared norm of multivector.*)
tPGA3/:SquaredNorm[a_tPGA3]:=Reverse[a]**a
tPieceValue/:SquaredNorm[a_tPieceValue]:=MapAt[SquaredNorm,a,{All,1}]
SquaredNorm/:SquaredNorm[scalar:Except[_tPGA3|_tPieceValue|_SquaredNorm]]:=SquaredNorm[ScalarPGA3[scalar]]

(*Returns the norm of multivector.*)
tPGA3/:Norm[a_tPGA3]:=StudySqrt[SquaredNorm[a]]
tPieceValue/:Norm[a_tPieceValue]:=Simplify[MapAt[Norm,a,{All,1}]]

(*Returns the normalized multivector such that Normalized[a]**Norm[a]==a and Norm[Normalized[a]]==1.*)
tPGA3/:Normalized[a_tPGA3]:=a**StudyInverse[Norm[a]]
tPieceValue/:Normalized[a_tPieceValue]:=Simplify[MapAt[Normalized,a,{All,1}]]

(*Returns the ideal norm of multivector.*)
tPGA3/:INorm[a_tPGA3]:=Norm[Dual[a]]
tPieceValue/:INorm[a_tPieceValue]:=Simplify[MapAt[INorm,a,{All,1}]]
INorm/:INorm[scalar:Except[_tPGA3|_tPieceValue|_INorm]]:=INorm[ScalarPGA3[scalar]]

(*Returns the cross product of 2 multivectors.*)
tPGA3/:Cross[a_tPGA3,b_tPGA3]:=Meet[a,b]**Inverse[e123]
tPieceValue/:Cross[a_tPieceValue,b_]:=Simplify[MapAt[Function[{x},Cross[x,b]],a,{All,1}]]
tPieceValue/:Cross[a_,b_tPieceValue]:=Simplify[MapAt[Function[{x},Cross[a,x]],b,{All,1}]]
tPGA3/:Cross[a_tPGA3,scalar:Except[_tPGA3|_tPieceValue|_Cross]]:=Cross[a,ScalarPGA3[scalar]]
tPGA3/:Cross[scalar:Except[_tPGA3|_tPieceValue|_Cross],b_tPGA3]:=Cross[ScalarPGA3[scalar],b]

(*Tests if 2 multivectors are equal.*)
tPGA3/:a_tPGA3==b_tPGA3:=
Simplify[a[[1]]==b[[1]]&&a[[2]]==b[[2]]&&a[[3]]==b[[3]]&&a[[4]]==b[[4]]&&a[[5]]==b[[5]]&&a[[6]]==b[[6]]&&a[[7]]==b[[7]]&&a[[8]]==b[[8]]&&a[[9]]==b[[9]]&&a[[10]]==b[[10]]&&a[[11]]==b[[11]]&&a[[12]]==b[[12]]&&a[[13]]==b[[13]]&&a[[14]]==b[[14]]&&a[[15]]==b[[15]]&&a[[16]]==b[[16]]]
tPieceValue/:Equal[a_tPieceValue,b_]:=Apply[tPieceBoolean,MapAt[Function[{x},Equal[x,b]],a,{All,1}]]
tPieceValue/:Equal[a_,b_tPieceValue]:=Apply[tPieceBoolean,MapAt[Function[{x},Equal[a,x]],b,{All,1}]]
tPGA3/:Equal[a_tPGA3,scalar:Except[_tPGA3|_tPieceValue|_tPieceBoolean|_Equal|_Unequal|_Not|_And|_Or]]:=Equal[a,ScalarPGA3[scalar]]
tPGA3/:Equal[scalar:Except[_tPGA3|_tPieceValue|_tPieceBoolean|_Equal|_Unequal|_Not|_And|_Or],b_tPGA3]:=Equal[ScalarPGA3[scalar],b]

(*Tests if 2 multivectors are not equal.*)
tPGA3/:a_tPGA3!=b_tPGA3:=
Simplify[a[[1]]!=b[[1]]||a[[2]]!=b[[2]]||a[[3]]!=b[[3]]||a[[4]]!=b[[4]]||a[[5]]!=b[[5]]||a[[6]]!=b[[6]]||a[[7]]!=b[[7]]||a[[8]]!=b[[8]]||a[[9]]!=b[[9]]||a[[10]]!=b[[10]]||a[[11]]!=b[[11]]||a[[12]]!=b[[12]]||a[[13]]!=b[[13]]||a[[14]]!=b[[14]]||a[[15]]!=b[[15]]||a[[16]]!=b[[16]]]
tPieceValue/:Unequal[a_tPieceValue,b_]:=Apply[tPieceBoolean,MapAt[Function[{x},Unequal[x,b]],a,{All,1}]]
tPieceValue/:Unequal[a_,b_tPieceValue]:=Apply[tPieceBoolean,MapAt[Function[{x},Unequal[a,x]],b,{All,1}]]
tPGA3/:Unequal[a_tPGA3,scalar:Except[_tPGA3|_tPieceValue|_tPieceBoolean|_Equal|_Unequal|_Not|_And|_Or]]:=Unequal[a,ScalarPGA3[scalar]]
tPGA3/:Unequal[scalar:Except[_tPGA3|_tPieceValue|_tPieceBoolean|_Equal|_Unequal|_Not|_And|_Or],b_tPGA3]:=Unequal[ScalarPGA3[scalar],b]

(*Tests if a < b.*)
tPGA3/:a_tPGA3<b_tPGA3:=Apply[tPieceBoolean,
COND[IsScalar[a]&&IsScalar[b],
a[[1]]<b[[1]]]]
tPieceValue/:Less[a_tPieceValue,b_]:=Apply[tPieceBoolean,MapAt[Function[{x},Less[x,b]],a,{All,1}]]
tPieceValue/:Less[a_,b_tPieceValue]:=Apply[tPieceBoolean,MapAt[Function[{x},Less[a,x]],b,{All,1}]]
tPGA3/:Less[a_tPGA3,scalar:Except[_tPGA3|_tPieceValue|_tPieceBoolean|_Less|_Equal|_Unequal|_Not|_And|_Or]]:=Less[a,ScalarPGA3[scalar]]
tPGA3/:Less[scalar:Except[_tPGA3|_tPieceValue|_tPieceBoolean|_Less|_Equal|_Unequal|_Not|_And|_Or],b_tPGA3]:=Less[ScalarPGA3[scalar],b]

(*Tests if a <= b.*)
tPGA3/:a_tPGA3<=b_tPGA3:=Apply[tPieceBoolean,
COND[IsScalar[a]&&IsScalar[b],
a[[1]]<=b[[1]]]]
tPieceValue/:LessEqual[a_tPieceValue,b_]:=Apply[tPieceBoolean,MapAt[Function[{x},LessEqual[x,b]],a,{All,1}]]
tPieceValue/:LessEqual[a_,b_tPieceValue]:=Apply[tPieceBoolean,MapAt[Function[{x},LessEqual[a,x]],b,{All,1}]]
tPGA3/:LessEqual[a_tPGA3,scalar:Except[_tPGA3|_tPieceValue|_tPieceBoolean|_LessEqual|_Equal|_Unequal|_Not|_And|_Or]]:=LessEqual[a,ScalarPGA3[scalar]]
tPGA3/:LessEqual[scalar:Except[_tPGA3|_tPieceValue|_tPieceBoolean|_LessEqual|_Equal|_Unequal|_Not|_And|_Or],b_tPGA3]:=LessEqual[ScalarPGA3[scalar],b]

(*Tests if the multivector is normalized.*)
tPGA3/:IsNormalized[a_tPGA3]:=(SquaredNorm[a]==1)
tPieceValue/:IsNormalized[a_tPieceValue]:=Apply[tPieceBoolean,MapAt[IsNormalized,a,{All,1}]]
IsNormalized/:IsNormalized[scalar:Except[_tPGA3|_tPieceValue|_IsNormalized]]:=IsNormalized[ScalarPGA3[scalar]]

(*Tests if the multivector is a blade.*)
tPGA3/:IsBlade[a_tPGA3]:=(Meet[a,a]==0)
tPieceValue/:IsBlade[a_tPieceValue]:=Apply[tPieceBoolean,MapAt[IsBlade,a,{All,1}]]
IsBlade/:IsBlade[scalar:Except[_tPGA3|_tPieceValue|_IsBlade]]:=IsBlade[ScalarPGA3[scalar]]

(*Tests if the multivector is a scalar.*)
tPGA3/:IsScalar[a_tPGA3]:=(Grade[0,a]==a)
tPieceValue/:IsScalar[a_tPieceValue]:=Apply[tPieceBoolean,MapAt[IsScalar,a,{All,1}]]
IsScalar/:IsScalar[scalar:Except[_tPGA3|_tPieceValue|_IsScalar]]:=IsScalar[ScalarPGA3[scalar]]

(*Tests if the multivector is a plane.*)
tPGA3/:IsPlane[a_tPGA3]:=(Grade[1,a]==a)
tPieceValue/:IsPlane[a_tPieceValue]:=Apply[tPieceBoolean,MapAt[IsPlane,a,{All,1}]]
IsPlane/:IsPlane[scalar:Except[_tPGA3|_tPieceValue|_IsPlane]]:=IsPlane[ScalarPGA3[scalar]]

(*Tests if the multivector is a bivector.*)
tPGA3/:IsBivector[a_tPGA3]:=(Grade[2,a]==a)
tPieceValue/:IsBivector[a_tPieceValue]:=Apply[tPieceBoolean,MapAt[IsBivector,a,{All,1}]]
IsBivector/:IsBivector[scalar:Except[_tPGA3|_tPieceValue|_IsBivector]]:=IsBivector[ScalarPGA3[scalar]]

(*Tests if the multivector is a point.*)
tPGA3/:IsPoint[a_tPGA3]:=(Grade[3,a]==a)
tPieceValue/:IsPoint[a_tPieceValue]:=Apply[tPieceBoolean,MapAt[IsPoint,a,{All,1}]]
IsPoint/:IsPoint[scalar:Except[_tPGA3|_tPieceValue|_IsPoint]]:=IsPoint[ScalarPGA3[scalar]]

(*Tests if the multivector is a pseudoscalar.*)
tPGA3/:IsPseudoscalar[a_tPGA3]:=(Grade[4,a]==a)
tPieceValue/:IsPseudoscalar[a_tPieceValue]:=Apply[tPieceBoolean,MapAt[IsPseudoscalar,a,{All,1}]]
IsPseudoscalar/:IsPseudoscalar[scalar:Except[_tPGA3|_tPieceValue|_IsPseudoscalar]]:=IsPseudoscalar[ScalarPGA3[scalar]]

(*Tests if the multivector is a motor.*)
tPGA3/:IsMotor[a_tPGA3]:=(MotorPart[a]==a)
tPieceValue/:IsMotor[a_tPieceValue]:=Apply[tPieceBoolean,MapAt[IsMotor,a,{All,1}]]
IsMotor/:IsMotor[scalar:Except[_tPGA3|_tPieceValue|_IsMotor]]:=IsMotor[ScalarPGA3[scalar]]

(*Tests if the multivector is a study number.*)
tPGA3/:IsStudy[a_tPGA3]:=(StudyPart[a]==a)
tPieceValue/:IsStudy[a_tPieceValue]:=Apply[tPieceBoolean,MapAt[IsStudy,a,{All,1}]]
IsStudy/:IsStudy[scalar:Except[_tPGA3|_tPieceValue|_IsStudy]]:=IsStudy[ScalarPGA3[scalar]]

(*Returns the full multivector with specified coefficient names.*)
full=Function[{name},tPGA3[
Symbol[SymbolName[name]<>"1"],
Symbol[SymbolName[name]<>"2"],
Symbol[SymbolName[name]<>"3"],
Symbol[SymbolName[name]<>"4"],
Symbol[SymbolName[name]<>"5"],
Symbol[SymbolName[name]<>"6"],
Symbol[SymbolName[name]<>"7"],
Symbol[SymbolName[name]<>"8"],
Symbol[SymbolName[name]<>"9"],
Symbol[SymbolName[name]<>"10"],
Symbol[SymbolName[name]<>"11"],
Symbol[SymbolName[name]<>"12"],
Symbol[SymbolName[name]<>"13"],
Symbol[SymbolName[name]<>"14"],
Symbol[SymbolName[name]<>"15"],
Symbol[SymbolName[name]<>"16"]]];

(*Returns the plane with specified coefficient names.*)
plane=Function[{name},Grade[1,full[name]]];

(*Returns the normalized plane with specified coefficient names.*)
planeN=Function[{name},With[{v=plane[name]},COND[IsNormalized[v],v]]];

(*Returns the plane through the origin with specified coefficient names.*)
planeO=Function[{name},OriginPart[plane[name]]];

(*Returns the normalized plane through the origin with specified coefficient names.*)
planeON=Function[{name},With[{v=planeO[name]},COND[IsNormalized[v],v]]];

(*Returns the bivector with specified coefficient names.*)
bivector=Function[{name},Grade[2,full[name]]];

(*Returns the line with specified coefficient names.*)
line=Function[{name},With[{v=bivector[name]},COND[IsBlade[v],v]]];

(*Returns the normalized line with specified coefficient names.*)
lineN=Function[{name},With[{v=bivector[name]},COND[IsNormalized[v],v]]];

(*Returns the ideal line with specified coefficient names.*)
lineI=Function[{name},IdealPart[bivector[name]]];

(*Returns the ideal line with unit ideal norm and specified coefficient names.*)
lineIN=Function[{name},With[{v=lineI[name]},COND[IsNormalized[Dual[v]],v]]];

(*Returns the line through the origin with specified coefficient names.*)
lineO=Function[{name},OriginPart[bivector[name]]];

(*Returns the normalized line through the origin with specified coefficient names.*)
lineON=Function[{name},With[{v=lineO[name]},COND[IsNormalized[v],v]]];

(*Returns the point with specified coefficient names.*)
point=Function[{name},Grade[3,full[name]]];

(*Returns the normalized point with specified coefficient names.*)
pointN=Function[{name},IdealPart[point[name]]+e123];

(*Returns the ideal point with specified coefficient names.*)
pointI=Function[{name},IdealPart[point[name]]];

(*Returns the ideal point with unit ideal norm and specified coefficient names.*)
pointIN=Function[{name},With[{v=pointI[name]},COND[IsNormalized[Dual[v]],v]]];

(*Returns the motor with specified coefficient names.*)
motor=Function[{name},MotorPart[full[name]]];

(*Returns the normalized motor with specified coefficient names.*)
motorN=Function[{name},With[{v=motor[name]},COND[IsNormalized[v],v]]];

(*Returns the study number with specified coefficient names.*)
study=Function[{name},StudyPart[full[name]]];

(*Returns simplified multivector over real numbers.*)
SimplifyReals/:SimplifyReals[v_tPieceValue]:=With[{processed=MapAt[SimplifyReals,v,{All,2}]},
If[TrueQ[Head[processed]==tPieceValue],
MapAt[Function[{x},{Assuming[x[[2]],SimplifyReals[x[[1]]]],x[[2]]}],processed,All],
SimplifyReals[processed]]]
SimplifyReals/:SimplifyReals[v_List]:=MapAt[SimplifyReals,v,All]
SimplifyReals/:SimplifyReals[v_]:=Assuming[Element[_Symbol,Reals],
With[{expanded=ExpandAll[LogicalExpand[v]]},
ReplaceAll[{Greater[a_,b_]:>Less[b,a],GreaterEqual[a_,b_]:>LessEqual[b,a]}][FullSimplify[
If[TrueQ[Head[expanded]==Or],
MapAt[FullSimplify,expanded,All],
expanded]]]]]

(*Shows the multivector in table form.*)
tPGA3/:MatrixForm[a_tPGA3]:=MatrixForm[Transpose[{{"1","e0","e1","e2","e3","e01","e02","e03","e12","e31","e23","e021","e013","e032","e123","e0123"},Apply[List,a]}]]
tPieceValue/:MatrixForm[a_tPieceValue]:=MatrixForm[MapAt[MatrixForm,Apply[List,a],{All,1}]]
tPieceBoolean/:MatrixForm[a_tPieceBoolean]:=MatrixForm[Apply[List,a]]
ScalarPGA3=Function[{scalar},tPGA3[scalar,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0]];
Unprotect[Real];
Real/:Equal[a_Real,b_]:=(Abs[a-b]<0.0000001)
Real/:Equal[a_,b_Real]:=(Abs[a-b]<0.0000001)
Protect[Real];
SetAttributes[tPieceValue,Orderless]
tPieceValue/:Simplify[v_tPieceValue]:=MapAt[Function[{x},{Assuming[x[[2]],Simplify[x[[1]]]],x[[2]]}],v,All]
tPGA3/:Greater[a_tPGA3,b_]:=Less[b,a]
tPGA3/:Greater[a_,b_tPGA3]:=Less[b,a]
tPieceValue/:Greater[a_tPieceValue,b_]:=Less[b,a]
tPieceValue/:Greater[a_,b_tPieceValue]:=Less[b,a]
tPGA3/:GreaterEqual[a_tPGA3,b_]:=LessEqual[b,a]
tPGA3/:GreaterEqual[a_,b_tPGA3]:=LessEqual[b,a]
tPieceValue/:GreaterEqual[a_tPieceValue,b_]:=LessEqual[b,a]
tPieceValue/:GreaterEqual[a_,b_tPieceValue]:=LessEqual[b,a]
tPieceValue/:tPieceValue[{v_,True}]:=v
tPieceValue/:tPieceValue[{_,False},t___]:=tPieceValue[t]
tPieceValue/:tPieceValue[{v_,c1_},{v_,c2_},t___]:=tPieceValue[{v,Simplify[c1||c2]},t]
tPieceValue/:tPieceValue[{v_tPieceValue,c_},t___]:=tPieceValue[Apply[Sequence,MapAt[Function[{x},Simplify[x&&c]],v,{All,2}]],t]
tPieceValue/:tPieceValue[{v_Piecewise,c_},t___]:=With[{orderless=PiecewiseExpand[v,Method->{"OrderlessConditions"->True}]},
tPieceValue[{Apply[tPieceValue][Append[orderless[[1]],{orderless[[2]],Simplify[Apply[And][MapAt[Function[{x},Not[x[[2]]]],orderless[[1]],All]]]}]],c},t]]
SetAttributes[tPieceBoolean,Orderless]
tPieceBoolean/:Not[a_tPieceBoolean]:=MapAt[Function[{x},Not[x]],a,{All,1}]
tPieceBoolean/:And[a_tPieceBoolean,b_]:=MapAt[Function[{x},Simplify[And[x,b]]],a,{All,1}]
tPieceBoolean/:And[a_,b_tPieceBoolean]:=MapAt[Function[{x},Simplify[And[a,x]]],b,{All,1}]
tPieceBoolean/:Or[a_tPieceBoolean,b_]:=MapAt[Function[{x},Simplify[Or[x,b]]],a,{All,1}]
tPieceBoolean/:Or[a_,b_tPieceBoolean]:=MapAt[Function[{x},Simplify[Or[a,x]]],b,{All,1}]
tPieceBoolean/:tPieceBoolean[{v_tPieceBoolean,c_},t___]:=tPieceBoolean[Apply[Sequence,MapAt[Function[{x},Simplify[x&&c]],v,{All,2}]],t]
Bake=Function[{a},
If[TrueQ[Head[a]==tPieceBoolean],
Simplify[Apply[Or,MapAt[Apply[And],a,All]]],
Simplify[a]]];
IF=Function[{c,t,f},tPieceValue[{t,Bake[c]},{f,Bake[!c]}]];
COND=Function[{c,t},tPieceValue[{t,Bake[c]}]];
1 Like