X-Git-Url: http://gitweb.michael.orlitzky.com/?p=sage.d.git;a=blobdiff_plain;f=mjo%2Feja%2Feja_utils.py;h=a8abeff6be073b2b0aea3dd4f5af933c251666a5;hp=a4328610e5e41db689455828bd0d8988225e745b;hb=HEAD;hpb=e8599960ef47e5a5af8aca360a041d30584f6c3f diff --git a/mjo/eja/eja_utils.py b/mjo/eja/eja_utils.py index a432861..a8abeff 100644 --- a/mjo/eja/eja_utils.py +++ b/mjo/eja/eja_utils.py @@ -1,43 +1,4 @@ -from sage.functions.other import sqrt from sage.structure.element import is_Matrix -from sage.matrix.constructor import matrix -from sage.modules.free_module_element import vector - -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""" @@ -171,13 +132,6 @@ def _all2list(x): return sum( map(_all2list, xl) , []) - -def _mat2vec(m): - return vector(m.base_ring(), m.list()) - -def _vec2mat(v): - return matrix(v.base_ring(), sqrt(v.degree()), v.list()) - def gram_schmidt(v, inner_product=None): """ Perform Gram-Schmidt on the list ``v`` which are assumed to be @@ -290,55 +244,39 @@ def gram_schmidt(v, inner_product=None): sage: v = [v1,v2,v3] sage: len(gram_schmidt(v)) == 2 True - """ - 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() + V = v[0].parent() - # Our "zero" needs to belong to the right space for sum() to work. - zero = v[0].parent().zero() + if inner_product is None: + inner_product = lambda x,y: x.inner_product(y) 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): + # Don't extend the given field with the necessary + # square roots. This will probably throw weird + # errors about the symbolic ring if you e.g. try + # to use it on a set of rational vectors that isn't + # already orthonormalized. + return sc(x, ~inner_product(x,x).sqrt()) - # 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