1 | // Copyright (c) 2019-2023 Alexander Medvednikov. All rights reserved. |
2 | // Use of this source code is governed by an MIT license |
3 | // that can be found in the LICENSE file. |
4 | module math |
5 | |
6 | const ( |
7 | uvnan = u64(0x7FF8000000000001) |
8 | uvinf = u64(0x7FF0000000000000) |
9 | uvneginf = u64(0xFFF0000000000000) |
10 | uvone = u64(0x3FF0000000000000) |
11 | mask = 0x7FF |
12 | shift = 64 - 11 - 1 |
13 | bias = 1023 |
14 | normalize_smallest_mask = u64(u64(1) << 52) |
15 | sign_mask = u64(0x8000000000000000) // (u64(1) << 63) |
16 | frac_mask = u64((u64(1) << u64(shift)) - u64(1)) |
17 | ) |
18 | |
19 | // inf returns positive infinity if sign >= 0, negative infinity if sign < 0. |
20 | pub fn inf(sign int) f64 { |
21 | v := if sign >= 0 { math.uvinf } else { math.uvneginf } |
22 | return f64_from_bits(v) |
23 | } |
24 | |
25 | // nan returns an IEEE 754 ``not-a-number'' value. |
26 | pub fn nan() f64 { |
27 | return f64_from_bits(math.uvnan) |
28 | } |
29 | |
30 | // is_nan reports whether f is an IEEE 754 ``not-a-number'' value. |
31 | pub fn is_nan(f f64) bool { |
32 | $if fast_math { |
33 | if f64_bits(f) == math.uvnan { |
34 | return true |
35 | } |
36 | } |
37 | // IEEE 754 says that only NaNs satisfy f != f. |
38 | // To avoid the floating-point hardware, could use: |
39 | // x := f64_bits(f); |
40 | // return u32(x>>shift)&mask == mask && x != uvinf && x != uvneginf |
41 | return f != f |
42 | } |
43 | |
44 | // is_inf reports whether f is an infinity, according to sign. |
45 | // If sign > 0, is_inf reports whether f is positive infinity. |
46 | // If sign < 0, is_inf reports whether f is negative infinity. |
47 | // If sign == 0, is_inf reports whether f is either infinity. |
48 | pub fn is_inf(f f64, sign int) bool { |
49 | // Test for infinity by comparing against maximum float. |
50 | // To avoid the floating-point hardware, could use: |
51 | // x := f64_bits(f); |
52 | // return sign >= 0 && x == uvinf || sign <= 0 && x == uvneginf; |
53 | return (sign >= 0 && f > max_f64) || (sign <= 0 && f < -max_f64) |
54 | } |
55 | |
56 | pub fn is_finite(f f64) bool { |
57 | return !is_nan(f) && !is_inf(f, 0) |
58 | } |
59 | |
60 | // normalize returns a normal number y and exponent exp |
61 | // satisfying x == y × 2**exp. It assumes x is finite and non-zero. |
62 | pub fn normalize(x f64) (f64, int) { |
63 | smallest_normal := 2.2250738585072014e-308 // 2**-1022 |
64 | if abs(x) < smallest_normal { |
65 | return x * math.normalize_smallest_mask, -52 |
66 | } |
67 | return x, 0 |
68 | } |