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