n = len(basis)
vector_basis = basis
- from sage.matrix.matrix import is_Matrix
+ from sage.structure.element import is_Matrix
basis_is_matrices = False
degree = 0
V = VectorSpace(field, degree)
- self._deorthonormalization_matrix = matrix.identity(field,n)
+ # Compute this from "Q" (obtained from Gram-Schmidt) below as
+ # R = Q.solve_right(A), where the rows of "Q" are the
+ # orthonormalized vector_basis and and the rows of "A" are the
+ # original vector_basis.
+ self._deorthonormalization_matrix = None
+
if orthonormalize:
+ from mjo.eja.eja_utils import gram_schmidt
A = matrix(field, vector_basis)
- # uh oh, this is only the "usual" inner product
- Q,R = A.gram_schmidt(orthonormal=True)
- self._deorthonormalization_matrix = R.inverse().transpose()
- vector_basis = Q.rows()
+ 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) )
+ W = V.span_of_basis( vector_basis )
+
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
if basis_is_matrices:
for m in basis:
m.set_immutable()
+ else:
+ basis = tuple( x.column() for x in basis )
super().__init__(field,
mult_table,
prefix,
category,
- basis,
+ basis, # matrix basis
check_field,
check_axioms)
return cls(n, field, **kwargs)
-class HadamardEJA(RationalBasisEuclideanJordanAlgebra,
+class HadamardEJA(RationalBasisEuclideanJordanAlgebraNg,
ConcreteEuclideanJordanAlgebra):
"""
Return the Euclidean Jordan Algebra corresponding to the set
"""
def __init__(self, n, field=AA, **kwargs):
V = VectorSpace(field, n)
- mult_table = [ [ V.gen(i)*(i == j) for j in range(n) ]
- for i in range(n) ]
+ basis = V.basis()
- # 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 = matrix.identity(field,n)
- self._inner_product_matrix = ip_table
+ def jordan_product(x,y):
+ return V([ xi*yi for (xi,yi) in zip(x,y) ])
+ def inner_product(x,y):
+ return x.inner_product(y)
super(HadamardEJA, self).__init__(field,
- mult_table,
- check_axioms=False,
+ basis,
+ jordan_product,
+ inner_product,
**kwargs)
self.rank.set_cache(n)
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)``
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