@@ -1189,10 +1189,6 @@ function modf(x::T) where T<:IEEEFloat
11891189 return (rx, ix)
11901190end
11911191
1192- @inline function use_power_by_squaring (n:: Integer )
1193- - 2 ^ 12 <= n <= 3 * 2 ^ 13
1194- end
1195-
11961192# @constprop aggressive to help the compiler see the switch between the integer and float
11971193# variants for callers with constant `y`
11981194@constprop :aggressive function ^ (x:: Float64 , y:: Float64 )
@@ -1205,33 +1201,24 @@ end
12051201 y = sign (y)* 0x1 .8 p62
12061202 end
12071203 yint = unsafe_trunc (Int64, y) # This is actually safe since julia freezes the result
1208- yisint = y == yint
1209- if yisint
1210- yint == 0 && return 1.0
1211- use_power_by_squaring (yint) && return @noinline pow_body (x, yint)
1212- end
1213- 2 * xu== 0 && return abs (y)* Inf * (! (y> 0 )) # if x === +0.0 or -0.0 (Inf * false === 0.0)
1214- s = 1
1215- if x < 0
1216- ! yisint && throw_exp_domainerror (x) # y isn't an integer
1217- s = ifelse (isodd (yint), - 1 , 1 )
1218- end
1219- ! isfinite (x) && return copysign (x,s)* (y> 0 || isnan (x)) # x is inf or NaN
1220- return copysign (pow_body (abs (x), y), s)
1221- end
1222-
1223- @assume_effects :foldable @noinline function pow_body (x:: Float64 , y:: Float64 )
1224- xu = reinterpret (UInt64, x)
1204+ y == yint && return @noinline x^ yint
1205+ 2 * xu== 0 && return abs (y)* Inf * (! (y> 0 )) # if x==0
1206+ x< 0 && throw_exp_domainerror (x) # |y| is small enough that y isn't an integer
1207+ ! isfinite (x) && return x* (y> 0 || isnan (x)) # x is inf or NaN
12251208 if xu < (UInt64 (1 )<< 52 ) # x is subnormal
12261209 xu = reinterpret (UInt64, x * 0x1 p52) # normalize x
12271210 xu &= ~ sign_mask (Float64)
12281211 xu -= UInt64 (52 ) << 52 # mess with the exponent
12291212 end
1230- logxhi,logxlo = _log_ext (xu)
1213+ return pow_body (xu, y)
1214+ end
1215+
1216+ @inline function pow_body (xu:: UInt64 , y:: Float64 )
1217+ logxhi,logxlo = Base. Math. _log_ext (xu)
12311218 xyhi, xylo = two_mul (logxhi,y)
12321219 xylo = muladd (logxlo, y, xylo)
12331220 hi = xyhi+ xylo
1234- return @inline Base. Math. exp_impl (hi, xylo- (hi- xyhi), Val (:ℯ ))
1221+ return Base. Math. exp_impl (hi, xylo- (hi- xyhi), Val (:ℯ ))
12351222end
12361223
12371224@constprop :aggressive function ^ (x:: T , y:: T ) where T <: Union{Float16, Float32}
@@ -1255,29 +1242,12 @@ end
12551242 return T (exp2 (log2 (abs (widen (x))) * y))
12561243end
12571244
1245+ # compensated power by squaring
12581246@constprop :aggressive @inline function ^ (x:: Float64 , n:: Integer )
1259- x^ clamp (n, Int64)
1260- end
1261- @constprop :aggressive @inline function ^ (x:: Float64 , n:: Int64 )
12621247 n == 0 && return one (x)
1263- if use_power_by_squaring (n)
1264- return pow_body (x, n)
1265- else
1266- s = ifelse (x < 0 && isodd (n), - 1.0 , 1.0 )
1267- x = abs (x)
1268- y = float (n)
1269- if y == n
1270- return copysign (pow_body (x, y), s)
1271- else
1272- n2 = n % 1024
1273- y = float (n - n2)
1274- return pow_body (x, y) * copysign (pow_body (x, n2), s)
1275- end
1276- end
1248+ return pow_body (x, n)
12771249end
12781250
1279- # compensated power by squaring
1280- # this method is only reliable for -2^20 < n < 2^20 (cf. #53881 #53886)
12811251@assume_effects :terminates_locally @noinline function pow_body (x:: Float64 , n:: Integer )
12821252 y = 1.0
12831253 xnlo = ynlo = 0.0
0 commit comments