float16 cbrt, 50% faster (#39441)

Co-authored-by: Kristoffer Carlsson <kcarlsson89@gmail.com>

Author: oscardssmith

base/special/cbrt.jl

@@ -146,3 +146,20 @@ function cbrt(x::Union{Float32,Float64})
     t = _approx_cbrt(x)
     return _improve_cbrt(x, t)
 end
+
+function cbrt(a::Float16)
+    if !isfinite(a) || iszero(a)
+        return a
+    end
+    x = Float32(a)
+
+    # 5 bit approximation. Simpler than _approx_cbrt since subnormals can not appear
+    u = highword(x) & 0x7fff_ffff
+    v = div(u, UInt32(3)) + 0x2a5119f2
+    t = copysign(fromhighword(Float32, v), x)
+
+    # 2 newton iterations
+    t = 0.33333334f0 * (2f0*t + x/(t*t))
+    t = 0.33333334f0 * (2f0*t + x/(t*t))
+    return Float16(t)
+end