Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
3 changes: 1 addition & 2 deletions src/Octonions.jl
Original file line number Diff line number Diff line change
Expand Up @@ -4,14 +4,13 @@ import Base: +, -, *, /, ^, ==
import Base: abs, abs2, conj, exp, inv, isreal, isfinite, isinf, iszero, isnan, log, real, sqrt
import Base: promote_rule, float
import Base: rand, randn
import LinearAlgebra: normalize
using Random

Base.@irrational INV_SQRT_EIGHT 0.3535533905932737622004 sqrt(big(0.125))

include("octonion.jl")

export Octonion, OctonionF16, OctonionF32, OctonionF64
export imag_part, normalize, normalizea, octo, octorand
export imag_part, octo, octorand

end # module
52 changes: 16 additions & 36 deletions src/octonion.jl
Original file line number Diff line number Diff line change
Expand Up @@ -7,16 +7,15 @@ struct Octonion{T<:Real} <: Number
v5::T
v6::T
v7::T
norm::Bool
end

Octonion{T}(x::Real) where {T<:Real} = Octonion(convert(T, x))
Octonion{T}(o::Octonion) where {T<:Real} =
Octonion{T}(o.s, o.v1, o.v2, o.v3, o.v4, o.v5, o.v6, o.v7, o.norm)
Octonion{T}(o.s, o.v1, o.v2, o.v3, o.v4, o.v5, o.v6, o.v7)

Octonion(s::Real, v1::Real, v2::Real, v3::Real, v4::Real, v5::Real, v6::Real, v7::Real, n::Bool = false) =
Octonion(promote(s, v1, v2, v3, v4, v5, v6, v7)..., n)
Octonion(x::Real) = Octonion(x, zero(x), zero(x), zero(x), zero(x), zero(x), zero(x), zero(x), abs(x) == one(x))
Octonion(s::Real, v1::Real, v2::Real, v3::Real, v4::Real, v5::Real, v6::Real, v7::Real) =
Octonion(promote(s, v1, v2, v3, v4, v5, v6, v7)...)
Octonion(x::Real) = Octonion(x, zero(x), zero(x), zero(x), zero(x), zero(x), zero(x), zero(x))
Octonion(s::Real, a::Vector) = Octonion(s, a[1], a[2], a[3], a[4], a[5], a[6], a[7])
Octonion(a::Vector) = Octonion(0, a[1], a[2], a[3], a[4], a[5], a[6], a[7])

Expand All @@ -28,7 +27,6 @@ promote_rule(::Type{Octonion{T}}, ::Type{S}) where {T <: Real, S <: Real} = Octo
promote_rule(::Type{Octonion{T}}, ::Type{Octonion{S}}) where {T <: Real, S <: Real} = Octonion{promote_type(T, S)}

octo(p, v1, v2, v3, v4, v5, v6, v7) = Octonion(p, v1, v2, v3, v4, v5, v6, v7)
octo(p, v1, v2, v3, v4, v5, v6, v7, n) = Octonion(p, v1, v2, v3, v4, v5, v6, v7, n)
octo(x) = Octonion(x)
octo(s, a) = Octonion(s, a)

Expand All @@ -40,37 +38,20 @@ imag_part(o::Octonion) = (o.v1, o.v2, o.v3, o.v4, o.v5, o.v6, o.v7)
(*)(o::Octonion, x::Real) = Octonion(o.s * x, o.v1 * x, o.v2 * x, o.v3 * x, o.v4 * x, o.v5 * x, o.v6 * x, o.v7 * x)
(*)(x::Real, o::Octonion) = o * x

