]> gitweb.michael.orlitzky.com - sage.d.git/blobdiff - mjo/eja/eja_algebra.py
eja: fix a copy/paste error in the docs.
[sage.d.git] / mjo / eja / eja_algebra.py
index 2bad32c2f500193e4126b7c5e209c0acb3116ede..afd0f55c103c885315e99bc90afe647262aadfbd 100644 (file)
@@ -230,7 +230,6 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
     We should compute that an element subalgebra is associative even
     if we circumvent the element method::
 
-        sage: set_random_seed()
         sage: J = random_eja(field=QQ,orthonormalize=False)
         sage: x = J.random_element()
         sage: A = x.subalgebra_generated_by(orthonormalize=False)
@@ -367,7 +366,7 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
 
         if orthonormalize:
             # Now "self._matrix_span" is the vector space of our
-            # algebra coordinates. The variables "X1", "X2",...  refer
+            # algebra coordinates. The variables "X0", "X1",...  refer
             # to the entries of vectors in self._matrix_span. Thus to
             # convert back and forth between the orthonormal
             # coordinates and the given ones, we need to stick the
@@ -432,7 +431,6 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
 
         TESTS::
 
-            sage: set_random_seed()
             sage: J = random_eja()
             sage: J(1)
             Traceback (most recent call last):
@@ -457,7 +455,6 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
 
         TESTS::
 
-            sage: set_random_seed()
             sage: J = random_eja()
             sage: n = J.dimension()
             sage: bi = J.zero()
@@ -499,7 +496,6 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
         Our inner product is "associative," which means the following for
         a symmetric bilinear form::
 
-            sage: set_random_seed()
             sage: J = random_eja()
             sage: x,y,z = J.random_elements(3)
             sage: (x*y).inner_product(z) == y.inner_product(x*z)
@@ -510,7 +506,6 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
         Ensure that this is the usual inner product for the algebras
         over `R^n`::
 
-            sage: set_random_seed()
             sage: J = HadamardEJA.random_instance()
             sage: x,y = J.random_elements(2)
             sage: actual = x.inner_product(y)
@@ -523,7 +518,6 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
         one). This is in Faraut and Koranyi, and also my "On the
         symmetry..." paper::
 
-            sage: set_random_seed()
             sage: J = BilinearFormEJA.random_instance()
             sage: n = J.dimension()
             sage: x = J.random_element()
@@ -636,7 +630,6 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
         The values we've presupplied to the constructors agree with
         the computation::
 
-            sage: set_random_seed()
             sage: J = random_eja()
             sage: J.is_associative() == J._jordan_product_is_associative()
             True
@@ -758,7 +751,6 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
         Ensure that we can convert any element back and forth
         faithfully between its matrix and algebra representations::
 
-            sage: set_random_seed()
             sage: J = random_eja()
             sage: x = J.random_element()
             sage: J(x.to_matrix()) == x
@@ -871,7 +863,7 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
 
             sage: J = JordanSpinEJA(3)
             sage: p = J.characteristic_polynomial_of(); p
-            X1^2 - X2^2 - X3^2 + (-2*t)*X1 + t^2
+            X0^2 - X1^2 - X2^2 + (-2*t)*X0 + t^2
             sage: xvec = J.one().to_vector()
             sage: p(*xvec)
             t^2 - 2*t + 1
@@ -920,13 +912,13 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
 
             sage: J = HadamardEJA(2)
             sage: J.coordinate_polynomial_ring()
-            Multivariate Polynomial Ring in X1, X2...
+            Multivariate Polynomial Ring in X0, X1...
             sage: J = RealSymmetricEJA(3,field=QQ,orthonormalize=False)
             sage: J.coordinate_polynomial_ring()
