eja: get the quaternions working.

[sage.d.git] / mjo / eja / euclidean_jordan_algebra.py

diff --git a/mjo/eja/euclidean_jordan_algebra.py b/mjo/eja/euclidean_jordan_algebra.py
index ca344edc709352dea2d1940a7646e058066d4b8d..b8c7a912b1ded8149f7aefdf14b90b99d1056b2b 100644 (file)
--- a/mjo/eja/euclidean_jordan_algebra.py
+++ b/mjo/eja/euclidean_jordan_algebra.py
@@ -94,6 +94,20 @@ class FiniteDimensionalEuclideanJordanAlgebra(FiniteDimensionalAlgebra):
         The inner product associated with this Euclidean Jordan algebra.
 
         Will default to the trace inner product if nothing else.
         The inner product associated with this Euclidean Jordan algebra.
 
         Will default to the trace inner product if nothing else.

+
+        EXAMPLES:
+
+        The inner product must satisfy its axiom for this algebra to truly
+        be a Euclidean Jordan Algebra::
+
+            sage: set_random_seed()
+            sage: J = random_eja()
+            sage: x = J.random_element()
+            sage: y = J.random_element()
+            sage: z = J.random_element()
+            sage: (x*y).inner_product(z) == y.inner_product(x*z)
+            True
+

         """
         if (not x in self) or (not y in self):
             raise TypeError("arguments must live in this algebra")
         """
         if (not x in self) or (not y in self):
             raise TypeError("arguments must live in this algebra")


@@ -160,6 +174,51 @@ class FiniteDimensionalEuclideanJordanAlgebra(FiniteDimensionalAlgebra):
         An element of a Euclidean Jordan algebra.
         """
 
         An element of a Euclidean Jordan algebra.
         """
 

+        def __init__(self, A, elt=None):
+            """
+            EXAMPLES:
+
+            The identity in `S^n` is converted to the identity in the EJA::
+
+                sage: J = RealSymmetricSimpleEJA(3)
+                sage: I = identity_matrix(QQ,3)
+                sage: J(I) == J.one()
+                True
+
+            This skew-symmetric matrix can't be represented in the EJA::
+
+                sage: J = RealSymmetricSimpleEJA(3)
+                sage: A = matrix(QQ,3, lambda i,j: i-j)
+                sage: J(A)
+                Traceback (most recent call last):
+                ...
+                ArithmeticError: vector is not in free module
+
+            """
+            # Goal: if we're given a matrix, and if it lives in our
+            # parent algebra's "natural ambient space," convert it
+            # into an algebra element.
+            #
+            # The catch is, we make a recursive call after converting
+            # the given matrix into a vector that lives in the algebra.
+            # This we need to try the parent class initializer first,
+            # to avoid recursing forever if we're given something that
+            # already fits into the algebra, but also happens to live
+            # in the parent's "natural ambient space" (this happens with
+            # vectors in R^n).
+            try:
+                FiniteDimensionalAlgebraElement.__init__(self, A, elt)
+            except ValueError:
+                natural_basis = A.natural_basis()
+                if elt in natural_basis[0].matrix_space():
+                    # Thanks for nothing! Matrix spaces aren't vector
+                    # spaces in Sage, so we have to figure out its
+                    # natural-basis coordinates ourselves.
+                    V = VectorSpace(elt.base_ring(), elt.nrows()**2)
+                    W = V.span( _mat2vec(s) for s in natural_basis )
+                    coords =  W.coordinates(_mat2vec(elt))
+                    FiniteDimensionalAlgebraElement.__init__(self, A, coords)
+

         def __pow__(self, n):
             """
             Return ``self`` raised to the power ``n``.
         def __pow__(self, n):
             """
             Return ``self`` raised to the power ``n``.


@@ -244,6 +303,30 @@ class FiniteDimensionalEuclideanJordanAlgebra(FiniteDimensionalAlgebra):
                 sage: J(x).inner_product(J(y))
                 32
 
                 sage: J(x).inner_product(J(y))
                 32
 

