1 // Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions.
2 // See https://llvm.org/LICENSE.txt for license information.
3 // SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception
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 #include "DD.h"
10 
11 long double __gcc_qmul(long double x, long double y) {
12   static const uint32_t infinityHi = UINT32_C(0x7ff00000);
13   DD dst = {.ld = x}, src = {.ld = y};
14 
15   register double A = dst.s.hi, a = dst.s.lo, B = src.s.hi, b = src.s.lo;
16 
17   double aHi, aLo, bHi, bLo;
18   double ab, tmp, tau;
19 
20   ab = A * B;
21 
22   // Detect special cases
23   if (ab == 0.0) {
24     dst.s.hi = ab;
25     dst.s.lo = 0.0;
26     return dst.ld;
27   }
28 
29   const doublebits abBits = {.d = ab};
30   if (((uint32_t)(abBits.x >> 32) & infinityHi) == infinityHi) {
31     dst.s.hi = ab;
32     dst.s.lo = 0.0;
33     return dst.ld;
34   }
35 
36   // Generic cases handled here.
37   aHi = high26bits(A);
38   bHi = high26bits(B);
39   aLo = A - aHi;
40   bLo = B - bHi;
41 
42   tmp = LOWORDER(ab, aHi, aLo, bHi, bLo);
43   tmp += (A * b + a * B);
44   tau = ab + tmp;
45 
46   dst.s.lo = (ab - tau) + tmp;
47   dst.s.hi = tau;
48 
49   return dst.ld;
50 }
51