-            Multivariate Polynomial Ring in X1, X2, X3, X4, X5, X6...
+            Multivariate Polynomial Ring in X0, X1, X2, X3, X4, X5...
 
         """
-        var_names = tuple( "X%d" % z for z in range(1, self.dimension()+1) )
+        var_names = tuple( "X%d" % z for z in range(self.dimension()) )
         return PolynomialRing(self.base_ring(), var_names)
 
     def inner_product(self, x, y):
@@ -948,7 +940,6 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
         Our inner product is "associative," which means the following for
         a symmetric bilinear form::
 
-            sage: set_random_seed()
             sage: J = random_eja()
             sage: x,y,z = J.random_elements(3)
             sage: (x*y).inner_product(z) == y.inner_product(x*z)
@@ -959,7 +950,6 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
         Ensure that this is the usual inner product for the algebras
         over `R^n`::
 
-            sage: set_random_seed()
             sage: J = HadamardEJA.random_instance()
             sage: x,y = J.random_elements(2)
             sage: actual = x.inner_product(y)
@@ -972,7 +962,6 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
         one). This is in Faraut and Koranyi, and also my "On the
         symmetry..." paper::
 
-            sage: set_random_seed()
             sage: J = BilinearFormEJA.random_instance()
             sage: n = J.dimension()
             sage: x = J.random_element()
@@ -1200,7 +1189,6 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
         The identity element acts like the identity, regardless of
         whether or not we orthonormalize::
 
-            sage: set_random_seed()
             sage: J = random_eja()
             sage: x = J.random_element()
             sage: J.one()*x == x and x*J.one() == x
@@ -1212,7 +1200,6 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
 
         ::
 
-            sage: set_random_seed()
             sage: J = random_eja(field=QQ, orthonormalize=False)
             sage: x = J.random_element()
             sage: J.one()*x == x and x*J.one() == x
@@ -1226,7 +1213,6 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
         regardless of the base field and whether or not we
         orthonormalize::
 
-            sage: set_random_seed()
             sage: J = random_eja()
             sage: actual = J.one().operator().matrix()
             sage: expected = matrix.identity(J.base_ring(), J.dimension())
@@ -1241,7 +1227,6 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
 
         ::
 
-            sage: set_random_seed()
             sage: J = random_eja(field=QQ, orthonormalize=False)
             sage: actual = J.one().operator().matrix()
             sage: expected = matrix.identity(J.base_ring(), J.dimension())
@@ -1257,7 +1242,6 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
         Ensure that the cached unit element (often precomputed by
         hand) agrees with the computed one::
 
-            sage: set_random_seed()
             sage: J = random_eja()
             sage: cached = J.one()
             sage: J.one.clear_cache()
@@ -1266,7 +1250,6 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
 
         ::
 
-            sage: set_random_seed()
             sage: J = random_eja(field=QQ, orthonormalize=False)
             sage: cached = J.one()
             sage: J.one.clear_cache()
@@ -1379,7 +1362,6 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
         Every algebra decomposes trivially with respect to its identity
         element::
 
-            sage: set_random_seed()
             sage: J = random_eja()
             sage: J0,J5,J1 = J.peirce_decomposition(J.one())
             sage: J0.dimension() == 0 and J5.dimension() == 0
@@ -1392,7 +1374,6 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
         elements in the two subalgebras are the projections onto their
         respective subspaces of the superalgebra's identity element::
 
-            sage: set_random_seed()
             sage: J = random_eja()
             sage: x = J.random_element()
             sage: if not J.is_trivial():
@@ -1518,6 +1499,64 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
                       for idx in range(count) )
 
 
+    def operator_polynomial_matrix(self):
+        r"""
+        Return the matrix of polynomials (over this algebra's
+        :meth:`coordinate_polynomial_ring`) that, when evaluated at
+        the basis coordinates of an element `x`, produces the basis
+        representation of `L_{x}`.
+
+        SETUP::
+
+            sage: from mjo.eja.eja_algebra import (HadamardEJA,
+            ....:                                  JordanSpinEJA)
+
+        EXAMPLES::
+
+            sage: J = HadamardEJA(4)
+            sage: L_x = J.operator_polynomial_matrix()
+            sage: L_x
+            [X0  0  0  0]
+            [ 0 X1  0  0]
+            [ 0  0 X2  0]
+            [ 0  0  0 X3]
+            sage: x = J.one()
+            sage: d = zip(J.coordinate_polynomial_ring().gens(), x.to_vector())
+            sage: L_x.subs(dict(d))
+            [1 0 0 0]
+            [0 1 0 0]
+            [0 0 1 0]
+            [0 0 0 1]
+
+        ::
+
+            sage: J = JordanSpinEJA(4)
+            sage: L_x = J.operator_polynomial_matrix()
+            sage: L_x
+            [X0 X1 X2 X3]
+            [X1 X0  0  0]
+            [X2  0 X0  0]
+            [X3  0  0 X0]
+            sage: x = J.one()
+            sage: d = zip(J.coordinate_polynomial_ring().gens(), x.to_vector())
+            sage: L_x.subs(dict(d))
+            [1 0 0 0]
+            [0 1 0 0]
+            [0 0 1 0]
+            [0 0 0 1]
+
+        """
+        R = self.coordinate_polynomial_ring()
+
+        def L_x_i_j(i,j):
+            # From a result in my book, these are the entries of the
+            # basis representation of L_x.
+            return sum( v*self.monomial(k).operator().matrix()[i,j]
+                        for (k,v) in enumerate(R.gens()) )
+
+        n = self.dimension()
+        return matrix(R, n, n, L_x_i_j)
+
     @cached_method
     def _charpoly_coefficients(self):
         r"""
