From c9db5fa5d2f30833e4da85d32418048d34483e43 Mon Sep 17 00:00:00 2001 From: Michael Orlitzky Date: Fri, 27 Nov 2020 10:38:18 -0500 Subject: [PATCH] eja: use new constructor for BilinearFormEJA. --- mjo/eja/eja_algebra.py | 75 ++++++++++++++++++++++++------------------ 1 file changed, 43 insertions(+), 32 deletions(-) diff --git a/mjo/eja/eja_algebra.py b/mjo/eja/eja_algebra.py index 4609ca2..9eba669 100644 --- a/mjo/eja/eja_algebra.py +++ b/mjo/eja/eja_algebra.py @@ -1062,8 +1062,15 @@ class RationalBasisEuclideanJordanAlgebraNg(FiniteDimensionalEuclideanJordanAlge if orthonormalize: from mjo.eja.eja_utils import gram_schmidt + A = matrix(field, vector_basis) vector_basis = gram_schmidt(vector_basis, inner_product) W = V.span_of_basis( vector_basis ) + Q = matrix(field, vector_basis) + # A = QR <==> A.T == R.T*Q.T + # So, Q.solve_right() is equivalent to the Q.T.solve_left() + # that we want. + self._deorthonormalization_matrix = Q.solve_right(A) + if basis_is_matrices: from mjo.eja.eja_utils import _vec2mat basis = tuple( map(_vec2mat,vector_basis) ) @@ -1073,15 +1080,30 @@ class RationalBasisEuclideanJordanAlgebraNg(FiniteDimensionalEuclideanJordanAlge mult_table = [ [0 for i in range(n)] for j in range(n) ] ip_table = [ [0 for i in range(n)] for j in range(n) ] + # Note: the Jordan and inner- products are defined in terms + # of the ambient basis. It's important that their arguments + # are in ambient coordinates as well. for i in range(n): for j in range(i+1): - # do another mat2vec because the multiplication - # table is in terms of vectors - elt = _mat2vec(jordan_product(basis[i],basis[j])) + # ortho basis w.r.t. ambient coords + q_i = vector_basis[i] + q_j = vector_basis[j] + + if basis_is_matrices: + q_i = _vec2mat(q_i) + q_j = _vec2mat(q_j) + + elt = jordan_product(q_i, q_j) + ip = inner_product(q_i, q_j) + + if basis_is_matrices: + # do another mat2vec because the multiplication + # table is in terms of vectors + elt = _mat2vec(elt) + elt = W.coordinate_vector(elt) mult_table[i][j] = elt mult_table[j][i] = elt - ip = inner_product(basis[i],basis[j]) ip_table[i][j] = ip ip_table[j][i] = ip @@ -2193,7 +2215,7 @@ class HadamardEJA(RationalBasisEuclideanJordanAlgebraNg, return cls(n, field, **kwargs) -class BilinearFormEJA(RationalBasisEuclideanJordanAlgebra, +class BilinearFormEJA(RationalBasisEuclideanJordanAlgebraNg, ConcreteEuclideanJordanAlgebra): r""" The rank-2 simple EJA consisting of real vectors ``x=(x0, x_bar)`` @@ -2274,39 +2296,28 @@ class BilinearFormEJA(RationalBasisEuclideanJordanAlgebra, True """ def __init__(self, B, field=AA, **kwargs): - n = B.nrows() - if not B.is_positive_definite(): raise ValueError("bilinear form is not positive-definite") + n = B.nrows() V = VectorSpace(field, n) - mult_table = [[V.zero() for j in range(n)] for i in range(n)] - for i in range(n): - for j in range(n): - x = V.gen(i) - y = V.gen(j) - x0 = x[0] - xbar = x[1:] - y0 = y[0] - ybar = y[1:] - z0 = (B*x).inner_product(y) - zbar = y0*xbar + x0*ybar - z = V([z0] + zbar.list()) - mult_table[i][j] = z - - # Inner products are real numbers and not algebra - # elements, so once we turn the algebra element - # into a vector in inner_product(), we never go - # back. As a result -- contrary to what we do with - # self._multiplication_table -- we store the inner - # product table as a plain old matrix and not as - # an algebra operator. - ip_table = B - self._inner_product_matrix = ip_table + + def inner_product(x,y): + return (B*x).inner_product(y) + + def jordan_product(x,y): + x0 = x[0] + xbar = x[1:] + y0 = y[0] + ybar = y[1:] + z0 = inner_product(x,y) + zbar = y0*xbar + x0*ybar + return V([z0] + zbar.list()) super(BilinearFormEJA, self).__init__(field, - mult_table, - check_axioms=False, + V.basis(), + jordan_product, + inner_product, **kwargs) # The rank of this algebra is two, unless we're in a -- 2.44.2