]>
gitweb.michael.orlitzky.com - sage.d.git/blob - eja_utils.py
a4328610e5e41db689455828bd0d8988225e745b
1 from sage
.functions
.other
import sqrt
2 from sage
.structure
.element
import is_Matrix
3 from sage
.matrix
.constructor
import matrix
4 from sage
.modules
.free_module_element
import vector
6 def _charpoly_sage_input(s
):
8 Helper function that you can use on the string output from sage
9 to convert a charpoly coefficient into the corresponding input
14 sage: from mjo.eja.eja_algebra import JordanSpinEJA
15 sage: from mjo.eja.eja_utils import _charpoly_sage_input
19 sage: J = JordanSpinEJA(4,QQ)
20 sage: a = J._charpoly_coefficients()
22 X1^2 - X2^2 - X3^2 - X4^2
23 sage: _charpoly_sage_input(str(a[0]))
24 'X[0]**2 - X[1]**2 - X[2]**2 - X[3]**2'
32 digit_out
= r
"X([0-9]+)"
36 return "X[" + str(int(m
.group(1)) - 1) + "]"
38 s
= re
.sub(exponent_out
, exponent_in
, s
)
39 return re
.sub(digit_out
, replace_digit
, s
)
44 Scale the vector, matrix, or cartesian-product-of-those-things
47 This works around the inability to scale certain elements of
48 Cartesian product spaces, as reported in
50 https://trac.sagemath.org/ticket/31435
54 This will do the wrong thing if you feed it a tuple or list.
58 sage: from mjo.eja.eja_utils import _scale
62 sage: v = vector(QQ, (1,2,3))
65 sage: m = matrix(QQ, [[1,2],[3,4]])
66 sage: M = cartesian_product([m.parent(), m.parent()])
67 sage: _scale(M((m,m)), 2)
73 if hasattr(x
, 'cartesian_factors'):
75 return P(tuple( _scale(x_i
, alpha
)
76 for x_i
in x
.cartesian_factors() ))
83 Flatten a vector, matrix, or cartesian product of those things
86 If the entries of the matrix themselves belong to a real vector
87 space (such as the complex numbers which can be thought of as
88 pairs of real numbers), they will also be expanded in vector form
89 and flattened into the list.
93 sage: from mjo.eja.eja_utils import _all2list
94 sage: from mjo.hurwitz import (QuaternionMatrixAlgebra,
96 ....: OctonionMatrixAlgebra)
100 sage: _all2list([[1]])
105 sage: V1 = VectorSpace(QQ,2)
106 sage: V2 = MatrixSpace(QQ,2)
108 sage: x2 = V1([1,-1])
110 sage: y2 = V2([0,1,1,0])
111 sage: _all2list((x1,y1))
113 sage: _all2list((x2,y2))
115 sage: M = cartesian_product([V1,V2])
116 sage: _all2list(M((x1,y1)))
118 sage: _all2list(M((x2,y2)))
123 sage: _all2list(Octonions().one())
124 [1, 0, 0, 0, 0, 0, 0, 0]
125 sage: _all2list(OctonionMatrixAlgebra(1).one())
126 [1, 0, 0, 0, 0, 0, 0, 0]
130 sage: _all2list(QuaternionAlgebra(QQ, -1, -1).one())
132 sage: _all2list(QuaternionMatrixAlgebra(1).one())
137 sage: V1 = VectorSpace(QQ,2)
138 sage: V2 = OctonionMatrixAlgebra(1,field=QQ)
139 sage: C = cartesian_product([V1,V2])
142 sage: _all2list(C( (x1,y1) ))
143 [3, 4, 1, 0, 0, 0, 0, 0, 0, 0]
146 if hasattr(x
, 'to_vector'):
147 # This works on matrices of e.g. octonions directly, without
148 # first needing to convert them to a list of octonions and
149 # then recursing down into the list. It also avoids the wonky
150 # list(x) when x is an element of a CFM. I don't know what it
151 # returns but it aint the coordinates. We don't recurse
152 # because vectors can only contain ring elements as entries.
153 return x
.to_vector().list()
156 # This sucks, but for performance reasons we don't want to
157 # call _all2list recursively on the contents of a matrix
158 # when we don't have to (they only contain ring elements
164 except TypeError: # x is not iterable
168 # Avoid the retardation of list(QQ(1)) == [1].
171 return sum( map(_all2list
, xl
) , [])
176 return vector(m
.base_ring(), m
.list())
179 return matrix(v
.base_ring(), sqrt(v
.degree()), v
.list())
181 def gram_schmidt(v
, inner_product
=None):
183 Perform Gram-Schmidt on the list ``v`` which are assumed to be
184 vectors over the same base ring. Returns a list of orthonormalized
185 vectors over the same base ring, which means that your base ring
186 needs to contain the appropriate roots.
190 sage: from mjo.eja.eja_utils import gram_schmidt
194 If you start with an orthonormal set, you get it back. We can use
195 the rationals here because we don't need any square roots::
197 sage: v1 = vector(QQ, (1,0,0))
198 sage: v2 = vector(QQ, (0,1,0))
199 sage: v3 = vector(QQ, (0,0,1))
201 sage: gram_schmidt(v) == v
204 The usual inner-product and norm are default::
206 sage: v1 = vector(AA,(1,2,3))
207 sage: v2 = vector(AA,(1,-1,6))
208 sage: v3 = vector(AA,(2,1,-1))
210 sage: u = gram_schmidt(v)
211 sage: all( u_i.inner_product(u_i).sqrt() == 1 for u_i in u )
213 sage: bool(u[0].inner_product(u[1]) == 0)
215 sage: bool(u[0].inner_product(u[2]) == 0)
217 sage: bool(u[1].inner_product(u[2]) == 0)
221 But if you supply a custom inner product, the result is
222 orthonormal with respect to that (and not the usual inner
225 sage: v1 = vector(AA,(1,2,3))
226 sage: v2 = vector(AA,(1,-1,6))
227 sage: v3 = vector(AA,(2,1,-1))
229 sage: B = matrix(AA, [ [6, 4, 2],
232 sage: ip = lambda x,y: (B*x).inner_product(y)
233 sage: norm = lambda x: ip(x,x)
234 sage: u = gram_schmidt(v,ip)
235 sage: all( norm(u_i) == 1 for u_i in u )
237 sage: ip(u[0],u[1]).is_zero()
239 sage: ip(u[0],u[2]).is_zero()
241 sage: ip(u[1],u[2]).is_zero()
244 This Gram-Schmidt routine can be used on matrices as well, so long
245 as an appropriate inner-product is provided::
247 sage: E11 = matrix(AA, [ [1,0],
249 sage: E12 = matrix(AA, [ [0,1],
251 sage: E22 = matrix(AA, [ [0,0],
253 sage: I = matrix.identity(AA,2)
254 sage: trace_ip = lambda X,Y: (X*Y).trace()
255 sage: gram_schmidt([E11,E12,I,E22], inner_product=trace_ip)
257 [1 0] [ 0 0.7071067811865475?] [0 0]
258 [0 0], [0.7071067811865475? 0], [0 1]
261 It even works on Cartesian product spaces whose factors are vector
264 sage: V1 = VectorSpace(AA,2)
265 sage: V2 = MatrixSpace(AA,2)
266 sage: M = cartesian_product([V1,V2])
268 sage: x2 = V1([1,-1])
270 sage: y2 = V2([0,1,1,0])
271 sage: z1 = M((x1,y1))
272 sage: z2 = M((x2,y2))
274 ....: return a[0].inner_product(b[0]) + (a[1]*b[1]).trace()
275 sage: U = gram_schmidt([z1,z2], inner_product=ip)
285 Ensure that zero vectors don't get in the way::
287 sage: v1 = vector(AA,(1,2,3))
288 sage: v2 = vector(AA,(1,-1,6))
289 sage: v3 = vector(AA,(0,0,0))
291 sage: len(gram_schmidt(v)) == 2
294 if inner_product
is None:
295 inner_product
= lambda x
,y
: x
.inner_product(y
)
297 ip
= inner_product(x
,x
)
298 # Don't expand the given field; the inner-product's codomain
299 # is already correct. For example QQ(2).sqrt() returns sqrt(2)
300 # in SR, and that will give you weird errors about symbolics
301 # when what's really going wrong is that you're trying to
302 # orthonormalize in QQ.
303 return ip
.parent()(ip
.sqrt())
305 v
= list(v
) # make a copy, don't clobber the input
307 # Drop all zero vectors before we start.
308 v
= [ v_i
for v_i
in v
if not v_i
.is_zero() ]
316 # Our "zero" needs to belong to the right space for sum() to work.
317 zero
= v
[0].parent().zero()
320 if hasattr(v
[0], 'cartesian_factors'):
321 # Only use the slow implementation if necessary.
325 return sc(x
, (inner_product(x
,y
)/inner_product(x
,x
)))
327 # First orthogonalize...
328 for i
in range(1,len(v
)):
329 # Earlier vectors can be made into zero so we have to ignore them.
330 v
[i
] -= sum( (proj(v
[j
],v
[i
])
332 if not v
[j
].is_zero() ),
335 # And now drop all zero vectors again if they were "orthogonalized out."
336 v
= [ v_i
for v_i
in v
if not v_i
.is_zero() ]
338 # Just normalize. If the algebra is missing the roots, we can't add
339 # them here because then our subalgebra would have a bigger field
340 # than the superalgebra.
341 for i
in range(len(v
)):
342 v
[i
] = sc(v
[i
], ~
norm(v
[i
]))