@@ -1533,7 +1572,6 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
         The theory shows that these are all homogeneous polynomials of
         a known degree::
 
-            sage: set_random_seed()
             sage: J = random_eja()
             sage: all(p.is_homogeneous() for p in J._charpoly_coefficients())
             True
@@ -1541,16 +1579,9 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
         """
         n = self.dimension()
         R = self.coordinate_polynomial_ring()
-        vars = R.gens()
         F = R.fraction_field()
 
-        def L_x_i_j(i,j):
-            # From a result in my book, these are the entries of the
-            # basis representation of L_x.
-            return sum( vars[k]*self.monomial(k).operator().matrix()[i,j]
-                        for k in range(n) )
-
-        L_x = matrix(F, n, n, L_x_i_j)
+        L_x = self.operator_polynomial_matrix()
 
         r = None
         if self.rank.is_in_cache():
@@ -1631,7 +1662,6 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
         positive integer rank, unless the algebra is trivial in
         which case its rank will be zero::
 
-            sage: set_random_seed()
             sage: J = random_eja()
             sage: r = J.rank()
             sage: r in ZZ
@@ -1642,7 +1672,6 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
         Ensure that computing the rank actually works, since the ranks
         of all simple algebras are known and will be cached by default::
 
-            sage: set_random_seed()    # long time
             sage: J = random_eja()     # long time
             sage: cached = J.rank()    # long time
             sage: J.rank.clear_cache() # long time
@@ -1771,7 +1800,7 @@ class RationalBasisEJA(FiniteDimensionalEJA):
 
             sage: J = JordanSpinEJA(3)
             sage: J._charpoly_coefficients()
-            (X1^2 - X2^2 - X3^2, -2*X1)
+            (X0^2 - X1^2 - X2^2, -2*X0)
             sage: a0 = J._charpoly_coefficients()[0]
             sage: J.base_ring()
             Algebraic Real Field
@@ -1817,7 +1846,6 @@ class ConcreteEJA(FiniteDimensionalEJA):
     Our basis is normalized with respect to the algebra's inner
     product, unless we specify otherwise::
 
-        sage: set_random_seed()
         sage: J = ConcreteEJA.random_instance()
         sage: all( b.norm() == 1 for b in J.gens() )
         True
@@ -1828,7 +1856,6 @@ class ConcreteEJA(FiniteDimensionalEJA):
     natural->EJA basis representation is an isometry and within the
     EJA the operator is self-adjoint by the Jordan axiom::
 
-        sage: set_random_seed()
         sage: J = ConcreteEJA.random_instance()
         sage: x = J.random_element()
         sage: x.operator().is_self_adjoint()
@@ -1902,11 +1929,11 @@ class ConcreteEJA(FiniteDimensionalEJA):
         return eja_class.random_instance(max_dimension, *args, **kwargs)
 
 
-class MatrixEJA(FiniteDimensionalEJA):
+class HermitianMatrixEJA(FiniteDimensionalEJA):
     @staticmethod
     def _denormalized_basis(A):
         """
