]>
Commit | Line | Data |
---|---|---|
777b1eea | 1 | /* logbf(). PowerPC/POWER7 version. |
d4697bc9 | 2 | Copyright (C) 2012-2014 Free Software Foundation, Inc. |
777b1eea AZ |
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 | #include "math_private.h" | |
20 | ||
21 | /* This implementation avoids FP to INT conversions by using VSX | |
22 | bitwise instructions over FP values. */ | |
23 | ||
24 | static const double two1div52 = 2.220446049250313e-16; /* 1/2**52 */ | |
25 | static const double two10m1 = -1023.0; /* -2**10 + 1 */ | |
26 | static const double two7m1 = -127.0; /* -2**7 + 1 */ | |
27 | ||
28 | /* FP mask to extract the exponent. */ | |
29 | static const union { | |
30 | unsigned long long mask; | |
31 | double d; | |
32 | } mask = { 0x7ff0000000000000ULL }; | |
33 | ||
34 | float | |
35 | __logbf (float x) | |
36 | { | |
37 | /* VSX operation are all done internally as double. */ | |
38 | double ret; | |
39 | ||
40 | if (__builtin_expect (x == 0.0, 0)) | |
41 | /* Raise FE_DIVBYZERO and return -HUGE_VAL[LF]. */ | |
42 | return -1.0 / __builtin_fabsf (x); | |
43 | ||
44 | /* ret = x & 0x7f800000; */ | |
45 | asm ( | |
46 | "xxland %x0,%x1,%x2\n" | |
47 | "fcfid %0,%0" | |
48 | : "=f"(ret) | |
49 | : "f" (x), "f" (mask.d)); | |
50 | /* ret = (ret >> 52) - 1023.0, since ret is double. */ | |
51 | ret = (ret * two1div52) + two10m1; | |
52 | if (__builtin_expect (ret > -two7m1, 0)) | |
53 | /* Multiplication is used to set logb (+-INF) = INF. */ | |
54 | return (x * x); | |
55 | /* Since operations are done with double we don't need | |
56 | additional tests for subnormal numbers. | |
57 | The test is to avoid logb_downward (0.0) == -0.0. */ | |
58 | return ret == -0.0 ? 0.0 : ret; | |
59 | } | |
60 | weak_alias (__logbf, logbf) |