| 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 | #include "../int_math.h" |
| 6 | #include "DD.h" |
| 7 | |
| 8 | #define makeFinite(x) \ |
| 9 | { \ |
| 10 | (x).s.hi = crt_copysign(crt_isinf((x).s.hi) ? 1.0 : 0.0, (x).s.hi); \ |
| 11 | (x).s.lo = 0.0; \ |
| 12 | } |
| 13 | |
| 14 | #define zeroNaN(x) \ |
| 15 | { \ |
| 16 | if (crt_isnan((x).s.hi)) { \ |
| 17 | (x).s.hi = crt_copysign(0.0, (x).s.hi); \ |
| 18 | (x).s.lo = 0.0; \ |
| 19 | } \ |
| 20 | } |
| 21 | |
| 22 | long double _Complex __multc3(long double a, long double b, long double c, |
| 23 | long double d) { |
| 24 | long double ac = __gcc_qmul(a, c); |
| 25 | long double bd = __gcc_qmul(b, d); |
| 26 | long double ad = __gcc_qmul(a, d); |
| 27 | long double bc = __gcc_qmul(b, c); |
| 28 | |
| 29 | DD real = {.ld = __gcc_qsub(ac, bd)}; |
| 30 | DD imag = {.ld = __gcc_qadd(ad, bc)}; |
| 31 | |
| 32 | if (crt_isnan(real.s.hi) && crt_isnan(imag.s.hi)) { |
| 33 | int recalc = 0; |
| 34 | |
| 35 | DD aDD = {.ld = a}; |
| 36 | DD bDD = {.ld = b}; |
| 37 | DD cDD = {.ld = c}; |
| 38 | DD dDD = {.ld = d}; |
| 39 | |
| 40 | if (crt_isinf(aDD.s.hi) || crt_isinf(bDD.s.hi)) { |
| 41 | makeFinite(aDD); |
| 42 | makeFinite(bDD); |
| 43 | zeroNaN(cDD); |
| 44 | zeroNaN(dDD); |
| 45 | recalc = 1; |
| 46 | } |
| 47 | |
| 48 | if (crt_isinf(cDD.s.hi) || crt_isinf(dDD.s.hi)) { |
| 49 | makeFinite(cDD); |
| 50 | makeFinite(dDD); |
| 51 | zeroNaN(aDD); |
| 52 | zeroNaN(bDD); |
| 53 | recalc = 1; |
| 54 | } |
| 55 | |
| 56 | if (!recalc) { |
| 57 | DD acDD = {.ld = ac}; |
| 58 | DD bdDD = {.ld = bd}; |
| 59 | DD adDD = {.ld = ad}; |
| 60 | DD bcDD = {.ld = bc}; |
| 61 | |
| 62 | if (crt_isinf(acDD.s.hi) || crt_isinf(bdDD.s.hi) || |
| 63 | crt_isinf(adDD.s.hi) || crt_isinf(bcDD.s.hi)) { |
| 64 | zeroNaN(aDD); |
| 65 | zeroNaN(bDD); |
| 66 | zeroNaN(cDD); |
| 67 | zeroNaN(dDD); |
| 68 | recalc = 1; |
| 69 | } |
| 70 | } |
| 71 | |
| 72 | if (recalc) { |
| 73 | real.s.hi = CRT_INFINITY * (aDD.s.hi * cDD.s.hi - bDD.s.hi * dDD.s.hi); |
| 74 | real.s.lo = 0.0; |
| 75 | imag.s.hi = CRT_INFINITY * (aDD.s.hi * dDD.s.hi + bDD.s.hi * cDD.s.hi); |
| 76 | imag.s.lo = 0.0; |
| 77 | } |
| 78 | } |
| 79 | |
| 80 | long double _Complex z; |
| 81 | __real__ z = real.ld; |
| 82 | __imag__ z = imag.ld; |
| 83 | |
| 84 | return z; |
| 85 | } |
| 86 | |