-        Returns a basis for the space of complex Hermitian n-by-n matrices.
+        Returns a basis for the given Hermitian matrix space.
 
         Why do we embed these? Basically, because all of numerical linear
         algebra assumes that you're working with vectors consisting of `n`
@@ -1919,41 +1946,37 @@ class MatrixEJA(FiniteDimensionalEJA):
             sage: from mjo.hurwitz import (ComplexMatrixAlgebra,
             ....:                          QuaternionMatrixAlgebra,
             ....:                          OctonionMatrixAlgebra)
-            sage: from mjo.eja.eja_algebra import MatrixEJA
+            sage: from mjo.eja.eja_algebra import HermitianMatrixEJA
 
         TESTS::
 
-            sage: set_random_seed()
             sage: n = ZZ.random_element(1,5)
             sage: A = MatrixSpace(QQ, n)
-            sage: B = MatrixEJA._denormalized_basis(A)
+            sage: B = HermitianMatrixEJA._denormalized_basis(A)
             sage: all( M.is_hermitian() for M in  B)
             True
 
         ::
 
-            sage: set_random_seed()
             sage: n = ZZ.random_element(1,5)
             sage: A = ComplexMatrixAlgebra(n, scalars=QQ)
-            sage: B = MatrixEJA._denormalized_basis(A)
+            sage: B = HermitianMatrixEJA._denormalized_basis(A)
             sage: all( M.is_hermitian() for M in  B)
             True
 
         ::
 
-            sage: set_random_seed()
             sage: n = ZZ.random_element(1,5)
             sage: A = QuaternionMatrixAlgebra(n, scalars=QQ)
-            sage: B = MatrixEJA._denormalized_basis(A)
+            sage: B = HermitianMatrixEJA._denormalized_basis(A)
             sage: all( M.is_hermitian() for M in B )
             True
 
         ::
 
-            sage: set_random_seed()
             sage: n = ZZ.random_element(1,5)
             sage: A = OctonionMatrixAlgebra(n, scalars=QQ)
-            sage: B = MatrixEJA._denormalized_basis(A)
+            sage: B = HermitianMatrixEJA._denormalized_basis(A)
             sage: all( M.is_hermitian() for M in B )
             True
 
@@ -2058,7 +2081,7 @@ class MatrixEJA(FiniteDimensionalEJA):
         self.rank.set_cache(matrix_space.nrows())
         self.one.set_cache( self(matrix_space.one()) )
 
-class RealSymmetricEJA(MatrixEJA, RationalBasisEJA, ConcreteEJA):
+class RealSymmetricEJA(HermitianMatrixEJA, RationalBasisEJA, ConcreteEJA):
     """
     The rank-n simple EJA consisting of real symmetric n-by-n
     matrices, the usual symmetric Jordan product, and the trace inner
@@ -2091,7 +2114,6 @@ class RealSymmetricEJA(MatrixEJA, RationalBasisEJA, ConcreteEJA):
 
     The dimension of this algebra is `(n^2 + n) / 2`::
 
-        sage: set_random_seed()
         sage: d = RealSymmetricEJA._max_random_instance_dimension()
         sage: n = RealSymmetricEJA._max_random_instance_size(d)
         sage: J = RealSymmetricEJA(n)
@@ -2100,7 +2122,6 @@ class RealSymmetricEJA(MatrixEJA, RationalBasisEJA, ConcreteEJA):
 
     The Jordan multiplication is what we think it is::
 
-        sage: set_random_seed()
         sage: J = RealSymmetricEJA.random_instance()
         sage: x,y = J.random_elements(2)
         sage: actual = (x*y).to_matrix()
@@ -2152,7 +2173,7 @@ class RealSymmetricEJA(MatrixEJA, RationalBasisEJA, ConcreteEJA):
 
 
 
-class ComplexHermitianEJA(MatrixEJA, RationalBasisEJA, ConcreteEJA):
+class ComplexHermitianEJA(HermitianMatrixEJA, RationalBasisEJA, ConcreteEJA):
     """
     The rank-n simple EJA consisting of complex Hermitian n-by-n
     matrices over the real numbers, the usual symmetric Jordan product,
