1 /* This file is distributed under the University of Illinois Open Source 2 * License. See LICENSE.TXT for details. 3 */ 4 5 /* long double __gcc_qmul(long double x, long double y); 6 * This file implements the PowerPC 128-bit double-double multiply operation. 7 * This implementation is shamelessly cribbed from Apple's DDRT, circa 1993(!) 8 */ 9 10 #include "DD.h" 11 12 long double __gcc_qmul(long double x, long double y) 13 { 14 static const uint32_t infinityHi = UINT32_C(0x7ff00000); 15 DD dst = { .ld = x }, src = { .ld = y }; 16 17 register double A = dst.s.hi, a = dst.s.lo, 18 B = src.s.hi, b = src.s.lo; 19 20 double aHi, aLo, bHi, bLo; 21 double ab, tmp, tau; 22 23 ab = A * B; 24 25 /* Detect special cases */ 26 if (ab == 0.0) { 27 dst.s.hi = ab; 28 dst.s.lo = 0.0; 29 return dst.ld; 30 } 31 32 const doublebits abBits = { .d = ab }; 33 if (((uint32_t)(abBits.x >> 32) & infinityHi) == infinityHi) { 34 dst.s.hi = ab; 35 dst.s.lo = 0.0; 36 return dst.ld; 37 } 38 39 /* Generic cases handled here. */ 40 aHi = high26bits(A); 41 bHi = high26bits(B); 42 aLo = A - aHi; 43 bLo = B - bHi; 44 45 tmp = LOWORDER(ab, aHi, aLo, bHi, bLo); 46 tmp += (A * b + a * B); 47 tau = ab + tmp; 48 49 dst.s.lo = (ab - tau) + tmp; 50 dst.s.hi = tau; 51 52 return dst.ld; 53 } 54