]>
Commit | Line | Data |
---|---|---|
bc9f6000 | 1 | /* Compute remainder and a congruent to the quotient. |
2b778ceb | 2 | Copyright (C) 1997-2021 Free Software Foundation, Inc. |
bc9f6000 UD |
3 | This file is part of the GNU C Library. |
4 | Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. | |
5 | ||
6 | The GNU C Library is free software; you can redistribute it and/or | |
41bdb6e2 AJ |
7 | modify it under the terms of the GNU Lesser General Public |
8 | License as published by the Free Software Foundation; either | |
9 | version 2.1 of the License, or (at your option) any later version. | |
bc9f6000 UD |
10 | |
11 | The GNU C Library is distributed in the hope that it will be useful, | |
12 | but WITHOUT ANY WARRANTY; without even the implied warranty of | |
13 | MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU | |
41bdb6e2 | 14 | Lesser General Public License for more details. |
bc9f6000 | 15 | |
41bdb6e2 | 16 | You should have received a copy of the GNU Lesser General Public |
59ba27a6 | 17 | License along with the GNU C Library; if not, see |
5a82c748 | 18 | <https://www.gnu.org/licenses/>. */ |
bc9f6000 UD |
19 | |
20 | #include <math.h> | |
21 | ||
1ed0291c | 22 | #include <math_private.h> |
86f9568a | 23 | #include <libm-alias-ldouble.h> |
bc9f6000 UD |
24 | |
25 | ||
63551311 | 26 | static const long double zero = 0.0; |
bc9f6000 UD |
27 | |
28 | ||
29 | long double | |
63551311 | 30 | __remquol (long double x, long double p, int *quo) |
bc9f6000 UD |
31 | { |
32 | int32_t ex,ep,hx,hp; | |
24ab7723 | 33 | uint32_t sx,lx,lp; |
63551311 | 34 | int cquo,qs; |
bc9f6000 UD |
35 | |
36 | GET_LDOUBLE_WORDS (ex, hx, lx, x); | |
37 | GET_LDOUBLE_WORDS (ep, hp, lp, p); | |
38 | sx = ex & 0x8000; | |
39 | qs = (sx ^ (ep & 0x8000)) >> 15; | |
40 | ep &= 0x7fff; | |
41 | ex &= 0x7fff; | |
42 | ||
43 | /* Purge off exception values. */ | |
44 | if ((ep | hp | lp) == 0) | |
45 | return (x * p) / (x * p); /* p = 0 */ | |
46 | if ((ex == 0x7fff) /* x not finite */ | |
47 | || ((ep == 0x7fff) /* p is NaN */ | |
3ce2232e | 48 | && (((hp & 0x7fffffff) | lp) != 0))) |
bc9f6000 UD |
49 | return (x * p) / (x * p); |
50 | ||
51 | if (ep <= 0x7ffb) | |
63551311 | 52 | x = __ieee754_fmodl (x, 8 * p); /* now x < 8p */ |
bc9f6000 UD |
53 | |
54 | if (((ex - ep) | (hx - hp) | (lx - lp)) == 0) | |
55 | { | |
56 | *quo = qs ? -1 : 1; | |
57 | return zero * x; | |
58 | } | |
59 | ||
60 | x = fabsl (x); | |
61 | p = fabsl (p); | |
62 | cquo = 0; | |
63 | ||
a820f9b3 | 64 | if (ep <= 0x7ffc && x >= 4 * p) |
bc9f6000 UD |
65 | { |
66 | x -= 4 * p; | |
63551311 UD |
67 | cquo += 4; |
68 | } | |
a820f9b3 | 69 | if (ep <= 0x7ffd && x >= 2 * p) |
63551311 UD |
70 | { |
71 | x -= 2 * p; | |
bc9f6000 UD |
72 | cquo += 2; |
73 | } | |
bc9f6000 UD |
74 | |
75 | if (ep < 0x0002) | |
76 | { | |
77 | if (x + x > p) | |
78 | { | |
79 | x -= p; | |
d705269e | 80 | ++cquo; |
bc9f6000 | 81 | if (x + x >= p) |
63551311 UD |
82 | { |
83 | x -= p; | |
84 | ++cquo; | |
85 | } | |
bc9f6000 UD |
86 | } |
87 | } | |
88 | else | |
89 | { | |
90 | long double p_half = 0.5 * p; | |
63551311 | 91 | if (x > p_half) |
bc9f6000 UD |
92 | { |
93 | x -= p; | |
d705269e | 94 | ++cquo; |
bc9f6000 | 95 | if (x >= p_half) |
63551311 UD |
96 | { |
97 | x -= p; | |
98 | ++cquo; | |
99 | } | |
bc9f6000 UD |
100 | } |
101 | } | |
102 | ||
103 | *quo = qs ? -cquo : cquo; | |
104 | ||
ce8fc784 JM |
105 | /* Ensure correct sign of zero result in round-downward mode. */ |
106 | if (x == 0.0L) | |
107 | x = 0.0L; | |
63551311 UD |
108 | if (sx) |
109 | x = -x; | |
bc9f6000 UD |
110 | return x; |
111 | } | |
86f9568a | 112 | libm_alias_ldouble (__remquo, remquo) |