]>
gitweb.michael.orlitzky.com - sage.d.git/blob - eja_utils.py
79d8ecfce61c555375deffd40501b2fb100c0379
1 from sage
.structure
.element
import is_Matrix
3 def _charpoly_sage_input(s
):
5 Helper function that you can use on the string output from sage
6 to convert a charpoly coefficient into the corresponding input
11 sage: from mjo.eja.eja_algebra import JordanSpinEJA
12 sage: from mjo.eja.eja_utils import _charpoly_sage_input
16 sage: J = JordanSpinEJA(4,QQ)
17 sage: a = J._charpoly_coefficients()
19 X1^2 - X2^2 - X3^2 - X4^2
20 sage: _charpoly_sage_input(str(a[0]))
21 'X[0]**2 - X[1]**2 - X[2]**2 - X[3]**2'
29 digit_out
= r
"X([0-9]+)"
33 return "X[" + str(int(m
.group(1)) - 1) + "]"
35 s
= re
.sub(exponent_out
, exponent_in
, s
)
36 return re
.sub(digit_out
, replace_digit
, s
)
41 Scale the vector, matrix, or cartesian-product-of-those-things
44 This works around the inability to scale certain elements of
45 Cartesian product spaces, as reported in
47 https://trac.sagemath.org/ticket/31435
51 This will do the wrong thing if you feed it a tuple or list.
55 sage: from mjo.eja.eja_utils import _scale
59 sage: v = vector(QQ, (1,2,3))
62 sage: m = matrix(QQ, [[1,2],[3,4]])
63 sage: M = cartesian_product([m.parent(), m.parent()])
64 sage: _scale(M((m,m)), 2)
70 if hasattr(x
, 'cartesian_factors'):
72 return P(tuple( _scale(x_i
, alpha
)
73 for x_i
in x
.cartesian_factors() ))
80 Flatten a vector, matrix, or cartesian product of those things
83 If the entries of the matrix themselves belong to a real vector
84 space (such as the complex numbers which can be thought of as
85 pairs of real numbers), they will also be expanded in vector form
86 and flattened into the list.
90 sage: from mjo.eja.eja_utils import _all2list
91 sage: from mjo.hurwitz import (QuaternionMatrixAlgebra,
93 ....: OctonionMatrixAlgebra)
97 sage: _all2list([[1]])
102 sage: V1 = VectorSpace(QQ,2)
103 sage: V2 = MatrixSpace(QQ,2)
105 sage: x2 = V1([1,-1])
107 sage: y2 = V2([0,1,1,0])
108 sage: _all2list((x1,y1))
110 sage: _all2list((x2,y2))
112 sage: M = cartesian_product([V1,V2])
113 sage: _all2list(M((x1,y1)))
115 sage: _all2list(M((x2,y2)))
120 sage: _all2list(Octonions().one())
121 [1, 0, 0, 0, 0, 0, 0, 0]
122 sage: _all2list(OctonionMatrixAlgebra(1).one())
123 [1, 0, 0, 0, 0, 0, 0, 0]
127 sage: _all2list(QuaternionAlgebra(QQ, -1, -1).one())
129 sage: _all2list(QuaternionMatrixAlgebra(1).one())
134 sage: V1 = VectorSpace(QQ,2)
135 sage: V2 = OctonionMatrixAlgebra(1,field=QQ)
136 sage: C = cartesian_product([V1,V2])
139 sage: _all2list(C( (x1,y1) ))
140 [3, 4, 1, 0, 0, 0, 0, 0, 0, 0]
143 if hasattr(x
, 'to_vector'):
144 # This works on matrices of e.g. octonions directly, without
145 # first needing to convert them to a list of octonions and
146 # then recursing down into the list. It also avoids the wonky
147 # list(x) when x is an element of a CFM. I don't know what it
148 # returns but it aint the coordinates. We don't recurse
149 # because vectors can only contain ring elements as entries.
150 return x
.to_vector().list()
153 # This sucks, but for performance reasons we don't want to
154 # call _all2list recursively on the contents of a matrix
155 # when we don't have to (they only contain ring elements
161 except TypeError: # x is not iterable
165 # Avoid the retardation of list(QQ(1)) == [1].
168 return sum( map(_all2list
, xl
) , [])
171 def gram_schmidt(v
, inner_product
=None):
173 Perform Gram-Schmidt on the list ``v`` which are assumed to be
174 vectors over the same base ring. Returns a list of orthonormalized
175 vectors over the same base ring, which means that your base ring
176 needs to contain the appropriate roots.
180 sage: from mjo.eja.eja_utils import gram_schmidt
184 If you start with an orthonormal set, you get it back. We can use
185 the rationals here because we don't need any square roots::
187 sage: v1 = vector(QQ, (1,0,0))
188 sage: v2 = vector(QQ, (0,1,0))
189 sage: v3 = vector(QQ, (0,0,1))
191 sage: gram_schmidt(v) == v
194 The usual inner-product and norm are default::
196 sage: v1 = vector(AA,(1,2,3))
197 sage: v2 = vector(AA,(1,-1,6))
198 sage: v3 = vector(AA,(2,1,-1))
200 sage: u = gram_schmidt(v)
201 sage: all( u_i.inner_product(u_i).sqrt() == 1 for u_i in u )
203 sage: bool(u[0].inner_product(u[1]) == 0)
205 sage: bool(u[0].inner_product(u[2]) == 0)
207 sage: bool(u[1].inner_product(u[2]) == 0)
211 But if you supply a custom inner product, the result is
212 orthonormal with respect to that (and not the usual inner
215 sage: v1 = vector(AA,(1,2,3))
216 sage: v2 = vector(AA,(1,-1,6))
217 sage: v3 = vector(AA,(2,1,-1))
219 sage: B = matrix(AA, [ [6, 4, 2],
222 sage: ip = lambda x,y: (B*x).inner_product(y)
223 sage: norm = lambda x: ip(x,x)
224 sage: u = gram_schmidt(v,ip)
225 sage: all( norm(u_i) == 1 for u_i in u )
227 sage: ip(u[0],u[1]).is_zero()
229 sage: ip(u[0],u[2]).is_zero()
231 sage: ip(u[1],u[2]).is_zero()
234 This Gram-Schmidt routine can be used on matrices as well, so long
235 as an appropriate inner-product is provided::
237 sage: E11 = matrix(AA, [ [1,0],
239 sage: E12 = matrix(AA, [ [0,1],
241 sage: E22 = matrix(AA, [ [0,0],
243 sage: I = matrix.identity(AA,2)
244 sage: trace_ip = lambda X,Y: (X*Y).trace()
245 sage: gram_schmidt([E11,E12,I,E22], inner_product=trace_ip)
247 [1 0] [ 0 0.7071067811865475?] [0 0]
248 [0 0], [0.7071067811865475? 0], [0 1]
251 It even works on Cartesian product spaces whose factors are vector
254 sage: V1 = VectorSpace(AA,2)
255 sage: V2 = MatrixSpace(AA,2)
256 sage: M = cartesian_product([V1,V2])
258 sage: x2 = V1([1,-1])
260 sage: y2 = V2([0,1,1,0])
261 sage: z1 = M((x1,y1))
262 sage: z2 = M((x2,y2))
264 ....: return a[0].inner_product(b[0]) + (a[1]*b[1]).trace()
265 sage: U = gram_schmidt([z1,z2], inner_product=ip)
275 Ensure that zero vectors don't get in the way::
277 sage: v1 = vector(AA,(1,2,3))
278 sage: v2 = vector(AA,(1,-1,6))
279 sage: v3 = vector(AA,(0,0,0))
281 sage: len(gram_schmidt(v)) == 2
284 if inner_product
is None:
285 inner_product
= lambda x
,y
: x
.inner_product(y
)
287 ip
= inner_product(x
,x
)
288 # Don't expand the given field; the inner-product's codomain
289 # is already correct. For example QQ(2).sqrt() returns sqrt(2)
290 # in SR, and that will give you weird errors about symbolics
291 # when what's really going wrong is that you're trying to
292 # orthonormalize in QQ.
293 return ip
.parent()(ip
.sqrt())
295 v
= list(v
) # make a copy, don't clobber the input
297 # Drop all zero vectors before we start.
298 v
= [ v_i
for v_i
in v
if not v_i
.is_zero() ]
306 # Our "zero" needs to belong to the right space for sum() to work.
307 zero
= v
[0].parent().zero()
310 if hasattr(v
[0], 'cartesian_factors'):
311 # Only use the slow implementation if necessary.
315 return sc(x
, (inner_product(x
,y
)/inner_product(x
,x
)))
317 # First orthogonalize...
318 for i
in range(1,len(v
)):
319 # Earlier vectors can be made into zero so we have to ignore them.
320 v
[i
] -= sum( (proj(v
[j
],v
[i
])
322 if not v
[j
].is_zero() ),
325 # And now drop all zero vectors again if they were "orthogonalized out."
326 v
= [ v_i
for v_i
in v
if not v_i
.is_zero() ]
328 # Just normalize. If the algebra is missing the roots, we can't add
329 # them here because then our subalgebra would have a bigger field
330 # than the superalgebra.
331 for i
in range(len(v
)):
332 v
[i
] = sc(v
[i
], ~
norm(v
[i
]))