@@ -2192,7 +2213,6 @@ class ComplexHermitianEJA(MatrixEJA, RationalBasisEJA, ConcreteEJA):
 
     The dimension of this algebra is `n^2`::
 
-        sage: set_random_seed()
         sage: d = ComplexHermitianEJA._max_random_instance_dimension()
         sage: n = ComplexHermitianEJA._max_random_instance_size(d)
         sage: J = ComplexHermitianEJA(n)
@@ -2201,7 +2221,6 @@ class ComplexHermitianEJA(MatrixEJA, RationalBasisEJA, ConcreteEJA):
 
     The Jordan multiplication is what we think it is::
 
-        sage: set_random_seed()
         sage: J = ComplexHermitianEJA.random_instance()
         sage: x,y = J.random_elements(2)
         sage: actual = (x*y).to_matrix()
@@ -2253,7 +2272,7 @@ class ComplexHermitianEJA(MatrixEJA, RationalBasisEJA, ConcreteEJA):
         return cls(n, **kwargs)
 
 
-class QuaternionHermitianEJA(MatrixEJA, RationalBasisEJA, ConcreteEJA):
+class QuaternionHermitianEJA(HermitianMatrixEJA, RationalBasisEJA, ConcreteEJA):
     r"""
     The rank-n simple EJA consisting of self-adjoint n-by-n quaternion
     matrices, the usual symmetric Jordan product, and the
@@ -2278,7 +2297,6 @@ class QuaternionHermitianEJA(MatrixEJA, RationalBasisEJA, ConcreteEJA):
 
     The dimension of this algebra is `2*n^2 - n`::
 
-        sage: set_random_seed()
         sage: d = QuaternionHermitianEJA._max_random_instance_dimension()
         sage: n = QuaternionHermitianEJA._max_random_instance_size(d)
         sage: J = QuaternionHermitianEJA(n)
@@ -2287,7 +2305,6 @@ class QuaternionHermitianEJA(MatrixEJA, RationalBasisEJA, ConcreteEJA):
 
     The Jordan multiplication is what we think it is::
 
-        sage: set_random_seed()
         sage: J = QuaternionHermitianEJA.random_instance()
         sage: x,y = J.random_elements(2)
         sage: actual = (x*y).to_matrix()
@@ -2343,7 +2360,7 @@ class QuaternionHermitianEJA(MatrixEJA, RationalBasisEJA, ConcreteEJA):
         n = ZZ.random_element(max_size + 1)
         return cls(n, **kwargs)
 
-class OctonionHermitianEJA(MatrixEJA, RationalBasisEJA, ConcreteEJA):
+class OctonionHermitianEJA(HermitianMatrixEJA, RationalBasisEJA, ConcreteEJA):
     r"""
     SETUP::
 
@@ -2433,7 +2450,7 @@ class OctonionHermitianEJA(MatrixEJA, RationalBasisEJA, ConcreteEJA):
     @staticmethod
     def _max_random_instance_size(max_dimension):
         r"""
-        The maximum rank of a random QuaternionHermitianEJA.
+        The maximum rank of a random OctonionHermitianEJA.
         """
         # There's certainly a formula for this, but with only four
         # cases to worry about, I'm not that motivated to derive it.
@@ -2660,7 +2677,6 @@ class BilinearFormEJA(RationalBasisEJA, ConcreteEJA):
     matrix.  We opt not to orthonormalize the basis, because if we
     did, we would have to normalize the `s_{i}` in a similar manner::
 
-        sage: set_random_seed()
         sage: n = ZZ.random_element(5)
         sage: M = matrix.random(QQ, max(0,n-1), algorithm='unimodular')
         sage: B11 = matrix.identity(QQ,1)
@@ -2822,7 +2838,6 @@ class JordanSpinEJA(BilinearFormEJA):
 
         Ensure that we have the usual inner product on `R^n`::
 
