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