test_mv.py
1 from packaging.version import Version 2 3 import pytest 4 import sympy 5 from sympy import symbols, S 6 from galgebra.ga import Ga 7 from galgebra.mv import proj, undual, g_invol, exp, norm, norm2, mag, mag2, ccon, rev, scalar, qform, sp, inv, shirokov_inverse, hitzer_inverse 8 9 10 class TestMv: 11 12 def test_deprecations(self): 13 ga, e_1, e_2, e_3 = Ga.build('e*1|2|3') 14 with pytest.warns(DeprecationWarning): 15 ga.mv_x 16 with pytest.warns(DeprecationWarning): 17 ga.mv_I 18 19 def test_is_base(self): 20 """ 21 Various tests on several multivectors. 22 """ 23 _g3d, e_1, e_2, e_3 = Ga.build('e*1|2|3') 24 25 assert (e_1).is_base() 26 assert (e_2).is_base() 27 assert (e_3).is_base() 28 assert (e_1 ^ e_2).is_base() 29 assert (e_2 ^ e_3).is_base() 30 assert (e_1 ^ e_3).is_base() 31 assert (e_1 ^ e_2 ^ e_3).is_base() 32 33 assert not (2 * e_1).is_base() 34 assert not (e_1 + e_2).is_base() 35 assert not (e_3 * 4).is_base() 36 assert not ((3 * e_1) ^ e_2).is_base() 37 assert not (2 * (e_2 ^ e_3)).is_base() 38 assert not (e_3 ^ e_1).is_base() 39 assert not (e_2 ^ e_1 ^ e_3).is_base() 40 41 def test_get_coefs(self): 42 _g3d, e_1, e_2, e_3 = Ga.build('e*1|2|3') 43 44 assert (e_1 * 3 + e_3).get_coefs(1) == [3, 0, 1] 45 46 # can always get coefficients of 0 47 assert (0*e_1).get_coefs(0) == [0] 48 assert (0*e_1).get_coefs(1) == [0, 0, 0] 49 assert (0*e_1).get_coefs(2) == [0, 0, 0] 50 assert (0*e_1).get_coefs(3) == [0] 51 52 # grade is wrong 53 with pytest.raises(ValueError): 54 (e_1 ^ e_2).get_coefs(1) 55 with pytest.raises(ValueError): 56 (e_1 ^ e_2).get_coefs(3) 57 58 def test_blade_coefs(self): 59 """ 60 Various tests on several multivectors. 61 """ 62 _g3d, e_1, e_2, e_3 = Ga.build('e*1|2|3') 63 64 m0 = 2 * e_1 + e_2 - e_3 + 3 * (e_1 ^ e_3) + (e_1 ^ e_3) + (e_2 ^ (3 * e_3)) 65 assert m0.blade_coefs([e_1]) == [2] 66 assert m0.blade_coefs([e_2]) == [1] 67 assert m0.blade_coefs([e_1, e_2]) == [2, 1] 68 assert m0.blade_coefs([e_1 ^ e_3]) == [4] 69 assert m0.blade_coefs([e_1 ^ e_3, e_2 ^ e_3]) == [4, 3] 70 assert m0.blade_coefs([e_2 ^ e_3, e_1 ^ e_3]) == [3, 4] 71 assert m0.blade_coefs([e_1, e_2 ^ e_3]) == [2, 3] 72 73 a = sympy.Symbol('a') 74 b = sympy.Symbol('b') 75 m1 = a * e_1 + e_2 - e_3 + b * (e_1 ^ e_2) 76 assert m1.blade_coefs([e_1]) == [a] 77 assert m1.blade_coefs([e_2]) == [1] 78 assert m1.blade_coefs([e_3]) == [-1] 79 assert m1.blade_coefs([e_1, e_2, e_3]) == [a, 1, -1] 80 assert m1.list() == [a, 1, -1] # alias 81 82 assert m1.blade_coefs([e_1 ^ e_2]) == [b] 83 assert m1.blade_coefs([e_2 ^ e_3]) == [0] 84 assert m1.blade_coefs([e_1 ^ e_3]) == [0] 85 assert m1.blade_coefs([e_1 ^ e_2 ^ e_3]) == [0] 86 87 # Invalid parameters 88 pytest.raises(ValueError, lambda: m1.blade_coefs([e_1 + e_2])) 89 pytest.raises(ValueError, lambda: m1.blade_coefs([e_2 ^ e_1])) 90 pytest.raises(ValueError, lambda: m1.blade_coefs([e_1, e_2 ^ e_1])) 91 pytest.raises(ValueError, lambda: m1.blade_coefs([a * e_1])) 92 pytest.raises(ValueError, lambda: m1.blade_coefs([3 * e_3])) 93 94 def test_rep_switching(self): 95 # this ga has a non-diagonal metric 96 _g3d, e_1, e_2, e_3 = Ga.build('e*1|2|3') 97 98 m0 = 2 * e_1 + e_2 - e_3 + 3 * (e_1 ^ e_3) + (e_1 ^ e_3) + (e_2 ^ (3 * e_3)) 99 m1 = (-4*(e_1 | e_3)-3*(e_2 | e_3))+2*e_1+e_2-e_3+4*e_1*e_3+3*e_2*e_3 100 # m1 was chosen to make this true 101 assert m0 == m1 102 103 # all objects start off in blade rep 104 assert m0.is_blade_rep 105 106 # convert to base rep 107 m0_base = m0.base_rep() 108 assert m0.is_blade_rep # original should not change 109 assert not m0_base.is_blade_rep 110 assert m0 == m0_base 111 112 # convert back 113 m0_base_blade = m0_base.blade_rep() 114 assert not m0_base.is_blade_rep # original should not change 115 assert m0_base_blade.is_blade_rep 116 assert m0 == m0_base_blade 117 118 def test_construction(self): 119 ga, e_1, e_2, e_3 = Ga.build('e*1|2|3') 120 121 def check(x, expected_grades): 122 assert x.grades == expected_grades 123 assert x != 0 124 125 # non-function symbol construction 126 check(ga.mv('A', 'scalar'), [0]) 127 check(ga.mv('A', 'grade', 0), [0]) 128 check(ga.mv('A', 0), [0]) 129 check(ga.mv('A', 'vector'), [1]) 130 check(ga.mv('A', 'grade', 1), [1]) 131 check(ga.mv('A', 1), [1]) 132 check(ga.mv('A', 'bivector'), [2]) 133 check(ga.mv('A', 'grade2'), [2]) 134 check(ga.mv('A', 2), [2]) 135 check(ga.mv('A', 'pseudo'), [3]) 136 check(ga.mv('A', 'spinor'), [0, 2]) 137 check(ga.mv('A', 'even'), [0, 2]) 138 check(ga.mv('A', 'odd'), [1, 3]) 139 check(ga.mv('A', 'mv'), [0, 1, 2, 3]) 140 141 # value construction 142 check(ga.mv([1, 2, 3], 'vector'), [1]) 143 144 # illegal arguments 145 with pytest.raises(TypeError): 146 ga.mv('A', 'vector', "too many arguments") 147 with pytest.raises(TypeError): 148 ga.mv('A', 'grade') # too few arguments 149 with pytest.raises(TypeError): 150 ga.mv('A', 'grade', not_an_argument=True) # invalid kwarg 151 with pytest.raises(TypeError): 152 ga.mv([1, 2, 3], 'vector', f=True) # can't pass f with coefficients 153 with pytest.raises(TypeError): 154 ga.mv(e_1, 'even') # Must be a string 155 156 def test_abs(self): 157 ga, e_1, e_2, e_3 = Ga.build('e*1|2|3', g=[1, 1, 1]) 158 B = ga.mv('B', 'bivector') 159 assert abs(B*B) == -(B*B).scalar() 160 161 def test_hashable(self): 162 ga, e_1, e_2, e_3 = Ga.build('e*1|2|3') 163 164 d = {} 165 d[e_1] = 1 166 d[e_2] = 2 167 assert d[e_1 + 0] == 1 168 d[10] = 3 # note: not a multivector key! 169 assert d[e_1 * 0 + 10] == 3 170 171 def test_subs(self): 172 ga, e_1, e_2, e_3 = Ga.build('e*1|2|3', g=[1, 1, 1]) 173 B = ga.mv('B', 'bivector') 174 B_inv = B.inv() 175 176 B_mag = sympy.Symbol('|B|') 177 178 # both of the sympy subs syntaxes work: 179 assert (-B / B_mag**2).subs(B_mag, abs(B)) == B_inv 180 assert (-B / B_mag**2).subs({B_mag: abs(B)}) == B_inv 181 182 def test_sympify(self): 183 ga, e_1, e_2, e_3 = Ga.build('e*1|2|3', g=[1, 1, 1]) 184 185 # Letting this succeed silently and not return an Mv instance would be 186 # dangerous. 187 with pytest.raises(sympy.SympifyError): 188 sympy.sympify(e_1) 189 190 # this is fine 191 sympy.sympify(e_1.obj) 192 193 def test_arithmetic(self): 194 ga, e_1, e_2, e_3 = Ga.build('e*1|2|3', g=[1, 1, 1]) 195 one = ga.mv(sympy.S.One) 196 197 # test that scalars are promoted to Mvs correctly 198 assert e_1 + 1 == e_1 + one 199 assert 1 + e_1 == one + e_1 200 assert e_1 - 1 == e_1 - one 201 assert 1 - e_1 == one - e_1 202 203 @pytest.mark.parametrize('make_one', [ 204 lambda ga: 1, 205 lambda ga: ga.mv(sympy.S.One), 206 pytest.param(lambda ga: sympy.S.One, marks=pytest.mark.skipif( 207 Version(sympy.__version__) < Version("1.6"), 208 # until sympy/sympy@bec42df53cf2486d485065ddad1c31011a48bf3b 209 reason="Cannot override < and > on sympy.Expr" 210 )) 211 ]) 212 def test_contraction(self, make_one): 213 ga, e_1, e_2 = Ga.build('e*1|2', g=[1, 1]) 214 e12 = e_1 ^ e_2 215 one = make_one(ga) 216 217 assert (one < e12) == e12 218 assert (e12 > one) == e12 219 assert (e12 < one) == 0 220 assert (one > e12) == 0 221 222 def test_proj(self): 223 g3coords = symbols('x y z', real=True) 224 V = Ga('e', g=[1, 1, 1], coords=g3coords) 225 226 u = V.mv("u", "vector") 227 v = V.mv("v", "vector") 228 w = V.mv("w", "vector") 229 B = V.mv("B", "mv") 230 231 assert proj(u, v) == v.project_in_blade(u) 232 assert proj(w, v) + proj(w, u) == proj(w, u + v) 233 assert proj(u, v) == (v | u) / u 234 235 Vr = u ^ v 236 assert proj(Vr, B) == ((B < Vr) * Vr.inv()) 237 assert proj(Vr, B) == ((B < Vr) < Vr.inv()) 238 239 def test_norm2(self): 240 g3coords = symbols('x y z', real=True) 241 g3 = Ga('e', g=[1, 1, 1], coords=g3coords) 242 A = g3.mv('A', 'mv') 243 244 assert A.norm2() == A.rev().sp(A) 245 assert A.norm2('+') == A.rev().sp(A) 246 assert A.norm2('-') == A.rev().sp(A) 247 248 assert norm2(A) == norm(A) * norm(A) 249 assert norm2(A, '+') == norm(A) * norm(A) 250 assert norm2(A, '-') == norm(A) * norm(A) 251 252 def test_norm_nonneg(self): 253 """Test that norm always returns a nonneg expression (issue 522).""" 254 from sympy import Abs 255 g3 = Ga('e', g=[1, 1, 1], coords=symbols('x y z', real=True)) 256 257 # scalar norm should include Abs and be nonneg 258 s = g3.mv('s', 'scalar') 259 s_norm = s.norm() 260 assert isinstance(s_norm, Abs) 261 assert s_norm.is_nonnegative 262 263 # vector norm is already nonneg (sum of squares under sqrt) 264 v = g3.mv('v', 'vector') 265 assert v.norm().is_nonnegative 266 267 # even-grade (bivector) norm is nonneg 268 B = g3.mv('B', 'bivector') 269 assert B.norm().is_nonnegative 270 271 # non-Euclidean: symbolic vector where norm may not be provably nonneg 272 g2mn = Ga('e', g=[1, -1], coords=symbols('x t', real=True)) 273 v = g2mn.mv('v', 'vector') 274 assert v.norm(hint='+').is_nonnegative 275 276 def test_mag2(self): 277 g3coords = symbols('x y z', real=True) 278 g3 = Ga('e', g=[1, 1, 1], coords=g3coords) 279 A = g3.mv('A', 'mv') 280 281 assert mag2(A) == mag(A) * mag(A) 282 283 def test_undual(self): 284 # A not a multivector in undual(A). 285 with pytest.raises(ValueError): 286 undual(1) 287 288 def test_g_invol(self): 289 g3coords = symbols('x y z', real=True) 290 g3 = Ga('e', g=[1, 1, 1], coords=g3coords) 291 A = g3.mv('A', 'mv') 292 293 assert g_invol(A.even()) == A.even() 294 assert g_invol(A.odd()) == -A.odd() 295 296 with pytest.raises(ValueError): 297 g_invol(1) 298 299 def test_exp(self): 300 g3coords = (x, y, z) = symbols('x y z', real=True) 301 g3 = Ga('e', g=[1, 1, 1], coords=g3coords) 302 303 u = g3.mv("u", "vector") 304 v = g3.mv("v", "vector") 305 306 assert exp(u ^ v) == exp(-v ^ u) 307 assert exp(x + y + z) == exp(z + y + x) 308 309 def test_ccon(self): 310 g3coords = symbols('x y z', real=True) 311 g3 = Ga('e', g=[1, 1, 1], coords=g3coords) 312 313 A = g3.mv('A', 'mv') 314 315 assert ccon(ccon(A)) == A 316 assert ccon(A) == g_invol(rev(A)) 317 assert ccon(A) == rev(g_invol(A)) 318 319 # not a multivector in ccon(A) 320 with pytest.raises(ValueError): 321 ccon(1) 322 323 def test_scalar(self): 324 g3coords = symbols('x y z', real=True) 325 g3 = Ga('e', g=[1, 1, 1], coords=g3coords) 326 327 A = g3.mv('A', 'mv') 328 u = g3.mv("u", "vector") 329 330 assert scalar(A) == A.scalar() 331 assert scalar(u) == 0 332 assert scalar(u * u) == qform(u) 333 334 # not a multivector in scalar(A) 335 with pytest.raises(ValueError): 336 scalar(1) 337 338 def test_sp(self): 339 g3coords = symbols('x y z', real=True) 340 g3 = Ga('e', g=[1, 1, 1], coords=g3coords) 341 342 A = g3.mv('A', 'mv') 343 B = g3.mv('B', 'mv') 344 345 assert sp(A, B) == A.sp(B) 346 assert sp(A, B) == (A*B).scalar() 347 348 # not a multivector in sp(A, B) 349 with pytest.raises(ValueError): 350 sp(1, 1) 351 352 @pytest.mark.parametrize('inv_func', [inv, shirokov_inverse, hitzer_inverse]) 353 @pytest.mark.parametrize('n', range(2, 5)) # tests are slow for larger n 354 def test_inv(self, inv_func, n): 355 gncoords = symbols(" ".join(f"x{i}" for i in range(n)), real=True) 356 gn = Ga('e', g=[1] * n, coords=gncoords) 357 358 # invert versors 359 A = gn.mv('A', 'vector') 360 Ainv = inv_func(A) 361 assert A * Ainv == 1 + 0 * A 362 363 # invert generic multivectors 364 if inv_func in [shirokov_inverse, hitzer_inverse] and n<4: 365 B = gn.mv('A', 'mv') 366 Binv = inv_func(B) 367 assert B * Binv == 1 + 0 * B 368 369 def test_rtruediv(self): 370 """Test that scalar/Mv works (issue 512).""" 371 ga, e_1, e_2, e_3 = Ga.build('e*1|2|3') 372 I = ga.I() 373 374 # 1/I should be the inverse of I 375 result = 1 / I 376 assert result * I == ga.mv(1) 377 378 # scalar / vector 379 v = e_1 + e_2 380 result = 2 / v 381 assert result * v == ga.mv(2) 382 383 # sympy scalar / Mv 384 from sympy import S, symbols 385 result = S(3) / I 386 assert result * I == ga.mv(3) 387 388 # symbolic numerator 389 a = symbols('a') 390 result = a / I 391 assert result * I == ga.mv(a) 392 393 # non-Euclidean metric (Minkowski-like) 394 ga_m, e_t, e_x = Ga.build('e*t|x', g=[1, -1]) 395 v = e_t + e_x 396 # v*v = 1 - 1 = 0 for null vector — skip (not invertible) 397 v2 = e_t + 2 * e_x # v2*v2 = 1 - 4 = -3 (invertible) 398 result = 1 / v2 399 assert result * v2 == ga_m.mv(1) 400 401 # reproduce #537: is_blade() must handle null vectors and null blades 402 403 def test_is_blade_zero_and_scalar(self): 404 """Zero mv is not a blade; grade-0 scalars are.""" 405 ga, e1, e2 = Ga.build('e*1|2', g=[1, 1]) 406 assert ga.mv(S.Zero).is_blade() is False 407 assert ga.mv(S.One).is_blade() is True 408 409 def test_is_blade_vectors(self): 410 """Grade-1: non-null and null vectors are both blades.""" 411 ga, e0, e1 = Ga.build('e*0|1', g=[1, -1]) 412 # non-null vectors 413 assert e0.is_blade() is True 414 assert e1.is_blade() is True 415 # null vector b = e0+e1 (b*b = 0): is_versor() fails but is_blade() must pass 416 b = e0 + e1 417 assert b.is_zero() is False 418 assert (b * b).is_zero() # confirm null 419 assert b.is_versor() is False # no inverse 420 assert b.is_blade() is True # still a 1-blade 421 422 def test_is_blade_bivectors(self): 423 """Grade-2: non-null and null bivectors.""" 424 # non-null bivector via is_versor() 425 ga, e1, e2, e3 = Ga.build('e*1|2|3', g=[1, 1, 1]) 426 assert (e1 ^ e2).is_blade() is True 427 # null 2-blade: outer-product squaring test must detect it 428 ga3, f0, f1, f2 = Ga.build('f*0|1|2', g=[1, -1, -1]) 429 B = (f0 + f1) ^ f2 430 assert B.is_zero() is False 431 assert B.is_versor() is False # null constituent 432 assert B.is_blade() is True # but it IS a 2-blade 433 434 def test_is_blade_trivector(self): 435 """Grade-3: pseudoscalar of R^3 is a blade.""" 436 ga, e1, e2, e3 = Ga.build('e*1|2|3', g=[1, 1, 1]) 437 assert (e1 ^ e2 ^ e3).is_blade() is True 438 439 @pytest.mark.slow 440 def test_is_blade_grade3_known_limitation(self): 441 """Grade >= 3: B^B=0 is not sufficient; is_blade() may give false positives. 442 443 Dorst/Fontijne/Mann counterexample in R^6: the 3-vector 444 X = e1^e2^e5 + e1^e3^e6 + e2^e4^e6 - e3^e4^e5 445 satisfies X^X = 0 but is NOT a blade. Our algorithm incorrectly 446 returns True; see is_blade() docstring. 447 """ 448 ga, e1, e2, e3, e4, e5, e6 = Ga.build('e*1|2|3|4|5|6', g=[1] * 6) 449 X = (e1 ^ e2 ^ e5) + (e1 ^ e3 ^ e6) + (e2 ^ e4 ^ e6) - (e3 ^ e4 ^ e5) 450 assert (X ^ X).is_zero() # necessary condition satisfied → triggers false positive 451 assert X.is_blade() is True # known limitation: algorithm returns True 452 453 def test_is_blade_non_simple_bivector(self): 454 """Non-simple bivector: B^B != 0, so is_blade() correctly returns False.""" 455 ga, e1, e2, e3, e4 = Ga.build('e*1|2|3|4', g=[1, 1, 1, 1]) 456 B = (e1 ^ e2) + (e3 ^ e4) # non-simple; B^B = 2*e1^e2^e3^e4 != 0 457 assert not (B ^ B).is_zero() 458 assert B.is_blade() is False 459 460 def test_is_blade_non_homogeneous(self): 461 """Non-grade-homogeneous mv is not a blade.""" 462 ga, e0, e1 = Ga.build('e*0|1', g=[1, -1]) 463 assert (e0 + (e0 ^ e1)).is_blade() is False 464 465 def test_is_blade_result_cached(self): 466 """is_blade() caches its result in blade_flg.""" 467 ga, e1, e2 = Ga.build('e*1|2', g=[1, 1]) 468 assert e1.is_blade() is True 469 assert e1.blade_flg is True # cache populated 470 assert e1.is_blade() is True # second call hits cached branch 471 472 def test_reflect_in_null_blade_raises(self): 473 """reflect_in_blade() must raise ValueError for null blades (#537).""" 474 ga, f0, f1, f2 = Ga.build('f*0|1|2', g=[1, -1, -1]) 475 with pytest.raises(ValueError, match='null blade'): 476 f0.reflect_in_blade((f0 + f1) ^ f2) 477 478 def test_project_in_null_blade_raises(self): 479 """project_in_blade() must raise ValueError for null blades (#537).""" 480 ga, f0, f1, f2 = Ga.build('f*0|1|2', g=[1, -1, -1]) 481 with pytest.raises(ValueError, match='null blade'): 482 f0.project_in_blade((f0 + f1) ^ f2) 483