-            sage: set_random_seed()
             sage: J = JordanSpinEJA.random_instance()
             sage: x,y = J.random_elements(2)
             sage: actual = x.inner_product(y)
@@ -2943,7 +2958,6 @@ class CartesianProductEJA(FiniteDimensionalEJA):
     The Jordan product is inherited from our factors and implemented by
     our CombinatorialFreeModule Cartesian product superclass::
 
-        sage: set_random_seed()
         sage: J1 = HadamardEJA(2)
         sage: J2 = RealSymmetricEJA(2)
         sage: J = cartesian_product([J1,J2])
@@ -3080,7 +3094,6 @@ class CartesianProductEJA(FiniteDimensionalEJA):
 
     The cached unit element is the same one that would be computed::
 
-        sage: set_random_seed()              # long time
         sage: J1 = random_eja()              # long time
         sage: J2 = random_eja()              # long time
         sage: J = cartesian_product([J1,J2]) # long time
@@ -3299,7 +3312,6 @@ class CartesianProductEJA(FiniteDimensionalEJA):
 
         The answer never changes::
 
-            sage: set_random_seed()
             sage: J1 = random_eja()
             sage: J2 = random_eja()
             sage: J = cartesian_product([J1,J2])
@@ -3389,7 +3401,6 @@ class CartesianProductEJA(FiniteDimensionalEJA):
 
         The answer never changes::
 
-            sage: set_random_seed()
             sage: J1 = random_eja()
             sage: J2 = random_eja()
             sage: J = cartesian_product([J1,J2])
@@ -3402,7 +3413,6 @@ class CartesianProductEJA(FiniteDimensionalEJA):
         produce the identity map, and mismatching them should produce
         the zero map::
 
-            sage: set_random_seed()
             sage: J1 = random_eja()
             sage: J2 = random_eja()
             sage: J = cartesian_product([J1,J2])
@@ -3507,7 +3517,6 @@ def random_eja(max_dimension=None, *args, **kwargs):
 
     TESTS::
 
-        sage: set_random_seed()
         sage: n = ZZ.random_element(1,5)
         sage: J = random_eja(max_dimension=n, field=QQ, orthonormalize=False)
         sage: J.dimension() <= n
@@ -3533,3 +3542,234 @@ def random_eja(max_dimension=None, *args, **kwargs):
         # if the sub-call also Decides on a cartesian product.
         J2 = random_eja(new_max_dimension, *args, **kwargs)
         return cartesian_product([J1,J2])
