Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Add isunit and change norm field to property #75

Closed
wants to merge 10 commits into from
Closed

Add isunit and change norm field to property #75

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
10 commits
Select commit Hold shift + click to select a range
e0d1aa8
Define isunit
sethaxen Feb 25, 2022
8abf205
Remove norm field
sethaxen Feb 25, 2022
62dc247
Define isunit for our types
sethaxen Feb 25, 2022
05747b5
Drop all calls to .norm
sethaxen Feb 25, 2022
2c3f48b
Use isunit for external type checks as well
sethaxen Feb 25, 2022
abd6b99
Use more reliable isunit implementation
sethaxen Feb 25, 2022
f61dac0
Don't use norm argument anywhere
sethaxen Feb 25, 2022
138c746
Merge branch 'master' into isunit
sethaxen Mar 7, 2022
8182877
Add missing deprecations
sethaxen Mar 7, 2022
76b15e4
Merge remote-tracking branch 'upstream/master' into isunit
sethaxen Mar 20, 2022
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
89 changes: 54 additions & 35 deletions src/DualQuaternion.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -1,39 +1,53 @@
using DualNumbers

struct DualQuaternion{T<:Real} <: Number
q0::Quaternion{T}
qe::Quaternion{T}
norm::Bool
end

DualQuaternion{T}(dq::DualQuaternion) where {T<:Real} = DualQuaternion{T}(dq.q0, dq.qe, dq.norm)
function DualQuaternion{T}(d1::Dual, d2::Dual, d3::Dual, d4::Dual, n::Bool=false) where {T<:Real}
Base.@deprecate(
DualQuaternion{T}(q0::Quaternion{T}, qe::Quaternion{T}, norm::Bool) where {T<:Real},
DualQuaternion{T}(q0, qe)
)
Base.@deprecate(
DualQuaternion{T}(d1::Dual, d2::Dual, d3::Dual, d4::Dual, norm::Bool) where {T<:Real},
DualQuaternion{T}(d1, d2, d3, d4)
)
Base.@deprecate(
DualQuaternion(q0::Quaternion, qe::Quaternion, norm::Bool),
DualQuaternion(q0, qe)
)
Base.@deprecate(
DualQuaternion(d1::Dual, d2::Dual, d3::Dual, d4::Dual, norm::Bool),
DualQuaternion(d1, d2, d3, d4)
)
Base.@deprecate dualquat(q0, qe, n) dualquat(q0, qe)
Base.@deprecate dualquat(d1, d2, d3, d4, n) dualquat(d1, d2, d3, d4)

DualQuaternion{T}(dq::DualQuaternion) where {T<:Real} = DualQuaternion{T}(dq.q0, dq.qe)
function DualQuaternion{T}(d1::Dual, d2::Dual, d3::Dual, d4::Dual) where {T<:Real}
return DualQuaternion{T}(
Quaternion(DualNumbers.value(d1), DualNumbers.value(d2), DualNumbers.value(d3), DualNumbers.value(d4), n),
Quaternion(DualNumbers.value(d1), DualNumbers.value(d2), DualNumbers.value(d3), DualNumbers.value(d4)),
Quaternion(DualNumbers.epsilon(d1), DualNumbers.epsilon(d2), DualNumbers.epsilon(d3), DualNumbers.epsilon(d4)),
n,
)
end
function DualQuaternion{T}(q0::Quaternion) where {T<:Real}
return DualQuaternion{T}(convert(Quaternion{T}, q0), zero(Quaternion{T}), q0.norm)
return DualQuaternion{T}(convert(Quaternion{T}, q0), zero(Quaternion{T}))
end
function DualQuaternion{T}(d::Dual) where {T<:Real}
return DualQuaternion(
Quaternion{T}(DualNumbers.value(d)),
Quaternion{T}(DualNumbers.epsilon(d)),
(DualNumbers.value(d)==one(DualNumbers.value(d))) & iszero(DualNumbers.epsilon(d)))
Quaternion{T}(DualNumbers.epsilon(d)))
end
function DualQuaternion{T}(x::Real) where {T<:Real}
return DualQuaternion(convert(Quaternion{T}, x), zero(Quaternion{T}), abs(x) == one(x))
return DualQuaternion(convert(Quaternion{T}, x), zero(Quaternion{T}))
end