+            The inner product on `S^n` is `<X,Y> = trace(X*Y)`, where
+            multiplication is the usual matrix multiplication in `S^n`,
+            so the inner product of the identity matrix with itself
+            should be the `n`::
+
+                sage: J = RealSymmetricSimpleEJA(3)
+                sage: J.one().inner_product(J.one())
+                3
+
+            Likewise, the inner product on `C^n` is `<X,Y> =
+            Re(trace(X*Y))`, where we must necessarily take the real
+            part because the product of Hermitian matrices may not be
+            Hermitian::
+
+                sage: J = ComplexHermitianSimpleEJA(3)
+                sage: J.one().inner_product(J.one())
+                3
+
+            Ditto for the quaternions::
+
+                sage: J = QuaternionHermitianSimpleEJA(3)
+                sage: J.one().inner_product(J.one())
+                3
+

             TESTS:
 
             Ensure that we can always compute an inner product, and that
             TESTS:
 
             Ensure that we can always compute an inner product, and that


@@ -597,6 +680,25 @@ class FiniteDimensionalEuclideanJordanAlgebra(FiniteDimensionalAlgebra):
                 [0 0 0 0 1 0]
                 [0 0 0 0 0 1]
 
                 [0 0 0 0 1 0]
                 [0 0 0 0 0 1]
 

+            ::
+
+                sage: J = QuaternionHermitianSimpleEJA(3)
+                sage: J.one()
+                e0 + e9 + e14
+                sage: J.one().natural_representation()
+                [1 0 0 0 0 0 0 0 0 0 0 0]
+                [0 1 0 0 0 0 0 0 0 0 0 0]
+                [0 0 1 0 0 0 0 0 0 0 0 0]
+                [0 0 0 1 0 0 0 0 0 0 0 0]
+                [0 0 0 0 1 0 0 0 0 0 0 0]
+                [0 0 0 0 0 1 0 0 0 0 0 0]
+                [0 0 0 0 0 0 1 0 0 0 0 0]
+                [0 0 0 0 0 0 0 1 0 0 0 0]
+                [0 0 0 0 0 0 0 0 1 0 0 0]
+                [0 0 0 0 0 0 0 0 0 1 0 0]
+                [0 0 0 0 0 0 0 0 0 0 1 0]
+                [0 0 0 0 0 0 0 0 0 0 0 1]
+

             """
             B = self.parent().natural_basis()
             W = B[0].matrix_space()
             """
             B = self.parent().natural_basis()
             W = B[0].matrix_space()


@@ -941,13 +1043,10 @@ def eja_rn(dimension, field=QQ):
     Qs = [ matrix(field, dimension, dimension, lambda k,j: 1*(k == j == i))
            for i in xrange(dimension) ]
 
     Qs = [ matrix(field, dimension, dimension, lambda k,j: 1*(k == j == i))
            for i in xrange(dimension) ]
 

-    # The usual inner product on R^n.
-    ip = lambda x, y: x.vector().inner_product(y.vector())
-

     return FiniteDimensionalEuclideanJordanAlgebra(field,
                                                    Qs,
                                                    rank=dimension,
     return FiniteDimensionalEuclideanJordanAlgebra(field,
                                                    Qs,
                                                    rank=dimension,

-                                                   inner_product=ip)
+                                                   inner_product=_usual_ip)

 
 
 
 
 
 


@@ -968,6 +1067,12 @@ def random_eja():
       * The ``n``-by-``n`` rational symmetric matrices with the symmetric
         product.
 
       * The ``n``-by-``n`` rational symmetric matrices with the symmetric
         product.
 

+      * The ``n``-by-``n`` complex-rational Hermitian matrices embedded
+        in the space of ``2n``-by-``2n`` real symmetric matrices.
+
+      * The ``n``-by-``n`` quaternion-rational Hermitian matrices embedded
+        in the space of ``4n``-by-``4n`` real symmetric matrices.
+

     Later this might be extended to return Cartesian products of the
     EJAs above.
 
     Later this might be extended to return Cartesian products of the
     EJAs above.
 