+
+
+class ComplexSkewSymmetricEJA(RationalBasisEJA, ConcreteEJA):
+    r"""
+    The skew-symmetric EJA of order `2n` described in Faraut and
+    Koranyi's Exercise III.1.b. It has dimension `2n^2 - n`.
+
+    It is (not obviously) isomorphic to the QuaternionHermitianEJA of
+    order `n`, as can be inferred by comparing rank/dimension or
+    explicitly from their "characteristic polynomial of" functions,
+    which just so happen to align nicely.
+
+    SETUP::
+
+        sage: from mjo.eja.eja_algebra import (ComplexSkewSymmetricEJA,
+        ....:                                  QuaternionHermitianEJA)
+        sage: from mjo.eja.eja_operator import FiniteDimensionalEJAOperator
+
+    EXAMPLES:
+
+    This EJA is isomorphic to the quaternions::
+
+        sage: J = ComplexSkewSymmetricEJA(2, field=QQ, orthonormalize=False)
+        sage: K = QuaternionHermitianEJA(2, field=QQ, orthonormalize=False)
+        sage: jordan_isom_matrix = matrix.diagonal(QQ,[-1,1,1,1,1,-1])
+        sage: phi = FiniteDimensionalEJAOperator(J,K,jordan_isom_matrix)
+        sage: all( phi(x*y) == phi(x)*phi(y)
+        ....:      for x in J.gens()
+        ....:      for y in J.gens() )
+        True
+        sage: x,y = J.random_elements(2)
+        sage: phi(x*y) == phi(x)*phi(y)
+        True
+
+    TESTS:
+
+    Random elements should satisfy the same conditions that the basis
+    elements do::
+
+        sage: K = ComplexSkewSymmetricEJA.random_instance(field=QQ,
+        ....:                                             orthonormalize=False)
+        sage: x,y = K.random_elements(2)
+        sage: z = x*y
+        sage: x = x.to_matrix()
+        sage: y = y.to_matrix()
+        sage: z = z.to_matrix()
+        sage: all( e.is_skew_symmetric() for e in (x,y,z) )
+        True
+        sage: J = -K.one().to_matrix()
+        sage: all( e*J == J*e.conjugate() for e in (x,y,z) )
+        True
+
+    The power law in Faraut & Koranyi's II.7.a is satisfied.
+    We're in a subalgebra of theirs, but powers are still
+    defined the same::
+
+        sage: K = ComplexSkewSymmetricEJA.random_instance(field=QQ,
+        ....:                                             orthonormalize=False)
+        sage: x = K.random_element()
+        sage: k = ZZ.random_element(5)
+        sage: actual = x^k
+        sage: J = -K.one().to_matrix()
+        sage: expected = K(-J*(J*x.to_matrix())^k)
+        sage: actual == expected
+        True
+
+    """
+    @staticmethod
+    def _max_random_instance_size(max_dimension):
+        # Obtained by solving d = 2n^2 - n, which comes from noticing
+        # that, in 2x2 block form, any element of this algebra has a
+        # free skew-symmetric top-left block, a Hermitian top-right
+        # block, and two bottom blocks that are determined by the top.
+        # The ZZ-int-ZZ thing is just "floor."
+        return ZZ(int(ZZ(8*max_dimension + 1).sqrt()/4 + 1/4))
+
+    @classmethod
+    def random_instance(cls, max_dimension=None, *args, **kwargs):
+        """
+        Return a random instance of this type of algebra.
+        """
+        class_max_d = cls._max_random_instance_dimension()
+        if (max_dimension is None or max_dimension > class_max_d):
+            max_dimension = class_max_d
+        max_size = cls._max_random_instance_size(max_dimension)
+        n = ZZ.random_element(max_size + 1)
+        return cls(n, **kwargs)
+
+    @staticmethod
+    def _denormalized_basis(A):
+        """
+        SETUP::
+
+            sage: from mjo.hurwitz import ComplexMatrixAlgebra
+            sage: from mjo.eja.eja_algebra import ComplexSkewSymmetricEJA
+
+        TESTS:
+
+        The basis elements are all skew-Hermitian::
+
+            sage: d_max = ComplexSkewSymmetricEJA._max_random_instance_dimension()
+            sage: n_max = ComplexSkewSymmetricEJA._max_random_instance_size(d_max)
+            sage: n = ZZ.random_element(n_max + 1)
+            sage: A = ComplexMatrixAlgebra(2*n, scalars=QQ)
+            sage: B = ComplexSkewSymmetricEJA._denormalized_basis(A)
+            sage: all( M.is_skew_symmetric() for M in  B)
+            True
+
+        The basis elements ``b`` all satisfy ``b*J == J*b.conjugate()``,
+        as in the definition of the algebra::
+
+            sage: d_max = ComplexSkewSymmetricEJA._max_random_instance_dimension()
+            sage: n_max = ComplexSkewSymmetricEJA._max_random_instance_size(d_max)
+            sage: n = ZZ.random_element(n_max + 1)
+            sage: A = ComplexMatrixAlgebra(2*n, scalars=QQ)
+            sage: I_n = matrix.identity(ZZ, n)
+            sage: J = matrix.block(ZZ, 2, 2, (0, I_n, -I_n, 0), subdivide=False)
+            sage: J = A.from_list(J.rows())
+            sage: B = ComplexSkewSymmetricEJA._denormalized_basis(A)
+            sage: all( b*J == J*b.conjugate()  for b in B )
+            True
+
+        """
+        es = A.entry_algebra_gens()
+        gen = lambda A,m: A.monomial(m)
+
+        basis = []
+
+        # The size of the blocks. We're going to treat these thing as
+        # 2x2 block matrices,
+        #
+        #   [  x1        x2      ]
+        #   [ -x2-conj   x1-conj ]
+        #
+        # where x1 is skew-symmetric and x2 is Hermitian.
+        #
+        m = A.nrows()/2
+
+        # We only loop through the top half of the matrix, because the
+        # bottom can be constructed from the top.
+        for i in range(m):
+            # First do the top-left block, which is skew-symmetric.
+            # We can compute the bottom-right block in the process.
+            for j in range(i+1):
+                if i != j:
+                    # Skew-symmetry implies zeros for (i == j).
+                    for e in es:
+                        # Top-left block's entry.
+                        E_ij  = gen(A, (i,j,e))
+                        E_ij -= gen(A, (j,i,e))
+
+                        # Bottom-right block's entry.
+                        F_ij  = gen(A, (i+m,j+m,e)).conjugate()
+                        F_ij -= gen(A, (j+m,i+m,e)).conjugate()
+
+                        basis.append(E_ij + F_ij)
+
+            # Now do the top-right block, which is Hermitian, and compute
+            # the bottom-left block along the way.
+            for j in range(m,i+m+1):
+                if (i+m) == j:
+                    # Hermitian matrices have real diagonal entries.
+                    # Top-right block's entry.
+                    E_ii = gen(A, (i,j,es[0]))
+
+                    # Bottom-left block's entry. Don't conjugate
+                    # 'cause it's real.
+                    E_ii -= gen(A, (i+m,j-m,es[0]))
+                    basis.append(E_ii)
+                else:
+                    for e in es:
+                        # Top-right block's entry. BEWARE! We're not
+                        # reflecting across the main diagonal as in
+                        # (i,j)~(j,i). We're only reflecting across
+                        # the diagonal for the top-right block.
+                        E_ij  = gen(A, (i,j,e))
+
+                        # Shift it back to non-offset coords, transpose,
+                        # conjugate, and put it back:
+                        #
+                        # (i,j) -> (i,j-m) -> (j-m, i) -> (j-m, i+m)
+                        E_ij += gen(A, (j-m,i+m,e)).conjugate()
+
+                        # Bottom-left's block's below-diagonal entry.
+                        # Just shift the top-right coords down m and
+                        # left m.
+                        F_ij  = -gen(A, (i+m,j-m,e)).conjugate()
+                        F_ij += -gen(A, (j,i,e)) # double-conjugate cancels
+
+                        basis.append(E_ij + F_ij)
+
+        return tuple( basis )
+
+    @staticmethod
+    @cached_method
+    def _J_matrix(matrix_space):
+        n = matrix_space.nrows() // 2
+        F = matrix_space.base_ring()
+        I_n = matrix.identity(F, n)
+        J = matrix.block(F, 2, 2, (0, I_n, -I_n, 0), subdivide=False)
+        return matrix_space.from_list(J.rows())
+
+    def J_matrix(self):
+        return ComplexSkewSymmetricEJA._J_matrix(self.matrix_space())
+
+    def __init__(self, n, field=AA, **kwargs):
+        # New code; always check the axioms.
+        #if "check_axioms" not in kwargs: kwargs["check_axioms"] = False
+
+        from mjo.hurwitz import ComplexMatrixAlgebra
+        A = ComplexMatrixAlgebra(2*n, scalars=field)
+        J = ComplexSkewSymmetricEJA._J_matrix(A)
+
+        def jordan_product(X,Y):
+            return (X*J*Y + Y*J*X)/2
+
+        def inner_product(X,Y):
+            return (X.conjugate_transpose()*Y).trace().real()
+
+        super().__init__(self._denormalized_basis(A),
+                         jordan_product,
+                         inner_product,
+                         field=field,
+                         matrix_space=A,
+                         **kwargs)
+
+        # This algebra is conjectured (by me) to be isomorphic to
+        # the quaternion Hermitian EJA of size n, and the rank
+        # would follow from that.
+        #self.rank.set_cache(n)
+        self.one.set_cache( self(-J) )