1 /* Copyright (C) 2005-2022 Free Software Foundation, Inc.
2 This file is part of the GNU C Library.
3
4 The GNU C Library is free software; you can redistribute it and/or
5 modify it under the terms of the GNU Lesser General Public
6 License as published by the Free Software Foundation; either
7 version 2.1 of the License, or (at your option) any later version.
8
9 The GNU C Library is distributed in the hope that it will be useful,
10 but WITHOUT ANY WARRANTY; without even the implied warranty of
11 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 Lesser General Public License for more details.
13
14 You should have received a copy of the GNU Lesser General Public
15 License along with the GNU C Library; if not, see
16 <https://www.gnu.org/licenses/>. */
17
18 #include <stdbool.h>
19 #include <math.h>
20 #include <complex.h>
21
22 attribute_hidden
23 long double _Complex
__multc3(long double a,long double b,long double c,long double d)24 __multc3 (long double a, long double b, long double c, long double d)
25 {
26 long double ac, bd, ad, bc, x, y;
27
28 ac = a * c;
29 bd = b * d;
30 ad = a * d;
31 bc = b * c;
32
33 x = ac - bd;
34 y = ad + bc;
35
36 if (isnan (x) && isnan (y))
37 {
38 /* Recover infinities that computed as NaN + iNaN. */
39 bool recalc = 0;
40 if (isinf (a) || isinf (b))
41 {
42 /* z is infinite. "Box" the infinity and change NaNs in
43 the other factor to 0. */
44 a = copysignl (isinf (a) ? 1 : 0, a);
45 b = copysignl (isinf (b) ? 1 : 0, b);
46 if (isnan (c)) c = copysignl (0, c);
47 if (isnan (d)) d = copysignl (0, d);
48 recalc = 1;
49 }
50 if (isinf (c) || isinf (d))
51 {
52 /* w is infinite. "Box" the infinity and change NaNs in
53 the other factor to 0. */
54 c = copysignl (isinf (c) ? 1 : 0, c);
55 d = copysignl (isinf (d) ? 1 : 0, d);
56 if (isnan (a)) a = copysignl (0, a);
57 if (isnan (b)) b = copysignl (0, b);
58 recalc = 1;
59 }
60 if (!recalc
61 && (isinf (ac) || isinf (bd)
62 || isinf (ad) || isinf (bc)))
63 {
64 /* Recover infinities from overflow by changing NaNs to 0. */
65 if (isnan (a)) a = copysignl (0, a);
66 if (isnan (b)) b = copysignl (0, b);
67 if (isnan (c)) c = copysignl (0, c);
68 if (isnan (d)) d = copysignl (0, d);
69 recalc = 1;
70 }
71 if (recalc)
72 {
73 x = INFINITY * (a * c - b * d);
74 y = INFINITY * (a * d + b * c);
75 }
76 }
77
78 return x + I * y;
79 }
80