conj(o::Octonion) = Octonion(o.s, -o.v1, -o.v2, -o.v3, -o.v4, -o.v5, -o.v6, -o.v7, o.norm)
conj(o::Octonion) = Octonion(o.s, -o.v1, -o.v2, -o.v3, -o.v4, -o.v5, -o.v6, -o.v7)
abs(o::Octonion) = sqrt(abs2(o))
float(q::Octonion{T}) where T = convert(Octonion{float(T)}, q)
abs_imag(o::Octonion) = sqrt((o.v4 * o.v4 + (o.v2 * o.v2 + o.v6 * o.v6)) + ((o.v1 * o.v1 + o.v5 * o.v5) + (o.v3 * o.v3 + o.v7 * o.v7))) # ordered to match abs2
abs2(o::Octonion) = ((o.s * o.s + o.v4 * o.v4) + (o.v2 * o.v2 + o.v6 * o.v6)) + ((o.v1 * o.v1 + o.v5 * o.v5) + (o.v3 * o.v3 + o.v7 * o.v7))
inv(o::Octonion) = o.norm ? conj(o) : conj(o) / abs2(o)
inv(o::Octonion) = conj(o) / abs2(o)

isreal(o::Octonion) = iszero(o.v1) & iszero(o.v2) & iszero(o.v3) & iszero(o.v4) & iszero(o.v5) & iszero(o.v6) & iszero(o.v7)
isfinite(o::Octonion) = o.norm | (isfinite(real(o)) & isfinite(o.v1) & isfinite(o.v2) & isfinite(o.v3) & isfinite(o.v4) & isfinite(o.v5) & isfinite(o.v6) & isfinite(o.v7))
iszero(o::Octonion) = ~o.norm & iszero(real(o)) & iszero(o.v1) & iszero(o.v2) & iszero(o.v3) & iszero(o.v4) & iszero(o.v5) & iszero(o.v6) & iszero(o.v7)
isfinite(o::Octonion) = isfinite(real(o)) & isfinite(o.v1) & isfinite(o.v2) & isfinite(o.v3) & isfinite(o.v4) & isfinite(o.v5) & isfinite(o.v6) & isfinite(o.v7)
iszero(o::Octonion) = iszero(real(o)) & iszero(o.v1) & iszero(o.v2) & iszero(o.v3) & iszero(o.v4) & iszero(o.v5) & iszero(o.v6) & iszero(o.v7)
isnan(o::Octonion) = isnan(real(o)) | isnan(o.v1) | isnan(o.v2) | isnan(o.v3) | isnan(o.v4) | isnan(o.v5) | isnan(o.v6) | isnan(o.v7)
isinf(o::Octonion) = ~o.norm & (isinf(real(o)) | isinf(o.v1) | isinf(o.v2) | isinf(o.v3) | isinf(o.v4) | isinf(o.v5) | isinf(o.v6) | isinf(o.v7))
isinf(o::Octonion) = isinf(real(o)) | isinf(o.v1) | isinf(o.v2) | isinf(o.v3) | isinf(o.v4) | isinf(o.v5) | isinf(o.v6) | isinf(o.v7)

function normalize(o::Octonion)
if (o.norm)
return o
end
o = o / abs(o)
Octonion(o.s, o.v1, o.v2, o.v3, o.v4, o.v5, o.v6, o.v7, true)
end

function normalizea(o::Octonion)
if (o.norm)
return (o, one(o.s))
end
a = abs(o)
o = o / a
(Octonion(o.s, o.v1, o.v2, o.v3, o.v4, o.v5, o.v6, o.v7, true), a)
end

(-)(o::Octonion) = Octonion(-o.s, -o.v1, -o.v2, -o.v3, -o.v4, -o.v5, -o.v6, -o.v7, o.norm)
(-)(o::Octonion) = Octonion(-o.s, -o.v1, -o.v2, -o.v3, -o.v4, -o.v5, -o.v6, -o.v7)

