mjo/eja/eja_element.py

   1 from sage.matrix.constructor import matrix
   2 from sage.modules.free_module import VectorSpace
   3 from sage.modules.with_basis.indexed_element import IndexedFreeModuleElement
   4
   5 from mjo.eja.eja_operator import FiniteDimensionalEuclideanJordanAlgebraOperator
   6 from mjo.eja.eja_utils import _mat2vec
   7
   8 class FiniteDimensionalEuclideanJordanAlgebraElement(IndexedFreeModuleElement):
   9     """
  10     An element of a Euclidean Jordan algebra.
  11     """
  12
  13     def __dir__(self):
  14         """
  15         Oh man, I should not be doing this. This hides the "disabled"
  16         methods ``left_matrix`` and ``matrix`` from introspection;
  17         in particular it removes them from tab-completion.
  18         """
  19         return filter(lambda s: s not in ['left_matrix', 'matrix'],
  20                       dir(self.__class__) )
  21
  22
  23
  24
  25     def __pow__(self, n):
  26         """
  27         Return ``self`` raised to the power ``n``.
  28
  29         Jordan algebras are always power-associative; see for
  30         example Faraut and Korányi, Proposition II.1.2 (ii).
  31
  32         We have to override this because our superclass uses row
  33         vectors instead of column vectors! We, on the other hand,
  34         assume column vectors everywhere.
  35
  36         SETUP::
  37
  38             sage: from mjo.eja.eja_algebra import random_eja
  39
  40         TESTS:
  41
  42         The definition of `x^2` is the unambiguous `x*x`::
  43
  44             sage: set_random_seed()
  45             sage: x = random_eja().random_element()
  46             sage: x*x == (x^2)
  47             True
  48
  49         A few examples of power-associativity::
  50
  51             sage: set_random_seed()
  52             sage: x = random_eja().random_element()
  53             sage: x*(x*x)*(x*x) == x^5
  54             True
  55             sage: (x*x)*(x*x*x) == x^5
  56             True
  57
  58         We also know that powers operator-commute (Koecher, Chapter
  59         III, Corollary 1)::
  60
  61             sage: set_random_seed()
  62             sage: x = random_eja().random_element()
  63             sage: m = ZZ.random_element(0,10)
  64             sage: n = ZZ.random_element(0,10)
  65             sage: Lxm = (x^m).operator()
  66             sage: Lxn = (x^n).operator()
  67             sage: Lxm*Lxn == Lxn*Lxm
  68             True
  69
  70         """
  71         if n == 0:
  72             return self.parent().one()
  73         elif n == 1:
  74             return self
  75         else:
  76             return (self**(n-1))*self
  77
  78
  79     def apply_univariate_polynomial(self, p):
  80         """
  81         Apply the univariate polynomial ``p`` to this element.
  82
  83         A priori, SageMath won't allow us to apply a univariate
  84         polynomial to an element of an EJA, because we don't know
  85         that EJAs are rings (they are usually not associative). Of
  86         course, we know that EJAs are power-associative, so the
  87         operation is ultimately kosher. This function sidesteps
  88         the CAS to get the answer we want and expect.
  89
  90         SETUP::
  91
  92             sage: from mjo.eja.eja_algebra import (HadamardEJA,
  93             ....:                                  random_eja)
  94
  95         EXAMPLES::
  96
  97             sage: R = PolynomialRing(QQ, 't')
  98             sage: t = R.gen(0)
  99             sage: p = t^4 - t^3 + 5*t - 2
 100             sage: J = HadamardEJA(5)
 101             sage: J.one().apply_univariate_polynomial(p) == 3*J.one()
 102             True
 103
 104         TESTS:
 105
 106         We should always get back an element of the algebra::
 107
 108             sage: set_random_seed()
 109             sage: p = PolynomialRing(AA, 't').random_element()
 110             sage: J = random_eja()
 111             sage: x = J.random_element()
 112             sage: x.apply_univariate_polynomial(p) in J
 113             True
 114
 115         """
 116         if len(p.variables()) > 1:
 117             raise ValueError("not a univariate polynomial")
 118         P = self.parent()
 119         R = P.base_ring()
 120         # Convert the coeficcients to the parent's base ring,
 121         # because a priori they might live in an (unnecessarily)
 122         # larger ring for which P.sum() would fail below.
 123         cs = [ R(c) for c in p.coefficients(sparse=False) ]
 124         return P.sum( cs[k]*(self**k) for k in range(len(cs)) )
 125
 126
 127     def characteristic_polynomial(self):
 128         """
 129         Return the characteristic polynomial of this element.
 130
 131         SETUP::
 132
 133             sage: from mjo.eja.eja_algebra import HadamardEJA
 134
 135         EXAMPLES:
 136
 137         The rank of `R^3` is three, and the minimal polynomial of
 138         the identity element is `(t-1)` from which it follows that
 139         the characteristic polynomial should be `(t-1)^3`::
 140
 141             sage: J = HadamardEJA(3)
 142             sage: J.one().characteristic_polynomial()
 143             t^3 - 3*t^2 + 3*t - 1
 144
 145         Likewise, the characteristic of the zero element in the
 146         rank-three algebra `R^{n}` should be `t^{3}`::
 147
 148             sage: J = HadamardEJA(3)
 149             sage: J.zero().characteristic_polynomial()
 150             t^3
 151
 152         TESTS:
 153
 154         The characteristic polynomial of an element should evaluate
 155         to zero on that element::
 156
 157             sage: set_random_seed()
 158             sage: x = HadamardEJA(3).random_element()
 159             sage: p = x.characteristic_polynomial()
 160             sage: x.apply_univariate_polynomial(p)
 161             0
 162
 163         The characteristic polynomials of the zero and unit elements
 164         should be what we think they are in a subalgebra, too::
 165
 166             sage: J = HadamardEJA(3)
 167             sage: p1 = J.one().characteristic_polynomial()
 168             sage: q1 = J.zero().characteristic_polynomial()
 169             sage: e0,e1,e2 = J.gens()
 170             sage: A = (e0 + 2*e1 + 3*e2).subalgebra_generated_by() # dim 3
 171             sage: p2 = A.one().characteristic_polynomial()
 172             sage: q2 = A.zero().characteristic_polynomial()
 173             sage: p1 == p2
 174             True
 175             sage: q1 == q2
 176             True
 177
 178         """
 179         p = self.parent().characteristic_polynomial_of()
 180         return p(*self.to_vector())
 181
 182
 183     def inner_product(self, other):
 184         """
 185         Return the parent algebra's inner product of myself and ``other``.
 186
 187         SETUP::
 188
 189             sage: from mjo.eja.eja_algebra import (
 190             ....:   ComplexHermitianEJA,
 191             ....:   JordanSpinEJA,
 192             ....:   QuaternionHermitianEJA,
 193             ....:   RealSymmetricEJA,
 194             ....:   random_eja)
 195
 196         EXAMPLES:
 197
 198         The inner product in the Jordan spin algebra is the usual
 199         inner product on `R^n` (this example only works because the
 200         basis for the Jordan algebra is the standard basis in `R^n`)::
 201
 202             sage: J = JordanSpinEJA(3)
 203             sage: x = vector(QQ,[1,2,3])
 204             sage: y = vector(QQ,[4,5,6])
 205             sage: x.inner_product(y)
 206             32
 207             sage: J.from_vector(x).inner_product(J.from_vector(y))
 208             32
 209
 210         The inner product on `S^n` is `<X,Y> = trace(X*Y)`, where
 211         multiplication is the usual matrix multiplication in `S^n`,
 212         so the inner product of the identity matrix with itself
 213         should be the `n`::
 214
 215             sage: J = RealSymmetricEJA(3)
 216             sage: J.one().inner_product(J.one())
 217             3
 218
 219         Likewise, the inner product on `C^n` is `<X,Y> =
 220         Re(trace(X*Y))`, where we must necessarily take the real
 221         part because the product of Hermitian matrices may not be
 222         Hermitian::
 223
 224             sage: J = ComplexHermitianEJA(3)
 225             sage: J.one().inner_product(J.one())
 226             3
 227
 228         Ditto for the quaternions::
 229
 230             sage: J = QuaternionHermitianEJA(3)
 231             sage: J.one().inner_product(J.one())
 232             3
 233
 234         TESTS:
 235
 236         Ensure that we can always compute an inner product, and that
 237         it gives us back a real number::
 238
 239             sage: set_random_seed()
 240             sage: J = random_eja()
 241             sage: x,y = J.random_elements(2)
 242             sage: x.inner_product(y) in RLF
 243             True
 244
 245         """
 246         P = self.parent()
 247         if not other in P:
 248             raise TypeError("'other' must live in the same algebra")
 249
 250         return P.inner_product(self, other)
 251
 252
 253     def operator_commutes_with(self, other):
 254         """
 255         Return whether or not this element operator-commutes
 256         with ``other``.
 257
 258         SETUP::
 259
 260             sage: from mjo.eja.eja_algebra import random_eja
 261
 262         EXAMPLES:
 263
 264         The definition of a Jordan algebra says that any element
 265         operator-commutes with its square::
 266
 267             sage: set_random_seed()
 268             sage: x = random_eja().random_element()
 269             sage: x.operator_commutes_with(x^2)
 270             True
 271
 272         TESTS:
 273
 274         Test Lemma 1 from Chapter III of Koecher::
 275
 276             sage: set_random_seed()
 277             sage: u,v = random_eja().random_elements(2)
 278             sage: lhs = u.operator_commutes_with(u*v)
 279             sage: rhs = v.operator_commutes_with(u^2)
 280             sage: lhs == rhs
 281             True
 282
 283         Test the first polarization identity from my notes, Koecher
 284         Chapter III, or from Baes (2.3)::
 285
 286             sage: set_random_seed()
 287             sage: x,y = random_eja().random_elements(2)
 288             sage: Lx = x.operator()
 289             sage: Ly = y.operator()
 290             sage: Lxx = (x*x).operator()
 291             sage: Lxy = (x*y).operator()
 292             sage: bool(2*Lx*Lxy + Ly*Lxx == 2*Lxy*Lx + Lxx*Ly)
 293             True
 294
 295         Test the second polarization identity from my notes or from
 296         Baes (2.4)::
 297
 298             sage: set_random_seed()
 299             sage: x,y,z = random_eja().random_elements(3)
 300             sage: Lx = x.operator()
 301             sage: Ly = y.operator()
 302             sage: Lz = z.operator()
 303             sage: Lzy = (z*y).operator()
 304             sage: Lxy = (x*y).operator()
 305             sage: Lxz = (x*z).operator()
 306             sage: bool(Lx*Lzy + Lz*Lxy + Ly*Lxz == Lzy*Lx + Lxy*Lz + Lxz*Ly)
 307             True
 308
 309         Test the third polarization identity from my notes or from
 310         Baes (2.5)::
 311
 312             sage: set_random_seed()
 313             sage: u,y,z = random_eja().random_elements(3)
 314             sage: Lu = u.operator()
 315             sage: Ly = y.operator()
 316             sage: Lz = z.operator()
 317             sage: Lzy = (z*y).operator()
 318             sage: Luy = (u*y).operator()
 319             sage: Luz = (u*z).operator()
 320             sage: Luyz = (u*(y*z)).operator()
 321             sage: lhs = Lu*Lzy + Lz*Luy + Ly*Luz
 322             sage: rhs = Luyz + Ly*Lu*Lz + Lz*Lu*Ly
 323             sage: bool(lhs == rhs)
 324             True
 325
 326         """
 327         if not other in self.parent():
 328             raise TypeError("'other' must live in the same algebra")
 329
 330         A = self.operator()
 331         B = other.operator()
 332         return (A*B == B*A)
 333
 334
 335     def det(self):
 336         """
 337         Return my determinant, the product of my eigenvalues.
 338
 339         SETUP::
 340
 341             sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
 342             ....:                                  TrivialEJA,
 343             ....:                                  random_eja)
 344
 345         EXAMPLES::
 346
 347             sage: J = JordanSpinEJA(2)
 348             sage: e0,e1 = J.gens()
 349             sage: x = sum( J.gens() )
 350             sage: x.det()
 351             0
 352
 353         ::
 354
 355             sage: J = JordanSpinEJA(3)
 356             sage: e0,e1,e2 = J.gens()
 357             sage: x = sum( J.gens() )
 358             sage: x.det()
 359             -1
 360
 361         The determinant of the sole element in the rank-zero trivial
 362         algebra is ``1``, by three paths of reasoning. First, its
 363         characteristic polynomial is a constant ``1``, so the constant
 364         term in that polynomial is ``1``. Second, the characteristic
 365         polynomial evaluated at zero is again ``1``. And finally, the
 366         (empty) product of its eigenvalues is likewise just unity::
 367
 368             sage: J = TrivialEJA()
 369             sage: J.zero().det()
 370             1
 371
 372         TESTS:
 373
 374         An element is invertible if and only if its determinant is
 375         non-zero::
 376
 377             sage: set_random_seed()
 378             sage: x = random_eja().random_element()
 379             sage: x.is_invertible() == (x.det() != 0)
 380             True
 381
 382         Ensure that the determinant is multiplicative on an associative
 383         subalgebra as in Faraut and Korányi's Proposition II.2.2::
 384
 385             sage: set_random_seed()
 386             sage: J = random_eja().random_element().subalgebra_generated_by()
 387             sage: x,y = J.random_elements(2)
 388             sage: (x*y).det() == x.det()*y.det()
 389             True
 390         """
 391         P = self.parent()
 392         r = P.rank()
 393
 394         if r == 0:
 395             # Special case, since we don't get the a0=1
 396             # coefficient when the rank of the algebra
 397             # is zero.
 398             return P.base_ring().one()
 399
 400         p = P._charpoly_coefficients()[0]
 401         # The _charpoly_coeff function already adds the factor of -1
 402         # to ensure that _charpoly_coefficients()[0] is really what
 403         # appears in front of t^{0} in the charpoly. However, we want
 404         # (-1)^r times THAT for the determinant.
 405         return ((-1)**r)*p(*self.to_vector())
 406
 407
 408     def inverse(self):
 409         """
 410         Return the Jordan-multiplicative inverse of this element.
 411
 412         ALGORITHM:
 413
 414         We appeal to the quadratic representation as in Koecher's
 415         Theorem 12 in Chapter III, Section 5.
 416
 417         SETUP::
 418
 419             sage: from mjo.eja.eja_algebra import (ComplexHermitianEJA,
 420             ....:                                  JordanSpinEJA,
 421             ....:                                  random_eja)
 422
 423         EXAMPLES:
 424
 425         The inverse in the spin factor algebra is given in Alizadeh's
 426         Example 11.11::
 427
 428             sage: set_random_seed()
 429             sage: J = JordanSpinEJA.random_instance()
 430             sage: x = J.random_element()
 431             sage: while not x.is_invertible():
 432             ....:     x = J.random_element()
 433             sage: x_vec = x.to_vector()
 434             sage: x0 = x_vec[:1]
 435             sage: x_bar = x_vec[1:]
 436             sage: coeff = x0.inner_product(x0) - x_bar.inner_product(x_bar)
 437             sage: x_inverse = x_vec.parent()(x0.list() + (-x_bar).list())
 438             sage: if not coeff.is_zero(): x_inverse = x_inverse/coeff
 439             sage: x.inverse() == J.from_vector(x_inverse)
 440             True
 441
 442         Trying to invert a non-invertible element throws an error:
 443
 444             sage: JordanSpinEJA(3).zero().inverse()
 445             Traceback (most recent call last):
 446             ...
 447             ValueError: element is not invertible
 448
 449         TESTS:
 450
 451         The identity element is its own inverse::
 452
 453             sage: set_random_seed()
 454             sage: J = random_eja()
 455             sage: J.one().inverse() == J.one()
 456             True
 457
 458         If an element has an inverse, it acts like one::
 459
 460             sage: set_random_seed()
 461             sage: J = random_eja()
 462             sage: x = J.random_element()
 463             sage: (not x.is_invertible()) or (x.inverse()*x == J.one())
 464             True
 465
 466         The inverse of the inverse is what we started with::
 467
 468             sage: set_random_seed()
 469             sage: J = random_eja()
 470             sage: x = J.random_element()
 471             sage: (not x.is_invertible()) or (x.inverse().inverse() == x)
 472             True
 473
 474         Proposition II.2.3 in Faraut and Korányi says that the inverse
 475         of an element is the inverse of its left-multiplication operator
 476         applied to the algebra's identity, when that inverse exists::
 477
 478             sage: set_random_seed()
 479             sage: J = random_eja()
 480             sage: x = J.random_element()
 481             sage: (not x.operator().is_invertible()) or (
 482             ....:    x.operator().inverse()(J.one()) == x.inverse() )
 483             True
 484
 485         Proposition II.2.4 in Faraut and Korányi gives a formula for
 486         the inverse based on the characteristic polynomial and the
 487         Cayley-Hamilton theorem for Euclidean Jordan algebras::
 488
 489             sage: set_random_seed()
 490             sage: J = ComplexHermitianEJA(3)
 491             sage: x = J.random_element()
 492             sage: while not x.is_invertible():
 493             ....:     x = J.random_element()
 494             sage: r = J.rank()
 495             sage: a = x.characteristic_polynomial().coefficients(sparse=False)
 496             sage: expected  = (-1)^(r+1)/x.det()
 497             sage: expected *= sum( a[i+1]*x^i for i in range(r) )
 498             sage: x.inverse() == expected
 499             True
 500
 501         """
 502         if not self.is_invertible():
 503             raise ValueError("element is not invertible")
 504
 505         if self.parent()._charpoly_coefficients.is_in_cache():
 506             # We can invert using our charpoly if it will be fast to
 507             # compute. If the coefficients are cached, our rank had
 508             # better be too!
 509             r = self.parent().rank()
 510             a = self.characteristic_polynomial().coefficients(sparse=False)
 511             return (-1)**(r+1)*sum(a[i+1]*self**i for i in range(r))/self.det()
 512
 513         return (~self.quadratic_representation())(self)
 514
 515
 516     def is_invertible(self):
 517         """
 518         Return whether or not this element is invertible.
 519
 520         ALGORITHM:
 521
 522         The usual way to do this is to check if the determinant is
 523         zero, but we need the characteristic polynomial for the
 524         determinant. The minimal polynomial is a lot easier to get,
 525         so we use Corollary 2 in Chapter V of Koecher to check
 526         whether or not the paren't algebra's zero element is a root
 527         of this element's minimal polynomial.
 528
 529         That is... unless the coefficients of our algebra's
 530         "characteristic polynomial of" function are already cached!
 531         In that case, we just use the determinant (which will be fast
 532         as a result).
 533
 534         Beware that we can't use the superclass method, because it
 535         relies on the algebra being associative.
 536
 537         SETUP::
 538
 539             sage: from mjo.eja.eja_algebra import random_eja
 540
 541         TESTS:
 542
 543         The identity element is always invertible::
 544
 545             sage: set_random_seed()
 546             sage: J = random_eja()
 547             sage: J.one().is_invertible()
 548             True
 549
 550         The zero element is never invertible in a non-trivial algebra::
 551
 552             sage: set_random_seed()
 553             sage: J = random_eja()
 554             sage: (not J.is_trivial()) and J.zero().is_invertible()
 555             False
 556
 557         """
 558         if self.is_zero():
 559             if self.parent().is_trivial():
 560                 return True
 561             else:
 562                 return False
 563
 564         if self.parent()._charpoly_coefficients.is_in_cache():
 565             # The determinant will be quicker than computing the minimal
 566             # polynomial from scratch, most likely.
 567             return (not self.det().is_zero())
 568
 569         # In fact, we only need to know if the constant term is non-zero,
 570         # so we can pass in the field's zero element instead.
 571         zero = self.base_ring().zero()
 572         p = self.minimal_polynomial()
 573         return not (p(zero) == zero)
 574
 575
 576     def is_primitive_idempotent(self):
 577         """
 578         Return whether or not this element is a primitive (or minimal)
 579         idempotent.
 580
 581         A primitive idempotent is a non-zero idempotent that is not
 582         the sum of two other non-zero idempotents. Remark 2.7.15 in
 583         Baes shows that this is what he refers to as a "minimal
 584         idempotent."
 585
 586         An element of a Euclidean Jordan algebra is a minimal idempotent
 587         if it :meth:`is_idempotent` and if its Peirce subalgebra
 588         corresponding to the eigenvalue ``1`` has dimension ``1`` (Baes,
 589         Proposition 2.7.17).
 590
 591         SETUP::
 592
 593             sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
 594             ....:                                  RealSymmetricEJA,
 595             ....:                                  TrivialEJA,
 596             ....:                                  random_eja)
 597
 598         WARNING::
 599
 600         This method is sloooooow.
 601
 602         EXAMPLES:
 603
 604         The spectral decomposition of a non-regular element should always
 605         contain at least one non-minimal idempotent::
 606
 607             sage: J = RealSymmetricEJA(3)
 608             sage: x = sum(J.gens())
 609             sage: x.is_regular()
 610             False
 611             sage: [ c.is_primitive_idempotent()
 612             ....:   for (l,c) in x.spectral_decomposition() ]
 613             [False, True]
 614
 615         On the other hand, the spectral decomposition of a regular
 616         element should always be in terms of minimal idempotents::
 617
 618             sage: J = JordanSpinEJA(4)
 619             sage: x = sum( i*J.gens()[i] for i in range(len(J.gens())) )
 620             sage: x.is_regular()
 621             True
 622             sage: [ c.is_primitive_idempotent()
 623             ....:   for (l,c) in x.spectral_decomposition() ]
 624             [True, True]
 625
 626         TESTS:
 627
 628         The identity element is minimal only in an EJA of rank one::
 629
 630             sage: set_random_seed()
 631             sage: J = random_eja()
 632             sage: J.rank() == 1 or not J.one().is_primitive_idempotent()
 633             True
 634
 635         A non-idempotent cannot be a minimal idempotent::
 636
 637             sage: set_random_seed()
 638             sage: J = JordanSpinEJA(4)
 639             sage: x = J.random_element()
 640             sage: (not x.is_idempotent()) and x.is_primitive_idempotent()
 641             False
 642
 643         Proposition 2.7.19 in Baes says that an element is a minimal
 644         idempotent if and only if it's idempotent with trace equal to
 645         unity::
 646
 647             sage: set_random_seed()
 648             sage: J = JordanSpinEJA(4)
 649             sage: x = J.random_element()
 650             sage: expected = (x.is_idempotent() and x.trace() == 1)
 651             sage: actual = x.is_primitive_idempotent()
 652             sage: actual == expected
 653             True
 654
 655         Primitive idempotents must be non-zero::
 656
 657             sage: set_random_seed()
 658             sage: J = random_eja()
 659             sage: J.zero().is_idempotent()
 660             True
 661             sage: J.zero().is_primitive_idempotent()
 662             False
 663
 664         As a consequence of the fact that primitive idempotents must
 665         be non-zero, there are no primitive idempotents in a trivial
 666         Euclidean Jordan algebra::
 667
 668             sage: J = TrivialEJA()
 669             sage: J.one().is_idempotent()
 670             True
 671             sage: J.one().is_primitive_idempotent()
 672             False
 673
 674         """
 675         if not self.is_idempotent():
 676             return False
 677
 678         if self.is_zero():
 679             return False
 680
 681         (_,_,J1) = self.parent().peirce_decomposition(self)
 682         return (J1.dimension() == 1)
 683
 684
 685     def is_nilpotent(self):
 686         """
 687         Return whether or not some power of this element is zero.
 688
 689         ALGORITHM:
 690
 691         We use Theorem 5 in Chapter III of Koecher, which says that
 692         an element ``x`` is nilpotent if and only if ``x.operator()``
 693         is nilpotent. And it is a basic fact of linear algebra that
 694         an operator on an `n`-dimensional space is nilpotent if and
 695         only if, when raised to the `n`th power, it equals the zero
 696         operator (for example, see Axler Corollary 8.8).
 697
 698         SETUP::
 699
 700             sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
 701             ....:                                  random_eja)
 702
 703         EXAMPLES::
 704
 705             sage: J = JordanSpinEJA(3)
 706             sage: x = sum(J.gens())
 707             sage: x.is_nilpotent()
 708             False
 709
 710         TESTS:
 711
 712         The identity element is never nilpotent, except in a trivial EJA::
 713
 714             sage: set_random_seed()
 715             sage: J = random_eja()
 716             sage: J.one().is_nilpotent() and not J.is_trivial()
 717             False
 718
 719         The additive identity is always nilpotent::
 720
 721             sage: set_random_seed()
 722             sage: random_eja().zero().is_nilpotent()
 723             True
 724
 725         """
 726         P = self.parent()
 727         zero_operator = P.zero().operator()
 728         return self.operator()**P.dimension() == zero_operator
 729
 730
 731     def is_regular(self):
 732         """
 733         Return whether or not this is a regular element.
 734
 735         SETUP::
 736
 737             sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
 738             ....:                                  random_eja)
 739
 740         EXAMPLES:
 741
 742         The identity element always has degree one, but any element
 743         linearly-independent from it is regular::
 744
 745             sage: J = JordanSpinEJA(5)
 746             sage: J.one().is_regular()
 747             False
 748             sage: e0, e1, e2, e3, e4 = J.gens() # e0 is the identity
 749             sage: for x in J.gens():
 750             ....:     (J.one() + x).is_regular()
 751             False
 752             True
 753             True
 754             True
 755             True
 756
 757         TESTS:
 758
 759         The zero element should never be regular, unless the parent
 760         algebra has dimension less than or equal to one::
 761
 762             sage: set_random_seed()
 763             sage: J = random_eja()
 764             sage: J.dimension() <= 1 or not J.zero().is_regular()
 765             True
 766
 767         The unit element isn't regular unless the algebra happens to
 768         consist of only its scalar multiples::
 769
 770             sage: set_random_seed()
 771             sage: J = random_eja()
 772             sage: J.dimension() <= 1 or not J.one().is_regular()
 773             True
 774
 775         """
 776         return self.degree() == self.parent().rank()
 777
 778
 779     def degree(self):
 780         """
 781         Return the degree of this element, which is defined to be
 782         the degree of its minimal polynomial.
 783
 784         ALGORITHM:
 785
 786         For now, we skip the messy minimal polynomial computation
 787         and instead return the dimension of the vector space spanned
 788         by the powers of this element. The latter is a bit more
 789         straightforward to compute.
 790
 791         SETUP::
 792
 793             sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
 794             ....:                                  random_eja)
 795
 796         EXAMPLES::
 797
 798             sage: J = JordanSpinEJA(4)
 799             sage: J.one().degree()
 800             1
 801             sage: e0,e1,e2,e3 = J.gens()
 802             sage: (e0 - e1).degree()
 803             2
 804
 805         In the spin factor algebra (of rank two), all elements that
 806         aren't multiples of the identity are regular::
 807
 808             sage: set_random_seed()
 809             sage: J = JordanSpinEJA.random_instance()
 810             sage: n = J.dimension()
 811             sage: x = J.random_element()
 812             sage: x.degree() == min(n,2) or (x == x.coefficient(0)*J.one())
 813             True
 814
 815         TESTS:
 816
 817         The zero and unit elements are both of degree one in nontrivial
 818         algebras::
 819
 820             sage: set_random_seed()
 821             sage: J = random_eja()
 822             sage: d = J.zero().degree()
 823             sage: (J.is_trivial() and d == 0) or d == 1
 824             True
 825             sage: d = J.one().degree()
 826             sage: (J.is_trivial() and d == 0) or d == 1
 827             True
 828
 829         Our implementation agrees with the definition::
 830
 831             sage: set_random_seed()
 832             sage: x = random_eja().random_element()
 833             sage: x.degree() == x.minimal_polynomial().degree()
 834             True
 835
 836         """
 837         if self.is_zero() and not self.parent().is_trivial():
 838             # The minimal polynomial of zero in a nontrivial algebra
 839             # is "t"; in a trivial algebra it's "1" by convention
 840             # (it's an empty product).
 841             return 1
 842         return self.subalgebra_generated_by().dimension()
 843
 844
 845     def left_matrix(self):
 846         """
 847         Our parent class defines ``left_matrix`` and ``matrix``
 848         methods whose names are misleading. We don't want them.
 849         """
 850         raise NotImplementedError("use operator().matrix() instead")
 851
 852     matrix = left_matrix
 853
 854
 855     def minimal_polynomial(self):
 856         """
 857         Return the minimal polynomial of this element,
 858         as a function of the variable `t`.
 859
 860         ALGORITHM:
 861
 862         We restrict ourselves to the associative subalgebra
 863         generated by this element, and then return the minimal
 864         polynomial of this element's operator matrix (in that
 865         subalgebra). This works by Baes Proposition 2.3.16.
 866
 867         SETUP::
 868
 869             sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
 870             ....:                                  RealSymmetricEJA,
 871             ....:                                  TrivialEJA,
 872             ....:                                  random_eja)
 873
 874         EXAMPLES:
 875
 876         Keeping in mind that the polynomial ``1`` evaluates the identity
 877         element (also the zero element) of the trivial algebra, it is clear
 878         that the polynomial ``1`` is the minimal polynomial of the only
 879         element in a trivial algebra::
 880
 881             sage: J = TrivialEJA()
 882             sage: J.one().minimal_polynomial()
 883             1
 884             sage: J.zero().minimal_polynomial()
 885             1
 886
 887         TESTS:
 888
 889         The minimal polynomial of the identity and zero elements are
 890         always the same, except in trivial algebras where the minimal
 891         polynomial of the unit/zero element is ``1``::
 892
 893             sage: set_random_seed()
 894             sage: J = random_eja()
 895             sage: mu = J.one().minimal_polynomial()
 896             sage: t = mu.parent().gen()
 897             sage: mu + int(J.is_trivial())*(t-2)
 898             t - 1
 899             sage: mu = J.zero().minimal_polynomial()
 900             sage: t = mu.parent().gen()
 901             sage: mu + int(J.is_trivial())*(t-1)
 902             t
 903
 904         The degree of an element is (by one definition) the degree
 905         of its minimal polynomial::
 906
 907             sage: set_random_seed()
 908             sage: x = random_eja().random_element()
 909             sage: x.degree() == x.minimal_polynomial().degree()
 910             True
 911
 912         The minimal polynomial and the characteristic polynomial coincide
 913         and are known (see Alizadeh, Example 11.11) for all elements of
 914         the spin factor algebra that aren't scalar multiples of the
 915         identity. We require the dimension of the algebra to be at least
 916         two here so that said elements actually exist::
 917
 918             sage: set_random_seed()
 919             sage: n_max = max(2, JordanSpinEJA._max_random_instance_size())
 920             sage: n = ZZ.random_element(2, n_max)
 921             sage: J = JordanSpinEJA(n)
 922             sage: y = J.random_element()
 923             sage: while y == y.coefficient(0)*J.one():
 924             ....:     y = J.random_element()
 925             sage: y0 = y.to_vector()[0]
 926             sage: y_bar = y.to_vector()[1:]
 927             sage: actual = y.minimal_polynomial()
 928             sage: t = PolynomialRing(J.base_ring(),'t').gen(0)
 929             sage: expected = t^2 - 2*y0*t + (y0^2 - norm(y_bar)^2)
 930             sage: bool(actual == expected)
 931             True
 932
 933         The minimal polynomial should always kill its element::
 934
 935             sage: set_random_seed()
 936             sage: x = random_eja().random_element()
 937             sage: p = x.minimal_polynomial()
 938             sage: x.apply_univariate_polynomial(p)
 939             0
 940
 941         The minimal polynomial is invariant under a change of basis,
 942         and in particular, a re-scaling of the basis::
 943
 944             sage: set_random_seed()
 945             sage: n_max = RealSymmetricEJA._max_random_instance_size()
 946             sage: n = ZZ.random_element(1, n_max)
 947             sage: J1 = RealSymmetricEJA(n)
 948             sage: J2 = RealSymmetricEJA(n,normalize_basis=False)
 949             sage: X = random_matrix(AA,n)
 950             sage: X = X*X.transpose()
 951             sage: x1 = J1(X)
 952             sage: x2 = J2(X)
 953             sage: x1.minimal_polynomial() == x2.minimal_polynomial()
 954             True
 955
 956         """
 957         if self.is_zero():
 958             # We would generate a zero-dimensional subalgebra
 959             # where the minimal polynomial would be constant.
 960             # That might be correct, but only if *this* algebra
 961             # is trivial too.
 962             if not self.parent().is_trivial():
 963                 # Pretty sure we know what the minimal polynomial of
 964                 # the zero operator is going to be. This ensures
 965                 # consistency of e.g. the polynomial variable returned
 966                 # in the "normal" case without us having to think about it.
 967                 return self.operator().minimal_polynomial()
 968
 969         A = self.subalgebra_generated_by(orthonormalize_basis=False)
 970         return A(self).operator().minimal_polynomial()
 971
 972
 973
 974     def natural_representation(self):
 975         """
 976         Return a more-natural representation of this element.
 977
 978         Every finite-dimensional Euclidean Jordan Algebra is a
 979         direct sum of five simple algebras, four of which comprise
 980         Hermitian matrices. This method returns the original
 981         "natural" representation of this element as a Hermitian
 982         matrix, if it has one. If not, you get the usual representation.
 983
 984         SETUP::
 985
 986             sage: from mjo.eja.eja_algebra import (ComplexHermitianEJA,
 987             ....:                                  QuaternionHermitianEJA)
 988
 989         EXAMPLES::
 990
 991             sage: J = ComplexHermitianEJA(3)
 992             sage: J.one()
 993             e0 + e3 + e8
 994             sage: J.one().natural_representation()
 995             [1 0 0 0 0 0]
 996             [0 1 0 0 0 0]
 997             [0 0 1 0 0 0]
 998             [0 0 0 1 0 0]
 999             [0 0 0 0 1 0]
1000             [0 0 0 0 0 1]
1001
1002         ::
1003
1004             sage: J = QuaternionHermitianEJA(3)
1005             sage: J.one()
1006             e0 + e5 + e14
1007             sage: J.one().natural_representation()
1008             [1 0 0 0 0 0 0 0 0 0 0 0]
1009             [0 1 0 0 0 0 0 0 0 0 0 0]
1010             [0 0 1 0 0 0 0 0 0 0 0 0]
1011             [0 0 0 1 0 0 0 0 0 0 0 0]
1012             [0 0 0 0 1 0 0 0 0 0 0 0]
1013             [0 0 0 0 0 1 0 0 0 0 0 0]
1014             [0 0 0 0 0 0 1 0 0 0 0 0]
1015             [0 0 0 0 0 0 0 1 0 0 0 0]
1016             [0 0 0 0 0 0 0 0 1 0 0 0]
1017             [0 0 0 0 0 0 0 0 0 1 0 0]
1018             [0 0 0 0 0 0 0 0 0 0 1 0]
1019             [0 0 0 0 0 0 0 0 0 0 0 1]
1020
1021         """
1022         B = self.parent().natural_basis()
1023         W = self.parent().natural_basis_space()
1024
1025         # This is just a manual "from_vector()", but of course
1026         # matrix spaces aren't vector spaces in sage, so they
1027         # don't have a from_vector() method.
1028         return W.linear_combination(zip(B,self.to_vector()))
1029
1030
1031     def norm(self):
1032         """
1033         The norm of this element with respect to :meth:`inner_product`.
1034
1035         SETUP::
1036
1037             sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
1038             ....:                                  HadamardEJA)
1039
1040         EXAMPLES::
1041
1042             sage: J = HadamardEJA(2)
1043             sage: x = sum(J.gens())
1044             sage: x.norm()
1045             1.414213562373095?
1046
1047         ::
1048
1049             sage: J = JordanSpinEJA(4)
1050             sage: x = sum(J.gens())
1051             sage: x.norm()
1052             2
1053
1054         """
1055         return self.inner_product(self).sqrt()
1056
1057
1058     def operator(self):
1059         """
1060         Return the left-multiplication-by-this-element
1061         operator on the ambient algebra.
1062
1063         SETUP::
1064
1065             sage: from mjo.eja.eja_algebra import random_eja
1066
1067         TESTS::
1068
1069             sage: set_random_seed()
1070             sage: J = random_eja()
1071             sage: x,y = J.random_elements(2)
1072             sage: x.operator()(y) == x*y
1073             True
1074             sage: y.operator()(x) == x*y
1075             True
1076
1077         """
1078         P = self.parent()
1079         left_mult_by_self = lambda y: self*y
1080         L = P.module_morphism(function=left_mult_by_self, codomain=P)
1081         return FiniteDimensionalEuclideanJordanAlgebraOperator(
1082                  P,
1083                  P,
1084                  L.matrix() )
1085
1086
1087     def quadratic_representation(self, other=None):
1088         """
1089         Return the quadratic representation of this element.
1090
1091         SETUP::
1092
1093             sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
1094             ....:                                  random_eja)
1095
1096         EXAMPLES:
1097
1098         The explicit form in the spin factor algebra is given by
1099         Alizadeh's Example 11.12::
1100
1101             sage: set_random_seed()
1102             sage: x = JordanSpinEJA.random_instance().random_element()
1103             sage: x_vec = x.to_vector()
1104             sage: Q = matrix.identity(x.base_ring(), 0)
1105             sage: n = x_vec.degree()
1106             sage: if n > 0:
1107             ....:     x0 = x_vec[0]
1108             ....:     x_bar = x_vec[1:]
1109             ....:     A = matrix(x.base_ring(), 1, [x_vec.inner_product(x_vec)])
1110             ....:     B = 2*x0*x_bar.row()
1111             ....:     C = 2*x0*x_bar.column()
1112             ....:     D = matrix.identity(x.base_ring(), n-1)
1113             ....:     D = (x0^2 - x_bar.inner_product(x_bar))*D
1114             ....:     D = D + 2*x_bar.tensor_product(x_bar)
1115             ....:     Q = matrix.block(2,2,[A,B,C,D])
1116             sage: Q == x.quadratic_representation().matrix()
1117             True
1118
1119         Test all of the properties from Theorem 11.2 in Alizadeh::
1120
1121             sage: set_random_seed()
1122             sage: J = random_eja()
1123             sage: x,y = J.random_elements(2)
1124             sage: Lx = x.operator()
1125             sage: Lxx = (x*x).operator()
1126             sage: Qx = x.quadratic_representation()
1127             sage: Qy = y.quadratic_representation()
1128             sage: Qxy = x.quadratic_representation(y)
1129             sage: Qex = J.one().quadratic_representation(x)
1130             sage: n = ZZ.random_element(10)
1131             sage: Qxn = (x^n).quadratic_representation()
1132
1133         Property 1:
1134
1135             sage: 2*Qxy == (x+y).quadratic_representation() - Qx - Qy
1136             True
1137
1138         Property 2 (multiply on the right for :trac:`28272`):
1139
1140             sage: alpha = J.base_ring().random_element()
1141             sage: (alpha*x).quadratic_representation() == Qx*(alpha^2)
1142             True
1143
1144         Property 3:
1145
1146             sage: not x.is_invertible() or ( Qx(x.inverse()) == x )
1147             True
1148
1149             sage: not x.is_invertible() or (
1150             ....:   ~Qx
1151             ....:   ==
1152             ....:   x.inverse().quadratic_representation() )
1153             True
1154
1155             sage: Qxy(J.one()) == x*y
1156             True
1157
1158         Property 4:
1159
1160             sage: not x.is_invertible() or (
1161             ....:   x.quadratic_representation(x.inverse())*Qx
1162             ....:   == Qx*x.quadratic_representation(x.inverse()) )
1163             True
1164
1165             sage: not x.is_invertible() or (
1166             ....:   x.quadratic_representation(x.inverse())*Qx
1167             ....:   ==
1168             ....:   2*Lx*Qex - Qx )
1169             True
1170
1171             sage: 2*Lx*Qex - Qx == Lxx
1172             True
1173
1174         Property 5:
1175
1176             sage: Qy(x).quadratic_representation() == Qy*Qx*Qy
1177             True
1178
1179         Property 6:
1180
1181             sage: Qxn == (Qx)^n
1182             True
1183
1184         Property 7:
1185
1186             sage: not x.is_invertible() or (
1187             ....:   Qx*x.inverse().operator() == Lx )
1188             True
1189
1190         Property 8:
1191
1192             sage: not x.operator_commutes_with(y) or (
1193             ....:   Qx(y)^n == Qxn(y^n) )
1194             True
1195
1196         """
1197         if other is None:
1198             other=self
1199         elif not other in self.parent():
1200             raise TypeError("'other' must live in the same algebra")
1201
1202         L = self.operator()
1203         M = other.operator()
1204         return ( L*M + M*L - (self*other).operator() )
1205
1206
1207
1208     def spectral_decomposition(self):
1209         """
1210         Return the unique spectral decomposition of this element.
1211
1212         ALGORITHM:
1213
1214         Following Faraut and Korányi's Theorem III.1.1, we restrict this
1215         element's left-multiplication-by operator to the subalgebra it
1216         generates. We then compute the spectral decomposition of that
1217         operator, and the spectral projectors we get back must be the
1218         left-multiplication-by operators for the idempotents we
1219         seek. Thus applying them to the identity element gives us those
1220         idempotents.
1221
1222         Since the eigenvalues are required to be distinct, we take
1223         the spectral decomposition of the zero element to be zero
1224         times the identity element of the algebra (which is idempotent,
1225         obviously).
1226
1227         SETUP::
1228
1229             sage: from mjo.eja.eja_algebra import RealSymmetricEJA
1230
1231         EXAMPLES:
1232
1233         The spectral decomposition of the identity is ``1`` times itself,
1234         and the spectral decomposition of zero is ``0`` times the identity::
1235
1236             sage: J = RealSymmetricEJA(3)
1237             sage: J.one()
1238             e0 + e2 + e5
1239             sage: J.one().spectral_decomposition()
1240             [(1, e0 + e2 + e5)]
1241             sage: J.zero().spectral_decomposition()
1242             [(0, e0 + e2 + e5)]
1243
1244         TESTS::
1245
1246             sage: J = RealSymmetricEJA(4)
1247             sage: x = sum(J.gens())
1248             sage: sd = x.spectral_decomposition()
1249             sage: l0 = sd[0][0]
1250             sage: l1 = sd[1][0]
1251             sage: c0 = sd[0][1]
1252             sage: c1 = sd[1][1]
1253             sage: c0.inner_product(c1) == 0
1254             True
1255             sage: c0.is_idempotent()
1256             True
1257             sage: c1.is_idempotent()
1258             True
1259             sage: c0 + c1 == J.one()
1260             True
1261             sage: l0*c0 + l1*c1 == x
1262             True
1263
1264         The spectral decomposition should work in subalgebras, too::
1265
1266             sage: J = RealSymmetricEJA(4)
1267             sage: (e0, e1, e2, e3, e4, e5, e6, e7, e8, e9) = J.gens()
1268             sage: A = 2*e5 - 2*e8
1269             sage: (lambda1, c1) = A.spectral_decomposition()[1]
1270             sage: (J0, J5, J1) = J.peirce_decomposition(c1)
1271             sage: (f0, f1, f2) = J1.gens()
1272             sage: f0.spectral_decomposition()
1273             [(0, f2), (1, f0)]
1274
1275         """
1276         A = self.subalgebra_generated_by(orthonormalize_basis=True)
1277         result = []
1278         for (evalue, proj) in A(self).operator().spectral_decomposition():
1279             result.append( (evalue, proj(A.one()).superalgebra_element()) )
1280         return result
1281
1282     def subalgebra_generated_by(self, orthonormalize_basis=False):
1283         """
1284         Return the associative subalgebra of the parent EJA generated
1285         by this element.
1286
1287         Since our parent algebra is unital, we want "subalgebra" to mean
1288         "unital subalgebra" as well; thus the subalgebra that an element
1289         generates will itself be a Euclidean Jordan algebra after
1290         restricting the algebra operations appropriately. This is the
1291         subalgebra that Faraut and Korányi work with in section II.2, for
1292         example.
1293
1294         SETUP::
1295
1296             sage: from mjo.eja.eja_algebra import random_eja
1297
1298         TESTS:
1299
1300         This subalgebra, being composed of only powers, is associative::
1301
1302             sage: set_random_seed()
1303             sage: x0 = random_eja().random_element()
1304             sage: A = x0.subalgebra_generated_by()
1305             sage: x,y,z = A.random_elements(3)
1306             sage: (x*y)*z == x*(y*z)
1307             True
1308
1309         Squaring in the subalgebra should work the same as in
1310         the superalgebra::
1311
1312             sage: set_random_seed()
1313             sage: x = random_eja().random_element()
1314             sage: A = x.subalgebra_generated_by()
1315             sage: A(x^2) == A(x)*A(x)
1316             True
1317
1318         By definition, the subalgebra generated by the zero element is
1319         the one-dimensional algebra generated by the identity
1320         element... unless the original algebra was trivial, in which
1321         case the subalgebra is trivial too::
1322
1323             sage: set_random_seed()
1324             sage: A = random_eja().zero().subalgebra_generated_by()
1325             sage: (A.is_trivial() and A.dimension() == 0) or A.dimension() == 1
1326             True
1327
1328         """
1329         from mjo.eja.eja_element_subalgebra import FiniteDimensionalEuclideanJordanElementSubalgebra
1330         return FiniteDimensionalEuclideanJordanElementSubalgebra(self, orthonormalize_basis)
1331
1332
1333     def subalgebra_idempotent(self):
1334         """
1335         Find an idempotent in the associative subalgebra I generate
1336         using Proposition 2.3.5 in Baes.
1337
1338         SETUP::
1339
1340             sage: from mjo.eja.eja_algebra import random_eja
1341
1342         TESTS:
1343
1344         Ensure that we can find an idempotent in a non-trivial algebra
1345         where there are non-nilpotent elements, or that we get the dumb
1346         solution in the trivial algebra::
1347
1348             sage: set_random_seed()
1349             sage: J = random_eja()
1350             sage: x = J.random_element()
1351             sage: while x.is_nilpotent() and not J.is_trivial():
1352             ....:     x = J.random_element()
1353             sage: c = x.subalgebra_idempotent()
1354             sage: c^2 == c
1355             True
1356
1357         """
1358         if self.parent().is_trivial():
1359             return self
1360
1361         if self.is_nilpotent():
1362             raise ValueError("this only works with non-nilpotent elements!")
1363
1364         J = self.subalgebra_generated_by()
1365         u = J(self)
1366
1367         # The image of the matrix of left-u^m-multiplication
1368         # will be minimal for some natural number s...
1369         s = 0
1370         minimal_dim = J.dimension()
1371         for i in range(1, minimal_dim):
1372             this_dim = (u**i).operator().matrix().image().dimension()
1373             if this_dim < minimal_dim:
1374                 minimal_dim = this_dim
1375                 s = i
1376
1377         # Now minimal_matrix should correspond to the smallest
1378         # non-zero subspace in Baes's (or really, Koecher's)
1379         # proposition.
1380         #
1381         # However, we need to restrict the matrix to work on the
1382         # subspace... or do we? Can't we just solve, knowing that
1383         # A(c) = u^(s+1) should have a solution in the big space,
1384         # too?
1385         #
1386         # Beware, solve_right() means that we're using COLUMN vectors.
1387         # Our FiniteDimensionalAlgebraElement superclass uses rows.
1388         u_next = u**(s+1)
1389         A = u_next.operator().matrix()
1390         c = J.from_vector(A.solve_right(u_next.to_vector()))
1391
1392         # Now c is the idempotent we want, but it still lives in the subalgebra.
1393         return c.superalgebra_element()
1394
1395
1396     def trace(self):
1397         """
1398         Return my trace, the sum of my eigenvalues.
1399
1400         In a trivial algebra, however you want to look at it, the trace is
1401         an empty sum for which we declare the result to be zero.
1402
1403         SETUP::
1404
1405             sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
1406             ....:                                  HadamardEJA,
1407             ....:                                  TrivialEJA,
1408             ....:                                  random_eja)
1409
1410         EXAMPLES::
1411
1412             sage: J = TrivialEJA()
1413             sage: J.zero().trace()
1414             0
1415
1416         ::
1417             sage: J = JordanSpinEJA(3)
1418             sage: x = sum(J.gens())
1419             sage: x.trace()
1420             2
1421
1422         ::
1423
1424             sage: J = HadamardEJA(5)
1425             sage: J.one().trace()
1426             5
1427
1428         TESTS:
1429
1430         The trace of an element is a real number::
1431
1432             sage: set_random_seed()
1433             sage: J = random_eja()
1434             sage: J.random_element().trace() in RLF
1435             True
1436
1437         """
1438         P = self.parent()
1439         r = P.rank()
1440
1441         if r == 0:
1442             # Special case for the trivial algebra where
1443             # the trace is an empty sum.
1444             return P.base_ring().zero()
1445
1446         p = P._charpoly_coefficients()[r-1]
1447         # The _charpoly_coeff function already adds the factor of
1448         # -1 to ensure that _charpoly_coeff(r-1) is really what
1449         # appears in front of t^{r-1} in the charpoly. However,
1450         # we want the negative of THAT for the trace.
1451         return -p(*self.to_vector())
1452
1453
1454     def trace_inner_product(self, other):
1455         """
1456         Return the trace inner product of myself and ``other``.
1457
1458         SETUP::
1459
1460             sage: from mjo.eja.eja_algebra import random_eja
1461
1462         TESTS:
1463
1464         The trace inner product is commutative, bilinear, and associative::
1465
1466             sage: set_random_seed()
1467             sage: J = random_eja()
1468             sage: x,y,z = J.random_elements(3)
1469             sage: # commutative
1470             sage: x.trace_inner_product(y) == y.trace_inner_product(x)
1471             True
1472             sage: # bilinear
1473             sage: a = J.base_ring().random_element();
1474             sage: actual = (a*(x+z)).trace_inner_product(y)
1475             sage: expected = ( a*x.trace_inner_product(y) +
1476             ....:              a*z.trace_inner_product(y) )
1477             sage: actual == expected
1478             True
1479             sage: actual = x.trace_inner_product(a*(y+z))
1480             sage: expected = ( a*x.trace_inner_product(y) +
1481             ....:              a*x.trace_inner_product(z) )
1482             sage: actual == expected
1483             True
1484             sage: # associative
1485             sage: (x*y).trace_inner_product(z) == y.trace_inner_product(x*z)
1486             True
1487
1488         """
1489         if not other in self.parent():
1490             raise TypeError("'other' must live in the same algebra")
1491
1492         return (self*other).trace()
1493
1494
1495     def trace_norm(self):
1496         """
1497         The norm of this element with respect to :meth:`trace_inner_product`.
1498
1499         SETUP::
1500
1501             sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
1502             ....:                                  HadamardEJA)
1503
1504         EXAMPLES::
1505
1506             sage: J = HadamardEJA(2)
1507             sage: x = sum(J.gens())
1508             sage: x.trace_norm()
1509             1.414213562373095?
1510
1511         ::
1512
1513             sage: J = JordanSpinEJA(4)
1514             sage: x = sum(J.gens())
1515             sage: x.trace_norm()
1516             2.828427124746190?
1517
1518         """
1519         return self.trace_inner_product(self).sqrt()