DualQuaternion(q0::Quaternion, qe::Quaternion) = DualQuaternion(promote(q0, qe)..., false)
DualQuaternion(d1::Dual, d2::Dual, d3::Dual, d4::Dual, n::Bool=false) =
DualQuaternion(Quaternion(DualNumbers.value(d1), DualNumbers.value(d2), DualNumbers.value(d3), DualNumbers.value(d4), n),
Quaternion(DualNumbers.epsilon(d1), DualNumbers.epsilon(d2), DualNumbers.epsilon(d3), DualNumbers.epsilon(d4)), n)
DualQuaternion(x::Real) = DualQuaternion(Quaternion(x), Quaternion(zero(x)), abs(x) == one(x))
DualQuaternion(d::Dual) = DualQuaternion(Quaternion(DualNumbers.value(d)), Quaternion(DualNumbers.epsilon(d)), (DualNumbers.value(d)==one(DualNumbers.value(d))) & iszero(DualNumbers.epsilon(d)))
DualQuaternion(q::Quaternion) = DualQuaternion(q, zero(q), q.norm)
DualQuaternion(q0::Quaternion, qe::Quaternion) = DualQuaternion(promote(q0, qe)...)
DualQuaternion(d1::Dual, d2::Dual, d3::Dual, d4::Dual) =
DualQuaternion(Quaternion(DualNumbers.value(d1), DualNumbers.value(d2), DualNumbers.value(d3), DualNumbers.value(d4)),
Quaternion(DualNumbers.epsilon(d1), DualNumbers.epsilon(d2), DualNumbers.epsilon(d3), DualNumbers.epsilon(d4)))
DualQuaternion(x::Real) = DualQuaternion(Quaternion(x), Quaternion(zero(x)))
DualQuaternion(d::Dual) = DualQuaternion(Quaternion(DualNumbers.value(d)), Quaternion(DualNumbers.epsilon(d)))
DualQuaternion(q::Quaternion) = DualQuaternion(q, zero(q))
DualQuaternion(a::Vector) = DualQuaternion(zero(Quaternion{typeof(a[1])}), Quaternion(a))

const DualQuaternionF16 = DualQuaternion{Float16}
Expand All @@ -46,11 +60,17 @@ promote_rule(::Type{DualQuaternion{T}}, ::Type{Quaternion{S}}) where {T <: Real,
promote_rule(::Type{DualQuaternion{T}}, ::Type{DualQuaternion{S}}) where {T <: Real, S <: Real} = DualQuaternion{promote_type(T, S)}

dualquat(q1, q2) = DualQuaternion(q1, q2)
dualquat(q1, q2, n) = DualQuaternion(q1, q2, n)
dualquat(d1, d2, d3, d4) = DualQuaternion(d1, d2, d3, d4)
dualquat(d1, d2, d3, d4, n) = DualQuaternion(d1, d2, d3, d4, n)
dualquat(x) = DualQuaternion(x)

function Base.getproperty(dq::DualQuaternion, k::Symbol)
if k === :norm
Base.depwarn("`dq.norm` is deprecated. Please use `isunit(dq)` instead.", Symbol("dq.norm"))
return isunit(dq)
end
return getfield(dq, k)
end

function show(io::IO, dq::DualQuaternion)
show(io, dq.q0)
print(io, " + ( ")
Expand All @@ -69,58 +89,57 @@ Qe(dq::DualQuaternion) = dq.qe
dual(dq.q0.v2, dq.qe.v2) / d,
dual(dq.q0.v3, dq.qe.v3) / d)

abs2(dq::DualQuaternion) = dq.norm ? dual(one(dq.q0.s)) :
abs2(dq::DualQuaternion) = isunit(dq) ? dual(one(dq.q0.s)) :
dual(abs2(dq.q0),
2.0 * (dq.q0.s * dq.qe.s +
dq.q0.v1 * dq.qe.v1 +
dq.q0.v2 * dq.qe.v2 +
dq.q0.v3 * dq.qe.v3))

