|
16519
|
1 /* Hypotenuse of a right-angled triangle. |
|
17848
|
2 Copyright (C) 2012-2015 Free Software Foundation, Inc. |
|
16519
|
3 |
|
|
4 This program is free software: you can redistribute it and/or modify |
|
|
5 it under the terms of the GNU General Public License as published by |
|
|
6 the Free Software Foundation; either version 3 of the License, or |
|
|
7 (at your option) any later version. |
|
|
8 |
|
|
9 This program 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 |
|
|
12 GNU General Public License for more details. |
|
|
13 |
|
|
14 You should have received a copy of the GNU General Public License |
|
|
15 along with this program. If not, see <http://www.gnu.org/licenses/>. */ |
|
|
16 |
|
|
17 /* Written by Bruno Haible <bruno@clisp.org>, 2012. */ |
|
|
18 |
|
|
19 #include <config.h> |
|
|
20 |
|
|
21 /* Specification. */ |
|
|
22 #include <math.h> |
|
|
23 |
|
|
24 #if HAVE_SAME_LONG_DOUBLE_AS_DOUBLE |
|
|
25 |
|
|
26 long double |
|
|
27 hypotl (long double x, long double y) |
|
|
28 { |
|
|
29 return hypot (x, y); |
|
|
30 } |
|
|
31 |
|
|
32 #else |
|
|
33 |
|
|
34 long double |
|
|
35 hypotl (long double x, long double y) |
|
|
36 { |
|
|
37 if (isfinite (x) && isfinite (y)) |
|
|
38 { |
|
|
39 /* Determine absolute values. */ |
|
|
40 x = fabsl (x); |
|
|
41 y = fabsl (y); |
|
|
42 |
|
|
43 { |
|
|
44 /* Find the bigger and the smaller one. */ |
|
|
45 long double a; |
|
|
46 long double b; |
|
|
47 |
|
|
48 if (x >= y) |
|
|
49 { |
|
|
50 a = x; |
|
|
51 b = y; |
|
|
52 } |
|
|
53 else |
|
|
54 { |
|
|
55 a = y; |
|
|
56 b = x; |
|
|
57 } |
|
|
58 /* Now 0 <= b <= a. */ |
|
|
59 |
|
|
60 { |
|
|
61 int e; |
|
|
62 long double an; |
|
|
63 long double bn; |
|
|
64 |
|
|
65 /* Write a = an * 2^e, b = bn * 2^e with 0 <= bn <= an < 1. */ |
|
|
66 an = frexpl (a, &e); |
|
|
67 bn = ldexpl (b, - e); |
|
|
68 |
|
|
69 { |
|
|
70 long double cn; |
|
|
71 |
|
|
72 /* Through the normalization, no unneeded overflow or underflow |
|
|
73 will occur here. */ |
|
|
74 cn = sqrtl (an * an + bn * bn); |
|
|
75 return ldexpl (cn, e); |
|
|
76 } |
|
|
77 } |
|
|
78 } |
|
|
79 } |
|
|
80 else |
|
|
81 { |
|
|
82 if (isinf (x) || isinf (y)) |
|
|
83 /* x or y is infinite. Return +Infinity. */ |
|
|
84 return HUGE_VALL; |
|
|
85 else |
|
|
86 /* x or y is NaN. Return NaN. */ |
|
|
87 return x + y; |
|
|
88 } |
|
|
89 } |
|
|
90 |
|
|
91 #endif |
