╭ fmod.zig ─────────────────────────────────────────────────────╮ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ ╰───────────────────────────────────────────────────────────────╯
╭ find file ─ ─ ─ ─ ─ ─ ─ ─ ─ ─╮ ╭────╮ ╭ clone ─ ─ ─ ─ ─ ─ ─ ─ ─ ─ ─ ─ ╮ │ │ │http│ ssh │ │ ╰ ─ ─ ─ ─ ─ ─ ─ ─ ─ ─ ─ ─ ─ ─ ─╯ ╰────╯ ╰ ─ ─ ─ ─ ─ ─ ─ ─ ─ ─ ─ ─ ─ ─ ─ ╯
.. aarch64_outline_atomics.zig absv.zig absvdi2.zig absvdi2_test.zig absvsi2.zig absvsi2_test.zig absvti2.zig absvti2_test.zig addf3.zig addf3_test.zig addvdi3.zig addvsi3.zig arm.zig atomics.zig aulldiv.zig aullrem.zig bitreverse.zig bitreversedi2_test.zig bitreversesi2_test.zig bitreverseti2_test.zig bswap.zig bswapdi2_test.zig bswapsi2_test.zig bswapti2_test.zig clear_cache.zig clzdi2_test.zig clzsi2_test.zig clzti2_test.zig cmp.zig cmpdi2_test.zig cmpsi2_test.zig cmpti2_test.zig comparedf2_test.zig comparef.zig comparesf2_test.zig cos.zig count0bits.zig ctzdi2_test.zig ctzsi2_test.zig ctzti2_test.zig divc3.zig divc3_test.zig divdf3.zig divdf3_test.zig divmodei4.zig divsf3.zig divsf3_test.zig divtf3.zig divtf3_test.zig divti3_test.zig divxf3.zig divxf3_test.zig emutls.zig exp.zig exp2.zig exp_f128.zig extendf.zig extendf_test.zig fabs.zig ffsdi2_test.zig ffssi2_test.zig ffsti2_test.zig float_from_int.zig float_from_int_test.zig floor_ceil.zig fma.zig fmax.zig fmin.zig ╭─────────────────────────────────╮ │fmod.zig │ ╰─────────────────────────────────╯ fmodq_test.zig fmodx_test.zig hexagon.zig int.zig int_from_float.zig int_from_float_test.zig limb64.zig log.zig log10.zig log2.zig log_f128.zig memcpy.zig memmove.zig modti3_test.zig mulXi3.zig mulXi3_test.zig mulc3.zig mulc3_test.zig mulf3.zig mulf3_test.zig mulo.zig mulodi4_test.zig mulosi4_test.zig muloti4_test.zig mulvsi3.zig negdi2_test.zig negsi2_test.zig negti2_test.zig negv.zig negvdi2_test.zig negvsi2_test.zig negvti2_test.zig os_version_check.zig parity.zig paritydi2_test.zig paritysi2_test.zig parityti2_test.zig popcount.zig popcountdi2_test.zig popcountsi2_test.zig popcountti2_test.zig powiXf2.zig powiXf2_test.zig rem_pio2.zig rem_pio2_large.zig rem_pio2f.zig rem_pio2l.zig round.zig shift.zig shift_test.zig sin.zig sincos.zig sqrt.zig ssp.zig stack_probe.zig subvdi3.zig subvsi3.zig tan.zig trig.zig trunc.zig truncf.zig truncf_test.zig ucmpdi2_test.zig ucmpsi2_test.zig ucmpti2_test.zig udivmod.zig udivmoddi4_test.zig udivmodei4.zig udivmodsi4_test.zig udivmodti4_test.zig
1 const builtin = @import("builtin"); 2 const std = @import("std"); 3 const math = std.math; 4 const assert = std.debug.assert; 5 const arch = builtin.cpu.arch; 6 const compiler_rt = @import("../compiler_rt.zig"); 7 const symbol = compiler_rt.symbol; 8 const normalize = compiler_rt.normalize; 9 10 comptime { 11 symbol(&__fmodh, "__fmodh"); 12 symbol(&fmodf, "fmodf"); 13 symbol(&fmod, "fmod"); 14 symbol(&__fmodx, "__fmodx"); 15 symbol(&fmodq, "fmodf128"); 16 symbol(&fmodl, "fmodl"); 17 } 18 19 fn __fmodh(a: compiler_rt.f16.Abi, b: compiler_rt.f16.Abi) callconv(.c) compiler_rt.f16.Abi { 20 return compiler_rt.f16.toAbi(fmod_f16(compiler_rt.f16.fromAbi(a), compiler_rt.f16.fromAbi(b))); 21 } 22 pub fn fmod_f16(x: f16, y: f16) f16 { 23 // TODO: more efficient implementation 24 return @floatCast(fmod_f32(x, y)); 25 } 26 27 fn fmodf(a: compiler_rt.f32.Abi, b: compiler_rt.f32.Abi) callconv(.c) compiler_rt.f32.Abi { 28 return compiler_rt.f32.toAbi(fmod_f32(compiler_rt.f32.fromAbi(a), compiler_rt.f32.fromAbi(b))); 29 } 30 pub fn fmod_f32(x: f32, y: f32) f32 { 31 return generic_fmod(f32, x, y); 32 } 33 34 fn fmod(a: compiler_rt.f64.Abi, b: compiler_rt.f64.Abi) callconv(.c) compiler_rt.f64.Abi { 35 return compiler_rt.f64.toAbi(fmod_f64(compiler_rt.f64.fromAbi(a), compiler_rt.f64.fromAbi(b))); 36 } 37 pub fn fmod_f64(x: f64, y: f64) f64 { 38 return generic_fmod(f64, x, y); 39 } 40 41 fn __fmodx(a: compiler_rt.f80.Abi, b: compiler_rt.f80.Abi) callconv(.c) compiler_rt.f80.Abi { 42 return compiler_rt.f80.toAbi(fmod_f80(compiler_rt.f80.fromAbi(a), compiler_rt.f80.fromAbi(b))); 43 } 44 /// fmodx - floating modulo large, returns the remainder of division for f80 types 45 /// Logic and flow heavily inspired by MUSL fmodl for 113 mantissa digits 46 pub fn fmod_f80(a: f80, b: f80) f80 { 47 const T = f80; 48 const Z = @Int(.unsigned, @bitSizeOf(T)); 49 50 const significandBits = math.floatMantissaBits(T); 51 const fractionalBits = math.floatFractionalBits(T); 52 const exponentBits = math.floatExponentBits(T); 53 54 const signBit = (@as(Z, 1) << (significandBits + exponentBits)); 55 const maxExponent = ((1 << exponentBits) - 1); 56 57 var aRep: Z = @bitCast(a); 58 var bRep: Z = @bitCast(b); 59 60 const signA = aRep & signBit; 61 var expA: i32 = @intCast((@as(Z, @bitCast(a)) >> significandBits) & maxExponent); 62 var expB: i32 = @intCast((@as(Z, @bitCast(b)) >> significandBits) & maxExponent); 63 64 // There are 3 cases where the answer is undefined, check for: 65 // - fmodx(val, 0) 66 // - fmodx(val, NaN) 67 // - fmodx(inf, val) 68 // The sign on checked values does not matter. 69 // Doing (a * b) / (a * b) produces undefined results 70 // because the three cases always produce undefined calculations: 71 // - 0 / 0 72 // - val * NaN 73 // - inf / inf 74 if (b == 0 or math.isNan(b) or expA == maxExponent) { 75 return (a * b) / (a * b); 76 } 77 78 // Remove the sign from both 79 aRep &= ~signBit; 80 bRep &= ~signBit; 81 if (aRep <= bRep) { 82 if (aRep == bRep) { 83 return 0 * a; 84 } 85 return a; 86 } 87 88 if (expA == 0) expA = normalize(f80, &aRep); 89 if (expB == 0) expB = normalize(f80, &bRep); 90 91 var highA: u64 = 0; 92 const highB: u64 = 0; 93 var lowA: u64 = @truncate(aRep); 94 const lowB: u64 = @truncate(bRep); 95 96 while (expA > expB) : (expA -= 1) { 97 var high = highA -% highB; 98 const low = lowA -% lowB; 99 if (lowA < lowB) { 100 high -%= 1; 101 } 102 if (high >> 63 == 0) { 103 if ((high | low) == 0) { 104 return 0 * a; 105 } 106 highA = 2 *% high + (low >> 63); 107 lowA = 2 *% low; 108 } else { 109 highA = 2 *% highA + (lowA >> 63); 110 lowA = 2 *% lowA; 111 } 112 } 113 114 var high = highA -% highB; 115 const low = lowA -% lowB; 116 if (lowA < lowB) { 117 high -%= 1; 118 } 119 if (high >> 63 == 0) { 120 if ((high | low) == 0) { 121 return 0 * a; 122 } 123 highA = high; 124 lowA = low; 125 } 126 127 while ((lowA >> fractionalBits) == 0) { 128 lowA = 2 *% lowA; 129 expA = expA - 1; 130 } 131 132 // Combine the exponent with the sign and significand, normalize if happened to be denormalized 133 if (expA < -fractionalBits) { 134 return @bitCast(signA); 135 } else if (expA <= 0) { 136 return @bitCast((lowA >> @intCast(1 - expA)) | signA); 137 } else { 138 return @bitCast(lowA | (@as(Z, @as(u16, @intCast(expA))) << significandBits) | signA); 139 } 140 } 141 142 fn fmodq(a: compiler_rt.f128.Abi, b: compiler_rt.f128.Abi) callconv(.c) compiler_rt.f128.Abi { 143 return compiler_rt.f128.toAbi(fmod_f128(compiler_rt.f128.fromAbi(a), compiler_rt.f128.fromAbi(b))); 144 } 145 /// fmodq - floating modulo large, returns the remainder of division for f128 types 146 /// Logic and flow heavily inspired by MUSL fmodl for 113 mantissa digits 147 pub fn fmod_f128(a: f128, b: f128) f128 { 148 var amod = a; 149 var bmod = b; 150 const aPtr_u64: [*]u64 = @ptrCast(&amod); 151 const bPtr_u64: [*]u64 = @ptrCast(&bmod); 152 const aPtr_u16: [*]u16 = @ptrCast(&amod); 153 const bPtr_u16: [*]u16 = @ptrCast(&bmod); 154 155 const exp_and_sign_index = comptime switch (builtin.target.cpu.arch.endian()) { 156 .little => 7, 157 .big => 0, 158 }; 159 const low_index = comptime switch (builtin.target.cpu.arch.endian()) { 160 .little => 0, 161 .big => 1, 162 }; 163 const high_index = comptime switch (builtin.target.cpu.arch.endian()) { 164 .little => 1, 165 .big => 0, 166 }; 167 168 const signA = aPtr_u16[exp_and_sign_index] & 0x8000; 169 var expA: i32 = @intCast((aPtr_u16[exp_and_sign_index] & 0x7fff)); 170 var expB: i32 = @intCast((bPtr_u16[exp_and_sign_index] & 0x7fff)); 171 172 // There are 3 cases where the answer is undefined, check for: 173 // - fmodq(val, 0) 174 // - fmodq(val, NaN) 175 // - fmodq(inf, val) 176 // The sign on checked values does not matter. 177 // Doing (a * b) / (a * b) produces undefined results 178 // because the three cases always produce undefined calculations: 179 // - 0 / 0 180 // - val * NaN 181 // - inf / inf 182 if (b == 0 or std.math.isNan(b) or expA == 0x7fff) { 183 return (a * b) / (a * b); 184 } 185 186 // Remove the sign from both 187 aPtr_u16[exp_and_sign_index] = @bitCast(@as(i16, @intCast(expA))); 188 bPtr_u16[exp_and_sign_index] = @bitCast(@as(i16, @intCast(expB))); 189 if (amod <= bmod) { 190 if (amod == bmod) { 191 return 0 * a; 192 } 193 return a; 194 } 195 196 if (expA == 0) { 197 amod *= 0x1p120; 198 expA = @as(i32, aPtr_u16[exp_and_sign_index]) - 120; 199 } 200 201 if (expB == 0) { 202 bmod *= 0x1p120; 203 expB = @as(i32, bPtr_u16[exp_and_sign_index]) - 120; 204 } 205 206 // OR in extra non-stored mantissa digit 207 var highA: u64 = (aPtr_u64[high_index] & (std.math.maxInt(u64) >> 16)) | 1 << 48; 208 const highB: u64 = (bPtr_u64[high_index] & (std.math.maxInt(u64) >> 16)) | 1 << 48; 209 var lowA: u64 = aPtr_u64[low_index]; 210 const lowB: u64 = bPtr_u64[low_index]; 211 212 while (expA > expB) : (expA -= 1) { 213 var high = highA -% highB; 214 const low = lowA -% lowB; 215 if (lowA < lowB) { 216 high -%= 1; 217 } 218 if (high >> 63 == 0) { 219 if ((high | low) == 0) { 220 return 0 * a; 221 } 222 highA = 2 *% high + (low >> 63); 223 lowA = 2 *% low; 224 } else { 225 highA = 2 *% highA + (lowA >> 63); 226 lowA = 2 *% lowA; 227 } 228 } 229 230 var high = highA -% highB; 231 const low = lowA -% lowB; 232 if (lowA < lowB) { 233 high -= 1; 234 } 235 if (high >> 63 == 0) { 236 if ((high | low) == 0) { 237 return 0 * a; 238 } 239 highA = high; 240 lowA = low; 241 } 242 243 while (highA >> 48 == 0) { 244 highA = 2 *% highA + (lowA >> 63); 245 lowA = 2 *% lowA; 246 expA = expA - 1; 247 } 248 249 // Overwrite the current amod with the values in highA and lowA 250 aPtr_u64[high_index] = highA; 251 aPtr_u64[low_index] = lowA; 252 253 // Combine the exponent with the sign, normalize if happened to be denormalized 254 if (expA <= 0) { 255 aPtr_u16[exp_and_sign_index] = @as(u16, @truncate(@as(u32, @bitCast((expA +% 120))))) | signA; 256 amod *= 0x1p-120; 257 } else { 258 aPtr_u16[exp_and_sign_index] = @as(u16, @truncate(@as(u32, @bitCast(expA)))) | signA; 259 } 260 261 return amod; 262 } 263 264 pub fn fmodl(a: c_longdouble, b: c_longdouble) callconv(.c) c_longdouble { 265 switch (@typeInfo(c_longdouble).float.bits) { 266 64 => return fmod_f64(a, b), 267 80 => return fmod_f80(a, b), 268 128 => return fmod_f128(a, b), 269 else => comptime unreachable, 270 } 271 } 272 273 inline fn generic_fmod(comptime T: type, x: T, y: T) T { 274 const bits = @typeInfo(T).float.bits; 275 const uint = @Int(.unsigned, bits); 276 comptime assert(T == f32 or T == f64); 277 const digits = if (T == f32) 23 else 52; 278 const exp_bits = if (T == f32) 9 else 12; 279 const bits_minus_1 = bits - 1; 280 const mask = if (T == f32) 0xff else 0x7ff; 281 var ux: uint = @bitCast(x); 282 var uy: uint = @bitCast(y); 283 var ex: i32 = @intCast((ux >> digits) & mask); 284 var ey: i32 = @intCast((uy >> digits) & mask); 285 const sx = if (T == f32) @as(u32, @intCast(ux & 0x80000000)) else @as(i32, @intCast(ux >> bits_minus_1)); 286 var i: uint = undefined; 287 288 if (uy << 1 == 0 or math.isNan(@as(T, @bitCast(uy))) or ex == mask) 289 return (x * y) / (x * y); 290 291 if (ux << 1 <= uy << 1) { 292 if (ux << 1 == uy << 1) 293 return 0 * x; 294 return x; 295 } 296 297 // normalize x and y 298 if (ex == 0) { 299 i = ux << exp_bits; 300 while (i >> bits_minus_1 == 0) : ({ 301 ex -= 1; 302 i <<= 1; 303 }) {} 304 ux <<= @intCast(@as(u32, @bitCast(-ex + 1))); 305 } else { 306 ux &= math.maxInt(uint) >> exp_bits; 307 ux |= 1 << digits; 308 } 309 if (ey == 0) { 310 i = uy << exp_bits; 311 while (i >> bits_minus_1 == 0) : ({ 312 ey -= 1; 313 i <<= 1; 314 }) {} 315 uy <<= @intCast(@as(u32, @bitCast(-ey + 1))); 316 } else { 317 uy &= math.maxInt(uint) >> exp_bits; 318 uy |= 1 << digits; 319 } 320 321 // x mod y 322 while (ex > ey) : (ex -= 1) { 323 i = ux -% uy; 324 if (i >> bits_minus_1 == 0) { 325 if (i == 0) 326 return 0 * x; 327 ux = i; 328 } 329 ux <<= 1; 330 } 331 i = ux -% uy; 332 if (i >> bits_minus_1 == 0) { 333 if (i == 0) 334 return 0 * x; 335 ux = i; 336 } 337 while (ux >> digits == 0) : ({ 338 ux <<= 1; 339 ex -= 1; 340 }) {} 341 342 // scale result up 343 if (ex > 0) { 344 ux -%= 1 << digits; 345 ux |= @as(uint, @as(u32, @bitCast(ex))) << digits; 346 } else { 347 ux >>= @intCast(@as(u32, @bitCast(-ex + 1))); 348 } 349 if (T == f32) { 350 ux |= sx; 351 } else { 352 ux |= @as(uint, @intCast(sx)) << bits_minus_1; 353 } 354 return @bitCast(ux); 355 } 356 357 test fmod_f32 { 358 const nan_val = math.nan(f32); 359 const inf_val = math.inf(f32); 360 361 try std.testing.expect(math.isNan(fmod_f32(nan_val, 1.0))); 362 try std.testing.expect(math.isNan(fmod_f32(1.0, nan_val))); 363 try std.testing.expect(math.isNan(fmod_f32(inf_val, 1.0))); 364 try std.testing.expect(math.isNan(fmod_f32(0.0, 0.0))); 365 try std.testing.expect(math.isNan(fmod_f32(1.0, 0.0))); 366 367 try std.testing.expectEqual(@as(f32, 0.0), fmod_f32(0.0, 2.0)); 368 try std.testing.expectEqual(@as(f32, -0.0), fmod_f32(-0.0, 2.0)); 369 370 try std.testing.expectEqual(@as(f32, -2.0), fmod_f32(-32.0, 10.0)); 371 try std.testing.expectEqual(@as(f32, -2.0), fmod_f32(-32.0, -10.0)); 372 try std.testing.expectEqual(@as(f32, 2.0), fmod_f32(32.0, 10.0)); 373 try std.testing.expectEqual(@as(f32, 2.0), fmod_f32(32.0, -10.0)); 374 } 375 376 test fmod_f64 { 377 const nan_val = math.nan(f64); 378 const inf_val = math.inf(f64); 379 380 try std.testing.expect(math.isNan(fmod_f64(nan_val, 1.0))); 381 try std.testing.expect(math.isNan(fmod_f64(1.0, nan_val))); 382 try std.testing.expect(math.isNan(fmod_f64(inf_val, 1.0))); 383 try std.testing.expect(math.isNan(fmod_f64(0.0, 0.0))); 384 try std.testing.expect(math.isNan(fmod_f64(1.0, 0.0))); 385 386 try std.testing.expectEqual(@as(f64, 0.0), fmod_f64(0.0, 2.0)); 387 try std.testing.expectEqual(@as(f64, -0.0), fmod_f64(-0.0, 2.0)); 388 389 try std.testing.expectEqual(@as(f64, -2.0), fmod_f64(-32.0, 10.0)); 390 try std.testing.expectEqual(@as(f64, -2.0), fmod_f64(-32.0, -10.0)); 391 try std.testing.expectEqual(@as(f64, 2.0), fmod_f64(32.0, 10.0)); 392 try std.testing.expectEqual(@as(f64, 2.0), fmod_f64(32.0, -10.0)); 393 } 394 395 test { 396 _ = @import("fmodq_test.zig"); 397 _ = @import("fmodx_test.zig"); 398 } 399