X-Git-Url: http://gitweb.michael.orlitzky.com/?a=blobdiff_plain;f=mjo%2Feja%2Feja_utils.py;h=8334f516094fcbea657637da037f825c1cf4cbd2;hb=dc92538d3fc92d16c9b6432ad17c37cb0d6b2be9;hp=79d8ecfce61c555375deffd40501b2fb100c0379;hpb=40850626cb85d115363995ad6beaee8cb17d83da;p=sage.d.git diff --git a/mjo/eja/eja_utils.py b/mjo/eja/eja_utils.py index 79d8ecf..8334f51 100644 --- a/mjo/eja/eja_utils.py +++ b/mjo/eja/eja_utils.py @@ -1,41 +1,5 @@ from sage.structure.element import is_Matrix -def _charpoly_sage_input(s): - r""" - Helper function that you can use on the string output from sage - to convert a charpoly coefficient into the corresponding input - to be cached. - - SETUP:: - - sage: from mjo.eja.eja_algebra import JordanSpinEJA - sage: from mjo.eja.eja_utils import _charpoly_sage_input - - EXAMPLES:: - - sage: J = JordanSpinEJA(4,QQ) - sage: a = J._charpoly_coefficients() - sage: a[0] - X1^2 - X2^2 - X3^2 - X4^2 - sage: _charpoly_sage_input(str(a[0])) - 'X[0]**2 - X[1]**2 - X[2]**2 - X[3]**2' - - """ - import re - - exponent_out = r"\^" - exponent_in = r"**" - - digit_out = r"X([0-9]+)" - - def replace_digit(m): - # m is a match object - return "X[" + str(int(m.group(1)) - 1) + "]" - - s = re.sub(exponent_out, exponent_in, s) - return re.sub(digit_out, replace_digit, s) - - def _scale(x, alpha): r""" Scale the vector, matrix, or cartesian-product-of-those-things @@ -281,54 +245,40 @@ def gram_schmidt(v, inner_product=None): sage: len(gram_schmidt(v)) == 2 True """ + if len(v) == 0: + # cool + return v + + V = v[0].parent() + if inner_product is None: inner_product = lambda x,y: x.inner_product(y) + def norm(x): - ip = inner_product(x,x) # Don't expand the given field; the inner-product's codomain # is already correct. For example QQ(2).sqrt() returns sqrt(2) # in SR, and that will give you weird errors about symbolics # when what's really going wrong is that you're trying to # orthonormalize in QQ. - return ip.parent()(ip.sqrt()) - - v = list(v) # make a copy, don't clobber the input - - # Drop all zero vectors before we start. - v = [ v_i for v_i in v if not v_i.is_zero() ] - - if len(v) == 0: - # cool - return v - - R = v[0].base_ring() - - # Our "zero" needs to belong to the right space for sum() to work. - zero = v[0].parent().zero() + return V.base_ring()(inner_product(x,x).sqrt()) sc = lambda x,a: a*x - if hasattr(v[0], 'cartesian_factors'): + if hasattr(V, 'cartesian_factors'): # Only use the slow implementation if necessary. sc = _scale def proj(x,y): + # project y onto the span of {x} return sc(x, (inner_product(x,y)/inner_product(x,x))) - # First orthogonalize... - for i in range(1,len(v)): - # Earlier vectors can be made into zero so we have to ignore them. - v[i] -= sum( (proj(v[j],v[i]) - for j in range(i) - if not v[j].is_zero() ), - zero ) + def normalize(x): + return sc(x, ~norm(x)) - # And now drop all zero vectors again if they were "orthogonalized out." - v = [ v_i for v_i in v if not v_i.is_zero() ] + v_out = [] # make a copy, don't clobber the input - # Just normalize. If the algebra is missing the roots, we can't add - # them here because then our subalgebra would have a bigger field - # than the superalgebra. - for i in range(len(v)): - v[i] = sc(v[i], ~norm(v[i])) + for (i, v_i) in enumerate(v): + ortho_v_i = v_i - V.sum( proj(v_out[j],v_i) for j in range(i) ) + if not ortho_v_i.is_zero(): + v_out.append(normalize(ortho_v_i)) - return v + return v_out