mjo/eja/euclidean_jordan_algebra.py

   1 """
   2 Euclidean Jordan Algebras. These are formally-real Jordan Algebras;
   3 specifically those where u^2 + v^2 = 0 implies that u = v = 0. They
   4 are used in optimization, and have some additional nice methods beyond
   5 what can be supported in a general Jordan Algebra.
   6 """
   7
   8 from sage.algebras.finite_dimensional_algebras.finite_dimensional_algebra import FiniteDimensionalAlgebra
   9 from sage.algebras.finite_dimensional_algebras.finite_dimensional_algebra_element import FiniteDimensionalAlgebraElement
  10
  11 class FiniteDimensionalEuclideanJordanAlgebra(FiniteDimensionalAlgebra):
  12     @staticmethod
  13     def __classcall__(cls, field, mult_table, names='e', category=None):
  14         fda = super(FiniteDimensionalEuclideanJordanAlgebra, cls)
  15         return fda.__classcall_private__(cls,
  16                                          field,
  17                                          mult_table,
  18                                          names,
  19                                          category)
  20
  21     def __init__(self, field, mult_table, names='e', category=None):
  22         fda = super(FiniteDimensionalEuclideanJordanAlgebra, self)
  23         fda.__init__(field, mult_table, names, category)
  24
  25
  26     def _repr_(self):
  27         """
  28         Return a string representation of ``self``.
  29         """
  30         return "Euclidean Jordan algebra of degree {} over {}".format(self.degree(), self.base_ring())
  31
  32     def rank(self):
  33         """
  34         Return the rank of this EJA.
  35         """
  36         raise NotImplementedError
  37
  38
  39     class Element(FiniteDimensionalAlgebraElement):
  40         """
  41         An element of a Euclidean Jordan algebra.
  42
  43         Since EJAs are commutative, the "right multiplication" matrix is
  44         also the left multiplication matrix and must be symmetric::
  45
  46             sage: set_random_seed()
  47             sage: J = eja_ln(5)
  48             sage: J.random_element().matrix().is_symmetric()
  49             True
  50
  51         """
  52
  53         def __pow__(self, n):
  54             """
  55             Return ``self`` raised to the power ``n``.
  56
  57             Jordan algebras are always power-associative; see for
  58             example Faraut and Koranyi, Proposition II.1.2 (ii).
  59             """
  60             A = self.parent()
  61             if n == 0:
  62                 return A.one()
  63             elif n == 1:
  64                 return self
  65             else:
  66                 return A.element_class(A, self.vector()*(self.matrix()**(n-1)))
  67
  68
  69         def span_of_powers(self):
  70             """
  71             Return the vector space spanned by successive powers of
  72             this element.
  73             """
  74             # The dimension of the subalgebra can't be greater than
  75             # the big algebra, so just put everything into a list
  76             # and let span() get rid of the excess.
  77             V = self.vector().parent()
  78             return V.span( (self**d).vector() for d in xrange(V.dimension()) )
  79
  80
  81         def degree(self):
  82             """
  83             Compute the degree of this element the straightforward way
  84             according to the definition; by appending powers to a list
  85             and figuring out its dimension (that is, whether or not
  86             they're linearly dependent).
  87
  88             EXAMPLES::
  89
  90                 sage: J = eja_ln(4)
  91                 sage: J.one().degree()
  92                 1
  93                 sage: e0,e1,e2,e3 = J.gens()
  94                 sage: (e0 - e1).degree()
  95                 2
  96
  97             In the spin factor algebra (of rank two), all elements that
  98             aren't multiples of the identity are regular::
  99
 100                 sage: set_random_seed()
 101                 sage: n = ZZ.random_element(1,10).abs()
 102                 sage: J = eja_ln(n)
 103                 sage: x = J.random_element()
 104                 sage: x == x.coefficient(0)*J.one() or x.degree() == 2
 105                 True
 106
 107             """
 108             return self.span_of_powers().dimension()
 109
 110
 111         def subalgebra_generated_by(self):
 112             """
 113             Return the subalgebra of the parent EJA generated by this element.
 114             """
 115             # First get the subspace spanned by the powers of myself...
 116             V = self.span_of_powers()
 117             F = self.base_ring()
 118
 119             # Now figure out the entries of the right-multiplication
 120             # matrix for the successive basis elements b0, b1,... of
 121             # that subspace.
 122             mats = []
 123             for b_right in V.basis():
 124                 eja_b_right = self.parent()(b_right)
 125                 b_right_rows = []
 126                 # The first row of the right-multiplication matrix by
 127                 # b1 is what we get if we apply that matrix to b1. The
 128                 # second row of the right multiplication matrix by b1
 129                 # is what we get when we apply that matrix to b2...
 130                 for b_left in V.basis():
 131                     eja_b_left = self.parent()(b_left)
 132                     # Multiply in the original EJA, but then get the
 133                     # coordinates from the subalgebra in terms of its
 134                     # basis.
 135                     this_row = V.coordinates((eja_b_left*eja_b_right).vector())
 136                     b_right_rows.append(this_row)
 137                 b_right_matrix = matrix(F, b_right_rows)
 138                 mats.append(b_right_matrix)
 139
 140             return FiniteDimensionalEuclideanJordanAlgebra(F, mats)
 141
 142
 143         def minimal_polynomial(self):
 144             """
 145             EXAMPLES::
 146
 147                 sage: set_random_seed()
 148                 sage: n = ZZ.random_element(1,10).abs()
 149                 sage: J = eja_rn(n)
 150                 sage: x = J.random_element()
 151                 sage: x.degree() == x.minimal_polynomial().degree()
 152                 True
 153
 154             ::
 155
 156                 sage: set_random_seed()
 157                 sage: n = ZZ.random_element(1,10).abs()
 158                 sage: J = eja_ln(n)
 159                 sage: x = J.random_element()
 160                 sage: x.degree() == x.minimal_polynomial().degree()
 161                 True
 162
 163             The minimal polynomial and the characteristic polynomial coincide
 164             and are known (see Alizadeh, Example 11.11) for all elements of
 165             the spin factor algebra that aren't scalar multiples of the
 166             identity::
 167
 168                 sage: set_random_seed()
 169                 sage: n = ZZ.random_element(2,10).abs()
 170                 sage: J = eja_ln(n)
 171                 sage: y = J.random_element()
 172                 sage: while y == y.coefficient(0)*J.one():
 173                 ....:     y = J.random_element()
 174                 sage: y0 = y.vector()[0]
 175                 sage: y_bar = y.vector()[1:]
 176                 sage: actual = y.minimal_polynomial()
 177                 sage: x = SR.symbol('x', domain='real')
 178                 sage: expected = x^2 - 2*y0*x + (y0^2 - norm(y_bar)^2)
 179                 sage: bool(actual == expected)
 180                 True
 181
 182             """
 183             V = self.span_of_powers()
 184             assoc_subalg = self.subalgebra_generated_by()
 185             # Mis-design warning: the basis used for span_of_powers()
 186             # and subalgebra_generated_by() must be the same, and in
 187             # the same order!
 188             subalg_self = assoc_subalg(V.coordinates(self.vector()))
 189             return subalg_self.matrix().minimal_polynomial()
 190
 191
 192         def characteristic_polynomial(self):
 193             return self.matrix().characteristic_polynomial()
 194
 195
 196 def eja_rn(dimension, field=QQ):
 197     """
 198     Return the Euclidean Jordan Algebra corresponding to the set
 199     `R^n` under the Hadamard product.
 200
 201     EXAMPLES:
 202
 203     This multiplication table can be verified by hand::
 204
 205         sage: J = eja_rn(3)
 206         sage: e0,e1,e2 = J.gens()
 207         sage: e0*e0
 208         e0
 209         sage: e0*e1
 210         0
 211         sage: e0*e2
 212         0
 213         sage: e1*e1
 214         e1
 215         sage: e1*e2
 216         0
 217         sage: e2*e2
 218         e2
 219
 220     """
 221     # The FiniteDimensionalAlgebra constructor takes a list of
 222     # matrices, the ith representing right multiplication by the ith
 223     # basis element in the vector space. So if e_1 = (1,0,0), then
 224     # right (Hadamard) multiplication of x by e_1 picks out the first
 225     # component of x; and likewise for the ith basis element e_i.
 226     Qs = [ matrix(field, dimension, dimension, lambda k,j: 1*(k == j == i))
 227            for i in xrange(dimension) ]
 228
 229     return FiniteDimensionalEuclideanJordanAlgebra(field,Qs)
 230
 231
 232 def eja_ln(dimension, field=QQ):
 233     """
 234     Return the Jordan algebra corresponding to the Lorentz "ice cream"
 235     cone of the given ``dimension``.
 236
 237     EXAMPLES:
 238
 239     This multiplication table can be verified by hand::
 240
 241         sage: J = eja_ln(4)
 242         sage: e0,e1,e2,e3 = J.gens()
 243         sage: e0*e0
 244         e0
 245         sage: e0*e1
 246         e1
 247         sage: e0*e2
 248         e2
 249         sage: e0*e3
 250         e3
 251         sage: e1*e2
 252         0
 253         sage: e1*e3
 254         0
 255         sage: e2*e3
 256         0
 257
 258     In one dimension, this is the reals under multiplication::
 259
 260       sage: J1 = eja_ln(1)
 261       sage: J2 = eja_rn(1)
 262       sage: J1 == J2
 263       True
 264
 265     """
 266     Qs = []
 267     id_matrix = identity_matrix(field,dimension)
 268     for i in xrange(dimension):
 269         ei = id_matrix.column(i)
 270         Qi = zero_matrix(field,dimension)
 271         Qi.set_row(0, ei)
 272         Qi.set_column(0, ei)
 273         Qi += diagonal_matrix(dimension, [ei[0]]*dimension)
 274         # The addition of the diagonal matrix adds an extra ei[0] in the
 275         # upper-left corner of the matrix.
 276         Qi[0,0] = Qi[0,0] * ~field(2)
 277         Qs.append(Qi)
 278
 279     return FiniteDimensionalEuclideanJordanAlgebra(field,Qs)