One thing I should add in addition is that I think one’s best bet to get started with GeometricAlgebra in Julia without too much work is to just use matrices, and interpret them in terms of geometric algebra.
Geometric algebra true believers generally sneer at the use of matrices and claim that using a matrix is inefficient, and sometimes less accurate. If you use StaticArrays.jl, the performance claims don’t really hold up, though the accuracy concerns maybe carry a kernel of legitimacy, I think they’re a little contrived though.
On the upside,you’d benefit from having the full julia linear algebra ecosystem already implemented for you.
Now, you asked for help first doing this in terms of complex numbers first so that you may bootstrap yourself up to the full geometric algebra, but it’s actually easier for me if I show you how to do this with the 2D real geometric algebra. I’ll leave it to you to figure out how to do the 3D version, okay?
Lets start with the usual building blocks. We want an identity element 𝟙, and two linearly independent basis 'vector’s σ1 and σ2 which square to the identity and anticommute. We can write these as
using StaticArrays
𝟙 = SA[1 0
0 1]
σ1 = SA[0 1
1 0]
σ2 = SA[1 0
0 -1]
Here’s your proof that these have the correct behaviour:
julia> σ1^2 == 𝟙
true
julia> σ2^2 == 𝟙
true
julia> σ1 * σ2 == - σ2 * σ1
true
Now, 2x2 matrices span a 4D vector space, but we currently only have 3 basis elements. What’s the third? We construct it with products of the existing elements. For the 2d algebra, there’s only one such independent product: σ1 * σ2, which we’ll call
ℐ = σ1 * σ2
(that’s typed \scrI<TAB>) and this has all the desired properties of complex numbers, e.g.
julia> ℐ^2 == -𝟙
true
julia> exp(π * ℐ) ≈ -𝟙
true
julia> exp(π/2 * ℐ) ≈ ℐ
true
julia> exp(π/4 * ℐ)' * σ1 * exp(π/4 * ℐ) ≈ σ2
true
So you can just use these things for all your algebraic manipulations, but one problem you’ll encounter is that the way these things print isn’t very informative. i.e.
julia> 1𝟙 + 2σ1 + 3σ2 + 4ℐ
2×2 SMatrix{2, 2, Int64, 4} with indices SOneTo(2)×SOneTo(2):
4 -2
6 -2
What can we do about this? Well, once we have our basis constructed, it’s actually quite easy to just project any old matrix onto that basis and read off their elements:
function decomp(basis::NamedTuple{basisnames}, M::AbstractMatrix) where {basisnames}
N = size(M, 1)
@assert N == size(M, 2)
map(Tuple(basis)) do xi
elem = (1//N) * (M ⋅ xi)
end |> NamedTuple{basisnames}
end
and print them nicely
GAprint(basis) = M -> GAprint(basis, M)
function GAprint(basis::NamedTuple{basisnames}, M::AbstractMatrix) where {basisnames}
elems = map(decomp(basis, M) |> Tuple) do elem
if elem isa Rational && denominator(elem) == 1
numerator(elem)
else
elem
end
end
s = join(("$(elems[i]) $(basisnames[i])" for i in eachindex(elems) if elems[i] != 0), " + ")
if s == ""
println("𝟘")
else
println(s)
end
end
e.g.,
basis2d = (;𝟙, σ1, σ2, ℐ)
julia> SA[4 -2; 6 -2] |> GAprint(basis2d)
1 𝟙 + 2 σ1 + 3 σ2 + 4 ℐ
julia> exp(π/4 * ℐ)' * σ1 * exp(π/4 * ℐ) |> GAprint(basis2d)
3.3306690738754696e-16 σ1 + 1.0000000000000002 σ2
Is this at all helpful @banksiaboy?
So I can see how hard it might be to do the geometric algebra basis. Complex numbers were a first step for me.