@@ -981,7 +1086,8 @@ def random_eja():
     constructor = choice([eja_rn,
                           JordanSpinSimpleEJA,
                           RealSymmetricSimpleEJA,
     constructor = choice([eja_rn,
                           JordanSpinSimpleEJA,
                           RealSymmetricSimpleEJA,

-                          ComplexHermitianSimpleEJA])
+                          ComplexHermitianSimpleEJA,
+                          QuaternionHermitianSimpleEJA])

     return constructor(n, field=QQ)
 
 
     return constructor(n, field=QQ)
 
 


@@ -1042,6 +1148,54 @@ def _complex_hermitian_basis(n, field=QQ):
     return tuple(S)
 
 
     return tuple(S)
 
 

+def _quaternion_hermitian_basis(n, field=QQ):
+    """
+    Returns a basis for the space of quaternion Hermitian n-by-n matrices.
+
+    TESTS::
+
+        sage: set_random_seed()
+        sage: n = ZZ.random_element(1,5)
+        sage: all( M.is_symmetric() for M in _quaternion_hermitian_basis(n) )
+        True
+
+    """
+    Q = QuaternionAlgebra(QQ,-1,-1)
+    I,J,K = Q.gens()
+
+    # This is like the symmetric case, but we need to be careful:
+    #
+    #   * We want conjugate-symmetry, not just symmetry.
+    #   * The diagonal will (as a result) be real.
+    #
+    S = []
+    for i in xrange(n):
+        for j in xrange(i+1):
+            Eij = matrix(Q, n, lambda k,l: k==i and l==j)
+            if i == j:
+                Sij = _embed_quaternion_matrix(Eij)
+                S.append(Sij)
+            else:
+                # Beware, orthogonal but not normalized! The second,
+                # third, and fourth ones have a minus because they're
+                # conjugated.
+                Sij_real = _embed_quaternion_matrix(Eij + Eij.transpose())
+                S.append(Sij_real)
+                Sij_I = _embed_quaternion_matrix(I*Eij - I*Eij.transpose())
+                S.append(Sij_I)
+                Sij_J = _embed_quaternion_matrix(J*Eij - J*Eij.transpose())
+                S.append(Sij_J)
+                Sij_K = _embed_quaternion_matrix(K*Eij - K*Eij.transpose())
+                S.append(Sij_K)
+    return tuple(S)
+
+
+def _mat2vec(m):
+        return vector(m.base_ring(), m.list())
+
+def _vec2mat(v):
+        return matrix(v.base_ring(), sqrt(v.degree()), v.list())
+

 def _multiplication_table_from_matrix_basis(basis):
     """
     At least three of the five simple Euclidean Jordan algebras have the
 def _multiplication_table_from_matrix_basis(basis):
     """
     At least three of the five simple Euclidean Jordan algebras have the


@@ -1062,19 +1216,13 @@ def _multiplication_table_from_matrix_basis(basis):
     field = basis[0].base_ring()
     dimension = basis[0].nrows()
 
     field = basis[0].base_ring()
     dimension = basis[0].nrows()
 

-    def mat2vec(m):
-        return vector(field, m.list())
-
-    def vec2mat(v):
-        return matrix(field, dimension, v.list())
-

     V = VectorSpace(field, dimension**2)
     V = VectorSpace(field, dimension**2)

-    W = V.span( mat2vec(s) for s in basis )
+    W = V.span( _mat2vec(s) for s in basis )

 
     # Taking the span above reorders our basis (thanks, jerk!) so we
     # need to put our "matrix basis" in the same order as the
     # (reordered) vector basis.
 
     # Taking the span above reorders our basis (thanks, jerk!) so we
     # need to put our "matrix basis" in the same order as the
     # (reordered) vector basis.

-    S = tuple( vec2mat(b) for b in W.basis() )
+    S = tuple( _vec2mat(b) for b in W.basis() )

 
     Qs = []
     for s in S:
 
     Qs = []
     for s in S:


@@ -1087,7 +1235,7 @@ def _multiplication_table_from_matrix_basis(basis):
         # why we're computing rows here and not columns.
         Q_rows = []
         for t in S:
         # why we're computing rows here and not columns.
         Q_rows = []
         for t in S:

-            this_row = mat2vec((s*t + t*s)/2)
+            this_row = _mat2vec((s*t + t*s)/2)

             Q_rows.append(W.coordinates(this_row))
         Q = matrix(field, W.dimension(), Q_rows)
         Qs.append(Q)
             Q_rows.append(W.coordinates(this_row))
         Q = matrix(field, W.dimension(), Q_rows)
         Qs.append(Q)


@@ -1108,13 +1256,13 @@ def _embed_complex_matrix(M):
         sage: x2 = F(1 + 2*i)
         sage: x3 = F(-i)
         sage: x4 = F(6)
         sage: x2 = F(1 + 2*i)
         sage: x3 = F(-i)
         sage: x4 = F(6)

-        sage: M = matrix(F,2,[x1,x2,x3,x4])
+        sage: M = matrix(F,2,[[x1,x2],[x3,x4]])

         sage: _embed_complex_matrix(M)
         sage: _embed_complex_matrix(M)

-        [ 4  2| 1 -2]
-        [-2  4| 2  1]
+        [ 4 -2| 1  2]
+        [ 2  4|-2  1]

         [-----+-----]
         [-----+-----]

-        [ 0  1| 6  0]
-        [-1  0| 0  6]
+        [ 0 -1| 6  0]
+        [ 1  0| 0  6]

 
     """
     n = M.nrows()
 
     """
     n = M.nrows()


@@ -1125,7 +1273,7 @@ def _embed_complex_matrix(M):
     for z in M.list():
         a = z.real()
         b = z.imag()
     for z in M.list():
         a = z.real()
         b = z.imag()

-        blocks.append(matrix(field, 2, [[a,-b],[b,a]]))
+        blocks.append(matrix(field, 2, [[a,b],[-b,a]]))

 
     # We can drop the imaginaries here.
     return block_matrix(field.base_ring(), n, blocks)
 
     # We can drop the imaginaries here.
     return block_matrix(field.base_ring(), n, blocks)


@@ -1142,8 +1290,17 @@ def _unembed_complex_matrix(M):
         ....:                 [ 9,  10, 11, 12],
         ....:                 [-10, 9, -12, 11] ])
         sage: _unembed_complex_matrix(A)
         ....:                 [ 9,  10, 11, 12],
         ....:                 [-10, 9, -12, 11] ])
         sage: _unembed_complex_matrix(A)

-        [  -2*i + 1   -4*i + 3]
-        [ -10*i + 9 -12*i + 11]
+        [  2*i + 1   4*i + 3]
+        [ 10*i + 9 12*i + 11]
+
+    TESTS::
+
+        sage: set_random_seed()
+        sage: F = QuadraticField(-1, 'i')
+        sage: M = random_matrix(F, 3)
+        sage: _unembed_complex_matrix(_embed_complex_matrix(M)) == M
+        True
+

     """
     n = ZZ(M.nrows())
     if M.ncols() != n:
     """
     n = ZZ(M.nrows())
     if M.ncols() != n:


@@ -1161,15 +1318,122 @@ def _unembed_complex_matrix(M):
         for j in xrange(n/2):
             submat = M[2*k:2*k+2,2*j:2*j+2]
             if submat[0,0] != submat[1,1]:
         for j in xrange(n/2):
             submat = M[2*k:2*k+2,2*j:2*j+2]
             if submat[0,0] != submat[1,1]:

-                raise ArgumentError('bad real submatrix')
+                raise ValueError('bad on-diagonal submatrix')

             if submat[0,1] != -submat[1,0]:
             if submat[0,1] != -submat[1,0]:

-                raise ArgumentError('bad imag submatrix')
-            z = submat[0,0] + submat[1,0]*i
+                raise ValueError('bad off-diagonal submatrix')
+            z = submat[0,0] + submat[0,1]*i

             elements.append(z)
 
     return matrix(F, n/2, elements)
 
 
             elements.append(z)
 
     return matrix(F, n/2, elements)
 
 

+def _embed_quaternion_matrix(M):
+    """
+    Embed the n-by-n quaternion matrix ``M`` into the space of real
+    matrices of size 4n-by-4n by first sending each quaternion entry
+    `z = a + bi + cj + dk` to the block-complex matrix
+    ``[[a + bi, c+di],[-c + di, a-bi]]`, and then embedding those into
+    a real matrix.
+
+    EXAMPLES::
+
+        sage: Q = QuaternionAlgebra(QQ,-1,-1)
+        sage: i,j,k = Q.gens()
+        sage: x = 1 + 2*i + 3*j + 4*k
+        sage: M = matrix(Q, 1, [[x]])
+        sage: _embed_quaternion_matrix(M)
+        [ 1  2  3  4]
+        [-2  1 -4  3]
+        [-3  4  1 -2]
+        [-4 -3  2  1]
+
+    """
+    quaternions = M.base_ring()
+    n = M.nrows()
+    if M.ncols() != n:
+        raise ValueError("the matrix 'M' must be square")
+
+    F = QuadraticField(-1, 'i')
+    i = F.gen()
+
+    blocks = []
+    for z in M.list():
+        t = z.coefficient_tuple()
+        a = t[0]
+        b = t[1]
+        c = t[2]
+        d = t[3]
+        cplx_matrix = matrix(F, 2, [[ a + b*i, c + d*i],
+                                    [-c + d*i, a - b*i]])
+        blocks.append(_embed_complex_matrix(cplx_matrix))
+
+    # We should have real entries by now, so use the realest field
+    # we've got for the return value.
+    return block_matrix(quaternions.base_ring(), n, blocks)
+
+
+def _unembed_quaternion_matrix(M):
+    """
+    The inverse of _embed_quaternion_matrix().
+
+    EXAMPLES::
+
+        sage: M = matrix(QQ, [[ 1,  2,  3,  4],
+        ....:                 [-2,  1, -4,  3],
+        ....:                 [-3,  4,  1, -2],
+        ....:                 [-4, -3,  2,  1]])
+        sage: _unembed_quaternion_matrix(M)
+        [1 + 2*i + 3*j + 4*k]
+
+    TESTS::
+
+        sage: set_random_seed()
+        sage: Q = QuaternionAlgebra(QQ, -1, -1)
+        sage: M = random_matrix(Q, 3)
+        sage: _unembed_quaternion_matrix(_embed_quaternion_matrix(M)) == M
+        True
+
+    """
+    n = ZZ(M.nrows())
+    if M.ncols() != n:
+        raise ValueError("the matrix 'M' must be square")
+    if not n.mod(4).is_zero():
+        raise ValueError("the matrix 'M' must be a complex embedding")
+
+    Q = QuaternionAlgebra(QQ,-1,-1)
+    i,j,k = Q.gens()
+
+    # Go top-left to bottom-right (reading order), converting every
+    # 4-by-4 block we see to a 2-by-2 complex block, to a 1-by-1
+    # quaternion block.
+    elements = []
+    for l in xrange(n/4):
+        for m in xrange(n/4):
+            submat = _unembed_complex_matrix(M[4*l:4*l+4,4*m:4*m+4])
+            if submat[0,0] != submat[1,1].conjugate():
+                raise ValueError('bad on-diagonal submatrix')
+            if submat[0,1] != -submat[1,0].conjugate():
+                raise ValueError('bad off-diagonal submatrix')
+            z  = submat[0,0].real() + submat[0,0].imag()*i
+            z += submat[0,1].real()*j + submat[0,1].imag()*k
+            elements.append(z)
+
+    return matrix(Q, n/4, elements)
+
+
+# The usual inner product on R^n.
+def _usual_ip(x,y):
+    return x.vector().inner_product(y.vector())
+
+# The inner product used for the real symmetric simple EJA.
+# We keep it as a separate function because e.g. the complex
+# algebra uses the same inner product, except divided by 2.
+def _matrix_ip(X,Y):
+    X_mat = X.natural_representation()
+    Y_mat = Y.natural_representation()
+    return (X_mat*Y_mat).trace()
+
+

 def RealSymmetricSimpleEJA(n, field=QQ):
     """
     The rank-n simple EJA consisting of real symmetric n-by-n
 def RealSymmetricSimpleEJA(n, field=QQ):
     """
     The rank-n simple EJA consisting of real symmetric n-by-n


