@@ -1189,10 +1189,6 @@ function modf(x::T) where T<:IEEEFloat
1189
1189
return (rx, ix)
1190
1190
end
1191
1191
1192
- @inline function use_power_by_squaring (n:: Integer )
1193
- - 2 ^ 12 <= n <= 3 * 2 ^ 13
1194
- end
1195
-
1196
1192
# @constprop aggressive to help the compiler see the switch between the integer and float
1197
1193
# variants for callers with constant `y`
1198
1194
@constprop :aggressive function ^ (x:: Float64 , y:: Float64 )
@@ -1205,33 +1201,24 @@ end
1205
1201
y = sign (y)* 0x1 .8 p62
1206
1202
end
1207
1203
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
1225
1208
if xu < (UInt64 (1 )<< 52 ) # x is subnormal
1226
1209
xu = reinterpret (UInt64, x * 0x1 p52) # normalize x
1227
1210
xu &= ~ sign_mask (Float64)
1228
1211
xu -= UInt64 (52 ) << 52 # mess with the exponent
1229
1212
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)
1231
1218
xyhi, xylo = two_mul (logxhi,y)
1232
1219
xylo = muladd (logxlo, y, xylo)
1233
1220
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 (:ℯ ))
1235
1222
end
1236
1223
1237
1224
@constprop :aggressive function ^ (x:: T , y:: T ) where T <: Union{Float16, Float32}
@@ -1255,29 +1242,12 @@ end
1255
1242
return T (exp2 (log2 (abs (widen (x))) * y))
1256
1243
end
1257
1244
1245
+ # compensated power by squaring
1258
1246
@constprop :aggressive @inline function ^ (x:: Float64 , n:: Integer )
1259
- x^ clamp (n, Int64)
1260
- end
1261
- @constprop :aggressive @inline function ^ (x:: Float64 , n:: Int64 )
1262
1247
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)
1277
1249
end
1278
1250
1279
- # compensated power by squaring
1280
- # this method is only reliable for -2^20 < n < 2^20 (cf. #53881 #53886)
1281
1251
@assume_effects :terminates_locally @noinline function pow_body (x:: Float64 , n:: Integer )
1282
1252
y = 1.0
1283
1253
xnlo = ynlo = 0.0
0 commit comments