abs(dq::DualQuaternion) = dq.norm ? dual(one(dq.q0.s)) : sqrt(abs2(dq))
abs(dq::DualQuaternion) = isunit(dq) ? dual(one(dq.q0.s)) : sqrt(abs2(dq))
float(dq::DualQuaternion{T}) where T = convert(DualQuaternion{float(T)}, dq)

conj(dq::DualQuaternion) = DualQuaternion(conj(dq.q0), conj(dq.qe), dq.norm)
dconj(dq::DualQuaternion) = DualQuaternion(dq.q0, -dq.qe, dq.norm)
conj(dq::DualQuaternion) = DualQuaternion(conj(dq.q0), conj(dq.qe))
dconj(dq::DualQuaternion) = DualQuaternion(dq.q0, -dq.qe)

inv(dq::DualQuaternion) = dq.norm ? conj(dq) : conj(dq) / abs2(dq)
inv(dq::DualQuaternion) = isunit(dq) ? conj(dq) : conj(dq) / abs2(dq)

function normalize(dq::DualQuaternion)
if (dq.norm)
if (isunit(dq))
return dq
end
a = abs(dq)
if abs(a) > 0
qa = dq / a
dualquat(qa.q0, qa.qe, true)
dualquat(qa.q0, qa.qe)
else
dq
end
end

function normalizea(dq::DualQuaternion{T}) where {T}
if (dq.norm)
if (isunit(dq))
return (dq, one(Dual{T}))
end
a = abs(dq)
if abs(a) > 0
qa = dq / a
dualquat(qa.q0, qa.qe, true), a
dualquat(qa.q0, qa.qe), a
else
dq, zero(Dual{T})
end
end

(-)(dq::DualQuaternion) = DualQuaternion(-dq.q0, -dq.qe, dq.norm)
(-)(dq::DualQuaternion) = DualQuaternion(-dq.q0, -dq.qe)

(+)(dq::DualQuaternion, dw::DualQuaternion) = DualQuaternion(dq.q0 + dw.q0, dq.qe + dw.qe)
(-)(dq::DualQuaternion, dw::DualQuaternion) = DualQuaternion(dq.q0 - dw.q0, dq.qe - dw.qe)
(*)(dq::DualQuaternion, dw::DualQuaternion) = DualQuaternion(dq.q0 * dw.q0,
dq.q0 * dw.qe + dq.qe * dw.q0,
dq.norm && dw.norm)
dq.q0 * dw.qe + dq.qe * dw.q0)
(*)(dq::DualQuaternion, d::Dual) = (*)(Base.promote(dq, d)...)
(*)(d::Dual, dq::DualQuaternion) = (*)(Base.promote(d, dq)...)
(/)(dq::DualQuaternion, dw::DualQuaternion) = dq * inv(dw)
(==)(q::DualQuaternion, w::DualQuaternion) = (q.q0 == w.q0) & (q.qe == w.qe) # ignore .norm field
(==)(q::DualQuaternion, w::DualQuaternion) = (q.q0 == w.q0) & (q.qe == w.qe)

function angleaxis(dq::DualQuaternion)
tq = dq.qe * conj(dq.q0)
Expand Down Expand Up @@ -156,7 +175,7 @@ function exp(dq::DualQuaternion)
se = exp(se)
dq = dualquat(quat(0.0, imag_part(dq.q0)...), quat(0.0, imag_part(dq.qe)...))
dq, th = normalizea(dq)
if dq.norm
if isunit(dq)
dualquat(se) * (dualquat(cos(th)) + dq * dualquat(sin(th)))
else
dualquat(se)
Expand All @@ -180,5 +199,5 @@ dualquatrand() = dualquat(quatrand(), quatrand())
ndualquatrand() = normalize(dualquatrand())

