1  """
  2  Regression tests validating galgebra against Russell Goyder's geometric
  3  algebra cheat sheet (https://russellgoyder.github.io/geometric-algebra-cheat-sheet/).
  4  
  5  Covers all 7 sections in G(3,0): geometric product, vector-multivector
  6  products, bivector identities, commutator product, pseudoscalar,
  7  outermorphism/rotation, and adjoint.
  8  """
  9  import pytest
 10  from sympy import symbols, cos, sin, S, simplify, trigsimp
 11  
 12  from galgebra.ga import Ga
 13  
 14  
 15  @pytest.fixture(scope='module')
 16  def g3():
 17      ga, e1, e2, e3 = Ga.build('e*1|2|3', g=[1, 1, 1])
 18      return ga, e1, e2, e3
 19  
 20  
 21  @pytest.fixture(scope='module')
 22  def vecs(g3):
 23      ga, e1, e2, e3 = g3
 24      a1, a2, a3 = symbols('a1 a2 a3', real=True)
 25      b1, b2, b3 = symbols('b1 b2 b3', real=True)
 26      c1, c2, c3 = symbols('c1 c2 c3', real=True)
 27      d1, d2, d3 = symbols('d1 d2 d3', real=True)
 28      a = a1*e1 + a2*e2 + a3*e3
 29      b = b1*e1 + b2*e2 + b3*e3
 30      c = c1*e1 + c2*e2 + c3*e3
 31      d = d1*e1 + d2*e2 + d3*e3
 32      return a, b, c, d
 33  
 34  
 35  def test_geometric_product_decomposition(g3, vecs):
 36      """Section 1: ab = a·b + a∧b, inner/outer as symmetric/antisymmetric parts."""
 37      ga, e1, e2, e3 = g3
 38      a, b, c, d = vecs
 39  
 40      assert simplify((a*b - (a|b) - (a^b)).obj) == 0
 41      assert simplify(((a|b) - ga.mv(S.Half)*(a*b + b*a)).obj) == 0
 42      assert simplify(((a^b) - ga.mv(S.Half)*(a*b - b*a)).obj) == 0
 43  
 44  
 45  def test_vector_multivector_products(g3, vecs):
 46      """Section 2: a*B = a<B + a∧B (grade raising/lowering)."""
 47      ga, e1, e2, e3 = g3
 48      a, b, c, d = vecs
 49  
 50      B = b ^ c
 51      assert simplify((a*B - (a < B) - (a ^ B)).obj) == 0
 52  
 53      T = b ^ c ^ d
 54      assert simplify((a*T - (a < T) - (a ^ T)).obj) == 0
 55  
 56  
 57  def test_bivector_identities(g3, vecs):
 58      """Section 3: a·(b∧c), (a∧b)·B associativity, (a∧b)·(c∧d) scalar formula."""
 59      ga, e1, e2, e3 = g3
 60      a, b, c, d = vecs
 61  
 62      # a·(b∧c) = (a·b)c - (a·c)b
 63      lhs1 = a < (b^c)
 64      rhs1 = (a|b)*c - (a|c)*b
 65      assert simplify((lhs1 - rhs1).obj) == 0
 66  
 67      # (a∧b)·B = a·(b·B)  (left contraction associates)
 68      lhs2 = (a^b) < (c^d)
 69      rhs2 = a < (b < (c^d))
 70      assert simplify((lhs2 - rhs2).obj) == 0
 71  
 72      # (a∧b)·(c∧d) = (b·c)(a·d) - (a·c)(b·d)
 73      lhs3 = (a^b) | (c^d)
 74      rhs3 = (b|c)*(a|d) - (a|c)*(b|d)
 75      assert simplify((lhs3 - rhs3).obj) == 0
 76  
 77  
 78  def test_commutator_product(g3, vecs):
 79      """Section 4: B×v = B·v, grade-preserving on bivectors, Jacobi identity."""
 80      ga, e1, e2, e3 = g3
 81      a, b, c, d = vecs
 82  
 83      def comm(X, Y):
 84          return ga.mv(S.Half) * (X*Y - Y*X)
 85  
 86      B2 = e1 ^ e2
 87      for v in [e1, e2, e3, a]:
 88          assert simplify((comm(B2, v) - (B2 | v)).obj) == 0
 89  
 90      B3 = e2 ^ e3
 91      assert (comm(B2, B3)).pure_grade() == 2
 92  
 93      B4 = e1 ^ e3
 94      jacobi = comm(B2, comm(B3, B4)) + comm(B3, comm(B4, B2)) + comm(B4, comm(B2, B3))
 95      assert simplify(jacobi.obj) == 0
 96  
 97  
 98  def test_pseudoscalar(g3):
 99      """Section 5: I²=-1, I central in G(3,0), dual via I⁻¹."""
100      ga, e1, e2, e3 = g3
101      I = ga.I()
102  
103      assert (I*I).scalar() == S.NegativeOne
104  
105      # (-1)^(r*(n-1)) = (-1)^(2r) = 1 for all r in 3D: I is central
106      for mv in [e1, e1^e2, I]:
107          assert simplify((I*mv - mv*I).obj) == 0
108  
109      # dual(A) = A * I^{-1}  (Iinv+ convention, not galgebra's built-in dual()
110      # which uses I+ i.e. dual(A) = I*A — a different sign convention)
111      assert simplify((I.inv() + I).obj) == 0   # I^{-1} = -I when I^2=-1
112      assert e1 * I.inv() == -(e2 ^ e3)
113  
114  
115  def test_outermorphism_rotation(g3):
116      """Section 6: rotor normalization, rotation sandwich, det(rotation)=1."""
117      ga, e1, e2, e3 = g3
118      theta = symbols('theta', real=True)
119      R = ga.mv(cos(theta/2)) - ga.mv(sin(theta/2)) * (e1^e2)
120  
121      assert trigsimp((R * ~R).scalar()) == 1
122  
123      rotated = ga.mv(trigsimp((R * e1 * ~R).obj))
124      assert trigsimp((rotated - ga.mv(cos(theta))*e1 - ga.mv(sin(theta))*e2).obj) == 0
125  
126      L = ga.lt(lambda x: R*x*~R)
127      assert trigsimp(L.det()) == 1
128  
129  
130  def test_adjoint(g3, vecs):
131      """Section 7: a·F̄(b) = F(a)·b and det(FG) = det(F)·det(G)."""
132      ga, e1, e2, e3 = g3
133      a, b, c, d = vecs
134  
135      L = ga.lt([e1, e1 + e2, e3])
136      assert simplify((a | L(b) - L.adj()(a) | b).obj) == 0
137  
138      theta, phi = symbols('theta phi', real=True)
139      R1 = ga.mv(cos(theta/2)) - ga.mv(sin(theta/2)) * (e1^e2)
140      R2 = ga.mv(cos(phi/2)) - ga.mv(sin(phi/2)) * (e2^e3)
141      L1 = ga.lt(lambda x: R1*x*~R1)
142      L2 = ga.lt(lambda x: R2*x*~R2)
143      L12 = ga.lt(lambda x: L1(L2(x)))
144      assert trigsimp(L12.det() - L1.det()*L2.det()) == 0