@@ -1197,6 +1461,22 @@ def RealSymmetricSimpleEJA(n, field=QQ):
         sage: J.degree() == (n^2 + n)/2
         True
 
         sage: J.degree() == (n^2 + n)/2
         True
 

+    The Jordan multiplication is what we think it is::
+
+        sage: set_random_seed()
+        sage: n = ZZ.random_element(1,5)
+        sage: J = RealSymmetricSimpleEJA(n)
+        sage: x = J.random_element()
+        sage: y = J.random_element()
+        sage: actual = (x*y).natural_representation()
+        sage: X = x.natural_representation()
+        sage: Y = y.natural_representation()
+        sage: expected = (X*Y + Y*X)/2
+        sage: actual == expected
+        True
+        sage: J(expected) == x*y
+        True
+

     """
     S = _real_symmetric_basis(n, field=field)
     (Qs, T) = _multiplication_table_from_matrix_basis(S)
     """
     S = _real_symmetric_basis(n, field=field)
     (Qs, T) = _multiplication_table_from_matrix_basis(S)


@@ -1204,7 +1484,8 @@ def RealSymmetricSimpleEJA(n, field=QQ):
     return FiniteDimensionalEuclideanJordanAlgebra(field,
                                                    Qs,
                                                    rank=n,
     return FiniteDimensionalEuclideanJordanAlgebra(field,
                                                    Qs,
                                                    rank=n,

-                                                   natural_basis=T)
+                                                   natural_basis=T,
+                                                   inner_product=_matrix_ip)

 
 
 def ComplexHermitianSimpleEJA(n, field=QQ):
 
 
 def ComplexHermitianSimpleEJA(n, field=QQ):


@@ -1224,23 +1505,96 @@ def ComplexHermitianSimpleEJA(n, field=QQ):
         sage: J.degree() == n^2
         True
 
         sage: J.degree() == n^2
         True
 

+    The Jordan multiplication is what we think it is::
+
+        sage: set_random_seed()
+        sage: n = ZZ.random_element(1,5)
+        sage: J = ComplexHermitianSimpleEJA(n)
+        sage: x = J.random_element()
+        sage: y = J.random_element()
+        sage: actual = (x*y).natural_representation()
+        sage: X = x.natural_representation()
+        sage: Y = y.natural_representation()
+        sage: expected = (X*Y + Y*X)/2
+        sage: actual == expected
+        True
+        sage: J(expected) == x*y
+        True
+

     """
     S = _complex_hermitian_basis(n)
     (Qs, T) = _multiplication_table_from_matrix_basis(S)
     """
     S = _complex_hermitian_basis(n)
     (Qs, T) = _multiplication_table_from_matrix_basis(S)

+
+    # Since a+bi on the diagonal is represented as
+    #
+    #   a + bi  = [  a  b  ]
+    #             [ -b  a  ],
+    #
+    # we'll double-count the "a" entries if we take the trace of
+    # the embedding.
+    ip = lambda X,Y: _matrix_ip(X,Y)/2
+

     return FiniteDimensionalEuclideanJordanAlgebra(field,
                                                    Qs,
                                                    rank=n,
     return FiniteDimensionalEuclideanJordanAlgebra(field,
                                                    Qs,
                                                    rank=n,

-                                                   natural_basis=T)
+                                                   natural_basis=T,
+                                                   inner_product=ip)

 
 
 
 

