summaryrefslogtreecommitdiff log msg author committer range
path: root/lib/math/scale2-impl.myr
blob: f75aac864e0c51714c9edaecf722b43603b8f948 (plain)
 ```1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 ``` ``````use std use "fpmath" use "util" /* The scaleB function recommended by the 2008 revision of IEEE 754. Since only integer exponents are considered, the naive approach works quite well. */ pkg math = const scale232 : (x : flt32, m : int32 -> flt32) const scale264 : (x : flt64, m : int64 -> flt64) ;; const scale232 = {x : flt32, m : int32 var n, e, s (n, e, s) = std.flt32explode(x) (n, e, s) = scale2gen(n, e, s, -126, 127, 24, m) -> std.flt32assem(n, e, s) } const scale264 = {x : flt64, m : int64 var n, e, s (n, e, s) = std.flt64explode(x) (n, e, s) = scale2gen(n, e, s, -1022, 1023, 53, m) -> std.flt64assem(n, e, s) } generic scale2gen = {n : bool, e : @i, s : @u, emin : @i, emax : @i, p : @u, m : @i :: numeric, integral @i, numeric, integral @u if e < emin && s == 0 -> (n, e, s) elif m == 0 -> (n, e, s) elif m < 0 var sprime = s var eprime = e if e < emin sprime = s >> (-m) if need_round_away(0, (s : uint64), (-m : int64)) sprime++ ;; eprime = emin - 1 elif e + m < emin - p - 2 sprime = 0 eprime = emin - 1 elif e + m < emin sprime = s >> ((emin - m - e) : @u) if need_round_away(0, (s : uint64), ((emin - m - e) : int64)) sprime++ ;; eprime = emin - 1 else eprime = e + m ;; -> (n, eprime, sprime) ;; if e < emin var first_1 = find_first1_64((s : uint64), (p : int64)) var shift = p - (first_1 : @u) - 1 if m >= (shift : @i) s = s << shift m -= (shift : @i) e = emin else -> (n, e, s << m) ;; ;; if e + m > emax && s != 0 -> (n, emax + 1, 1 << (p - 1)) ;; -> (n, e + m, s) } ``````