(+)(o::Octonion, w::Octonion) = Octonion(o.s + w.s,
o.v1 + w.v1,
Expand Down Expand Up @@ -99,13 +80,13 @@ function (*)(o::Octonion, w::Octonion)
v5 = ((o.s * w.v5 + o.v5 * w.s) + (o.v2 * w.v7 - o.v7 * w.v2)) + ((o.v1 * w.v6 - o.v6 * w.v1) + (o.v4 * w.v3 - o.v3 * w.v4))
v6 = ((o.s * w.v6 + o.v6 * w.s) + (o.v3 * w.v7 - o.v7 * w.v3)) + ((o.v2 * w.v4 - o.v4 * w.v2) + (o.v5 * w.v1 - o.v1 * w.v5))
v7 = ((o.s * w.v7 + o.v7 * w.s) + (o.v5 * w.v2 - o.v2 * w.v5)) + ((o.v4 * w.v1 - o.v1 * w.v4) + (o.v6 * w.v3 - o.v3 * w.v6))
return Octonion(s, v1, v2, v3, v4, v5, v6, v7, o.norm & w.norm)
return Octonion(s, v1, v2, v3, v4, v5, v6, v7)
end

(/)(o::Octonion, w::Octonion) = o * inv(w)

(==)(q::Octonion, w::Octonion) = (q.s == w.s) & (q.v1 == w.v1) & (q.v2 == w.v2) & (q.v3 == w.v3) &
(q.v4 == w.v4) & (q.v5 == w.v5) & (q.v6 == w.v6) & (q.v7 == w.v7) # ignore .norm field
(q.v4 == w.v4) & (q.v5 == w.v5) & (q.v6 == w.v6) & (q.v7 == w.v7)

function exp(o::Octonion)
s = o.s
Expand All @@ -122,12 +103,12 @@ function exp(o::Octonion)
scale * o.v4,
scale * o.v5,
scale * o.v6,
scale * o.v7,
abs(s) == 0)
scale * o.v7)
end

function log(o::Octonion)
o, a = normalizea(o)
a = abs(o)
o = o / a
s = o.s
M = abs_imag(o)
th = atan(M, s)
Expand Down Expand Up @@ -157,7 +138,7 @@ octorand() = octo(randn(), randn(), randn(), randn(), randn(), randn(), randn(),

function rand(rng::AbstractRNG, ::Random.SamplerType{Octonion{T}}) where {T<:Real}
Octonion{T}(rand(rng, T), rand(rng, T), rand(rng, T), rand(rng, T),
rand(rng, T), rand(rng, T), rand(rng, T), rand(rng, T), false)
rand(rng, T), rand(rng, T), rand(rng, T), rand(rng, T))
end

