1 /* Inline math functions for i387 and SSE.
2 Copyright (C) 1995-2012 Free Software Foundation, Inc.
3 This file is part of the GNU C Library.
4
5 The GNU C Library is free software; you can redistribute it and/or
6 modify it under the terms of the GNU Lesser General Public
7 License as published by the Free Software Foundation; either
8 version 2.1 of the License, or (at your option) any later version.
9
10 The GNU C Library is distributed in the hope that it will be useful,
11 but WITHOUT ANY WARRANTY; without even the implied warranty of
12 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
13 Lesser General Public License for more details.
14
15 You should have received a copy of the GNU Lesser General Public
16 License along with the GNU C Library; if not, see
17 <http://www.gnu.org/licenses/>. */
18
19 #ifndef _MATH_H
20 # error "Never use <bits/mathinline.h> directly; include <math.h> instead."
21 #endif
22
23 #ifndef __extern_always_inline
24 # define __MATH_INLINE __inline
25 #else
26 # define __MATH_INLINE __extern_always_inline
27 #endif
28
29
30 #if defined __USE_ISOC99 && defined __GNUC__ && __GNUC__ >= 2
31 /* GCC 2.97 and up have builtins that actually can be used. */
32 # if !__GNUC_PREREQ (2,97)
33 /* ISO C99 defines some macros to perform unordered comparisons. The
34 ix87 FPU supports this with special opcodes and we should use them.
35 These must not be inline functions since we have to be able to handle
36 all floating-point types. */
37 # undef isgreater
38 # undef isgreaterequal
39 # undef isless
40 # undef islessequal
41 # undef islessgreater
42 # undef isunordered
43 # ifdef __i686__
44 /* For the PentiumPro and more recent processors we can provide
45 better code. */
46 # define isgreater(x, y) \
47 ({ register char __result; \
48 __asm__ ("fucomip %%st(1), %%st; seta %%al" \
49 : "=a" (__result) : "u" (y), "t" (x) : "cc", "st"); \
50 __result; })
51 # define isgreaterequal(x, y) \
52 ({ register char __result; \
53 __asm__ ("fucomip %%st(1), %%st; setae %%al" \
54 : "=a" (__result) : "u" (y), "t" (x) : "cc", "st"); \
55 __result; })
56
57 # define isless(x, y) \
58 ({ register char __result; \
59 __asm__ ("fucomip %%st(1), %%st; seta %%al" \
60 : "=a" (__result) : "u" (x), "t" (y) : "cc", "st"); \
61 __result; })
62
63 # define islessequal(x, y) \
64 ({ register char __result; \
65 __asm__ ("fucomip %%st(1), %%st; setae %%al" \
66 : "=a" (__result) : "u" (x), "t" (y) : "cc", "st"); \
67 __result; })
68
69 # define islessgreater(x, y) \
70 ({ register char __result; \
71 __asm__ ("fucomip %%st(1), %%st; setne %%al" \
72 : "=a" (__result) : "u" (y), "t" (x) : "cc", "st"); \
73 __result; })
74
75 # define isunordered(x, y) \
76 ({ register char __result; \
77 __asm__ ("fucomip %%st(1), %%st; setp %%al" \
78 : "=a" (__result) : "u" (y), "t" (x) : "cc", "st"); \
79 __result; })
80 # else
81 /* This is the dumb, portable code for i386 and above. */
82 # define isgreater(x, y) \
83 ({ register char __result; \
84 __asm__ ("fucompp; fnstsw; testb $0x45, %%ah; setz %%al" \
85 : "=a" (__result) : "u" (y), "t" (x) : "cc", "st", "st(1)"); \
86 __result; })
87
88 # define isgreaterequal(x, y) \
89 ({ register char __result; \
90 __asm__ ("fucompp; fnstsw; testb $0x05, %%ah; setz %%al" \
91 : "=a" (__result) : "u" (y), "t" (x) : "cc", "st", "st(1)"); \
92 __result; })
93
94 # define isless(x, y) \
95 ({ register char __result; \
96 __asm__ ("fucompp; fnstsw; testb $0x45, %%ah; setz %%al" \
97 : "=a" (__result) : "u" (x), "t" (y) : "cc", "st", "st(1)"); \
98 __result; })
99
100 # define islessequal(x, y) \
101 ({ register char __result; \
102 __asm__ ("fucompp; fnstsw; testb $0x05, %%ah; setz %%al" \
103 : "=a" (__result) : "u" (x), "t" (y) : "cc", "st", "st(1)"); \
104 __result; })
105
106 # define islessgreater(x, y) \
107 ({ register char __result; \
108 __asm__ ("fucompp; fnstsw; testb $0x44, %%ah; setz %%al" \
109 : "=a" (__result) : "u" (y), "t" (x) : "cc", "st", "st(1)"); \
110 __result; })
111
112 # define isunordered(x, y) \
113 ({ register char __result; \
114 __asm__ ("fucompp; fnstsw; sahf; setp %%al" \
115 : "=a" (__result) : "u" (y), "t" (x) : "cc", "st", "st(1)"); \
116 __result; })
117 # endif /* __i686__ */
118 # endif /* GCC 2.97 */
119
120 /* The gcc, version 2.7 or below, has problems with all this inlining
121 code. So disable it for this version of the compiler. */
122 # if __GNUC_PREREQ (2, 8)
123 __BEGIN_NAMESPACE_C99
124
125 /* Test for negative number. Used in the signbit() macro. */
126 __MATH_INLINE int
__NTH(__signbitf (float __x))127 __NTH (__signbitf (float __x))
128 {
129 # ifdef __SSE2_MATH__
130 int __m;
131 __asm ("pmovmskb %1, %0" : "=r" (__m) : "x" (__x));
132 return (__m & 0x8) != 0;
133 # else
134 __extension__ union { float __f; int __i; } __u = { __f: __x };
135 return __u.__i < 0;
136 # endif
137 }
138 __MATH_INLINE int
__NTH(__signbit (double __x))139 __NTH (__signbit (double __x))
140 {
141 # ifdef __SSE2_MATH__
142 int __m;
143 __asm ("pmovmskb %1, %0" : "=r" (__m) : "x" (__x));
144 return (__m & 0x80) != 0;
145 # else
146 __extension__ union { double __d; int __i[2]; } __u = { __d: __x };
147 return __u.__i[1] < 0;
148 # endif
149 }
150 __MATH_INLINE int
__NTH(__signbitl (long double __x))151 __NTH (__signbitl (long double __x))
152 {
153 __extension__ union { long double __l; int __i[3]; } __u = { __l: __x };
154 return (__u.__i[2] & 0x8000) != 0;
155 }
156
157 __END_NAMESPACE_C99
158 # endif
159 #endif
160
161
162 /* The gcc, version 2.7 or below, has problems with all this inlining
163 code. So disable it for this version of the compiler. */
164 #if __GNUC_PREREQ (2, 8)
165 # if !__GNUC_PREREQ (3, 4) && !defined __NO_MATH_INLINES \
166 && defined __OPTIMIZE__
167 /* GCC 3.4 introduced builtins for all functions below, so
168 there's no need to define any of these inline functions. */
169
170 # ifdef __USE_ISOC99
171 __BEGIN_NAMESPACE_C99
172
173 /* Round to nearest integer. */
174 # ifdef __SSE_MATH__
175 __MATH_INLINE long int
__NTH(lrintf (float __x))176 __NTH (lrintf (float __x))
177 {
178 long int __res;
179 /* Mark as volatile since the result is dependent on the state of
180 the SSE control register (the rounding mode). Otherwise GCC might
181 remove these assembler instructions since it does not know about
182 the rounding mode change and cannot currently be told. */
183 __asm __volatile__ ("cvtss2si %1, %0" : "=r" (__res) : "xm" (__x));
184 return __res;
185 }
186 # endif
187 # ifdef __SSE2_MATH__
188 __MATH_INLINE long int
__NTH(lrint (double __x))189 __NTH (lrint (double __x))
190 {
191 long int __res;
192 /* Mark as volatile since the result is dependent on the state of
193 the SSE control register (the rounding mode). Otherwise GCC might
194 remove these assembler instructions since it does not know about
195 the rounding mode change and cannot currently be told. */
196 __asm __volatile__ ("cvtsd2si %1, %0" : "=r" (__res) : "xm" (__x));
197 return __res;
198 }
199 # endif
200 # ifdef __x86_64__
201 __MATH_INLINE long long int
__NTH(llrintf (float __x))202 __NTH (llrintf (float __x))
203 {
204 long long int __res;
205 /* Mark as volatile since the result is dependent on the state of
206 the SSE control register (the rounding mode). Otherwise GCC might
207 remove these assembler instructions since it does not know about
208 the rounding mode change and cannot currently be told. */
209 __asm __volatile__ ("cvtss2si %1, %0" : "=r" (__res) : "xm" (__x));
210 return __res;
211 }
212 __MATH_INLINE long long int
__NTH(llrint (double __x))213 __NTH (llrint (double __x))
214 {
215 long long int __res;
216 /* Mark as volatile since the result is dependent on the state of
217 the SSE control register (the rounding mode). Otherwise GCC might
218 remove these assembler instructions since it does not know about
219 the rounding mode change and cannot currently be told. */
220 __asm __volatile__ ("cvtsd2si %1, %0" : "=r" (__res) : "xm" (__x));
221 return __res;
222 }
223 # endif
224
225 # if defined __FINITE_MATH_ONLY__ && __FINITE_MATH_ONLY__ > 0 \
226 && defined __SSE2_MATH__
227 /* Determine maximum of two values. */
228 __MATH_INLINE float
__NTH(fmaxf (float __x,float __y))229 __NTH (fmaxf (float __x, float __y))
230 {
231 # ifdef __AVX__
232 float __res;
233 __asm ("vmaxss %2, %1, %0" : "=x" (__res) : "x" (x), "xm" (__y));
234 return __res;
235 # else
236 __asm ("maxss %1, %0" : "+x" (__x) : "xm" (__y));
237 return __x;
238 # endif
239 }
240 __MATH_INLINE double
__NTH(fmax (double __x,double __y))241 __NTH (fmax (double __x, double __y))
242 {
243 # ifdef __AVX__
244 float __res;
245 __asm ("vmaxsd %2, %1, %0" : "=x" (__res) : "x" (x), "xm" (__y));
246 return __res;
247 # else
248 __asm ("maxsd %1, %0" : "+x" (__x) : "xm" (__y));
249 return __x;
250 # endif
251 }
252
253 /* Determine minimum of two values. */
254 __MATH_INLINE float
__NTH(fminf (float __x,float __y))255 __NTH (fminf (float __x, float __y))
256 {
257 # ifdef __AVX__
258 float __res;
259 __asm ("vminss %2, %1, %0" : "=x" (__res) : "x" (x), "xm" (__y));
260 return __res;
261 # else
262 __asm ("minss %1, %0" : "+x" (__x) : "xm" (__y));
263 return __x;
264 # endif
265 }
266 __MATH_INLINE double
__NTH(fmin (double __x,double __y))267 __NTH (fmin (double __x, double __y))
268 {
269 # ifdef __AVX__
270 float __res;
271 __asm ("vminsd %2, %1, %0" : "=x" (__res) : "x" (x), "xm" (__y));
272 return __res;
273 # else
274 __asm ("minsd %1, %0" : "+x" (__x) : "xm" (__y));
275 return __x;
276 # endif
277 }
278 # endif
279
280 __END_NAMESPACE_C99
281 # endif
282
283 # if defined __SSE4_1__ && defined __SSE2_MATH__
284 # if defined __USE_MISC || defined __USE_XOPEN_EXTENDED || defined __USE_ISOC99
285 __BEGIN_NAMESPACE_C99
286
287 /* Round to nearest integer. */
288 __MATH_INLINE double
__NTH(rint (double __x))289 __NTH (rint (double __x))
290 {
291 double __res;
292 /* Mark as volatile since the result is dependent on the state of
293 the SSE control register (the rounding mode). Otherwise GCC might
294 remove these assembler instructions since it does not know about
295 the rounding mode change and cannot currently be told. */
296 __asm __volatile__ ("roundsd $4, %1, %0" : "=x" (__res) : "xm" (__x));
297 return __res;
298 }
299 __MATH_INLINE float
__NTH(rintf (float __x))300 __NTH (rintf (float __x))
301 {
302 float __res;
303 /* Mark as volatile since the result is dependent on the state of
304 the SSE control register (the rounding mode). Otherwise GCC might
305 remove these assembler instructions since it does not know about
306 the rounding mode change and cannot currently be told. */
307 __asm __volatile__ ("roundss $4, %1, %0" : "=x" (__res) : "xm" (__x));
308 return __res;
309 }
310
311 # ifdef __USE_ISOC99
312 /* Round to nearest integer without raising inexact exception. */
313 __MATH_INLINE double
__NTH(nearbyint (double __x))314 __NTH (nearbyint (double __x))
315 {
316 double __res;
317 /* Mark as volatile since the result is dependent on the state of
318 the SSE control register (the rounding mode). Otherwise GCC might
319 remove these assembler instructions since it does not know about
320 the rounding mode change and cannot currently be told. */
321 __asm __volatile__ ("roundsd $0xc, %1, %0" : "=x" (__res) : "xm" (__x));
322 return __res;
323 }
324 __MATH_INLINE float
__NTH(nearbyintf (float __x))325 __NTH (nearbyintf (float __x))
326 {
327 float __res;
328 /* Mark as volatile since the result is dependent on the state of
329 the SSE control register (the rounding mode). Otherwise GCC might
330 remove these assembler instructions since it does not know about
331 the rounding mode change and cannot currently be told. */
332 __asm __volatile__ ("roundss $0xc, %1, %0" : "=x" (__res) : "xm" (__x));
333 return __res;
334 }
335 # endif
336
337 __END_NAMESPACE_C99
338 # endif
339
340 __BEGIN_NAMESPACE_STD
341 /* Smallest integral value not less than X. */
342 __MATH_INLINE double
__NTH(ceil (double __x))343 __NTH (ceil (double __x))
344 {
345 double __res;
346 __asm ("roundsd $2, %1, %0" : "=x" (__res) : "xm" (__x));
347 return __res;
348 }
349 __END_NAMESPACE_STD
350
351 __BEGIN_NAMESPACE_C99
352 __MATH_INLINE float
__NTH(ceilf (float __x))353 __NTH (ceilf (float __x))
354 {
355 float __res;
356 __asm ("roundss $2, %1, %0" : "=x" (__res) : "xm" (__x));
357 return __res;
358 }
359 __END_NAMESPACE_C99
360
361 __BEGIN_NAMESPACE_STD
362 /* Largest integer not greater than X. */
363 __MATH_INLINE double
__NTH(floor (double __x))364 __NTH (floor (double __x))
365 {
366 double __res;
367 __asm ("roundsd $1, %1, %0" : "=x" (__res) : "xm" (__x));
368 return __res;
369 }
370 __END_NAMESPACE_STD
371
372 __BEGIN_NAMESPACE_C99
373 __MATH_INLINE float
__NTH(floorf (float __x))374 __NTH (floorf (float __x))
375 {
376 float __res;
377 __asm ("roundss $1, %1, %0" : "=x" (__res) : "xm" (__x));
378 return __res;
379 }
380 __END_NAMESPACE_C99
381 # endif
382 # endif
383 #endif
384
385 #ifndef __x86_64__
386 # if ((!defined __NO_MATH_INLINES || defined __LIBC_INTERNAL_MATH_INLINES) \
387 && defined __OPTIMIZE__)
388
389 /* The inline functions do not set errno or raise necessarily the
390 correct exceptions. */
391 # undef math_errhandling
392
393 /* A macro to define float, double, and long double versions of various
394 math functions for the ix87 FPU. FUNC is the function name (which will
395 be suffixed with f and l for the float and long double version,
396 respectively). OP is the name of the FPU operation.
397 We define two sets of macros. The set with the additional NP
398 doesn't add a prototype declaration. */
399
400 # if defined __USE_MISC || defined __USE_ISOC99
401 # define __inline_mathop(func, op) \
402 __inline_mathop_ (double, func, op) \
403 __inline_mathop_ (float, __CONCAT(func,f), op) \
404 __inline_mathop_ (long double, __CONCAT(func,l), op)
405 # define __inline_mathopNP(func, op) \
406 __inline_mathopNP_ (double, func, op) \
407 __inline_mathopNP_ (float, __CONCAT(func,f), op) \
408 __inline_mathopNP_ (long double, __CONCAT(func,l), op)
409 # else
410 # define __inline_mathop(func, op) \
411 __inline_mathop_ (double, func, op)
412 # define __inline_mathopNP(func, op) \
413 __inline_mathopNP_ (double, func, op)
414 # endif
415
416 # define __inline_mathop_(float_type, func, op) \
417 __inline_mathop_decl_ (float_type, func, op, "0" (__x))
418 # define __inline_mathopNP_(float_type, func, op) \
419 __inline_mathop_declNP_ (float_type, func, op, "0" (__x))
420
421
422 # if defined __USE_MISC || defined __USE_ISOC99
423 # define __inline_mathop_decl(func, op, params...) \
424 __inline_mathop_decl_ (double, func, op, params) \
425 __inline_mathop_decl_ (float, __CONCAT(func,f), op, params) \
426 __inline_mathop_decl_ (long double, __CONCAT(func,l), op, params)
427 # define __inline_mathop_declNP(func, op, params...) \
428 __inline_mathop_declNP_ (double, func, op, params) \
429 __inline_mathop_declNP_ (float, __CONCAT(func,f), op, params) \
430 __inline_mathop_declNP_ (long double, __CONCAT(func,l), op, params)
431 # else
432 # define __inline_mathop_decl(func, op, params...) \
433 __inline_mathop_decl_ (double, func, op, params)
434 # define __inline_mathop_declNP(func, op, params...) \
435 __inline_mathop_declNP_ (double, func, op, params)
436 # endif
437
438 # define __inline_mathop_decl_(float_type, func, op, params...) \
439 __MATH_INLINE float_type func (float_type) __THROW; \
440 __inline_mathop_declNP_ (float_type, func, op, params)
441
442 # define __inline_mathop_declNP_(float_type, func, op, params...) \
443 __MATH_INLINE float_type __NTH (func (float_type __x)) \
444 { \
445 register float_type __result; \
446 __asm __volatile__ (op : "=t" (__result) : params); \
447 return __result; \
448 }
449
450
451 # if defined __USE_MISC || defined __USE_ISOC99
452 # define __inline_mathcode(func, arg, code) \
453 __inline_mathcode_ (double, func, arg, code) \
454 __inline_mathcode_ (float, __CONCAT(func,f), arg, code) \
455 __inline_mathcode_ (long double, __CONCAT(func,l), arg, code)
456 # define __inline_mathcodeNP(func, arg, code) \
457 __inline_mathcodeNP_ (double, func, arg, code) \
458 __inline_mathcodeNP_ (float, __CONCAT(func,f), arg, code) \
459 __inline_mathcodeNP_ (long double, __CONCAT(func,l), arg, code)
460 # define __inline_mathcode2(func, arg1, arg2, code) \
461 __inline_mathcode2_ (double, func, arg1, arg2, code) \
462 __inline_mathcode2_ (float, __CONCAT(func,f), arg1, arg2, code) \
463 __inline_mathcode2_ (long double, __CONCAT(func,l), arg1, arg2, code)
464 # define __inline_mathcodeNP2(func, arg1, arg2, code) \
465 __inline_mathcodeNP2_ (double, func, arg1, arg2, code) \
466 __inline_mathcodeNP2_ (float, __CONCAT(func,f), arg1, arg2, code) \
467 __inline_mathcodeNP2_ (long double, __CONCAT(func,l), arg1, arg2, code)
468 # define __inline_mathcode3(func, arg1, arg2, arg3, code) \
469 __inline_mathcode3_ (double, func, arg1, arg2, arg3, code) \
470 __inline_mathcode3_ (float, __CONCAT(func,f), arg1, arg2, arg3, code) \
471 __inline_mathcode3_ (long double, __CONCAT(func,l), arg1, arg2, arg3, code)
472 # define __inline_mathcodeNP3(func, arg1, arg2, arg3, code) \
473 __inline_mathcodeNP3_ (double, func, arg1, arg2, arg3, code) \
474 __inline_mathcodeNP3_ (float, __CONCAT(func,f), arg1, arg2, arg3, code) \
475 __inline_mathcodeNP3_ (long double, __CONCAT(func,l), arg1, arg2, arg3, code)
476 # else
477 # define __inline_mathcode(func, arg, code) \
478 __inline_mathcode_ (double, func, (arg), code)
479 # define __inline_mathcodeNP(func, arg, code) \
480 __inline_mathcodeNP_ (double, func, (arg), code)
481 # define __inline_mathcode2(func, arg1, arg2, code) \
482 __inline_mathcode2_ (double, func, arg1, arg2, code)
483 # define __inline_mathcodeNP2(func, arg1, arg2, code) \
484 __inline_mathcodeNP2_ (double, func, arg1, arg2, code)
485 # define __inline_mathcode3(func, arg1, arg2, arg3, code) \
486 __inline_mathcode3_ (double, func, arg1, arg2, arg3, code)
487 # define __inline_mathcodeNP3(func, arg1, arg2, arg3, code) \
488 __inline_mathcodeNP3_ (double, func, arg1, arg2, arg3, code)
489 # endif
490
491 # define __inline_mathcode_(float_type, func, arg, code) \
492 __MATH_INLINE float_type func (float_type) __THROW; \
493 __inline_mathcodeNP_(float_type, func, arg, code)
494
495 # define __inline_mathcodeNP_(float_type, func, arg, code) \
496 __MATH_INLINE float_type __NTH (func (float_type arg)) \
497 { \
498 code; \
499 }
500
501
502 # define __inline_mathcode2_(float_type, func, arg1, arg2, code) \
503 __MATH_INLINE float_type func (float_type, float_type) __THROW; \
504 __inline_mathcodeNP2_ (float_type, func, arg1, arg2, code)
505
506 # define __inline_mathcodeNP2_(float_type, func, arg1, arg2, code) \
507 __MATH_INLINE float_type __NTH (func (float_type arg1, float_type arg2)) \
508 { \
509 code; \
510 }
511
512 # define __inline_mathcode3_(float_type, func, arg1, arg2, arg3, code) \
513 __MATH_INLINE float_type func (float_type, float_type, float_type) __THROW; \
514 __inline_mathcodeNP3_(float_type, func, arg1, arg2, arg3, code)
515
516 # define __inline_mathcodeNP3_(float_type, func, arg1, arg2, arg3, code) \
517 __MATH_INLINE float_type __NTH (func (float_type arg1, float_type arg2, \
518 float_type arg3)) \
519 { \
520 code; \
521 }
522 # endif
523
524
525 # if !defined __NO_MATH_INLINES && defined __OPTIMIZE__
526 /* Miscellaneous functions */
527
528 /* __FAST_MATH__ is defined by gcc -ffast-math. */
529 # ifdef __FAST_MATH__
530 # ifdef __USE_GNU
531 # define __sincos_code \
532 register long double __cosr; \
533 register long double __sinr; \
534 register unsigned int __swtmp; \
535 __asm __volatile__ \
536 ("fsincos\n\t" \
537 "fnstsw %w2\n\t" \
538 "testl $0x400, %2\n\t" \
539 "jz 1f\n\t" \
540 "fldpi\n\t" \
541 "fadd %%st(0)\n\t" \
542 "fxch %%st(1)\n\t" \
543 "2: fprem1\n\t" \
544 "fnstsw %w2\n\t" \
545 "testl $0x400, %2\n\t" \
546 "jnz 2b\n\t" \
547 "fstp %%st(1)\n\t" \
548 "fsincos\n\t" \
549 "1:" \
550 : "=t" (__cosr), "=u" (__sinr), "=a" (__swtmp) : "0" (__x)); \
551 *__sinx = __sinr; \
552 *__cosx = __cosr
553
554 __MATH_INLINE void
__NTH(__sincos (double __x,double * __sinx,double * __cosx))555 __NTH (__sincos (double __x, double *__sinx, double *__cosx))
556 {
557 __sincos_code;
558 }
559
560 __MATH_INLINE void
__NTH(__sincosf (float __x,float * __sinx,float * __cosx))561 __NTH (__sincosf (float __x, float *__sinx, float *__cosx))
562 {
563 __sincos_code;
564 }
565
566 __MATH_INLINE void
__NTH(__sincosl (long double __x,long double * __sinx,long double * __cosx))567 __NTH (__sincosl (long double __x, long double *__sinx, long double *__cosx))
568 {
569 __sincos_code;
570 }
571 # endif
572
573
574 /* Optimized inline implementation, sometimes with reduced precision
575 and/or argument range. */
576
577 # if __GNUC_PREREQ (3, 5)
578 # define __expm1_code \
579 register long double __temp; \
580 __temp = __builtin_expm1l (__x); \
581 return __temp ? __temp : __x
582 # else
583 # define __expm1_code \
584 register long double __value; \
585 register long double __exponent; \
586 register long double __temp; \
587 __asm __volatile__ \
588 ("fldl2e # e^x - 1 = 2^(x * log2(e)) - 1\n\t" \
589 "fmul %%st(1) # x * log2(e)\n\t" \
590 "fst %%st(1)\n\t" \
591 "frndint # int(x * log2(e))\n\t" \
592 "fxch\n\t" \
593 "fsub %%st(1) # fract(x * log2(e))\n\t" \
594 "f2xm1 # 2^(fract(x * log2(e))) - 1\n\t" \
595 "fscale # 2^(x * log2(e)) - 2^(int(x * log2(e)))\n\t" \
596 : "=t" (__value), "=u" (__exponent) : "0" (__x)); \
597 __asm __volatile__ \
598 ("fscale # 2^int(x * log2(e))\n\t" \
599 : "=t" (__temp) : "0" (1.0), "u" (__exponent)); \
600 __temp -= 1.0; \
601 __temp += __value; \
602 return __temp ? __temp : __x
603 # endif
__inline_mathcodeNP_(long double,__expm1l,__x,__expm1_code)604 __inline_mathcodeNP_ (long double, __expm1l, __x, __expm1_code)
605
606 # if __GNUC_PREREQ (3, 4)
607 __inline_mathcodeNP_ (long double, __expl, __x, return __builtin_expl (__x))
608 # else
609 # define __exp_code \
610 register long double __value; \
611 register long double __exponent; \
612 __asm __volatile__ \
613 ("fldl2e # e^x = 2^(x * log2(e))\n\t" \
614 "fmul %%st(1) # x * log2(e)\n\t" \
615 "fst %%st(1)\n\t" \
616 "frndint # int(x * log2(e))\n\t" \
617 "fxch\n\t" \
618 "fsub %%st(1) # fract(x * log2(e))\n\t" \
619 "f2xm1 # 2^(fract(x * log2(e))) - 1\n\t" \
620 : "=t" (__value), "=u" (__exponent) : "0" (__x)); \
621 __value += 1.0; \
622 __asm __volatile__ \
623 ("fscale" \
624 : "=t" (__value) : "0" (__value), "u" (__exponent)); \
625 return __value
626 __inline_mathcodeNP (exp, __x, __exp_code)
627 __inline_mathcodeNP_ (long double, __expl, __x, __exp_code)
628 # endif
629
630
631 # if !__GNUC_PREREQ (3, 5)
632 __inline_mathcodeNP (tan, __x, \
633 register long double __value; \
634 register long double __value2 __attribute__ ((__unused__)); \
635 __asm __volatile__ \
636 ("fptan" \
637 : "=t" (__value2), "=u" (__value) : "0" (__x)); \
638 return __value)
639 # endif
640 # endif /* __FAST_MATH__ */
641
642
643 # if __GNUC_PREREQ (3, 4)
644 __inline_mathcodeNP2_ (long double, __atan2l, __y, __x,
645 return __builtin_atan2l (__y, __x))
646 # else
647 # define __atan2_code \
648 register long double __value; \
649 __asm __volatile__ \
650 ("fpatan" \
651 : "=t" (__value) : "0" (__x), "u" (__y) : "st(1)"); \
652 return __value
653 # ifdef __FAST_MATH__
654 __inline_mathcodeNP2 (atan2, __y, __x, __atan2_code)
655 # endif
656 __inline_mathcodeNP2_ (long double, __atan2l, __y, __x, __atan2_code)
657 # endif
658
659
660 # if defined __FAST_MATH__ && !__GNUC_PREREQ (3, 5)
661 __inline_mathcodeNP2 (fmod, __x, __y, \
662 register long double __value; \
663 __asm __volatile__ \
664 ("1: fprem\n\t" \
665 "fnstsw %%ax\n\t" \
666 "sahf\n\t" \
667 "jp 1b" \
668 : "=t" (__value) : "0" (__x), "u" (__y) : "ax", "cc"); \
669 return __value)
670 # endif
671
672
673 # ifdef __FAST_MATH__
674 # if !__GNUC_PREREQ (3,3)
675 __inline_mathopNP (sqrt, "fsqrt")
676 __inline_mathopNP_ (long double, __sqrtl, "fsqrt")
677 # define __libc_sqrtl(n) __sqrtl (n)
678 # else
679 # define __libc_sqrtl(n) __builtin_sqrtl (n)
680 # endif
681 # endif
682
683 # if __GNUC_PREREQ (2, 8)
684 __inline_mathcodeNP_ (double, fabs, __x, return __builtin_fabs (__x))
685 # if defined __USE_MISC || defined __USE_ISOC99
686 __inline_mathcodeNP_ (float, fabsf, __x, return __builtin_fabsf (__x))
687 __inline_mathcodeNP_ (long double, fabsl, __x, return __builtin_fabsl (__x))
688 # endif
689 __inline_mathcodeNP_ (long double, __fabsl, __x, return __builtin_fabsl (__x))
690 # else
691 __inline_mathop (fabs, "fabs")
692 __inline_mathop_ (long double, __fabsl, "fabs")
693 # endif
694
695 # ifdef __FAST_MATH__
696 # if !__GNUC_PREREQ (3, 4)
697 /* The argument range of this inline version is reduced. */
698 __inline_mathopNP (sin, "fsin")
699 /* The argument range of this inline version is reduced. */
700 __inline_mathopNP (cos, "fcos")
701
702 __inline_mathop_declNP (log, "fldln2; fxch; fyl2x", "0" (__x) : "st(1)")
703 # endif
704
705 # if !__GNUC_PREREQ (3, 5)
706 __inline_mathop_declNP (log10, "fldlg2; fxch; fyl2x", "0" (__x) : "st(1)")
707
708 __inline_mathcodeNP (asin, __x, return __atan2l (__x, __libc_sqrtl (1.0 - __x * __x)))
709 __inline_mathcodeNP (acos, __x, return __atan2l (__libc_sqrtl (1.0 - __x * __x), __x))
710 # endif
711
712 # if !__GNUC_PREREQ (3, 4)
713 __inline_mathop_declNP (atan, "fld1; fpatan", "0" (__x) : "st(1)")
714 # endif
715 # endif /* __FAST_MATH__ */
716
717 __inline_mathcode_ (long double, __sgn1l, __x, \
718 __extension__ union { long double __xld; unsigned int __xi[3]; } __n = \
719 { __xld: __x }; \
720 __n.__xi[2] = (__n.__xi[2] & 0x8000) | 0x3fff; \
721 __n.__xi[1] = 0x80000000; \
722 __n.__xi[0] = 0; \
723 return __n.__xld)
724
725
726 # ifdef __FAST_MATH__
727 /* The argument range of the inline version of sinhl is slightly reduced. */
728 __inline_mathcodeNP (sinh, __x, \
729 register long double __exm1 = __expm1l (__fabsl (__x)); \
730 return 0.5 * (__exm1 / (__exm1 + 1.0) + __exm1) * __sgn1l (__x))
731
732 __inline_mathcodeNP (cosh, __x, \
733 register long double __ex = __expl (__x); \
734 return 0.5 * (__ex + 1.0 / __ex))
735
736 __inline_mathcodeNP (tanh, __x, \
737 register long double __exm1 = __expm1l (-__fabsl (__x + __x)); \
738 return __exm1 / (__exm1 + 2.0) * __sgn1l (-__x))
739 # endif
740
741 __inline_mathcodeNP (floor, __x, \
742 register long double __value; \
743 register int __ignore; \
744 unsigned short int __cw; \
745 unsigned short int __cwtmp; \
746 __asm __volatile ("fnstcw %3\n\t" \
747 "movzwl %3, %1\n\t" \
748 "andl $0xf3ff, %1\n\t" \
749 "orl $0x0400, %1\n\t" /* rounding down */ \
750 "movw %w1, %2\n\t" \
751 "fldcw %2\n\t" \
752 "frndint\n\t" \
753 "fldcw %3" \
754 : "=t" (__value), "=&q" (__ignore), "=m" (__cwtmp), \
755 "=m" (__cw) \
756 : "0" (__x)); \
757 return __value)
758
759 __inline_mathcodeNP (ceil, __x, \
760 register long double __value; \
761 register int __ignore; \
762 unsigned short int __cw; \
763 unsigned short int __cwtmp; \
764 __asm __volatile ("fnstcw %3\n\t" \
765 "movzwl %3, %1\n\t" \
766 "andl $0xf3ff, %1\n\t" \
767 "orl $0x0800, %1\n\t" /* rounding up */ \
768 "movw %w1, %2\n\t" \
769 "fldcw %2\n\t" \
770 "frndint\n\t" \
771 "fldcw %3" \
772 : "=t" (__value), "=&q" (__ignore), "=m" (__cwtmp), \
773 "=m" (__cw) \
774 : "0" (__x)); \
775 return __value)
776
777 # ifdef __FAST_MATH__
778 # define __ldexp_code \
779 register long double __value; \
780 __asm __volatile__ \
781 ("fscale" \
782 : "=t" (__value) : "0" (__x), "u" ((long double) __y)); \
783 return __value
784
785 __MATH_INLINE double
786 __NTH (ldexp (double __x, int __y))
787 {
788 __ldexp_code;
789 }
790 # endif
791
792
793 /* Optimized versions for some non-standardized functions. */
794 # if defined __USE_ISOC99 || defined __USE_MISC
795
796 # ifdef __FAST_MATH__
__inline_mathcodeNP(expm1,__x,__expm1_code)797 __inline_mathcodeNP (expm1, __x, __expm1_code)
798
799 /* We cannot rely on M_SQRT being defined. So we do it for ourself
800 here. */
801 # define __M_SQRT2 1.41421356237309504880L /* sqrt(2) */
802
803 # if !__GNUC_PREREQ (3, 5)
804 __inline_mathcodeNP (log1p, __x, \
805 register long double __value; \
806 if (__fabsl (__x) >= 1.0 - 0.5 * __M_SQRT2) \
807 __value = logl (1.0 + __x); \
808 else \
809 __asm __volatile__ \
810 ("fldln2\n\t" \
811 "fxch\n\t" \
812 "fyl2xp1" \
813 : "=t" (__value) : "0" (__x) : "st(1)"); \
814 return __value)
815 # endif
816
817
818 /* The argument range of the inline version of asinhl is slightly reduced. */
819 __inline_mathcodeNP (asinh, __x, \
820 register long double __y = __fabsl (__x); \
821 return (log1pl (__y * __y / (__libc_sqrtl (__y * __y + 1.0) + 1.0) + __y) \
822 * __sgn1l (__x)))
823
824 __inline_mathcodeNP (acosh, __x, \
825 return logl (__x + __libc_sqrtl (__x - 1.0) * __libc_sqrtl (__x + 1.0)))
826
827 __inline_mathcodeNP (atanh, __x, \
828 register long double __y = __fabsl (__x); \
829 return -0.5 * log1pl (-(__y + __y) / (1.0 + __y)) * __sgn1l (__x))
830
831 /* The argument range of the inline version of hypotl is slightly reduced. */
832 __inline_mathcodeNP2 (hypot, __x, __y,
833 return __libc_sqrtl (__x * __x + __y * __y))
834
835 # if !__GNUC_PREREQ (3, 5)
836 __inline_mathcodeNP(logb, __x, \
837 register long double __value; \
838 register long double __junk; \
839 __asm __volatile__ \
840 ("fxtract\n\t" \
841 : "=t" (__junk), "=u" (__value) : "0" (__x)); \
842 return __value)
843 # endif
844
845 # endif
846 # endif
847
848 # ifdef __USE_ISOC99
849 # ifdef __FAST_MATH__
850
851 # if !__GNUC_PREREQ (3, 5)
852 __inline_mathop_declNP (log2, "fld1; fxch; fyl2x", "0" (__x) : "st(1)")
853 # endif
854
855 __MATH_INLINE float
856 __NTH (ldexpf (float __x, int __y))
857 {
858 __ldexp_code;
859 }
860
861 __MATH_INLINE long double
__NTH(ldexpl (long double __x,int __y))862 __NTH (ldexpl (long double __x, int __y))
863 {
864 __ldexp_code;
865 }
866
867 __inline_mathopNP (rint, "frndint")
868 # endif /* __FAST_MATH__ */
869
870 # define __lrint_code \
871 long int __lrintres; \
872 __asm__ __volatile__ \
873 ("fistpl %0" \
874 : "=m" (__lrintres) : "t" (__x) : "st"); \
875 return __lrintres
876 __MATH_INLINE long int
__NTH(lrintf (float __x))877 __NTH (lrintf (float __x))
878 {
879 __lrint_code;
880 }
881 __MATH_INLINE long int
__NTH(lrint (double __x))882 __NTH (lrint (double __x))
883 {
884 __lrint_code;
885 }
886 __MATH_INLINE long int
__NTH(lrintl (long double __x))887 __NTH (lrintl (long double __x))
888 {
889 __lrint_code;
890 }
891 # undef __lrint_code
892
893 # define __llrint_code \
894 long long int __llrintres; \
895 __asm__ __volatile__ \
896 ("fistpll %0" \
897 : "=m" (__llrintres) : "t" (__x) : "st"); \
898 return __llrintres
899 __MATH_INLINE long long int
__NTH(llrintf (float __x))900 __NTH (llrintf (float __x))
901 {
902 __llrint_code;
903 }
904 __MATH_INLINE long long int
__NTH(llrint (double __x))905 __NTH (llrint (double __x))
906 {
907 __llrint_code;
908 }
909 __MATH_INLINE long long int
__NTH(llrintl (long double __x))910 __NTH (llrintl (long double __x))
911 {
912 __llrint_code;
913 }
914 # undef __llrint_code
915
916 # endif
917
918
919 # ifdef __USE_MISC
920
921 # if defined __FAST_MATH__ && !__GNUC_PREREQ (3, 5)
922 __inline_mathcodeNP2 (drem, __x, __y, \
923 register double __value; \
924 register int __clobbered; \
925 __asm __volatile__ \
926 ("1: fprem1\n\t" \
927 "fstsw %%ax\n\t" \
928 "sahf\n\t" \
929 "jp 1b" \
930 : "=t" (__value), "=&a" (__clobbered) : "0" (__x), "u" (__y) : "cc"); \
931 return __value)
932 # endif
933
934
935 /* This function is used in the `isfinite' macro. */
936 __MATH_INLINE int
__NTH(__finite (double __x))937 __NTH (__finite (double __x))
938 {
939 return (__extension__
940 (((((union { double __d; int __i[2]; }) {__d: __x}).__i[1]
941 | 0x800fffffu) + 1) >> 31));
942 }
943
944 # endif /* __USE_MISC */
945
946 /* Undefine some of the large macros which are not used anymore. */
947 # undef __atan2_code
948 # ifdef __FAST_MATH__
949 # undef __expm1_code
950 # undef __exp_code
951 # undef __sincos_code
952 # endif /* __FAST_MATH__ */
953
954 # endif /* __NO_MATH_INLINES */
955
956
957 /* This code is used internally in the GNU libc. */
958 # ifdef __LIBC_INTERNAL_MATH_INLINES
959 __inline_mathop (__ieee754_sqrt, "fsqrt")
960 __inline_mathcode2 (__ieee754_atan2, __y, __x,
961 register long double __value;
962 __asm __volatile__ ("fpatan\n\t"
963 : "=t" (__value)
964 : "0" (__x), "u" (__y) : "st(1)");
965 return __value;)
966 # endif
967
968 #endif /* !__x86_64__ */
969