`power_by_squaring(::Float16, ::Integer)` involves unnecessary conversion to `Float64`

[giordano]

On any target which has native support for Float16 (e.g. Apple Silicon) you can see that power_by_squaring is doing unnecessary promotion to higher precision (Float64, not just Float32):

julia> code_llvm((Float16,); debuginfo=:none) do (x) x ^ 2 end

; Function Signature: var"#9"(Float16)
define half @"julia_#9_2873"(half %"x::Float16") #0 {
top:
  %0 = fpext half %"x::Float16" to double
  %1 = call double @"j_#power_by_squaring#522_2882"(double %0, i16 signext 2)
  %2 = fptrunc double %1 to float
  %3 = fptrunc float %2 to half
  ret half %3
}

julia> code_llvm((Float16,); debuginfo=:none) do (x) abs2(x) end

; Function Signature: var"#11"(Float16)
define half @"julia_#11_2885"(half %"x::Float16") #0 {
top:
  %0 = fmul half %"x::Float16", %"x::Float16"
  ret half %0
}

Edit: @jakobnissen pointed out this is happening around

julia/base/math.jl

Lines 1343 to 1349 in 912460b

	function ^(x::Float32, n::Integer)
	n == -2 && return (i=inv(x); i*i)
	n == 3 && return x*x*x #keep compatibility with literal_pow
	n < 0 && return Float32(Base.power_by_squaring(inv(Float64(x)),-n))
	Float32(Base.power_by_squaring(Float64(x),n))
	end
	@inline ^(x::Float16, y::Integer) = Float16(Float32(x) ^ y)

Note that this doesn't happen for Float32

julia> code_llvm((Float32,); debuginfo=:none) do x x^2 end

; Function Signature: var"#59"(Float32)
define float @"julia_#59_13034"(float %"x::Float32") #0 {
top:
  %0 = fmul float %"x::Float32", %"x::Float32"
  ret float %0
}

@MasonProtter showed that's because

julia/base/intfuncs.jl

Line 370 in 912460b

@inline literal_pow(::typeof(^), x::HWNumber, ::Val{2}) = x*x

and Float16 isn't in

julia/base/intfuncs.jl

Lines 362 to 363 in 912460b

	const HWReal = Union{Int8,Int16,Int32,Int64,UInt8,UInt16,UInt32,UInt64,Float32,Float64}
	const HWNumber = Union{HWReal, Complex{<:HWReal}, Rational{<:HWReal}}

[giordano]

Now that llvm automatically handles Float16 demotion correctly (#55486), should we be able to extend the literal_pow method to Float16 for all targets?

[oscardssmith]

we could, but it's slightly tricky since for hardware that doesn't support Float16, doing the math in Float16 is probably a bad idea for speed reasons. otoh, it is reasonable to argue that if you do float16 computation without hardware support, it's your fault that it's slow...

[giordano]

To be clear, you aren't really doing math in Float16 on hardware without native support: with the demotion working correctly in LLVM all Float16s are automatically converted to Float32 and the result of each operation is converted back to Float16, and so on. Yes, all these Float16 -> Float32 -> Float16 conversions are going to be slow on hardware which doesn't have native support for Float16, but as you said that should be on the users. Float16 is useful only in limited cases anyway, and users should be aware that it's worth it only on dedicated hardware. The advantage for julia would be in not treating Float16 differently from Float32 and Float64 anymore.