]> gitweb.michael.orlitzky.com - sage.d.git/blob - mjo/eja/eja_element_subalgebra.py
mjo/ldlt.py: various efficiency improvements to ldlt_fast().
[sage.d.git] / mjo / eja / eja_element_subalgebra.py
1 from sage.matrix.constructor import matrix
2
3 from mjo.eja.eja_subalgebra import FiniteDimensionalEuclideanJordanSubalgebra
4
5
6 class FiniteDimensionalEuclideanJordanElementSubalgebra(FiniteDimensionalEuclideanJordanSubalgebra):
7 def __init__(self, elt, orthonormalize_basis):
8 self._superalgebra = elt.parent()
9 category = self._superalgebra.category().Associative()
10 V = self._superalgebra.vector_space()
11 field = self._superalgebra.base_ring()
12
13 # This list is guaranteed to contain all independent powers,
14 # because it's the maximal set of powers that could possibly
15 # be independent (by a dimension argument).
16 powers = [ elt**k for k in range(V.dimension()) ]
17 power_vectors = [ p.to_vector() for p in powers ]
18 P = matrix(field, power_vectors)
19
20 if orthonormalize_basis == False:
21 # In this case, we just need to figure out which elements
22 # of the "powers" list are redundant... First compute the
23 # vector subspace spanned by the powers of the given
24 # element.
25
26 # Figure out which powers form a linearly-independent set.
27 ind_rows = P.pivot_rows()
28
29 # Pick those out of the list of all powers.
30 superalgebra_basis = tuple(map(powers.__getitem__, ind_rows))
31
32 # If our superalgebra is a subalgebra of something else, then
33 # these vectors won't have the right coordinates for
34 # V.span_of_basis() unless we use V.from_vector() on them.
35 basis_vectors = map(power_vectors.__getitem__, ind_rows)
36 else:
37 # If we're going to orthonormalize the basis anyway, we
38 # might as well just do Gram-Schmidt on the whole list of
39 # powers. The redundant ones will get zero'd out. If this
40 # looks like a roundabout way to orthonormalize, it is.
41 # But converting everything from algebra elements to vectors
42 # to matrices and then back again turns out to be about
43 # as fast as reimplementing our own Gram-Schmidt that
44 # works in an EJA.
45 G,_ = P.gram_schmidt(orthonormal=True)
46 basis_vectors = [ g for g in G.rows() if not g.is_zero() ]
47 superalgebra_basis = [ self._superalgebra.from_vector(b)
48 for b in basis_vectors ]
49
50 W = V.span_of_basis( V.from_vector(v) for v in basis_vectors )
51
52 # The rank is the highest possible degree of a minimal
53 # polynomial, and is bounded above by the dimension. We know
54 # in this case that there's an element whose minimal
55 # polynomial has the same degree as the space's dimension
56 # (remember how we constructed the space?), so that must be
57 # its rank too.
58 rank = W.dimension()
59
60 fdeja = super(FiniteDimensionalEuclideanJordanElementSubalgebra, self)
61 return fdeja.__init__(self._superalgebra,
62 superalgebra_basis,
63 rank=rank,
64 category=category)
65
66
67 def _a_regular_element(self):
68 """
69 Override the superalgebra method to return the one
70 regular element that is sure to exist in this
71 subalgebra, namely the element that generated it.
72
73 SETUP::
74
75 sage: from mjo.eja.eja_algebra import random_eja
76
77 TESTS::
78
79 sage: set_random_seed()
80 sage: J = random_eja().random_element().subalgebra_generated_by()
81 sage: J._a_regular_element().is_regular()
82 True
83
84 """
85 if self.dimension() == 0:
86 return self.zero()
87 else:
88 return self.monomial(1)
89
90
91 def one(self):
92 """
93 Return the multiplicative identity element of this algebra.
94
95 The superclass method computes the identity element, which is
96 beyond overkill in this case: the superalgebra identity
97 restricted to this algebra is its identity. Note that we can't
98 count on the first basis element being the identity -- it migth
99 have been scaled if we orthonormalized the basis.
100
101 SETUP::
102
103 sage: from mjo.eja.eja_algebra import (HadamardEJA,
104 ....: random_eja)
105
106 EXAMPLES::
107
108 sage: J = HadamardEJA(5)
109 sage: J.one()
110 e0 + e1 + e2 + e3 + e4
111 sage: x = sum(J.gens())
112 sage: A = x.subalgebra_generated_by()
113 sage: A.one()
114 f0
115 sage: A.one().superalgebra_element()
116 e0 + e1 + e2 + e3 + e4
117
118 TESTS:
119
120 The identity element acts like the identity over the rationals::
121
122 sage: set_random_seed()
123 sage: x = random_eja().random_element()
124 sage: A = x.subalgebra_generated_by()
125 sage: x = A.random_element()
126 sage: A.one()*x == x and x*A.one() == x
127 True
128
129 The identity element acts like the identity over the algebraic
130 reals with an orthonormal basis::
131
132 sage: set_random_seed()
133 sage: x = random_eja(AA).random_element()
134 sage: A = x.subalgebra_generated_by(orthonormalize_basis=True)
135 sage: x = A.random_element()
136 sage: A.one()*x == x and x*A.one() == x
137 True
138
139 The matrix of the unit element's operator is the identity over
140 the rationals::
141
142 sage: set_random_seed()
143 sage: x = random_eja().random_element()
144 sage: A = x.subalgebra_generated_by()
145 sage: actual = A.one().operator().matrix()
146 sage: expected = matrix.identity(A.base_ring(), A.dimension())
147 sage: actual == expected
148 True
149
150 The matrix of the unit element's operator is the identity over
151 the algebraic reals with an orthonormal basis::
152
153 sage: set_random_seed()
154 sage: x = random_eja(AA).random_element()
155 sage: A = x.subalgebra_generated_by(orthonormalize_basis=True)
156 sage: actual = A.one().operator().matrix()
157 sage: expected = matrix.identity(A.base_ring(), A.dimension())
158 sage: actual == expected
159 True
160
161 """
162 if self.dimension() == 0:
163 return self.zero()
164 else:
165 sa_one = self.superalgebra().one().to_vector()
166 sa_coords = self.vector_space().coordinate_vector(sa_one)
167 return self.from_vector(sa_coords)