1 // Copyright 2010 The Go Authors. All rights reserved. 2 // Use of this source code is governed by a BSD-style 3 // license that can be found in the LICENSE file. 4 5 // Software IEEE754 64-bit floating point. 6 // Only referred to (and thus linked in) by arm port 7 // and by tests in this directory. 8 9 package runtime 10 11 const ( 12 mantbits64 uint = 52 13 expbits64 uint = 11 14 bias64 = -1<<(expbits64-1) + 1 15 16 nan64 uint64 = (1<<expbits64-1)<<mantbits64 + 1 17 inf64 uint64 = (1<<expbits64 - 1) << mantbits64 18 neg64 uint64 = 1 << (expbits64 + mantbits64) 19 20 mantbits32 uint = 23 21 expbits32 uint = 8 22 bias32 = -1<<(expbits32-1) + 1 23 24 nan32 uint32 = (1<<expbits32-1)<<mantbits32 + 1 25 inf32 uint32 = (1<<expbits32 - 1) << mantbits32 26 neg32 uint32 = 1 << (expbits32 + mantbits32) 27 ) 28 29 func funpack64(f uint64) (sign, mant uint64, exp int, inf, nan bool) { 30 sign = f & (1 << (mantbits64 + expbits64)) 31 mant = f & (1<<mantbits64 - 1) 32 exp = int(f>>mantbits64) & (1<<expbits64 - 1) 33 34 switch exp { 35 case 1<<expbits64 - 1: 36 if mant != 0 { 37 nan = true 38 return 39 } 40 inf = true 41 return 42 43 case 0: 44 // denormalized 45 if mant != 0 { 46 exp += bias64 + 1 47 for mant < 1<<mantbits64 { 48 mant <<= 1 49 exp-- 50 } 51 } 52 53 default: 54 // add implicit top bit 55 mant |= 1 << mantbits64 56 exp += bias64 57 } 58 return 59 } 60 61 func funpack32(f uint32) (sign, mant uint32, exp int, inf, nan bool) { 62 sign = f & (1 << (mantbits32 + expbits32)) 63 mant = f & (1<<mantbits32 - 1) 64 exp = int(f>>mantbits32) & (1<<expbits32 - 1) 65 66 switch exp { 67 case 1<<expbits32 - 1: 68 if mant != 0 { 69 nan = true 70 return 71 } 72 inf = true 73 return 74 75 case 0: 76 // denormalized 77 if mant != 0 { 78 exp += bias32 + 1 79 for mant < 1<<mantbits32 { 80 mant <<= 1 81 exp-- 82 } 83 } 84 85 default: 86 // add implicit top bit 87 mant |= 1 << mantbits32 88 exp += bias32 89 } 90 return 91 } 92 93 func fpack64(sign, mant uint64, exp int, trunc uint64) uint64 { 94 mant0, exp0, trunc0 := mant, exp, trunc 95 if mant == 0 { 96 return sign 97 } 98 for mant < 1<<mantbits64 { 99 mant <<= 1 100 exp-- 101 } 102 for mant >= 4<<mantbits64 { 103 trunc |= mant & 1 104 mant >>= 1 105 exp++ 106 } 107 if mant >= 2<<mantbits64 { 108 if mant&1 != 0 && (trunc != 0 || mant&2 != 0) { 109 mant++ 110 if mant >= 4<<mantbits64 { 111 mant >>= 1 112 exp++ 113 } 114 } 115 mant >>= 1 116 exp++ 117 } 118 if exp >= 1<<expbits64-1+bias64 { 119 return sign ^ inf64 120 } 121 if exp < bias64+1 { 122 if exp < bias64-int(mantbits64) { 123 return sign | 0 124 } 125 // repeat expecting denormal 126 mant, exp, trunc = mant0, exp0, trunc0 127 for exp < bias64 { 128 trunc |= mant & 1 129 mant >>= 1 130 exp++ 131 } 132 if mant&1 != 0 && (trunc != 0 || mant&2 != 0) { 133 mant++ 134 } 135 mant >>= 1 136 exp++ 137 if mant < 1<<mantbits64 { 138 return sign | mant 139 } 140 } 141 return sign | uint64(exp-bias64)<<mantbits64 | mant&(1<<mantbits64-1) 142 } 143 144 func fpack32(sign, mant uint32, exp int, trunc uint32) uint32 { 145 mant0, exp0, trunc0 := mant, exp, trunc 146 if mant == 0 { 147 return sign 148 } 149 for mant < 1<<mantbits32 { 150 mant <<= 1 151 exp-- 152 } 153 for mant >= 4<<mantbits32 { 154 trunc |= mant & 1 155 mant >>= 1 156 exp++ 157 } 158 if mant >= 2<<mantbits32 { 159 if mant&1 != 0 && (trunc != 0 || mant&2 != 0) { 160 mant++ 161 if mant >= 4<<mantbits32 { 162 mant >>= 1 163 exp++ 164 } 165 } 166 mant >>= 1 167 exp++ 168 } 169 if exp >= 1<<expbits32-1+bias32 { 170 return sign ^ inf32 171 } 172 if exp < bias32+1 { 173 if exp < bias32-int(mantbits32) { 174 return sign | 0 175 } 176 // repeat expecting denormal 177 mant, exp, trunc = mant0, exp0, trunc0 178 for exp < bias32 { 179 trunc |= mant & 1 180 mant >>= 1 181 exp++ 182 } 183 if mant&1 != 0 && (trunc != 0 || mant&2 != 0) { 184 mant++ 185 } 186 mant >>= 1 187 exp++ 188 if mant < 1<<mantbits32 { 189 return sign | mant 190 } 191 } 192 return sign | uint32(exp-bias32)<<mantbits32 | mant&(1<<mantbits32-1) 193 } 194 195 func fadd64(f, g uint64) uint64 { 196 fs, fm, fe, fi, fn := funpack64(f) 197 gs, gm, ge, gi, gn := funpack64(g) 198 199 // Special cases. 200 switch { 201 case fn || gn: // NaN + x or x + NaN = NaN 202 return nan64 203 204 case fi && gi && fs != gs: // +Inf + -Inf or -Inf + +Inf = NaN 205 return nan64 206 207 case fi: // Inf + g = Inf 208 return f 209 210 case gi: // f + Inf = Inf 211 return g 212 213 case fm == 0 && gm == 0 && fs != 0 && gs != 0: // -0 + -0 = -0 214 return f 215 216 case fm == 0: // 0 + g = g but 0 + -0 = +0 217 if gm == 0 { 218 g ^= gs 219 } 220 return g 221 222 case gm == 0: // f + 0 = f 223 return f 224 225 } 226 227 if fe < ge || fe == ge && fm < gm { 228 f, g, fs, fm, fe, gs, gm, ge = g, f, gs, gm, ge, fs, fm, fe 229 } 230 231 shift := uint(fe - ge) 232 fm <<= 2 233 gm <<= 2 234 trunc := gm & (1<<shift - 1) 235 gm >>= shift 236 if fs == gs { 237 fm += gm 238 } else { 239 fm -= gm 240 if trunc != 0 { 241 fm-- 242 } 243 } 244 if fm == 0 { 245 fs = 0 246 } 247 return fpack64(fs, fm, fe-2, trunc) 248 } 249 250 func fsub64(f, g uint64) uint64 { 251 return fadd64(f, fneg64(g)) 252 } 253 254 func fneg64(f uint64) uint64 { 255 return f ^ (1 << (mantbits64 + expbits64)) 256 } 257 258 func fmul64(f, g uint64) uint64 { 259 fs, fm, fe, fi, fn := funpack64(f) 260 gs, gm, ge, gi, gn := funpack64(g) 261 262 // Special cases. 263 switch { 264 case fn || gn: // NaN * g or f * NaN = NaN 265 return nan64 266 267 case fi && gi: // Inf * Inf = Inf (with sign adjusted) 268 return f ^ gs 269 270 case fi && gm == 0, fm == 0 && gi: // 0 * Inf = Inf * 0 = NaN 271 return nan64 272 273 case fm == 0: // 0 * x = 0 (with sign adjusted) 274 return f ^ gs 275 276 case gm == 0: // x * 0 = 0 (with sign adjusted) 277 return g ^ fs 278 } 279 280 // 53-bit * 53-bit = 107- or 108-bit 281 lo, hi := mullu(fm, gm) 282 shift := mantbits64 - 1 283 trunc := lo & (1<<shift - 1) 284 mant := hi<<(64-shift) | lo>>shift 285 return fpack64(fs^gs, mant, fe+ge-1, trunc) 286 } 287 288 func fdiv64(f, g uint64) uint64 { 289 fs, fm, fe, fi, fn := funpack64(f) 290 gs, gm, ge, gi, gn := funpack64(g) 291 292 // Special cases. 293 switch { 294 case fn || gn: // NaN / g = f / NaN = NaN 295 return nan64 296 297 case fi && gi: // Inf / Inf = NaN 298 return nan64 299 300 case !fi && !gi && fm == 0 && gm == 0: // 0 / 0 = NaN 301 return nan64 302 303 case fi, !gi && gm == 0: // Inf / g = f / 0 = Inf 304 return fs ^ gs ^ inf64 305 306 case gi, fm == 0: // f / Inf = 0 / g = Inf 307 return fs ^ gs ^ 0 308 } 309 _, _, _, _ = fi, fn, gi, gn 310 311 // 53-bit<<54 / 53-bit = 53- or 54-bit. 312 shift := mantbits64 + 2 313 q, r := divlu(fm>>(64-shift), fm<<shift, gm) 314 return fpack64(fs^gs, q, fe-ge-2, r) 315 } 316 317 func f64to32(f uint64) uint32 { 318 fs, fm, fe, fi, fn := funpack64(f) 319 if fn { 320 return nan32 321 } 322 fs32 := uint32(fs >> 32) 323 if fi { 324 return fs32 ^ inf32 325 } 326 const d = mantbits64 - mantbits32 - 1 327 return fpack32(fs32, uint32(fm>>d), fe-1, uint32(fm&(1<<d-1))) 328 } 329 330 func f32to64(f uint32) uint64 { 331 const d = mantbits64 - mantbits32 332 fs, fm, fe, fi, fn := funpack32(f) 333 if fn { 334 return nan64 335 } 336 fs64 := uint64(fs) << 32 337 if fi { 338 return fs64 ^ inf64 339 } 340 return fpack64(fs64, uint64(fm)<<d, fe, 0) 341 } 342 343 func fcmp64(f, g uint64) (cmp int32, isnan bool) { 344 fs, fm, _, fi, fn := funpack64(f) 345 gs, gm, _, gi, gn := funpack64(g) 346 347 switch { 348 case fn, gn: // flag NaN 349 return 0, true 350 351 case !fi && !gi && fm == 0 && gm == 0: // 0 == 0 352 return 0, false 353 354 case fs > gs: // f < 0, g > 0 355 return -1, false 356 357 case fs < gs: // f > 0, g < 0 358 return +1, false 359 360 // Same sign, not NaN. 361 // Can compare encodings directly now. 362 // Reverse for sign. 363 case fs == 0 && f < g, fs != 0 && f > g: 364 return -1, false 365 366 case fs == 0 && f > g, fs != 0 && f < g: 367 return +1, false 368 } 369 370 // f == g 371 return 0, false 372 } 373 374 func f64toint(f uint64) (val int64, ok bool) { 375 fs, fm, fe, fi, fn := funpack64(f) 376 377 switch { 378 case fi, fn: // NaN 379 return 0, false 380 381 case fe < -1: // f < 0.5 382 return 0, false 383 384 case fe > 63: // f >= 2^63 385 if fs != 0 && fm == 0 { // f == -2^63 386 return -1 << 63, true 387 } 388 if fs != 0 { 389 return 0, false 390 } 391 return 0, false 392 } 393 394 for fe > int(mantbits64) { 395 fe-- 396 fm <<= 1 397 } 398 for fe < int(mantbits64) { 399 fe++ 400 fm >>= 1 401 } 402 val = int64(fm) 403 if fs != 0 { 404 val = -val 405 } 406 return val, true 407 } 408 409 func fintto64(val int64) (f uint64) { 410 fs := uint64(val) & (1 << 63) 411 mant := uint64(val) 412 if fs != 0 { 413 mant = -mant 414 } 415 return fpack64(fs, mant, int(mantbits64), 0) 416 } 417 418 // 64x64 -> 128 multiply. 419 // adapted from hacker's delight. 420 func mullu(u, v uint64) (lo, hi uint64) { 421 const ( 422 s = 32 423 mask = 1<<s - 1 424 ) 425 u0 := u & mask 426 u1 := u >> s 427 v0 := v & mask 428 v1 := v >> s 429 w0 := u0 * v0 430 t := u1*v0 + w0>>s 431 w1 := t & mask 432 w2 := t >> s 433 w1 += u0 * v1 434 return u * v, u1*v1 + w2 + w1>>s 435 } 436 437 // 128/64 -> 64 quotient, 64 remainder. 438 // adapted from hacker's delight 439 func divlu(u1, u0, v uint64) (q, r uint64) { 440 const b = 1 << 32 441 442 if u1 >= v { 443 return 1<<64 - 1, 1<<64 - 1 444 } 445 446 // s = nlz(v); v <<= s 447 s := uint(0) 448 for v&(1<<63) == 0 { 449 s++ 450 v <<= 1 451 } 452 453 vn1 := v >> 32 454 vn0 := v & (1<<32 - 1) 455 un32 := u1<<s | u0>>(64-s) 456 un10 := u0 << s 457 un1 := un10 >> 32 458 un0 := un10 & (1<<32 - 1) 459 q1 := un32 / vn1 460 rhat := un32 - q1*vn1 461 462 again1: 463 if q1 >= b || q1*vn0 > b*rhat+un1 { 464 q1-- 465 rhat += vn1 466 if rhat < b { 467 goto again1 468 } 469 } 470 471 un21 := un32*b + un1 - q1*v 472 q0 := un21 / vn1 473 rhat = un21 - q0*vn1 474 475 again2: 476 if q0 >= b || q0*vn0 > b*rhat+un0 { 477 q0-- 478 rhat += vn1 479 if rhat < b { 480 goto again2 481 } 482 } 483 484 return q1*b + q0, (un21*b + un0 - q0*v) >> s 485 } 486