| 1 | // Ported from musl, which is licensed under the MIT license: |
| 2 | // https://git.musl-libc.org/cgit/musl/tree/COPYRIGHT |
| 3 | // |
| 4 | // https://git.musl-libc.org/cgit/musl/tree/src/math/atanf.c |
| 5 | // https://git.musl-libc.org/cgit/musl/tree/src/math/atan.c |
| 6 | // https://git.musl-libc.org/cgit/musl/tree/src/math/atanl.c |
| 7 | // |
| 8 | // Ported from ARM-software, which is licensed under the MIT license: |
| 9 | // https://github.com/ARM-software/optimized-routines/blob/master/LICENSE |
| 10 | // |
| 11 | // https://github.com/ARM-software/optimized-routines/blob/master/math/aarch64/advsimd/atanf.c |
| 12 | // https://github.com/ARM-software/optimized-routines/blob/master/math/aarch64/advsimd/atan.c |
| 13 | |
| 14 | const builtin = @import("builtin"); |
| 15 | const std = @import("../std.zig"); |
| 16 | const math = std.math; |
| 17 | const mem = std.mem; |
| 18 | const testing = std.testing; |
| 19 | |
| 20 | /// Returns the arc-tangent of x. |
| 21 | /// |
| 22 | /// Special Cases: |
| 23 | /// - atan(+-0) = +-0 |
| 24 | /// - atan(+-inf) = +-pi/2 |
| 25 | pub fn atan(x: anytype) @TypeOf(x) { |
| 26 | const T = @TypeOf(x); |
| 27 | switch (@typeInfo(T)) { |
| 28 | .float => |info| switch (info.bits) { |
| 29 | 16 => return atanBinary16(x), |
| 30 | 32 => return atanBinary32(x), |
| 31 | 64 => return atanBinary64(x), |
| 32 | 80 => return atanExtended80(x), |
| 33 | 128 => return atanBinary128(x), |
| 34 | else => comptime unreachable, |
| 35 | }, |
| 36 | .vector => |info| switch (info.child) { |
| 37 | f32 => return atanBinary32Vec(info.len, x), |
| 38 | f64 => return atanBinary64Vec(info.len, x), |
| 39 | else => @compileError("unimplemented"), |
| 40 | }, |
| 41 | else => comptime unreachable, |
| 42 | } |
| 43 | } |
| 44 | |
| 45 | fn atanBinary16(x: f16) f16 { |
| 46 | const atanhi: []const f32 = &.{ |
| 47 | 4.6364760399e-01, // atan(0.5)hi 0x3eed6338 |
| 48 | 7.8539812565e-01, // atan(1.0)hi 0x3f490fda |
| 49 | 9.8279368877e-01, // atan(1.5)hi 0x3f7b985e |
| 50 | 1.5707962513e+00, // atan(inf)hi 0x3fc90fda |
| 51 | }; |
| 52 | const aT: []const f32 = &.{ |
| 53 | 0x1.fffcccp-1, |
| 54 | -0x1.52e8ccp-2, |
| 55 | 0x1.522336p-3, |
| 56 | }; |
| 57 | |
| 58 | const hx: u16 = @bitCast(x); |
| 59 | const ix = hx & 0x7fff; |
| 60 | const sign = (hx >> 15) != 0; |
| 61 | // if |x| >= 2^11 |
| 62 | if (ix >= 0x6800) { |
| 63 | if (math.isNan(x)) { |
| 64 | return x; |
| 65 | } |
| 66 | const z = atanhi[3] + 0x1p-120; |
| 67 | return @floatCast(if (sign) -z else z); |
| 68 | } |
| 69 | const x_: f32, const id: ?usize = blk: { |
| 70 | // |x| < 0.4375 |
| 71 | if (ix < 0x3700) { |
| 72 | // |x| < 2^(-6) |
| 73 | if (ix < 0x2400) { |
| 74 | if (ix < 0x400) { |
| 75 | // raise underflow for subnormal x |
| 76 | mem.doNotOptimizeAway(x * x); |
| 77 | } |
| 78 | return x; |
| 79 | } |
| 80 | break :blk .{ @floatCast(x), null }; |
| 81 | } else { |
| 82 | const x_: f32 = @floatCast(@abs(x)); |
| 83 | // |x| < 1.1875 |
| 84 | if (ix < 0x3cc0) { |
| 85 | // 7/16 <= |x| < 11/16 |
| 86 | if (ix < 0x3980) { |
| 87 | break :blk .{ (2.0 * x_ - 1.0) / (2.0 + x_), 0 }; |
| 88 | } |
| 89 | // 11/16 <= |x| < 19/16 |
| 90 | else { |
| 91 | break :blk .{ (x_ - 1.0) / (x_ + 1.0), 1 }; |
| 92 | } |
| 93 | } else { |
| 94 | // |x| < 2.4375 |
| 95 | if (ix < 0x40e0) { |
| 96 | break :blk .{ (x_ - 1.5) / (1.0 + 1.5 * x_), 2 }; |
| 97 | } |
| 98 | // 2.4375 <= |x| < 2^11 |
| 99 | else { |
| 100 | break :blk .{ -1.0 / x_, 3 }; |
| 101 | } |
| 102 | } |
| 103 | } |
| 104 | }; |
| 105 | // end of argument reduction |
| 106 | const z = x_ * x_; |
| 107 | const s = aT[0] + z * (aT[1] + z * aT[2]); |
| 108 | if (id) |id_| { |
| 109 | const z_ = atanhi[id_] + x_ * s; |
| 110 | return @floatCast(if (sign) -z_ else z_); |
| 111 | } else { |
| 112 | return @floatCast(x_ * s); |
| 113 | } |
| 114 | } |
| 115 | |
| 116 | fn atanBinary32(x: f32) f32 { |
| 117 | const atanhi: []const f32 = &.{ |
| 118 | 4.6364760399e-01, // atan(0.5)hi 0x3eed6338 |
| 119 | 7.8539812565e-01, // atan(1.0)hi 0x3f490fda |
| 120 | 9.8279368877e-01, // atan(1.5)hi 0x3f7b985e |
| 121 | 1.5707962513e+00, // atan(inf)hi 0x3fc90fda |
| 122 | }; |
| 123 | const atanlo: []const f32 = &.{ |
| 124 | 5.0121582440e-09, // atan(0.5)lo 0x31ac3769 |
| 125 | 3.7748947079e-08, // atan(1.0)lo 0x33222168 |
| 126 | 3.4473217170e-08, // atan(1.5)lo 0x33140fb4 |
| 127 | 7.5497894159e-08, // atan(inf)lo 0x33a22168 |
| 128 | }; |
| 129 | const aT: []const f32 = &.{ |
| 130 | 3.3333328366e-01, |
| 131 | -1.9999158382e-01, |
| 132 | 1.4253635705e-01, |
| 133 | -1.0648017377e-01, |
| 134 | 6.1687607318e-02, |
| 135 | }; |
| 136 | |
| 137 | const hx: u32 = @bitCast(x); |
| 138 | const ix = hx & 0x7fff_ffff; |
| 139 | const sign = (hx >> 31) != 0; |
| 140 | // if |x| >= 2^26 |
| 141 | if (ix >= 0x4c80_0000) { |
| 142 | if (math.isNan(x)) { |
| 143 | return x; |
| 144 | } |
| 145 | const z = atanhi[3] + 0x1p-120; |
| 146 | return if (sign) -z else z; |
| 147 | } |
| 148 | const x_, const id: ?usize = blk: { |
| 149 | // |x| < 0.4375 |
| 150 | if (ix < 0x3ee00000) { |
| 151 | // |x| < 2^(-12) |
| 152 | if (ix < 0x39800000) { |
| 153 | if (ix < 0x00800000) { |
| 154 | // raise underflow for subnormal x |
| 155 | mem.doNotOptimizeAway(x * x); |
| 156 | } |
| 157 | return x; |
| 158 | } |
| 159 | break :blk .{ x, null }; |
| 160 | } else { |
| 161 | const x_ = @abs(x); |
| 162 | // |x| < 1.1875 |
| 163 | if (ix < 0x3f98_0000) { |
| 164 | // 7/16 <= |x| < 11/16 |
| 165 | if (ix < 0x3f30_0000) { |
| 166 | break :blk .{ (2.0 * x_ - 1.0) / (2.0 + x_), 0 }; |
| 167 | } |
| 168 | // 11/16 <= |x| < 19/16 |
| 169 | else { |
| 170 | break :blk .{ (x_ - 1.0) / (x_ + 1.0), 1 }; |
| 171 | } |
| 172 | } else { |
| 173 | // |x| < 2.4375 |
| 174 | if (ix < 0x401c_0000) { |
| 175 | break :blk .{ (x_ - 1.5) / (1.0 + 1.5 * x_), 2 }; |
| 176 | } |
| 177 | // 2.4375 <= |x| < 2^26 |
| 178 | else { |
| 179 | break :blk .{ -1.0 / x_, 3 }; |
| 180 | } |
| 181 | } |
| 182 | } |
| 183 | }; |
| 184 | // end of argument reduction |
| 185 | const z = x_ * x_; |
| 186 | const w = z * z; |
| 187 | // break sum from i=0 to 10 aT[i]z^(i+1) into odd and even poly |
| 188 | const s1 = z * (aT[0] + w * (aT[2] + w * aT[4])); |
| 189 | const s2 = w * (aT[1] + w * aT[3]); |
| 190 | if (id) |id_| { |
| 191 | const z_ = atanhi[id_] - ((x_ * (s1 + s2) - atanlo[id_]) - x_); |
| 192 | return if (sign) -z_ else z_; |
| 193 | } else { |
| 194 | return x_ - x_ * (s1 + s2); |
| 195 | } |
| 196 | } |
| 197 | |
| 198 | fn atanBinary64(x: f64) f64 { |
| 199 | const atanhi: []const f64 = &.{ |
| 200 | 4.63647609000806093515e-01, // atan(0.5)hi 0x3FDDAC67, 0x0561BB4F |
| 201 | 7.85398163397448278999e-01, // atan(1.0)hi 0x3FE921FB, 0x54442D18 |
| 202 | 9.82793723247329054082e-01, // atan(1.5)hi 0x3FEF730B, 0xD281F69B |
| 203 | 1.57079632679489655800e+00, // atan(inf)hi 0x3FF921FB, 0x54442D18 |
| 204 | }; |
| 205 | const atanlo: []const f64 = &.{ |
| 206 | 2.26987774529616870924e-17, // atan(0.5)lo 0x3C7A2B7F, 0x222F65E2 |
| 207 | 3.06161699786838301793e-17, // atan(1.0)lo 0x3C81A626, 0x33145C07 |
| 208 | 1.39033110312309984516e-17, // atan(1.5)lo 0x3C700788, 0x7AF0CBBD |
| 209 | 6.12323399573676603587e-17, // atan(inf)lo 0x3C91A626, 0x33145C07 |
| 210 | }; |
| 211 | const aT: []const f64 = &.{ |
| 212 | 3.33333333333329318027e-01, // 0x3FD55555, 0x5555550D |
| 213 | -1.99999999998764832476e-01, // 0xBFC99999, 0x9998EBC4 |
| 214 | 1.42857142725034663711e-01, // 0x3FC24924, 0x920083FF |
| 215 | -1.11111104054623557880e-01, // 0xBFBC71C6, 0xFE231671 |
| 216 | 9.09088713343650656196e-02, // 0x3FB745CD, 0xC54C206E |
| 217 | -7.69187620504482999495e-02, // 0xBFB3B0F2, 0xAF749A6D |
| 218 | 6.66107313738753120669e-02, // 0x3FB10D66, 0xA0D03D51 |
| 219 | -5.83357013379057348645e-02, // 0xBFADDE2D, 0x52DEFD9A |
| 220 | 4.97687799461593236017e-02, // 0x3FA97B4B, 0x24760DEB |
| 221 | -3.65315727442169155270e-02, // 0xBFA2B444, 0x2C6A6C2F |
| 222 | 1.62858201153657823623e-02, // 0x3F90AD3A, 0xE322DA11 |
| 223 | }; |
| 224 | |
| 225 | const hx: u64 = @bitCast(x); |
| 226 | const ix: u32 = @truncate((hx >> 32) & 0x7fffffff); |
| 227 | const sign = (hx >> 63) != 0; |
| 228 | // if |x| >= 2^66 |
| 229 | if (ix >= 0x44100000) { |
| 230 | if (math.isNan(x)) { |
| 231 | return x; |
| 232 | } |
| 233 | const z = atanhi[3] + 0x1p-120; |
| 234 | return if (sign) -z else z; |
| 235 | } |
| 236 | const x_, const id: ?usize = blk: { |
| 237 | // |x| < 0.4375 |
| 238 | if (ix < 0x3fdc_0000) { |
| 239 | // |x| < 2^(-27) |
| 240 | if (ix < 0x3e40_0000) { |
| 241 | if (ix < 0x0010_0000) { |
| 242 | // raise underflow for subnormal x |
| 243 | mem.doNotOptimizeAway(@as(f32, @floatCast(x))); |
| 244 | } |
| 245 | return x; |
| 246 | } |
| 247 | break :blk .{ x, null }; |
| 248 | } else { |
| 249 | const x_ = @abs(x); |
| 250 | // |x| < 1.1875 |
| 251 | if (ix < 0x3ff3_0000) { |
| 252 | // 7/16 <= |x| < 11/16 |
| 253 | if (ix < 0x3fe6_0000) { |
| 254 | break :blk .{ (2.0 * x_ - 1.0) / (2.0 + x_), 0 }; |
| 255 | } |
| 256 | // 11/16 <= |x| < 19/16 |
| 257 | else { |
| 258 | break :blk .{ (x_ - 1.0) / (x_ + 1.0), 1 }; |
| 259 | } |
| 260 | } else { |
| 261 | // |x| < 2.4375 |
| 262 | if (ix < 0x4003_8000) { |
| 263 | break :blk .{ (x_ - 1.5) / (1.0 + 1.5 * x_), 2 }; |
| 264 | } |
| 265 | // 2.4375 <= |x| < 2^66 |
| 266 | else { |
| 267 | break :blk .{ -1.0 / x_, 3 }; |
| 268 | } |
| 269 | } |
| 270 | } |
| 271 | }; |
| 272 | // end of argument reduction |
| 273 | const z = x_ * x_; |
| 274 | const w = z * z; |
| 275 | // break sum from i=0 to 10 aT[i]z^(i+1) into odd and even poly |
| 276 | const s1 = z * (aT[0] + w * (aT[2] + w * (aT[4] + w * (aT[6] + w * (aT[8] + w * aT[10]))))); |
| 277 | const s2 = w * (aT[1] + w * (aT[3] + w * (aT[5] + w * (aT[7] + w * aT[9])))); |
| 278 | if (id) |id_| { |
| 279 | const z_ = atanhi[id_] - (x_ * (s1 + s2) - atanlo[id_] - x_); |
| 280 | return if (sign) -z_ else z_; |
| 281 | } else { |
| 282 | return x_ - x_ * (s1 + s2); |
| 283 | } |
| 284 | } |
| 285 | |
| 286 | fn atanExtended80(x: f80) f80 { |
| 287 | const atanhi: []const f80 = &.{ |
| 288 | 4.63647609000806116202e-01, |
| 289 | 7.85398163397448309628e-01, |
| 290 | 9.82793723247329067960e-01, |
| 291 | 1.57079632679489661926e+00, |
| 292 | }; |
| 293 | const atanlo: []const f80 = &.{ |
| 294 | 1.18469937025062860669e-20, |
| 295 | -1.25413940316708300586e-20, |
| 296 | 2.55232234165405176172e-20, |
| 297 | -2.50827880633416601173e-20, |
| 298 | }; |
| 299 | const aT: []const f80 = &.{ |
| 300 | 3.33333333333333333017e-01, |
| 301 | -1.99999999999999632011e-01, |
| 302 | 1.42857142857046531280e-01, |
| 303 | -1.11111111100562372733e-01, |
| 304 | 9.09090902935647302252e-02, |
| 305 | -7.69230552476207730353e-02, |
| 306 | 6.66661718042406260546e-02, |
| 307 | -5.88158892835030888692e-02, |
| 308 | 5.25499891539726639379e-02, |
| 309 | -4.70119845393155721494e-02, |
| 310 | 4.03539201366454414072e-02, |
| 311 | -2.91303858419364158725e-02, |
| 312 | 1.24822046299269234080e-02, |
| 313 | }; |
| 314 | |
| 315 | const hx: u80 = @bitCast(x); |
| 316 | const se: u16 = @truncate(hx >> 64); |
| 317 | const e = se & 0x7fff; |
| 318 | const sign = se >> 15 != 0; |
| 319 | // if |x| is large, atan(x)~=pi/2 |
| 320 | if (e >= 0x3fff + math.floatMantissaBits(f80) + 1) { |
| 321 | if (math.isNan(x)) { |
| 322 | return x; |
| 323 | } |
| 324 | return if (sign) -atanhi[3] else atanhi[3]; |
| 325 | } |
| 326 | // Extract the exponent and the first few bits of the mantissa. |
| 327 | const m: u64 = @truncate(hx & 0x0000_ffff_ffff_ffff_ffff); |
| 328 | const expman = ((@as(u32, @intCast(se)) & 0x7fff) << 8) | (@as(u32, @truncate(m >> 55)) & 0xff); |
| 329 | const x_, const id: ?usize = blk: { |
| 330 | // |x| < 0.4375 |
| 331 | if (expman < ((0x3fff - 2) << 8) + 0xc0) { |
| 332 | // if |x| is small, atanl(x)~=x |
| 333 | if (e < 0x3fff - (math.floatMantissaBits(f80) + 1) / 2) { |
| 334 | // raise underflow if subnormal |
| 335 | if (e == 0) { |
| 336 | std.mem.doNotOptimizeAway(@as(f32, @floatCast(x))); |
| 337 | } |
| 338 | return x; |
| 339 | } |
| 340 | break :blk .{ x, null }; |
| 341 | } else { |
| 342 | const x_ = @abs(x); |
| 343 | // |x| < 1.1875 |
| 344 | if (expman < (0x3fff << 8) + 0x30) { |
| 345 | // 7/16 <= |x| < 11/16 |
| 346 | if (expman < ((0x3fff - 1) << 8) + 0x60) { |
| 347 | break :blk .{ (2.0 * x_ - 1.0) / (2.0 + x_), 0 }; |
| 348 | } |
| 349 | // 11/16 <= |x| < 19/16 |
| 350 | else { |
| 351 | break :blk .{ (x_ - 1.0) / (x_ + 1.0), 1 }; |
| 352 | } |
| 353 | } else { |
| 354 | // |x| < 2.4375 |
| 355 | if (expman < ((0x3fff + 1) << 8) + 0x38) { |
| 356 | break :blk .{ (x_ - 1.5) / (1.0 + 1.5 * x_), 2 }; |
| 357 | } |
| 358 | // 2.4375 <= |x| |
| 359 | else { |
| 360 | break :blk .{ -1.0 / x_, 3 }; |
| 361 | } |
| 362 | } |
| 363 | } |
| 364 | }; |
| 365 | // end of argument reduction |
| 366 | const z = x_ * x_; |
| 367 | const w = z * z; |
| 368 | // break sum aT[i]z^(i+1) into odd and even poly |
| 369 | const s1 = z * (aT[0] + w * (aT[2] + w * (aT[4] + w * (aT[6] + w * (aT[8] + w * (aT[10] + w * aT[12])))))); |
| 370 | const s2 = w * (aT[1] + w * (aT[3] + w * (aT[5] + w * (aT[7] + w * (aT[9] + w * aT[11]))))); |
| 371 | if (id) |id_| { |
| 372 | const z_ = atanhi[id_] - ((x_ * (s1 + s2) - atanlo[id_]) - x_); |
| 373 | return if (sign) -z_ else z_; |
| 374 | } else { |
| 375 | return x_ - x_ * (s1 + s2); |
| 376 | } |
| 377 | } |
| 378 | |
| 379 | fn atanBinary128(x: f128) f128 { |
| 380 | const atanhi: []const f128 = &.{ |
| 381 | 4.63647609000806116214256231461214397e-01, |
| 382 | 7.85398163397448309615660845819875699e-01, |
| 383 | 9.82793723247329067985710611014666038e-01, |
| 384 | 1.57079632679489661923132169163975140e+00, |
| 385 | }; |
| 386 | const atanlo: []const f128 = &.{ |
| 387 | 4.89509642257333492668618435220297706e-36, |
| 388 | 2.16795253253094525619926100651083806e-35, |
| 389 | -2.31288434538183565909319952098066272e-35, |
| 390 | 4.33590506506189051239852201302167613e-35, |
| 391 | }; |
| 392 | const aT: []const f128 = &.{ |
| 393 | 3.33333333333333333333333333333333125e-01, |
| 394 | -1.99999999999999999999999999999180430e-01, |
| 395 | 1.42857142857142857142857142125269827e-01, |
| 396 | -1.11111111111111111111110834490810169e-01, |
| 397 | 9.09090909090909090908522355708623681e-02, |
| 398 | -7.69230769230769230696553844935357021e-02, |
| 399 | 6.66666666666666660390096773046256096e-02, |
| 400 | -5.88235294117646671706582985209643694e-02, |
| 401 | 5.26315789473666478515847092020327506e-02, |
| 402 | -4.76190476189855517021024424991436144e-02, |
| 403 | 4.34782608678695085948531993458097026e-02, |
| 404 | -3.99999999632663469330634215991142368e-02, |
| 405 | 3.70370363987423702891250829918659723e-02, |
| 406 | -3.44827496515048090726669907612335954e-02, |
| 407 | 3.22579620681420149871973710852268528e-02, |
| 408 | -3.03020767654269261041647570626778067e-02, |
| 409 | 2.85641979882534783223403715930946138e-02, |
| 410 | -2.69824879726738568189929461383741323e-02, |
| 411 | 2.54194698498808542954187110873675769e-02, |
| 412 | -2.35083879708189059926183138130183215e-02, |
| 413 | 2.04832358998165364349957325067131428e-02, |
| 414 | -1.54489555488544397858507248612362957e-02, |
| 415 | 8.64492360989278761493037861575248038e-03, |
| 416 | -2.58521121597609872727919154569765469e-03, |
| 417 | }; |
| 418 | |
| 419 | const hx: u128 = @bitCast(x); |
| 420 | const se: u16 = @truncate(hx >> 112); |
| 421 | const e = se & 0x7fff; |
| 422 | const sign = se >> 15 != 0; |
| 423 | // if |x| is large, atan(x)~=pi/2 |
| 424 | if (e >= 0x3fff + math.floatMantissaBits(f128) + 2) { |
| 425 | if (math.isNan(x)) { |
| 426 | return x; |
| 427 | } |
| 428 | return if (sign) -atanhi[3] else atanhi[3]; |
| 429 | } |
| 430 | // Extract the exponent and the first few bits of the mantissa. |
| 431 | const top: u16 = @truncate((hx >> 96) & 0x0000_ffff); |
| 432 | const expman = ((@as(u32, @intCast(se)) & 0x7fff) << 8) | (@as(u32, @intCast(top)) >> 8); |
| 433 | const x_, const id: ?usize = blk: { |
| 434 | // |x| < 0.4375 |
| 435 | if (expman < ((0x3fff - 2) << 8) + 0xc0) { |
| 436 | // if |x| is small, atanl(x)~=x |
| 437 | if (e < 0x3fff - (math.floatMantissaBits(f128) + 2) / 2) { |
| 438 | // raise underflow if subnormal |
| 439 | if (e == 0) { |
| 440 | mem.doNotOptimizeAway(@as(f32, @floatCast(x))); |
| 441 | } |
| 442 | return x; |
| 443 | } |
| 444 | break :blk .{ x, null }; |
| 445 | } else { |
| 446 | const x_ = @abs(x); |
| 447 | // |x| < 1.1875 |
| 448 | if (expman < (0x3fff << 8) + 0x30) { |
| 449 | // 7/16 <= |x| < 11/16 |
| 450 | if (expman < ((0x3fff - 1) << 8) + 0x60) { |
| 451 | break :blk .{ (2.0 * x_ - 1.0) / (2.0 + x_), 0 }; |
| 452 | } |
| 453 | // 11/16 <= |x| < 19/16 |
| 454 | else { |
| 455 | break :blk .{ (x_ - 1.0) / (x_ + 1.0), 1 }; |
| 456 | } |
| 457 | } else { |
| 458 | // |x| < 2.4375 |
| 459 | if (expman < ((0x3fff + 1) << 8) + 0x38) { |
| 460 | break :blk .{ (x_ - 1.5) / (1.0 + 1.5 * x_), 2 }; |
| 461 | } |
| 462 | // 2.4375 <= |x| |
| 463 | else { |
| 464 | break :blk .{ -1.0 / x_, 3 }; |
| 465 | } |
| 466 | } |
| 467 | } |
| 468 | }; |
| 469 | // end of argument reduction |
| 470 | const z = x_ * x_; |
| 471 | const w = z * z; |
| 472 | // break sum aT[i]z^(i+1) into odd and even poly |
| 473 | const s1 = z * (aT[0] + w * (aT[2] + w * (aT[4] + w * (aT[6] + w * (aT[8] + w * (aT[10] + w * (aT[12] + w * (aT[14] + w * (aT[16] + w * (aT[18] + w * (aT[20] + w * aT[22]))))))))))); |
| 474 | const s2 = w * (aT[1] + w * (aT[3] + w * (aT[5] + w * (aT[7] + w * (aT[9] + w * (aT[11] + w * (aT[13] + w * (aT[15] + w * (aT[17] + w * (aT[19] + w * (aT[21] + w * aT[23]))))))))))); |
| 475 | if (id) |id_| { |
| 476 | const z_ = atanhi[id_] - ((x_ * (s1 + s2) - atanlo[id_]) - x_); |
| 477 | return if (sign) -z_ else z_; |
| 478 | } else { |
| 479 | return x_ - x_ * (s1 + s2); |
| 480 | } |
| 481 | } |
| 482 | |
| 483 | test "atanBinary16.special" { |
| 484 | try testing.expectEqual(0x0p+0, atanBinary16(0x0p+0)); |
| 485 | try testing.expectEqual(-0x0p+0, atanBinary16(-0x0p+0)); |
| 486 | try testing.expectApproxEqAbs(0x1.92p-1, atanBinary16(0x1p+0), math.floatEpsAt(f16, 0x1.92p-1)); |
| 487 | try testing.expectApproxEqAbs(-0x1.92p-1, atanBinary16(-0x1p+0), math.floatEpsAt(f16, -0x1.92p-1)); |
| 488 | try testing.expectApproxEqAbs(0x1.92p0, atanBinary16(math.inf(f16)), math.floatEpsAt(f16, 0x1.92p0)); |
| 489 | try testing.expectApproxEqAbs(-0x1.92p0, atanBinary16(-math.inf(f16)), math.floatEpsAt(f16, -0x1.92p0)); |
| 490 | try testing.expect(math.isNan(atanBinary16(math.nan(f16)))); |
| 491 | } |
| 492 | |
| 493 | test "atanBinary16" { |
| 494 | try testing.expectApproxEqAbs(-0x1.74cp-2, atanBinary16(-0x1.864p-2), math.floatEpsAt(f16, -0x1.74cp-2)); |
| 495 | try testing.expectApproxEqAbs(-0x1.374p0, atanBinary16(-0x1.59cp1), math.floatEpsAt(f16, -0x1.374p0)); |
| 496 | try testing.expectApproxEqAbs(-0x1.11cp0, atanBinary16(-0x1.d2cp0), math.floatEpsAt(f16, -0x1.11cp0)); |
| 497 | try testing.expectApproxEqAbs(-0x1.33cp-1, atanBinary16(-0x1.5f4p-1), math.floatEpsAt(f16, -0x1.33cp-1)); |
| 498 | try testing.expectApproxEqAbs(0x1.37p0, atanBinary16(0x1.588p1), math.floatEpsAt(f16, 0x1.37p0)); |
| 499 | try testing.expectApproxEqAbs(-0x1.99cp-2, atanBinary16(-0x1.b14p-2), math.floatEpsAt(f16, -0x1.99cp-2)); |
| 500 | try testing.expectApproxEqAbs(0x1.2fcp0, atanBinary16(0x1.3ccp1), math.floatEpsAt(f16, 0x1.2fcp0)); |
| 501 | try testing.expectApproxEqAbs(-0x1.08cp-2, atanBinary16(-0x1.0ecp-2), math.floatEpsAt(f16, -0x1.08cp-2)); |
| 502 | try testing.expectApproxEqAbs(0x1.2ap0, atanBinary16(0x1.298p1), math.floatEpsAt(f16, 0x1.2ap0)); |
| 503 | try testing.expectApproxEqAbs(-0x1.1c8p0, atanBinary16(-0x1.028p1), math.floatEpsAt(f16, -0x1.1c8p0)); |
| 504 | } |
| 505 | |
| 506 | test "atanBinary32.special" { |
| 507 | try testing.expectEqual(0x0p+0, atanBinary32(0x0p+0)); |
| 508 | try testing.expectEqual(-0x0p+0, atanBinary32(-0x0p+0)); |
| 509 | try testing.expectApproxEqAbs(0x1.921fb6p-1, atanBinary32(0x1p+0), math.floatEpsAt(f32, 0x1.921fb6p-1)); |
| 510 | try testing.expectApproxEqAbs(-0x1.921fb6p-1, atanBinary32(-0x1p+0), math.floatEpsAt(f32, -0x1.921fb6p-1)); |
| 511 | try testing.expectApproxEqAbs(0x1.921fb6p+0, atanBinary32(math.inf(f32)), math.floatEpsAt(f32, 0x1.921fb6p+0)); |
| 512 | try testing.expectApproxEqAbs(-0x1.921fb6p+0, atanBinary32(-math.inf(f32)), math.floatEpsAt(f32, -0x1.921fb6p+0)); |
| 513 | try testing.expect(math.isNan(atanBinary32(math.nan(f32)))); |
| 514 | } |
| 515 | |
| 516 | test "atanBinary32" { |
| 517 | try testing.expectApproxEqAbs(-0x1.74c62p-2, atanBinary32(-0x1.8629dp-2), math.floatEpsAt(f32, -0x1.74c62p-2)); |
| 518 | try testing.expectApproxEqAbs(-0x1.375fd8p0, atanBinary32(-0x1.59d42ep1), math.floatEpsAt(f32, -0x1.375fd8p0)); |
| 519 | try testing.expectApproxEqAbs(-0x1.11b8aep0, atanBinary32(-0x1.d2dbe2p0), math.floatEpsAt(f32, -0x1.11b8aep0)); |
| 520 | try testing.expectApproxEqAbs(-0x1.33d28cp-1, atanBinary32(-0x1.5f314ep-1), math.floatEpsAt(f32, -0x1.33d28cp-1)); |
| 521 | try testing.expectApproxEqAbs(0x1.37082ep0, atanBinary32(0x1.5869bp1), math.floatEpsAt(f32, 0x1.37082ep0)); |
| 522 | try testing.expectApproxEqAbs(-0x1.99d7cap-2, atanBinary32(-0x1.b13a06p-2), math.floatEpsAt(f32, -0x1.99d7cap-2)); |
| 523 | try testing.expectApproxEqAbs(0x1.2fcb12p0, atanBinary32(0x1.3cb0f2p1), math.floatEpsAt(f32, 0x1.2fcb12p0)); |
| 524 | try testing.expectApproxEqAbs(-0x1.08c71ap-2, atanBinary32(-0x1.0ed746p-2), math.floatEpsAt(f32, -0x1.08c71ap-2)); |
| 525 | try testing.expectApproxEqAbs(0x1.2a24e2p0, atanBinary32(0x1.299d54p1), math.floatEpsAt(f32, 0x1.2a24e2p0)); |
| 526 | try testing.expectApproxEqAbs(-0x1.1c6178p0, atanBinary32(-0x1.0264fcp1), math.floatEpsAt(f32, -0x1.1c6178p0)); |
| 527 | } |
| 528 | |
| 529 | test "atanBinary64.special" { |
| 530 | try testing.expectEqual(0x0p+0, atanBinary64(0x0p+0)); |
| 531 | try testing.expectEqual(-0x0p+0, atanBinary64(-0x0p+0)); |
| 532 | try testing.expectApproxEqAbs(0x1.921fb54442d18p-1, atanBinary64(0x1p+0), math.floatEpsAt(f64, 0x1.921fb54442d18p-1)); |
| 533 | try testing.expectApproxEqAbs(-0x1.921fb54442d18p-1, atanBinary64(-0x1p+0), math.floatEpsAt(f64, -0x1.921fb54442d18p-1)); |
| 534 | try testing.expectApproxEqAbs(0x1.921fb54442d18p+0, atanBinary64(math.inf(f64)), math.floatEpsAt(f64, 0x1.921fb54442d18p+0)); |
| 535 | try testing.expectApproxEqAbs(-0x1.921fb54442d18p+0, atanBinary64(-math.inf(f64)), math.floatEpsAt(f64, -0x1.921fb54442d18p+0)); |
| 536 | try testing.expect(math.isNan(atanBinary64(math.nan(f64)))); |
| 537 | } |
| 538 | |
| 539 | test "atanBinary64" { |
| 540 | try testing.expectApproxEqAbs(-0x1.74c61f4377016p-2, atanBinary64(-0x1.8629d0244cdccp-2), math.floatEpsAt(f64, -0x1.74c61f4377016p-2)); |
| 541 | try testing.expectApproxEqAbs(-0x1.375fd7987cc2p0, atanBinary64(-0x1.59d42d4659937p1), math.floatEpsAt(f64, -0x1.375fd7987cc2p0)); |
| 542 | try testing.expectApproxEqAbs(-0x1.11b8adeba5616p0, atanBinary64(-0x1.d2dbe23d04f06p0), math.floatEpsAt(f64, -0x1.11b8adeba5616p0)); |
| 543 | try testing.expectApproxEqAbs(-0x1.33d28ca762539p-1, atanBinary64(-0x1.5f314e72398e8p-1), math.floatEpsAt(f64, -0x1.33d28ca762539p-1)); |
| 544 | try testing.expectApproxEqAbs(0x1.37082ce2dd03p0, atanBinary64(0x1.5869af37b7d08p1), math.floatEpsAt(f64, 0x1.37082ce2dd03p0)); |
| 545 | try testing.expectApproxEqAbs(-0x1.99d7cac66dd44p-2, atanBinary64(-0x1.b13a05a662618p-2), math.floatEpsAt(f64, -0x1.99d7cac66dd44p-2)); |
| 546 | try testing.expectApproxEqAbs(0x1.2fcb120468e8ep0, atanBinary64(0x1.3cb0f12f39d8ap1), math.floatEpsAt(f64, 0x1.2fcb120468e8ep0)); |
| 547 | try testing.expectApproxEqAbs(-0x1.08c71aa0e509p-2, atanBinary64(-0x1.0ed746b39cbb7p-2), math.floatEpsAt(f64, -0x1.08c71aa0e509p-2)); |
| 548 | try testing.expectApproxEqAbs(0x1.2a24e22d861dfp0, atanBinary64(0x1.299d54ac7d6bp1), math.floatEpsAt(f64, 0x1.2a24e22d861dfp0)); |
| 549 | try testing.expectApproxEqAbs(-0x1.1c617825f9751p0, atanBinary64(-0x1.0264fb9f3d50ep1), math.floatEpsAt(f64, -0x1.1c617825f9751p0)); |
| 550 | } |
| 551 | |
| 552 | test "atanExtended80.special" { |
| 553 | try testing.expectEqual(0x0p+0, atanExtended80(0x0p+0)); |
| 554 | try testing.expectEqual(-0x0p+0, atanExtended80(-0x0p+0)); |
| 555 | try testing.expectApproxEqAbs(0x1.921fb54442d1846ap-1, atanExtended80(0x1p+0), math.floatEpsAt(f80, 0x1.921fb54442d1846ap-1)); |
| 556 | try testing.expectApproxEqAbs(-0x1.921fb54442d1846ap-1, atanExtended80(-0x1p+0), math.floatEpsAt(f80, -0x1.921fb54442d1846ap-1)); |
| 557 | try testing.expectApproxEqAbs(0x1.921fb54442d1846ap0, atanExtended80(math.inf(f80)), math.floatEpsAt(f80, 0x1.921fb54442d1846ap0)); |
| 558 | try testing.expectApproxEqAbs(-0x1.921fb54442d1846ap0, atanExtended80(-math.inf(f80)), math.floatEpsAt(f80, -0x1.921fb54442d1846ap0)); |
| 559 | try testing.expect(math.isNan(atanExtended80(math.nan(f80)))); |
| 560 | } |
| 561 | |
| 562 | test "atanExtended80" { |
| 563 | try testing.expectApproxEqAbs(-0x1.74c61f437701661p-2, atanExtended80(-0x1.8629d0244cdcbed8p-2), math.floatEpsAt(f80, -0x1.74c61f437701661p-2)); |
| 564 | try testing.expectApproxEqAbs(-0x1.375fd7987cc1fd02p0, atanExtended80(-0x1.59d42d4659936d9ep1), math.floatEpsAt(f80, -0x1.375fd7987cc1fd02p0)); |
| 565 | try testing.expectApproxEqAbs(-0x1.11b8adeba5615e04p0, atanExtended80(-0x1.d2dbe23d04f067b4p0), math.floatEpsAt(f80, -0x1.11b8adeba5615e04p0)); |
| 566 | try testing.expectApproxEqAbs(-0x1.33d28ca76253964cp-1, atanExtended80(-0x1.5f314e72398e7dbcp-1), math.floatEpsAt(f80, -0x1.33d28ca76253964cp-1)); |
| 567 | try testing.expectApproxEqAbs(0x1.37082ce2dd03010cp0, atanExtended80(0x1.5869af37b7d078cap1), math.floatEpsAt(f80, 0x1.37082ce2dd03010cp0)); |
| 568 | try testing.expectApproxEqAbs(-0x1.99d7cac66dd4438p-2, atanExtended80(-0x1.b13a05a66261821ap-2), math.floatEpsAt(f80, -0x1.99d7cac66dd4438p-2)); |
| 569 | try testing.expectApproxEqAbs(0x1.2fcb120468e8d9ecp0, atanExtended80(0x1.3cb0f12f39d899cp1), math.floatEpsAt(f80, 0x1.2fcb120468e8d9ecp0)); |
| 570 | try testing.expectApproxEqAbs(-0x1.08c71aa0e5090998p-2, atanExtended80(-0x1.0ed746b39cbb7614p-2), math.floatEpsAt(f80, -0x1.08c71aa0e5090998p-2)); |
| 571 | try testing.expectApproxEqAbs(0x1.2a24e22d861debfep0, atanExtended80(0x1.299d54ac7d6afc52p1), math.floatEpsAt(f80, 0x1.2a24e22d861debfep0)); |
| 572 | try testing.expectApproxEqAbs(-0x1.1c617825f97512b8p0, atanExtended80(-0x1.0264fb9f3d50e4fp1), math.floatEpsAt(f80, -0x1.1c617825f97512b8p0)); |
| 573 | } |
| 574 | |
| 575 | test "atanBinary128.special" { |
| 576 | try testing.expectEqual(0x0p+0, atanBinary128(0x0p+0)); |
| 577 | try testing.expectEqual(-0x0p+0, atanBinary128(-0x0p+0)); |
| 578 | try testing.expectApproxEqAbs(0x1.921fb54442d18469898cc51701b8p-1, atanBinary128(0x1p+0), math.floatEpsAt(f128, 0x1.921fb54442d18469898cc51701b8p-1)); |
| 579 | try testing.expectApproxEqAbs(-0x1.921fb54442d18469898cc51701b8p-1, atanBinary128(-0x1p+0), math.floatEpsAt(f128, -0x1.921fb54442d18469898cc51701b8p-1)); |
| 580 | try testing.expectApproxEqAbs(0x1.921fb54442d18469898cc51701b8p0, atanBinary128(math.inf(f128)), math.floatEpsAt(f128, 0x1.921fb54442d18469898cc51701b8p0)); |
| 581 | try testing.expectApproxEqAbs(-0x1.921fb54442d18469898cc51701b8p0, atanBinary128(-math.inf(f128)), math.floatEpsAt(f128, -0x1.921fb54442d18469898cc51701b8p0)); |
| 582 | try testing.expect(math.isNan(atanBinary128(math.nan(f128)))); |
| 583 | } |
| 584 | |
| 585 | test "atanBinary128" { |
| 586 | try testing.expectApproxEqAbs(-0x1.74c61f437701660ff76989d23707p-2, atanBinary128(-0x1.8629d0244cdcbed71792ccdec26dp-2), math.floatEpsAt(f128, -0x1.74c61f437701660ff76989d23707p-2)); |
| 587 | try testing.expectApproxEqAbs(-0x1.375fd7987cc1fd0119cf0cc5b708p0, atanBinary128(-0x1.59d42d4659936d9e22b5dea4faefp1), math.floatEpsAt(f128, -0x1.375fd7987cc1fd0119cf0cc5b708p0)); |
| 588 | try testing.expectApproxEqAbs(-0x1.11b8adeba5615e0370722b511231p0, atanBinary128(-0x1.d2dbe23d04f067b42da3f8efdf57p0), math.floatEpsAt(f128, -0x1.11b8adeba5615e0370722b511231p0)); |
| 589 | try testing.expectApproxEqAbs(-0x1.33d28ca76253964cb5d3581cdd88p-1, atanBinary128(-0x1.5f314e72398e7dbbe70fb072983ep-1), math.floatEpsAt(f128, -0x1.33d28ca76253964cb5d3581cdd88p-1)); |
| 590 | try testing.expectApproxEqAbs(0x1.37082ce2dd03010bbea814dc5882p0, atanBinary128(0x1.5869af37b7d078caa3456c44aecep1), math.floatEpsAt(f128, 0x1.37082ce2dd03010bbea814dc5882p0)); |
| 591 | try testing.expectApproxEqAbs(-0x1.99d7cac66dd4438077284b491a91p-2, atanBinary128(-0x1.b13a05a66261821a364ad8c6c999p-2), math.floatEpsAt(f128, -0x1.99d7cac66dd4438077284b491a91p-2)); |
| 592 | try testing.expectApproxEqAbs(0x1.2fcb120468e8d9ebdb74702314c8p0, atanBinary128(0x1.3cb0f12f39d899c0d963ac413297p1), math.floatEpsAt(f128, 0x1.2fcb120468e8d9ebdb74702314c8p0)); |
| 593 | try testing.expectApproxEqAbs(-0x1.08c71aa0e5090998206fbbe2090fp-2, atanBinary128(-0x1.0ed746b39cbb7614d8735e8315a8p-2), math.floatEpsAt(f128, -0x1.08c71aa0e5090998206fbbe2090fp-2)); |
| 594 | try testing.expectApproxEqAbs(0x1.2a24e22d861debfd6f974500567fp0, atanBinary128(0x1.299d54ac7d6afc5154643b601519p1), math.floatEpsAt(f128, 0x1.2a24e22d861debfd6f974500567fp0)); |
| 595 | try testing.expectApproxEqAbs(-0x1.1c617825f97512b7f38656ab12cdp0, atanBinary128(-0x1.0264fb9f3d50e4f0f966f0686064p1), math.floatEpsAt(f128, -0x1.1c617825f97512b7f38656ab12cdp0)); |
| 596 | } |
| 597 | |
| 598 | fn atanBinary32Vec(comptime vec_len: comptime_int, x: @Vector(vec_len, f32)) @TypeOf(x) { |
| 599 | const sign_mask: @Vector(vec_len, u32) = @splat(0x80000000); |
| 600 | const neg_one: @Vector(vec_len, f32) = @splat(-1.0); |
| 601 | const pi_over_2: @Vector(vec_len, u32) = @splat(0x3fc90fdb); |
| 602 | const zero: @Vector(vec_len, u32) = @splat(0); |
| 603 | const c0: @Vector(vec_len, f32) = @splat(-0x1.5554dcp-2); |
| 604 | const c1: @Vector(vec_len, f32) = @splat(0x1.9978ecp-3); |
| 605 | const c2: @Vector(vec_len, f32) = @splat(-0x1.230a94p-3); |
| 606 | const c3: @Vector(vec_len, f32) = @splat(0x1.b4debp-4); |
| 607 | const c4: @Vector(vec_len, f32) = @splat(-0x1.3550dap-4); |
| 608 | const c5: @Vector(vec_len, f32) = @splat(0x1.61eebp-5); |
| 609 | const c6: @Vector(vec_len, f32) = @splat(-0x1.0c17d4p-6); |
| 610 | const c7: @Vector(vec_len, f32) = @splat(0x1.7ea694p-9); |
| 611 | |
| 612 | const ix: @Vector(vec_len, u32) = @bitCast(x); |
| 613 | const sign = ix & sign_mask; |
| 614 | const pred = @abs(x) > @abs(neg_one); |
| 615 | const z = @select(f32, pred, neg_one / x, x); |
| 616 | const shift: @Vector(vec_len, f32) = @bitCast(@select(u32, pred, pi_over_2 ^ sign, zero)); |
| 617 | const z2 = z * z; |
| 618 | const z3 = z * z2; |
| 619 | const z4 = z2 * z2; |
| 620 | const z8 = z4 * z4; |
| 621 | const p0_1 = @mulAdd(@Vector(vec_len, f32), z2, c1, c0); |
| 622 | const p2_3 = @mulAdd(@Vector(vec_len, f32), z2, c3, c2); |
| 623 | const p4_5 = @mulAdd(@Vector(vec_len, f32), z2, c5, c4); |
| 624 | const p6_7 = @mulAdd(@Vector(vec_len, f32), z2, c7, c6); |
| 625 | const p0_3 = @mulAdd(@Vector(vec_len, f32), z4, p2_3, p0_1); |
| 626 | const p4_7 = @mulAdd(@Vector(vec_len, f32), z4, p6_7, p4_5); |
| 627 | const p0_7 = @mulAdd(@Vector(vec_len, f32), z8, p4_7, p0_3); |
| 628 | return @mulAdd(@Vector(vec_len, f32), z3, p0_7, shift + z); |
| 629 | } |
| 630 | |
| 631 | fn atanBinary64Vec(comptime vec_len: comptime_int, x: @Vector(vec_len, f64)) @TypeOf(x) { |
| 632 | const sign_mask: @Vector(vec_len, u64) = @splat(0x8000000000000000); |
| 633 | const neg_one: @Vector(vec_len, f64) = @splat(-1.0); |
| 634 | const pi_over_2: @Vector(vec_len, u64) = @splat(0x3ff921fb54442d18); |
| 635 | const zero: @Vector(vec_len, u64) = @splat(0); |
| 636 | const c0: @Vector(vec_len, f64) = @splat(-0x1.555555555552ap-2); |
| 637 | const c1: @Vector(vec_len, f64) = @splat(0x1.9999999995aebp-3); |
| 638 | const c2: @Vector(vec_len, f64) = @splat(-0x1.24924923923f6p-3); |
| 639 | const c3: @Vector(vec_len, f64) = @splat(0x1.c71c7184288a2p-4); |
| 640 | const c4: @Vector(vec_len, f64) = @splat(-0x1.745d11fb3d32bp-4); |
| 641 | const c5: @Vector(vec_len, f64) = @splat(0x1.3b136a18051b9p-4); |
| 642 | const c6: @Vector(vec_len, f64) = @splat(-0x1.110e6d985f496p-4); |
| 643 | const c7: @Vector(vec_len, f64) = @splat(0x1.e1bcf7f08801dp-5); |
| 644 | const c8: @Vector(vec_len, f64) = @splat(-0x1.ae644e28058c3p-5); |
| 645 | const c9: @Vector(vec_len, f64) = @splat(0x1.82eeb1fed85c6p-5); |
| 646 | const c10: @Vector(vec_len, f64) = @splat(-0x1.59d7f901566cbp-5); |
| 647 | const c11: @Vector(vec_len, f64) = @splat(0x1.2c982855ab069p-5); |
| 648 | const c12: @Vector(vec_len, f64) = @splat(-0x1.eb49592998177p-6); |
| 649 | const c13: @Vector(vec_len, f64) = @splat(0x1.69d8b396e3d38p-6); |
| 650 | const c14: @Vector(vec_len, f64) = @splat(-0x1.ca980345c4204p-7); |
| 651 | const c15: @Vector(vec_len, f64) = @splat(0x1.dc050eafde0b3p-8); |
| 652 | const c16: @Vector(vec_len, f64) = @splat(-0x1.7ea70755b8eccp-9); |
| 653 | const c17: @Vector(vec_len, f64) = @splat(0x1.ba3da3de903e8p-11); |
| 654 | const c18: @Vector(vec_len, f64) = @splat(-0x1.44a4b059b6f67p-13); |
| 655 | const c19: @Vector(vec_len, f64) = @splat(0x1.c4a45029e5a91p-17); |
| 656 | |
| 657 | const ix: @Vector(vec_len, u64) = @bitCast(x); |
| 658 | const sign = ix & sign_mask; |
| 659 | const pred = @abs(x) > @abs(neg_one); |
| 660 | const shift: @Vector(vec_len, f64) = @bitCast(@select(u64, pred, pi_over_2 ^ sign, zero)); |
| 661 | const z = @select(f64, pred, neg_one / x, x); |
| 662 | const z2 = z * z; |
| 663 | const z3 = z * z2; |
| 664 | const z4 = z2 * z2; |
| 665 | const z8 = z4 * z4; |
| 666 | const z16 = z8 * z8; |
| 667 | const p0_1 = @mulAdd(@Vector(vec_len, f64), z2, c1, c0); |
| 668 | const p2_3 = @mulAdd(@Vector(vec_len, f64), z2, c3, c2); |
| 669 | const p0_3 = @mulAdd(@Vector(vec_len, f64), z4, p2_3, p0_1); |
| 670 | const p4_5 = @mulAdd(@Vector(vec_len, f64), z2, c5, c4); |
| 671 | const p6_7 = @mulAdd(@Vector(vec_len, f64), z2, c7, c6); |
| 672 | const p4_7 = @mulAdd(@Vector(vec_len, f64), z4, p6_7, p4_5); |
| 673 | const p0_7 = @mulAdd(@Vector(vec_len, f64), z8, p4_7, p0_3); |
| 674 | const p8_9 = @mulAdd(@Vector(vec_len, f64), z2, c9, c8); |
| 675 | const p10_11 = @mulAdd(@Vector(vec_len, f64), z2, c11, c10); |
| 676 | const p8_11 = @mulAdd(@Vector(vec_len, f64), z4, p10_11, p8_9); |
| 677 | const p12_13 = @mulAdd(@Vector(vec_len, f64), z2, c13, c12); |
| 678 | const p14_15 = @mulAdd(@Vector(vec_len, f64), z2, c15, c14); |
| 679 | const p12_15 = @mulAdd(@Vector(vec_len, f64), z4, p14_15, p12_13); |
| 680 | const p16_17 = @mulAdd(@Vector(vec_len, f64), z2, c17, c16); |
| 681 | const p18_19 = @mulAdd(@Vector(vec_len, f64), z2, c19, c18); |
| 682 | const p16_19 = @mulAdd(@Vector(vec_len, f64), z4, p18_19, p16_17); |
| 683 | const p8_15 = @mulAdd(@Vector(vec_len, f64), z8, p12_15, p8_11); |
| 684 | const p8_19 = @mulAdd(@Vector(vec_len, f64), z16, p16_19, p8_15); |
| 685 | const p0_19 = @mulAdd(@Vector(vec_len, f64), p8_19, z16, p0_7); |
| 686 | return @mulAdd(@Vector(vec_len, f64), z3, p0_19, shift + z); |
| 687 | } |
| 688 | |
| 689 | test "atanBinary32Vec.special" { |
| 690 | const input: @Vector(7, f32) = .{ |
| 691 | 0x0p+0, |
| 692 | -0x0p+0, |
| 693 | 0x1p+0, |
| 694 | -0x1p+0, |
| 695 | math.inf(f32), |
| 696 | -math.inf(f32), |
| 697 | math.nan(f32), |
| 698 | }; |
| 699 | const output = atanBinary32Vec(7, input); |
| 700 | try testing.expectEqual(0x0p+0, output[0]); |
| 701 | try testing.expectEqual(-0x0p+0, output[1]); |
| 702 | try testing.expectApproxEqAbs(0x1.921fb6p-1, output[2], math.floatEpsAt(f32, 0x1.921fb6p-1)); |
| 703 | try testing.expectApproxEqAbs(-0x1.921fb6p-1, output[3], math.floatEpsAt(f32, -0x1.921fb6p-1)); |
| 704 | try testing.expectApproxEqAbs(0x1.921fb6p+0, output[4], math.floatEpsAt(f32, 0x1.921fb6p+0)); |
| 705 | try testing.expectApproxEqAbs(-0x1.921fb6p+0, output[5], math.floatEpsAt(f32, -0x1.921fb6p+0)); |
| 706 | try testing.expect(math.isNan(output[6])); |
| 707 | } |
| 708 | |
| 709 | test "atanBinary32Vec" { |
| 710 | const input: @Vector(10, f32) = .{ |
| 711 | -0x1.8629dp-2, |
| 712 | -0x1.59d42ep1, |
| 713 | -0x1.d2dbe2p0, |
| 714 | -0x1.5f314ep-1, |
| 715 | 0x1.5869bp1, |
| 716 | -0x1.b13a06p-2, |
| 717 | 0x1.3cb0f2p1, |
| 718 | -0x1.0ed746p-2, |
| 719 | 0x1.299d54p1, |
| 720 | -0x1.0264fcp1, |
| 721 | }; |
| 722 | const output = atanBinary32Vec(10, input); |
| 723 | try testing.expectApproxEqAbs(-0x1.74c62p-2, output[0], math.floatEpsAt(f32, -0x1.74c62p-2)); |
| 724 | try testing.expectApproxEqAbs(-0x1.375fd8p0, output[1], math.floatEpsAt(f32, -0x1.375fd8p0)); |
| 725 | try testing.expectApproxEqAbs(-0x1.11b8aep0, output[2], math.floatEpsAt(f32, -0x1.11b8aep0)); |
| 726 | try testing.expectApproxEqAbs(-0x1.33d28cp-1, output[3], math.floatEpsAt(f32, -0x1.33d28cp-1)); |
| 727 | try testing.expectApproxEqAbs(0x1.37082ep0, output[4], math.floatEpsAt(f32, 0x1.37082ep0)); |
| 728 | try testing.expectApproxEqAbs(-0x1.99d7cap-2, output[5], math.floatEpsAt(f32, -0x1.99d7cap-2)); |
| 729 | try testing.expectApproxEqAbs(0x1.2fcb12p0, output[6], math.floatEpsAt(f32, 0x1.2fcb12p0)); |
| 730 | try testing.expectApproxEqAbs(-0x1.08c71ap-2, output[7], math.floatEpsAt(f32, -0x1.08c71ap-2)); |
| 731 | try testing.expectApproxEqAbs(0x1.2a24e2p0, output[8], math.floatEpsAt(f32, 0x1.2a24e2p0)); |
| 732 | try testing.expectApproxEqAbs(-0x1.1c6178p0, output[9], math.floatEpsAt(f32, -0x1.1c6178p0)); |
| 733 | } |
| 734 | |
| 735 | test "atanBinary64Vec.special" { |
| 736 | const input: @Vector(7, f64) = .{ |
| 737 | 0x0p+0, |
| 738 | -0x0p+0, |
| 739 | 0x1p+0, |
| 740 | -0x1p+0, |
| 741 | math.inf(f64), |
| 742 | -math.inf(f64), |
| 743 | math.nan(f64), |
| 744 | }; |
| 745 | const output = atanBinary64Vec(7, input); |
| 746 | try testing.expectEqual(0x0p+0, output[0]); |
| 747 | try testing.expectEqual(-0x0p+0, output[1]); |
| 748 | try testing.expectApproxEqAbs(0x1.921fb54442d18p-1, output[2], math.floatEpsAt(f64, 0x1.921fb54442d18p-1)); |
| 749 | try testing.expectApproxEqAbs(-0x1.921fb54442d18p-1, output[3], math.floatEpsAt(f64, -0x1.921fb54442d18p-1)); |
| 750 | try testing.expectApproxEqAbs(0x1.921fb54442d18p+0, output[4], math.floatEpsAt(f64, 0x1.921fb54442d18p+0)); |
| 751 | try testing.expectApproxEqAbs(-0x1.921fb54442d18p+0, output[5], math.floatEpsAt(f64, -0x1.921fb54442d18p+0)); |
| 752 | try testing.expect(math.isNan(output[6])); |
| 753 | } |
| 754 | |
| 755 | test "atanBinary64Vec" { |
| 756 | const input: @Vector(10, f64) = .{ |
| 757 | -0x1.8629d0244cdccp-2, |
| 758 | -0x1.59d42d4659937p1, |
| 759 | -0x1.d2dbe23d04f06p0, |
| 760 | -0x1.5f314e72398e8p-1, |
| 761 | 0x1.5869af37b7d08p1, |
| 762 | -0x1.b13a05a662618p-2, |
| 763 | 0x1.3cb0f12f39d8ap1, |
| 764 | -0x1.0ed746b39cbb7p-2, |
| 765 | 0x1.299d54ac7d6bp1, |
| 766 | -0x1.0264fb9f3d50ep1, |
| 767 | }; |
| 768 | const output = atanBinary64Vec(10, input); |
| 769 | try testing.expectApproxEqAbs(-0x1.74c61f4377016p-2, output[0], math.floatEpsAt(f64, -0x1.74c61f4377016p-2)); |
| 770 | try testing.expectApproxEqAbs(-0x1.375fd7987cc2p0, output[1], math.floatEpsAt(f64, -0x1.375fd7987cc2p0)); |
| 771 | try testing.expectApproxEqAbs(-0x1.11b8adeba5616p0, output[2], math.floatEpsAt(f64, -0x1.11b8adeba5616p0)); |
| 772 | try testing.expectApproxEqAbs(-0x1.33d28ca762539p-1, output[3], math.floatEpsAt(f64, -0x1.33d28ca762539p-1)); |
| 773 | try testing.expectApproxEqAbs(0x1.37082ce2dd03p0, output[4], math.floatEpsAt(f64, 0x1.37082ce2dd03p0)); |
| 774 | try testing.expectApproxEqAbs(-0x1.99d7cac66dd44p-2, output[5], math.floatEpsAt(f64, -0x1.99d7cac66dd44p-2)); |
| 775 | try testing.expectApproxEqAbs(0x1.2fcb120468e8ep0, output[6], math.floatEpsAt(f64, 0x1.2fcb120468e8ep0)); |
| 776 | try testing.expectApproxEqAbs(-0x1.08c71aa0e509p-2, output[7], math.floatEpsAt(f64, -0x1.08c71aa0e509p-2)); |
| 777 | try testing.expectApproxEqAbs(0x1.2a24e22d861dfp0, output[8], math.floatEpsAt(f64, 0x1.2a24e22d861dfp0)); |
| 778 | try testing.expectApproxEqAbs(-0x1.1c617825f9751p0, output[9], math.floatEpsAt(f64, -0x1.1c617825f9751p0)); |
| 779 | } |