mjo/matrix_algebra.py

   1 from sage.misc.table import table
   2 from sage.categories.magmatic_algebras import MagmaticAlgebras
   3 from sage.misc.cachefunc import cached_method
   4 from sage.combinat.free_module import CombinatorialFreeModule
   5 from sage.modules.with_basis.indexed_element import IndexedFreeModuleElement
   6
   7 class MatrixAlgebraElement(IndexedFreeModuleElement):
   8     def nrows(self):
   9         return self.parent().nrows()
  10     ncols = nrows
  11
  12     @cached_method
  13     def rows(self):
  14         r"""
  15         SETUP::
  16
  17             sage: from mjo.matrix_algebra import MatrixAlgebra
  18
  19         EXAMPLES::
  20
  21             sage: M = MatrixAlgebra(2, QQbar,RDF)
  22             sage: A = M.monomial((0,0,1)) + 4*M.monomial((0,1,1))
  23             sage: A
  24             +-----+-----+
  25             | 1.0 | 4.0 |
  26             +-----+-----+
  27             | 0   | 0   |
  28             +-----+-----+
  29             sage: A.rows()
  30             [[1.0, 4.0], [0, 0]]
  31
  32         """
  33         zero = self.parent().entry_algebra().zero()
  34         l = [[zero for j in range(self.ncols())] for i in range(self.nrows())]
  35         for (k,v) in self.monomial_coefficients().items():
  36             (i,j,e) = k
  37             l[i][j] += v*e
  38         return l
  39
  40     def _repr_(self):
  41         r"""
  42         Display this matrix as a table.
  43
  44         The SageMath Matrix class representation is not easily reusable,
  45         but using a table fakes it.
  46
  47         SETUP::
  48
  49             sage: from mjo.matrix_algebra import MatrixAlgebra
  50
  51         EXAMPLES::
  52
  53             sage: MatrixAlgebra(2,ZZ,ZZ).zero()
  54             +---+---+
  55             | 0 | 0 |
  56             +---+---+
  57             | 0 | 0 |
  58             +---+---+
  59
  60         TESTS::
  61
  62             sage: MatrixAlgebra(0,ZZ,ZZ).zero()
  63             []
  64
  65         """
  66         if self.nrows() == 0 or self.ncols() == 0:
  67             # Otherwise we get a crash or a blank space, depending on
  68             # how hard we work for it. This is what MatrixSpace(...,
  69             # 0) returns.
  70             return "[]"
  71
  72         return table(self.rows(), frame=True)._repr_()
  73
  74
  75     def list(self):
  76         r"""
  77         Return one long list of this matrix's entries.
  78
  79         SETUP::
  80
  81             sage: from mjo.matrix_algebra import MatrixAlgebra
  82
  83         EXAMPLES::
  84
  85             sage: A = MatrixAlgebra(2,ZZ,ZZ)
  86             sage: A([[1,2],[3,4]]).list()
  87             [1, 2, 3, 4]
  88
  89         """
  90         return sum( self.rows(), [] )
  91
  92
  93     def __getitem__(self, indices):
  94         r"""
  95
  96         SETUP::
  97
  98             sage: from mjo.matrix_algebra import MatrixAlgebra
  99
 100         EXAMPLES::
 101
 102             sage: M = MatrixAlgebra(2,ZZ,ZZ)([[1,2],[3,4]])
 103             sage: M[0,0]
 104             1
 105             sage: M[0,1]
 106             2
 107             sage: M[1,0]
 108             3
 109             sage: M[1,1]
 110             4
 111
 112         """
 113         i,j = indices
 114         return self.rows()[i][j]
 115
 116     def trace(self):
 117         r"""
 118         Return the sum of this matrix's diagonal entries.
 119
 120         SETUP::
 121
 122             sage: from mjo.matrix_algebra import MatrixAlgebra
 123
 124         EXAMPLES:
 125
 126         The trace (being a sum of entries) belongs to the same algebra
 127         as those entries, and NOT the scalar ring::
 128
 129             sage: entries = MatrixSpace(ZZ,2)
 130             sage: scalars = ZZ
 131             sage: M = MatrixAlgebra(2, entries, scalars)
 132             sage: I = entries.one()
 133             sage: Z = entries.zero()
 134             sage: M([[I,Z],[Z,I]]).trace()
 135             [2 0]
 136             [0 2]
 137
 138         """
 139         zero = self.parent().entry_algebra().zero()
 140         return sum( (self[i,i] for i in range(self.nrows())), zero )
 141
 142     def matrix_space(self):
 143         r"""
 144
 145         SETUP::
 146
 147             sage: from mjo.matrix_algebra import MatrixAlgebra
 148
 149         TESTS::
 150
 151             sage: set_random_seed()
 152             sage: entries = QuaternionAlgebra(QQ,-1,-1)
 153             sage: M = MatrixAlgebra(3, entries, QQ)
 154             sage: M.random_element().matrix_space() == M
 155             True
 156
 157         """
 158         return self.parent()
 159
 160
 161 class MatrixAlgebra(CombinatorialFreeModule):
 162     r"""
 163     An algebra of ``n``-by-``n`` matrices over an arbitrary scalar
 164     ring whose entries come from a magmatic algebra that need not
 165     be the same as the scalars.
 166
 167     The usual matrix spaces in SageMath don't support separate spaces
 168     for the entries and the scalars; in particular they assume that
 169     the entries come from a commutative and associative ring. This
 170     is problematic in several interesting matrix algebras, like those
 171     where the entries are quaternions or octonions.
 172
 173     SETUP::
 174
 175         sage: from mjo.matrix_algebra import MatrixAlgebra
 176
 177     EXAMPLES::
 178
 179     The existence of a unit element is determined dynamically::
 180
 181         sage: MatrixAlgebra(2,ZZ,ZZ).one()
 182         +---+---+
 183         | 1 | 0 |
 184         +---+---+
 185         | 0 | 1 |
 186         +---+---+
 187
 188     """
 189     Element = MatrixAlgebraElement
 190
 191     def __init__(self, n, entry_algebra, scalars, prefix="A", **kwargs):
 192
 193         category = MagmaticAlgebras(scalars).FiniteDimensional()
 194         category = category.WithBasis()
 195
 196         if "Unital" in entry_algebra.category().axioms():
 197             category = category.Unital()
 198             entry_one = entry_algebra.one()
 199             self.one = lambda: self.sum( (self.monomial((i,i,entry_one))
 200                                      for i in range(self.nrows()) ) )
 201
 202         if "Associative" in entry_algebra.category().axioms():
 203             category = category.Associative()
 204
 205         self._nrows = n
 206
 207         # Since the scalar ring is real but the entries are not,
 208         # sticking a "1" in each position doesn't give us a basis for
 209         # the space. We actually need to stick each of e0, e1, ...  (a
 210         # basis for the entry algebra itself) into each position.
 211         self._entry_algebra = entry_algebra
 212
 213         # Needs to make the (overridden) method call when, for example,
 214         # the entry algebra is the complex numbers and its gens() method
 215         # lies to us.
 216         entry_basis = self.entry_algebra_gens()
 217
 218         basis_indices = [(i,j,e) for i in range(n)
 219                                  for j in range(n)
 220                                  for e in entry_basis]
 221
 222         super().__init__(scalars,
 223                          basis_indices,
 224                          category=category,
 225                          prefix=prefix,
 226                          bracket='(')
 227
 228     def _repr_(self):
 229         return ("Module of %d by %d matrices with entries in %s"
 230                 " over the scalar ring %s" %
 231                 (self.nrows(),
 232                  self.ncols(),
 233                  self.entry_algebra(),
 234                  self.base_ring()) )
 235
 236     def entry_algebra(self):
 237         r"""
 238         Return the algebra that our elements' entries come from.
 239         """
 240         return self._entry_algebra
 241
 242     def entry_algebra_gens(self):
 243         r"""
 244         Return a tuple of the generators of (that is, a basis for) the
 245         entries of this matrix algebra.
 246
 247         This can be overridden in subclasses to work around the
 248         inconsistency in the ``gens()`` methods of the various
 249         entry algebras.
 250         """
 251         return self.entry_algebra().gens()
 252
 253     def _entry_algebra_element_to_vector(self, entry):
 254         r"""
 255         Return a vector representation (of length equal to the cardinality
 256         of :meth:`entry_algebra_gens`) of the given ``entry``.
 257
 258         This can be overridden in subclasses to work around the fact that
 259         real numbers, complex numbers, quaternions, et cetera, all require
 260         different incantations to turn them into a vector.
 261
 262         It only makes sense to "guess" here in the superclass when no
 263         subclass that overrides :meth:`entry_algebra_gens` exists. So
 264         if you have a special subclass for your annoying entry algebra,
 265         override this with the correct implementation there instead of
 266         adding a bunch of awkward cases to this superclass method.
 267
 268         SETUP::
 269
 270             sage: from mjo.hurwitz import Octonions
 271             sage: from mjo.matrix_algebra import MatrixAlgebra
 272
 273         EXAMPLES:
 274
 275         Real numbers::
 276
 277             sage: A = MatrixAlgebra(1, AA, QQ)
 278             sage: A._entry_algebra_element_to_vector(AA(17))
 279             (17)
 280
 281         Octonions::
 282
 283             sage: A = MatrixAlgebra(1, Octonions(), QQ)
 284             sage: e = A.entry_algebra_gens()
 285             sage: A._entry_algebra_element_to_vector(e[0])
 286             (1, 0, 0, 0, 0, 0, 0, 0)
 287             sage: A._entry_algebra_element_to_vector(e[1])
 288             (0, 1, 0, 0, 0, 0, 0, 0)
 289             sage: A._entry_algebra_element_to_vector(e[2])
 290             (0, 0, 1, 0, 0, 0, 0, 0)
 291             sage: A._entry_algebra_element_to_vector(e[3])
 292             (0, 0, 0, 1, 0, 0, 0, 0)
 293             sage: A._entry_algebra_element_to_vector(e[4])
 294             (0, 0, 0, 0, 1, 0, 0, 0)
 295             sage: A._entry_algebra_element_to_vector(e[5])
 296             (0, 0, 0, 0, 0, 1, 0, 0)
 297             sage: A._entry_algebra_element_to_vector(e[6])
 298             (0, 0, 0, 0, 0, 0, 1, 0)
 299             sage: A._entry_algebra_element_to_vector(e[7])
 300             (0, 0, 0, 0, 0, 0, 0, 1)
 301
 302         Sage matrices::
 303
 304             sage: MS = MatrixSpace(QQ,2)
 305             sage: A = MatrixAlgebra(1, MS, QQ)
 306             sage: A._entry_algebra_element_to_vector(MS([[1,2],[3,4]]))
 307             (1, 2, 3, 4)
 308
 309         """
 310         if hasattr(entry, 'to_vector'):
 311             return entry.to_vector()
 312
 313         from sage.modules.free_module import FreeModule
 314         d = len(self.entry_algebra_gens())
 315         V = FreeModule(self.entry_algebra().base_ring(), d)
 316
 317         if hasattr(entry, 'list'):
 318             # sage matrices
 319             return V(entry.list())
 320
 321         # This works in AA, and will crash if it doesn't know what to
 322         # do, and that's fine because then I don't know what to do
 323         # either.
 324         return V((entry,))
 325
 326
 327
 328     def nrows(self):
 329         return self._nrows
 330     ncols = nrows
 331
 332     def product_on_basis(self, mon1, mon2):
 333         r"""
 334
 335         SETUP::
 336
 337             sage: from mjo.hurwitz import Octonions
 338             sage: from mjo.matrix_algebra import MatrixAlgebra
 339
 340         TESTS::
 341
 342             sage: O = Octonions(QQ)
 343             sage: e = O.gens()
 344             sage: e[2]*e[1]
 345             -e3
 346             sage: A = MatrixAlgebra(2,O,QQ)
 347             sage: A.product_on_basis( (0,0,e[2]), (0,0,e[1]) )
 348             +-----+---+
 349             | -e3 | 0 |
 350             +-----+---+
 351             | 0   | 0 |
 352             +-----+---+
 353
 354         """
 355         (i,j,e1) = mon1
 356         (k,l,e2) = mon2
 357         if j == k:
 358             # There's no reason to expect e1*e2 to itself be a monomial,
 359             # so we have to do some manual conversion to get one.
 360             p = self._entry_algebra_element_to_vector(e1*e2)
 361
 362             # We have to convert alpha_g because a priori it lives in the
 363             # base ring of the entry algebra.
 364             R = self.base_ring()
 365             return self.sum_of_terms( (((i,l,g), R(alpha_g))
 366                                       for (alpha_g, g)
 367                                       in zip(p, self.entry_algebra_gens()) ),
 368                                       distinct=True)
 369         else:
 370             return self.zero()
 371
 372     def from_list(self, entries):
 373         r"""
 374         Construct an element of this algebra from a list of lists of
 375         entries.
 376
 377         SETUP::
 378
 379             sage: from mjo.hurwitz import ComplexMatrixAlgebra
 380
 381         EXAMPLES::
 382
 383             sage: A = ComplexMatrixAlgebra(2, QQbar, ZZ)
 384             sage: M = A.from_list([[0,I],[-I,0]])
 385             sage: M
 386             +----+---+
 387             | 0  | I |
 388             +----+---+
 389             | -I | 0 |
 390             +----+---+
 391             sage: M.to_vector()
 392             (0, 0, 0, 1, 0, -1, 0, 0)
 393
 394         """
 395         nrows = len(entries)
 396         ncols = 0
 397         if nrows > 0:
 398             ncols = len(entries[0])
 399
 400         if (not all( len(r) == ncols for r in entries )) or (ncols != nrows):
 401             raise ValueError("list must be square")
 402
 403         def convert(e_ij):
 404             if e_ij in self.entry_algebra():
 405                 # Don't re-create an element if it already lives where
 406                 # it should!
 407                 return e_ij
 408
 409             try:
 410                 # This branch works with e.g. QQbar, where no
 411                 # to/from_vector() methods are available.
 412                 return self.entry_algebra()(e_ij)
 413             except TypeError:
 414                 # We have to pass through vectors to convert from the
 415                 # given entry algebra to ours. Otherwise we can fail to
 416                 # convert an element of (for example) Octonions(QQ) to
 417                 # Octonions(AA).
 418                 return self.entry_algebra().from_vector(e_ij.to_vector())
 419
 420         def entry_to_element(i,j,entry):
 421             # Convert an entry at i,j to a matrix whose only non-zero
 422             # entry is i,j and corresponds to the entry.
 423             p = self._entry_algebra_element_to_vector(entry)
 424
 425             # We have to convert alpha_g because a priori it lives in the
 426             # base ring of the entry algebra.
 427             R = self.base_ring()
 428             return self.sum_of_terms( (((i,j,g), R(alpha_g))
 429                                       for (alpha_g, g)
 430                                       in zip(p, self.entry_algebra_gens()) ),
 431                                       distinct=True)
 432
 433         return self.sum( entry_to_element(i,j,entries[i][j])
 434                          for j in range(ncols)
 435                          for i in range(nrows) )
 436
 437
 438     def _element_constructor_(self, elt):
 439         if elt in self:
 440             return self
 441         else:
 442             return self.from_list(elt)