function rand(rng::AbstractRNG, ::Random.SamplerType{DualQuaternion{T}}) where {T<:Real}
return DualQuaternion{T}(rand(rng, Quaternion{T}), rand(rng, Quaternion{T}), false)
return DualQuaternion{T}(rand(rng, Quaternion{T}), rand(rng, Quaternion{T}))
end
60 changes: 37 additions & 23 deletions src/Octonion.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -7,20 +7,19 @@ 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}(x::Complex) where {T<:Real} = Octonion(convert(Complex{T}, x))
Octonion{T}(q::Quaternion) where {T<:Real} = Octonion(convert(Quaternion{T}, q))
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(z::Complex) = Octonion(z.re, z.im, zero(z.re), zero(z.re), zero(z.re), zero(z.re), zero(z.re), zero(z.re), abs(z) == one(z.re))
Octonion(q::Quaternion) = Octonion(q.s, q.v1, q.v2, q.v3, zero(q.s), zero(q.s), zero(q.s), zero(q.s), q.norm)
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(z::Complex) = Octonion(z.re, z.im, zero(z.re), zero(z.re), zero(z.re), zero(z.re), zero(z.re), zero(z.re))
Octonion(q::Quaternion) = Octonion(q.s, q.v1, q.v2, q.v3, zero(q.s), zero(q.s), zero(q.s), zero(q.s))
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 @@ -34,10 +33,27 @@ promote_rule(::Type{Octonion{T}}, ::Type{Quaternion{S}}) where {T <: Real, S <:
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)

Base.@deprecate(
Octonion{T}(s::Real, v1::Real, v2::Real, v3::Real, v4::Real, v5::Real, v6::Real, v7::Real, norm::Bool) where {T<:Real},
Octonion{T}(s, v1, v2, v3, v4, v5, v6, v7)
)
Base.@deprecate(
Octonion(s::Real, v1::Real, v2::Real, v3::Real, v4::Real, v5::Real, v6::Real, v7::Real, norm::Bool),
Octonion(s, v1, v2, v3, v4, v5, v6, v7)
)
Base.@deprecate octo(p, v1, v2, v3, v4, v5, v6, v7, n) octo(p, v1, v2, v3, v4, v5, v6, v7)

function Base.getproperty(o::Octonion, k::Symbol)
if k === :norm
Base.depwarn("`o.norm` is deprecated. Please use `isunit(o)` instead.", Symbol("o.norm"))
return isunit(o)
end
return getfield(o, k)
end
hyrodium marked this conversation as resolved.
Show resolved Hide resolved

function show(io::IO, o::Octonion)
pm(x) = x < 0 ? " - $(-x)" : " + $x"
print(io, o.s, pm(o.v1), "im", pm(o.v2), "jm", pm(o.v3), "km", pm(o.v4), "ilm", pm(o.v5), "jlm", pm(o.v6), "klm", pm(o.v7), "lm")
Expand All @@ -49,37 +65,37 @@ 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)

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(o.s * o.s + o.v1 * o.v1 + o.v2 * o.v2 + o.v3 * o.v3 + o.v4 * o.v4 + o.v5 * o.v5 + o.v6 * o.v6 + o.v7 * o.v7)
float(q::Octonion{T}) where T = convert(Octonion{float(T)}, q)
abs_imag(o::Octonion) = sqrt(o.v1 * o.v1 + o.v2 * o.v2 + o.v3 * o.v3 + o.v4 * o.v4 + o.v5 * o.v5 + o.v6 * o.v6 + o.v7 * o.v7)
abs2(o::Octonion) = o.s * o.s + o.v1 * o.v1 + o.v2 * o.v2 + o.v3 * o.v3 + o.v4 * o.v4 + o.v5 * o.v5 + o.v6 * o.v6 + o.v7 * o.v7
inv(o::Octonion) = o.norm ? conj(o) : conj(o) / abs2(o)
inv(o::Octonion) = isunit(o) ? conj(o) : 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)
if (isunit(o))
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)
Octonion(o.s, o.v1, o.v2, o.v3, o.v4, o.v5, o.v6, o.v7)
end

function normalizea(o::Octonion)
if (o.norm)
if (isunit(o))
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)
(Octonion(o.s, o.v1, o.v2, o.v3, o.v4, o.v5, o.v6, o.v7), 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 @@ -121,7 +137,7 @@ 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 @@ -138,8 +154,7 @@ 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)
Expand Down Expand Up @@ -172,7 +187,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 @@ -185,6 +200,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
Loading