ftos_multiprecision.ha (7247B)
1 // SPDX-License-Identifier: MPL-2.0 2 // (c) Hare authors <https://harelang.org> 3 4 // Multiprecision float to string based on golang's strconv/decimal.go. 5 6 use ascii; 7 use math; 8 use strings; 9 10 type mp = struct { 11 // Numbers 0-9, not ascii. Length for small numbers is 12 // log10(mantissa * 5^-exp). Subnormal doubles have min exp -1074 and 13 // max mantissa 4e16, giving at most 767 digits. 14 buf: [768]u8, 15 // Number of valid digits in buf. May be 0 if the number rounds to 0. 16 nd: uint, 17 // Decimal point index, may be negative. 18 // -1 means 0.0ddd... 19 // 0 means 0.ddd... 20 // 1 means d.dd... 21 // and so on 22 dp: int, 23 }; 24 25 // These come from golang. The index into the table is amount of shift, up to 26 // 60. The number is the count of new digits that will be added by the shift, 27 // but one fewer if the number's prefix is smaller than the string prefix. 28 // 29 // For example, leftcheats[2] is (1, "25"). Any number left shifted by 2 will 30 // therefore be 1 digit longer, or zero digits longer if its first two digits 31 // are smaller than 25. 32 const leftcheats: [](size, str) = [ 33 (0, ""), 34 (1, "5"), 35 (1, "25"), 36 (1, "125"), 37 (2, "625"), 38 (2, "3125"), 39 (2, "15625"), 40 (3, "78125"), 41 (3, "390625"), 42 (3, "1953125"), 43 (4, "9765625"), 44 (4, "48828125"), 45 (4, "244140625"), 46 (4, "1220703125"), 47 (5, "6103515625"), 48 (5, "30517578125"), 49 (5, "152587890625"), 50 (6, "762939453125"), 51 (6, "3814697265625"), 52 (6, "19073486328125"), 53 (7, "95367431640625"), 54 (7, "476837158203125"), 55 (7, "2384185791015625"), 56 (7, "11920928955078125"), 57 (8, "59604644775390625"), 58 (8, "298023223876953125"), 59 (8, "1490116119384765625"), 60 (9, "7450580596923828125"), 61 (9, "37252902984619140625"), 62 (9, "186264514923095703125"), 63 (10, "931322574615478515625"), 64 (10, "4656612873077392578125"), 65 (10, "23283064365386962890625"), 66 (10, "116415321826934814453125"), 67 (11, "582076609134674072265625"), 68 (11, "2910383045673370361328125"), 69 (11, "14551915228366851806640625"), 70 (12, "72759576141834259033203125"), 71 (12, "363797880709171295166015625"), 72 (12, "1818989403545856475830078125"), 73 (13, "9094947017729282379150390625"), 74 (13, "45474735088646411895751953125"), 75 (13, "227373675443232059478759765625"), 76 (13, "1136868377216160297393798828125"), 77 (14, "5684341886080801486968994140625"), 78 (14, "28421709430404007434844970703125"), 79 (14, "142108547152020037174224853515625"), 80 (15, "710542735760100185871124267578125"), 81 (15, "3552713678800500929355621337890625"), 82 (15, "17763568394002504646778106689453125"), 83 (16, "88817841970012523233890533447265625"), 84 (16, "444089209850062616169452667236328125"), 85 (16, "2220446049250313080847263336181640625"), 86 (16, "11102230246251565404236316680908203125"), 87 (17, "55511151231257827021181583404541015625"), 88 (17, "277555756156289135105907917022705078125"), 89 (17, "1387778780781445675529539585113525390625"), 90 (18, "6938893903907228377647697925567626953125"), 91 (18, "34694469519536141888238489627838134765625"), 92 (18, "173472347597680709441192448139190673828125"), 93 (19, "867361737988403547205962240695953369140625"), 94 ]; 95 96 fn prefix_less_than_mp(m: *mp, s: str) bool = { 97 const u = strings::toutf8(s); 98 for (let i = 0z; i < len(s); i += 1) { 99 if (i >= m.nd) { 100 return true; 101 }; 102 if (m.buf[i] + '0': u8 != u[i]) { 103 return m.buf[i] + '0': u8 < u[i]; 104 }; 105 }; 106 return false; 107 }; 108 109 // Shift left by k. 110 fn shl_mp(m: *mp, k: u64) void = { 111 let delta = leftcheats[k].0; 112 if (prefix_less_than_mp(m, leftcheats[k].1)) 113 delta -= 1; 114 let r = (m.nd - 1): int; 115 let w = m.nd + delta; 116 let n = 0u64; 117 for (r >= 0; r -= 1) { 118 n += m.buf[r]: u64 << k; 119 const quo = n / 10; 120 const rem = n - 10 * quo; 121 w -= 1; 122 m.buf[w] = rem: u8; 123 n = quo; 124 }; 125 for (n > 0) { 126 const quo = n / 10; 127 const rem = n - 10 * quo; 128 w -= 1; 129 m.buf[w] = rem: u8; 130 n = quo; 131 }; 132 m.nd += delta: uint; 133 m.dp += delta: int; 134 }; 135 136 // Shift right by k. 137 fn shr_mp(m: *mp, k: u64) void = { 138 let r = 0z; 139 let w = 0z; 140 let n = 0u64; 141 const mask = (1 << k) - 1; 142 143 for (n >> k == 0; r += 1) { 144 if (r >= m.nd) { 145 for (n >> k == 0) { 146 n *= 10; 147 r += 1; 148 }; 149 break; 150 }; 151 n = 10 * n + m.buf[r]; 152 }; 153 m.dp -= r: int - 1; 154 155 for (r < m.nd; r += 1) { 156 const c = m.buf[r]; 157 const dig = n >> k; 158 n &= mask; 159 m.buf[w] = dig: u8; 160 w += 1; 161 n = n * 10 + c; 162 }; 163 164 for (n > 0; w += 1) { 165 const dig = n >> k; 166 n &= mask; 167 m.buf[w] = dig: u8; 168 n = n * 10; 169 }; 170 m.nd = w: uint; 171 }; 172 173 // Shift right (k < 0) or left (k > 0). We can only shift up to 60 at a time 174 // without losing bits, so break up big shifts. 175 fn shift_mp(m: *mp, k: int) void = { 176 if (k < 0) { 177 let nk = (-k): uint; 178 for (nk > 60) { 179 shr_mp(m, 60); 180 nk -= 60; 181 }; 182 shr_mp(m, nk); 183 } else if (k > 0) { 184 for (k > 60) { 185 shl_mp(m, 60); 186 k -= 60; 187 }; 188 shl_mp(m, k: uint); 189 }; 190 }; 191 192 fn init_mp(m: *mp, mantissa: u64, exponent: u32, eb: u64, mb: u64) void = { 193 let e2 = (eb + mb): i32; 194 let m2: u64 = 0; 195 if (exponent == 0) { 196 e2 = 1 - e2; 197 m2 = mantissa; 198 } else { 199 e2 = (exponent: i32) - e2; 200 m2 = (1u64 << mb) | mantissa; 201 }; 202 203 m.nd = declen(m2); 204 m.dp = m.nd: int; 205 for (let i = 0z; i < m.nd; i += 1) { 206 m.buf[m.nd - i - 1] = (m2 % 10): u8; 207 m2 /= 10; 208 }; 209 shift_mp(m, e2); 210 }; 211 212 fn init_mp_dec(m: *mp, mantissa: u64, exponent: i32) void = { 213 const dl = declen(mantissa); 214 for (let i = 0u; i < dl; i += 1) { 215 m.buf[dl - i - 1] = (mantissa % 10): u8; 216 mantissa /= 10; 217 }; 218 m.nd = dl; 219 m.dp = dl: i32 + exponent; 220 }; 221 222 fn round_up_mp(m: *mp) void = { 223 for (let i = 1z; i <= m.nd; i += 1) { 224 if (m.buf[m.nd - i] < 9) { 225 m.buf[m.nd - i] += 1; 226 return; 227 } else { 228 m.buf[m.nd - i] = 0; 229 }; 230 }; 231 // All high 232 m.buf[0] = 1; 233 m.nd = 1; 234 m.dp += 1; 235 }; 236 237 // Compute the number of figs to round to for the given arguments. 238 fn compute_round_mp(m: *mp, f: ffmt, prec: (void | uint), flag: fflags) uint = { 239 // nd is the number of sig figs that we want to end up with 240 let nd: int = match (prec) { 241 case void => 242 // we should only get here if Ryu did not extend past the 243 // decimal point 244 assert(ffpoint(flag)); 245 yield m.nd: int + (if (m.dp > 0) m.dp else 0); 246 case let u: uint => 247 yield switch (f) { 248 case ffmt::E => 249 yield u: int + 1; 250 case ffmt::F => 251 yield u: int + m.dp; 252 case ffmt::G => 253 yield if (u == 0) 1 else u: int; 254 }; 255 }; 256 const nde = if (nd < 2) 2 else nd; 257 const ndf = if (m.dp >= 0 && nd < m.dp + 1) m.dp + 1 else nd; 258 if (ffpoint(flag)) { 259 nd = switch (f) { 260 case ffmt::E => 261 // need at least two digits, d.de0. 262 yield nde; 263 case ffmt::F => 264 // need enough to clear the decimal point by one. 265 yield ndf; 266 case ffmt::G => 267 // XXX: dup'd with the condition in ftosf_handle 268 if (m.dp < -1 || m.dp - m.nd: int > 2) 269 yield nde; 270 yield ndf; 271 }; 272 }; 273 if (nd <= 0) { 274 nd = 0; 275 }; 276 return if (nd: uint > m.nd) m.nd else nd: uint; 277 }; 278 279 fn round_mp(m: *mp, nd: uint) void = { 280 assert(nd <= m.nd); 281 if (nd == m.nd) 282 return; 283 const oldnd = m.nd; 284 m.nd = nd; 285 if (m.buf[nd] > 5) { 286 round_up_mp(m); 287 } else if (m.buf[nd] == 5) { 288 let gt = false; 289 for (let i = m.nd + 1; i < oldnd; i += 1) { 290 if (m.buf[i] > 0) { 291 round_up_mp(m); 292 gt = true; 293 break; 294 }; 295 }; 296 if (!gt && nd > 0 && m.buf[nd - 1] & 1 > 0) { 297 round_up_mp(m); 298 }; 299 }; 300 }; 301 302 // Remove trailing zeros. 303 fn trim_mp(m: *mp) void = { 304 for (m.nd > 1 && m.buf[m.nd - 1] == 0) { 305 m.nd -= 1; 306 }; 307 };