function randn(rng::AbstractRNG, ::Type{Octonion{T}}) where {T<:AbstractFloat}
Expand All @@ -170,6 +151,5 @@ function randn(rng::AbstractRNG, ::Type{Octonion{T}}) where {T<:AbstractFloat}
randn(rng, T) * INV_SQRT_EIGHT,
randn(rng, T) * INV_SQRT_EIGHT,
randn(rng, T) * INV_SQRT_EIGHT,
false,
)
end
68 changes: 18 additions & 50 deletions test/octonion.jl
Original file line number Diff line number Diff line change
Expand Up @@ -24,32 +24,28 @@ end
@testset "Constructors" begin
@testset "from coefficients" begin
cs = [(1, 2.0, 3.0f0, 4//1, 5, 6, 7, 8), (1//1, 2.0f0, 3.0f0, 4, 5, 6, 7, 8)]
@testset for coef in cs, T in (Float32, Float64, Int), norm in (true, false)
q = @inferred Octonion{T}(coef..., norm)
@testset for coef in cs, T in (Float32, Float64, Int)
q = @inferred Octonion{T}(coef...)
@test q isa Octonion{T}
@test q.norm === norm
@test q === Octonion{T}(convert.(T, coef)..., norm)
q2 = @inferred Octonion(convert.(T, coef)..., norm)
@test Octonion(convert.(T, coef)..., norm) === q
if !norm
@test Octonion(convert.(T, coef)...) === q
end
@test q === Octonion{T}(convert.(T, coef)...)
q2 = @inferred Octonion(convert.(T, coef)...)
@test Octonion(convert.(T, coef)...) === q
end
end
@testset "from real" begin
@testset for x in (-1//1, 1.0, 2.0), T in (Float32, Float64, Int, Rational{Int})
coef = T.((x, zeros(7)...))
@test @inferred(Octonion{T}(x)) === Octonion{T}(coef..., isone(abs(x)))
@test @inferred(Octonion(T(x))) === Octonion{T}(coef..., isone(abs(x)))
@test @inferred(Octonion{T}(x)) === Octonion{T}(coef...)
@test @inferred(Octonion(T(x))) === Octonion{T}(coef...)
end
end
@testset "from octonion" begin
os = (
Octonion(1, 2, 3, 4, 5, 6, 7, 8), OctonionF64(0, 1, 0, 0, 0, 0, 0, 0, true)
Octonion(1, 2, 3, 4, 5, 6, 7, 8), OctonionF64(0, 1, 0, 0, 0, 0, 0, 0)
)
@testset for o in os, T in (Float32, Float64)
coef = T.((o.s, o.v1, o.v2, o.v3, o.v4, o.v5, o.v6, o.v7))
@test @inferred(Octonion{T}(o)) === Octonion{T}(coef..., o.norm)
@test @inferred(Octonion{T}(o)) === Octonion{T}(coef...)
@test @inferred(Octonion(o)) === o
end
end
Expand All @@ -64,16 +60,14 @@ end
@testset "==" begin
@test Octonion(1.0, 2, 3, 4, 5, 6, 7, 8) == Octonion(1, 2, 3, 4, 5, 6, 7, 8)
@test Octonion(1.0, 2, 3, 4, 5, 6, 7, 8) != Octonion(1, 2, 3, 4, 1, 2, 3, 4)
@test Octonion(1, 0, 0, 0, 0, 0, 0, 0, false) ==
Octonion(1, 0, 0, 0, 0, 0, 0, 0, true) # test that .norm field does not affect equality
end

@testset "convert" begin
@test convert(Octonion{Float64}, 1) === Octonion(1.0)
@test convert(Octonion{Float64}, Octonion(1, 2, 3, 4, 5, 6, 7, 8)) ===
Octonion(1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0)
@test convert(Octonion{Float64}, Octonion(0, 0, 0, 0, 1, 0, 0, 0, true)) ===
Octonion(0.0, 0.0, 0.0, 0.0, 1, 0, 0, 0, true)
@test convert(Octonion{Float64}, Octonion(0, 0, 0, 0, 1, 0, 0, 0)) ===
Octonion(0.0, 0.0, 0.0, 0.0, 1, 0, 0, 0)
end

@testset "promote" begin
Expand All @@ -86,13 +80,9 @@ end
end

@testset "shorthands" begin
@test octo(1) === Octonion(1) # checking the .norm field in particular
@test octo(1, 0, 0, 0, 0, 0, 0, 0) === Octonion(1, 0, 0, 0, 0, 0, 0, 0) # checking the .norm field in particular
@test octo(1) === Octonion(1)
@test octo(1, 2, 3, 4, 5, 6, 7, 8) === Octonion(1, 2, 3, 4, 5, 6, 7, 8)
@test octo(Octonion(1, 0, 0, 0, 0, 0, 0, 0)) === Octonion(1, 0, 0, 0, 0, 0, 0, 0) # checking the .norm field in particular
@test octo(Octonion(1, 2, 3, 4, 5, 6, 7, 8)) === Octonion(1, 2, 3, 4, 5, 6, 7, 8)
@test octo(1, 0, 0, 0, 0, 0, 0, 0, false).norm == false # respect the .norm input (even if wrong)
@test octo(1, 2, 3, 4, 5, 6, 7, 8, true).norm == true # respect the .norm input (even if wrong)
@test octo(1, collect(2:8)) === Octonion(1:8...)
@test octo(collect(2:8)) === Octonion(0, 2:8...)
end
Expand All @@ -101,14 +91,12 @@ end
@testset "octorand" begin
o = octorand()
@test o isa Octonion
@test !o.norm
end

@testset "rand($H)" for H in (OctonionF32, OctonionF64)
rng = Random.MersenneTwister(42)
o = rand(rng, H)
@test o isa H
@test !o.norm

os = rand(rng, H, 1000)
@test eltype(os) === H
Expand All @@ -126,7 +114,6 @@ end
rng = Random.MersenneTwister(42)
o = randn(rng, H)
@test o isa H
@test !o.norm

os = randn(rng, H, 10000)
@test eltype(os) === H
Expand All @@ -143,13 +130,13 @@ end

@testset "basic" begin
q = randn(OctonionF64)
qnorm = normalize(q)
qnorm = sign(q)
@test real(q) === q.s
@test_throws MethodError imag(q)
@test @test_deprecated(Octonions.imag(q)) == [q.v1, q.v2, q.v3, q.v4, q.v5, q.v6, q.v7]
@test imag_part(q) === (q.v1, q.v2, q.v3, q.v4, q.v5, q.v6, q.v7)
@test conj(q) ===
Octonion(q.s, -q.v1, -q.v2, -q.v3, -q.v4, -q.v5, -q.v6, -q.v7, q.norm)
Octonion(q.s, -q.v1, -q.v2, -q.v3, -q.v4, -q.v5, -q.v6, -q.v7)
@test conj(qnorm) === Octonion(
qnorm.s,
-qnorm.v1,
Expand All @@ -159,7 +146,6 @@ end
-qnorm.v5,
-qnorm.v6,
-qnorm.v7,
qnorm.norm,
)
@test conj(conj(q)) === q
@test conj(conj(qnorm)) === qnorm
Expand Down Expand Up @@ -212,7 +198,6 @@ end
@test !iszero(octo(0, 0, 0, 0, 0, 1, 0, 0))
@test !iszero(octo(0, 0, 0, 0, 0, 0, 1, 0))
@test !iszero(octo(0, 0, 0, 0, 0, 0, 0, 1))
@test !iszero(octo(0, 0, 0, 0, 0, 0, 0, 0, true))
end

@testset "isone" begin
Expand All @@ -239,7 +224,6 @@ end
@test !isfinite(octo(0, 0, 0, 0, 0, val, 0, 0))
@test !isfinite(octo(0, 0, 0, 0, 0, 0, val, 0))
@test !isfinite(octo(0, 0, 0, 0, 0, 0, 0, val))
@test isfinite(octo(fill(val, 8)..., true))
end
end

Expand All @@ -255,7 +239,6 @@ end
@test isinf(octo(0, 0, 0, 0, 0, inf, 0, 0))
@test isinf(octo(0, 0, 0, 0, 0, 0, inf, 0))
@test isinf(octo(0, 0, 0, 0, 0, 0, 0, inf))
@test !isinf(octo(fill(inf, 8)..., true))
end
end

Expand Down Expand Up @@ -456,29 +439,14 @@ end
end
end

@testset "normalize" begin
@testset "sign" begin
for _ in 1:100
q = randn(OctonionF64)
qnorm = @inferred normalize(q)
qnorm = @inferred sign(q)
@test abs(qnorm) ≈ 1
@test qnorm.norm
@test q ≈ abs(q) * qnorm
@test normalize(qnorm) === qnorm
@test sign(qnorm) qnorm
end
@test_broken @inferred(normalize(octo(1:8...)))
end

@testset "normalizea" begin
for _ in 1:100
q = randn(OctonionF64)
qnorm, a = @inferred normalizea(q)
@test abs(qnorm) ≈ 1
@test qnorm.norm
@test a isa Real
@test a ≈ abs(q)
@test q ≈ a * qnorm
@test normalizea(qnorm) === (qnorm, one(real(q)))
end
@test_broken @inferred(normalizea(octo(1:8...)))
@inferred(sign(octo(1:8...)))
end
end