-def QuaternionHermitianSimpleEJA(n):
+def QuaternionHermitianSimpleEJA(n, field=QQ):

     """
     The rank-n simple EJA consisting of self-adjoint n-by-n quaternion
     matrices, the usual symmetric Jordan product, and the
     real-part-of-trace inner product. It has dimension `2n^2 - n` over
     the reals.
     """
     The rank-n simple EJA consisting of self-adjoint n-by-n quaternion
     matrices, the usual symmetric Jordan product, and the
     real-part-of-trace inner product. It has dimension `2n^2 - n` over
     the reals.

+
+    TESTS:
+
+    The degree of this algebra is `n^2`::
+
+        sage: set_random_seed()
+        sage: n = ZZ.random_element(1,5)
+        sage: J = QuaternionHermitianSimpleEJA(n)
+        sage: J.degree() == 2*(n^2) - n
+        True
+
+    The Jordan multiplication is what we think it is::
+
+        sage: set_random_seed()
+        sage: n = ZZ.random_element(1,5)
+        sage: J = QuaternionHermitianSimpleEJA(n)
+        sage: x = J.random_element()
+        sage: y = J.random_element()
+        sage: actual = (x*y).natural_representation()
+        sage: X = x.natural_representation()
+        sage: Y = y.natural_representation()
+        sage: expected = (X*Y + Y*X)/2
+        sage: actual == expected
+        True
+        sage: J(expected) == x*y
+        True
+

     """
     """

-    pass
+    S = _quaternion_hermitian_basis(n)
+    (Qs, T) = _multiplication_table_from_matrix_basis(S)
+
+    # Since a+bi+cj+dk on the diagonal is represented as
+    #
+    #   a + bi +cj + dk = [  a  b  c  d]
+    #                     [ -b  a -d  c]
+    #                     [ -c  d  a -b]
+    #                     [ -d -c  b  a],
+    #
+    # we'll quadruple-count the "a" entries if we take the trace of
+    # the embedding.
+    ip = lambda X,Y: _matrix_ip(X,Y)/4
+
+    return FiniteDimensionalEuclideanJordanAlgebra(field,
+                                                   Qs,
+                                                   rank=n,
+                                                   natural_basis=T,
+                                                   inner_product=ip)
+

 
 def OctonionHermitianSimpleEJA(n):
     """
 
 def OctonionHermitianSimpleEJA(n):
     """


@@ -1277,6 +1631,13 @@ def JordanSpinSimpleEJA(n, field=QQ):
         sage: e2*e3
         0
 
         sage: e2*e3
         0
 

+    In one dimension, this is the reals under multiplication::
+
+        sage: J1 = JordanSpinSimpleEJA(1)
+        sage: J2 = eja_rn(1)
+        sage: J1 == J2
+        True
+

     """
     Qs = []
     id_matrix = identity_matrix(field, n)
     """
     Qs = []
     id_matrix = identity_matrix(field, n)


@@ -1291,13 +1652,10 @@ def JordanSpinSimpleEJA(n, field=QQ):
         Qi[0,0] = Qi[0,0] * ~field(2)
         Qs.append(Qi)
 
         Qi[0,0] = Qi[0,0] * ~field(2)
         Qs.append(Qi)
 

-    # The usual inner product on R^n.
-    ip = lambda x, y: x.vector().inner_product(y.vector())
-

     # The rank of the spin factor algebra is two, UNLESS we're in a
     # one-dimensional ambient space (the rank is bounded by the
     # ambient dimension).
     return FiniteDimensionalEuclideanJordanAlgebra(field,
                                                    Qs,
                                                    rank=min(n,2),
     # The rank of the spin factor algebra is two, UNLESS we're in a
     # one-dimensional ambient space (the rank is bounded by the
     # ambient dimension).
     return FiniteDimensionalEuclideanJordanAlgebra(field,
                                                    Qs,
                                                    rank=min(n,2),

-                                                   inner_product=ip)
+                                                   inner_product=_usual_ip)