Gentle Reader,

This file is a little rough. It needs improvement in a couple of ways:

Add clConstruction, clCoefficients, clFromCoefficients as in Clifford algebras.

Change multiplication formula from Baez's version to (a, b) (c, d) = (a c - d* b, a d + c* b); then check associativity, etc. This multiplication is supposed to make the isomorphism to the quaternions more natural (no minus sign on the k term).

The Cayley-Dickson Construction

The Cayley-Dickson construction produces the complex numbers from the real numbers, the quaternions from the complex numbers, the octonions from the quaternions, the sedenions from the octonions, etc. This section provides basic arithmetical routines for dealing with these and all other numbers produced by the Cayley-Dickson construction.

The real number x is represented here by the Mathematica list {0, x}. The complex number x + i y is represented by {1, x, y}, meaning {1, {0, x}, {0, y}}. That is, a complex number is of the form {1, a, b}, where a and b are real numbers. A quaternion is of the form {2, a, b}, where a and b are real or complex numbers. An octonion is of the form {3, a, b}, where a and b are real, complex, or quaternionic. In general, an nth-order Cayley-Dickson number is of the form {n, a, b}, where a and b are of order less than n. The functions in this section automatically coerce lower-order numbers to higher-order as needed.

cdZeroQ tests whether a number is 0. cdReduction returns the reduction of a number to the lowest possible order; for example, the quaternion 0 + i 3 + j 0 + k 0 reduces to the complex number 0 + i 3.

cdZeroQ[a_] := 
If[a[[1]] == 0, 
(a[[2]] == 0), 
(cdZeroQ[a[[2]]] && cdZeroQ[a[[3]]])];
cdReduction[a_] := 
If[a[[1]] == 0,
a,
If[cdZeroQ[a[[3]]],
cdReduction[a[[2]]],
{a[[1]], cdReduction[a[[2]]], cdReduction[a[[3]]]}]];

The following functions build real numbers x, complex numbers x + i y, and quaternions x + i y1 + j y2 + k y3, just for example. Notice that there is a sign change in the k-component of quaternions. I think this is right. Implicitly we are writing complex numbers as x + i y and quaternions as (x + i y1) + j (y2 + i y3), so that the last term is j i y3 = - k y3. If we define k = j i instead of the usual k = i j, then this sign change is unnecessary --- and that's how we do Clifford algebras below. But does this convention for quaternions agree with the famous dot/cross product formulas? Does it perhaps amount to a left-handed cross product?! I don't know.

cdRealNumber[x_] := 
{0, x};
cdComplexNumber[x_, y_] := 
{1, cdRealNumber[x], cdRealNumber[y]};
cdQuaternion[x_, y1_, y2_, y3_] := 
{2, cdComplexNumber[x, y1], cdComplexNumber[y2, -y3]};

Numbers come with a negation operation defined by -{0, x} = {0, -x} and -{n, a, b} = {n, -a, -b}. They also possess a conjugation operation defined by {0, x}* = {0, x} and {n, a, b}* = {n, a*, -b}.

cdNegation[a_] := 
If[a[[1]] == 0,
{0, -a[[2]]},
{a[[1]], cdNegation[a[[2]]], cdNegation[a[[3]]]}];
cdConjugate[a_] := 
If[a[[1]] == 0,
{0, a[[2]]},
{a[[1]], cdConjugate[a[[2]]], cdNegation[a[[3]]]}];

For two numbers of the same order, the sum operation is defined by {n, a, b} + {n, c, d} = {n, a + c, b + d} and the product operation is defined by {n, a, b}.{n, c, d} = {n, a c - d b*, a* d + cb}. The following functions perform these operations, with automatic coercion of lower-order numbers up to higher order.

cdSum[a_, b_] := 
If[a[[1]] < b[[1]],
{b[[1]], cdSum[a, b[[2]]], b[[3]]},
If[a[[1]] > b[[1]],
{a[[1]], cdSum[a[[2]], b], a[[3]]},
If[a[[1]] == 0,
{0, a[[2]] + b[[2]]},
{a[[1]], cdSum[a[[2]], b[[2]]], cdSum[a[[3]], b[[3]]]}]]];
cdProduct[a_, b_] := 
If[a[[1]] < b[[1]],
{b[[1]], cdProduct[a, b[[2]]], cdProduct[cdConjugate[a], b[[3]]]},
If[a[[1]] > b[[1]],
{a[[1]], cdProduct[a[[2]], b], cdProduct[b, a[[3]]]},
If[a[[1]] == 0,
{0, a[[2]] * b[[2]]},
{a[[1]], 
cdSum[cdProduct[a[[2]], b[[2]]], cdNegation[cdProduct[b[[3]], cdConjugate[a[[3]]]]]],
cdSum[cdProduct[cdConjugate[a[[2]]], b[[3]]], cdProduct[b[[2]], a[[3]]]]}]]];

The conjugation operation is "nicely normed", meaning that (a b)* = b* a* and a a* = a* a > 0 (hence real). The norm |a| of a is then the square root of a a*. The multiplicative inverse of a is a* / (a a*) = (1 / |a|^2) a*.

cdNormSquared[a_] := 
cdReduction[cdProduct[a, cdConjugate[a]]];
cdInverse[a_] := 
Module[{scalar}, 
If[a[[1]] == 0,
{0, 1/ a[[2]]},
scalar = cdNormSquared[a];
scalar = {0, 1 / scalar[[2]]};
cdProduct[scalar, cdConjugate[a]]]];

Here are some example calculations.

cdInverse[cdComplexNumber[x, y]]

quata = cdQuaternion[x, y1, y2, y3];
quatb = cdQuaternion[z, w1, w2, w3];
oct = {3, quata, quatb};
octNormSquared = cdProduct[oct, cdConjugate[oct]]
cdReduction[octNormSquared]

cdProduct[cdComplexNumber[3, 7], cdRealNumber[11]]

cdProduct[cdQuaternion[0, b1, b2, b3], cdQuaternion[0, d1, d2, d3]]