Skip to content

Type instability in axis #109

@mikmoore

Description

@mikmoore

Return type is always Vector{Float64} when angle(q) is zero or NaN, leading to a type instability.

julia> code_warntype(axis,(QuaternionF32,))
MethodInstance for Quaternions.axis(::QuaternionF32)
  from axis(q::Quaternion) in Quaternions at /home/mike/.julia/packages/Quaternions/GVtxM/src/Quaternion.jl:205
Arguments
  #self#::Core.Const(Quaternions.axis)
  q@_2::QuaternionF32
Locals
  s::Float32
  q@_4::QuaternionF32
Body::Union{Vector{Float32}, Vector{Float64}}
1 ─       (q@_4 = q@_2)
│         (q@_4 = Quaternions.normalize(q@_4))
│   %3  = Quaternions.angle(q@_4)::Float32%4  = (%3 / 2)::Float32
│         (s = Quaternions.sin(%4))
│   %6  = Quaternions.abs(s)::Float32%7  = (%6 > 0)::Bool
└──       goto #3 if not %7
2%9  = Base.getproperty(q@_4, :v1)::Float32%10 = Base.getproperty(q@_4, :v2)::Float32%11 = Base.getproperty(q@_4, :v3)::Float32%12 = Base.vect(%9, %10, %11)::Vector{Float32}%13 = (%12 / s)::Vector{Float32}
└──       return %13
3%15 = Base.vect(1.0, 0.0, 0.0)::Vector{Float64}
└──       return %15

I'd have made a PR but am not at a machine where I can do so right now. I think something like

function axis(q::Quaternion)
    s = sin(angle(q) / 2) * abs(q)
    iszero(s) ?
        [one(q.v1/s), zero(q.v2/s), zero(q.v3/s)] : 
        [q.v1/s, q.v2/s, q.v3/s]
end

will fix the instability and also allow NaN propagation. This still doesn't fix the case where only one imaginary component isinf, where we could give a reasonable answer (assuming we care to have axis defined for infinite quaternions), but that was broken under the old version too.

EDIT: Alternatively,

function axis(q::Quaternion)
    s = abs_imag(q)
    iszero(s) ?
        [one(q.v1/s), zero(q.v2/s), zero(q.v3/s)] : 
        [q.v1/s, q.v2/s, q.v3/s]
end

would be more accurate, propagate NaN on only the imaginary part, and would return a well-defined vector for the zero quaternion. If we don't want zero (or any not-convincingly unit quaternion?) to be a valid input, then I'd suggest a DomainError.

Metadata

Metadata

Assignees

No one assigned

    Labels

    No labels
    No labels

    Type

    No type

    Projects

    No projects

    Milestone

    No milestone

    Relationships

    None yet

    Development

    No branches or pull requests

    Issue actions