123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505 |
- /* Software floating-point emulation.
- Definitions for IEEE Extended Precision.
- Copyright (C) 1999-2014 Free Software Foundation, Inc.
- This file is part of the GNU C Library.
- Contributed by Jakub Jelinek (jj@ultra.linux.cz).
- The GNU C Library is free software; you can redistribute it and/or
- modify it under the terms of the GNU Lesser General Public
- License as published by the Free Software Foundation; either
- version 2.1 of the License, or (at your option) any later version.
- In addition to the permissions in the GNU Lesser General Public
- License, the Free Software Foundation gives you unlimited
- permission to link the compiled version of this file into
- combinations with other programs, and to distribute those
- combinations without any restriction coming from the use of this
- file. (The Lesser General Public License restrictions do apply in
- other respects; for example, they cover modification of the file,
- and distribution when not linked into a combine executable.)
- The GNU C Library is distributed in the hope that it will be useful,
- but WITHOUT ANY WARRANTY; without even the implied warranty of
- MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
- Lesser General Public License for more details.
- You should have received a copy of the GNU Lesser General Public
- License along with the GNU C Library; if not, see
- <http://www.gnu.org/licenses/>. */
- #if _FP_W_TYPE_SIZE < 32
- # error "Here's a nickel, kid. Go buy yourself a real computer."
- #endif
- #if _FP_W_TYPE_SIZE < 64
- # define _FP_FRACTBITS_E (4*_FP_W_TYPE_SIZE)
- # define _FP_FRACTBITS_DW_E (8*_FP_W_TYPE_SIZE)
- #else
- # define _FP_FRACTBITS_E (2*_FP_W_TYPE_SIZE)
- # define _FP_FRACTBITS_DW_E (4*_FP_W_TYPE_SIZE)
- #endif
- #define _FP_FRACBITS_E 64
- #define _FP_FRACXBITS_E (_FP_FRACTBITS_E - _FP_FRACBITS_E)
- #define _FP_WFRACBITS_E (_FP_WORKBITS + _FP_FRACBITS_E)
- #define _FP_WFRACXBITS_E (_FP_FRACTBITS_E - _FP_WFRACBITS_E)
- #define _FP_EXPBITS_E 15
- #define _FP_EXPBIAS_E 16383
- #define _FP_EXPMAX_E 32767
- #define _FP_QNANBIT_E \
- ((_FP_W_TYPE) 1 << (_FP_FRACBITS_E-2) % _FP_W_TYPE_SIZE)
- #define _FP_QNANBIT_SH_E \
- ((_FP_W_TYPE) 1 << (_FP_FRACBITS_E-2+_FP_WORKBITS) % _FP_W_TYPE_SIZE)
- #define _FP_IMPLBIT_E \
- ((_FP_W_TYPE) 1 << (_FP_FRACBITS_E-1) % _FP_W_TYPE_SIZE)
- #define _FP_IMPLBIT_SH_E \
- ((_FP_W_TYPE) 1 << (_FP_FRACBITS_E-1+_FP_WORKBITS) % _FP_W_TYPE_SIZE)
- #define _FP_OVERFLOW_E \
- ((_FP_W_TYPE) 1 << (_FP_WFRACBITS_E % _FP_W_TYPE_SIZE))
- #define _FP_WFRACBITS_DW_E (2 * _FP_WFRACBITS_E)
- #define _FP_WFRACXBITS_DW_E (_FP_FRACTBITS_DW_E - _FP_WFRACBITS_DW_E)
- #define _FP_HIGHBIT_DW_E \
- ((_FP_W_TYPE) 1 << (_FP_WFRACBITS_DW_E - 1) % _FP_W_TYPE_SIZE)
- typedef float XFtype __attribute__ ((mode (XF)));
- #if _FP_W_TYPE_SIZE < 64
- union _FP_UNION_E
- {
- XFtype flt;
- struct _FP_STRUCT_LAYOUT
- {
- # if __BYTE_ORDER == __BIG_ENDIAN
- unsigned long pad1 : _FP_W_TYPE_SIZE;
- unsigned long pad2 : (_FP_W_TYPE_SIZE - 1 - _FP_EXPBITS_E);
- unsigned long sign : 1;
- unsigned long exp : _FP_EXPBITS_E;
- unsigned long frac1 : _FP_W_TYPE_SIZE;
- unsigned long frac0 : _FP_W_TYPE_SIZE;
- # else
- unsigned long frac0 : _FP_W_TYPE_SIZE;
- unsigned long frac1 : _FP_W_TYPE_SIZE;
- unsigned exp : _FP_EXPBITS_E;
- unsigned sign : 1;
- # endif /* not bigendian */
- } bits __attribute__ ((packed));
- };
- # define FP_DECL_E(X) _FP_DECL (4, X)
- # define FP_UNPACK_RAW_E(X, val) \
- do \
- { \
- union _FP_UNION_E FP_UNPACK_RAW_E_flo; \
- FP_UNPACK_RAW_E_flo.flt = (val); \
- \
- X##_f[2] = 0; \
- X##_f[3] = 0; \
- X##_f[0] = FP_UNPACK_RAW_E_flo.bits.frac0; \
- X##_f[1] = FP_UNPACK_RAW_E_flo.bits.frac1; \
- X##_e = FP_UNPACK_RAW_E_flo.bits.exp; \
- X##_s = FP_UNPACK_RAW_E_flo.bits.sign; \
- } \
- while (0)
- # define FP_UNPACK_RAW_EP(X, val) \
- do \
- { \
- union _FP_UNION_E *FP_UNPACK_RAW_EP_flo \
- = (union _FP_UNION_E *) (val); \
- \
- X##_f[2] = 0; \
- X##_f[3] = 0; \
- X##_f[0] = FP_UNPACK_RAW_EP_flo->bits.frac0; \
- X##_f[1] = FP_UNPACK_RAW_EP_flo->bits.frac1; \
- X##_e = FP_UNPACK_RAW_EP_flo->bits.exp; \
- X##_s = FP_UNPACK_RAW_EP_flo->bits.sign; \
- } \
- while (0)
- # define FP_PACK_RAW_E(val, X) \
- do \
- { \
- union _FP_UNION_E FP_PACK_RAW_E_flo; \
- \
- if (X##_e) \
- X##_f[1] |= _FP_IMPLBIT_E; \
- else \
- X##_f[1] &= ~(_FP_IMPLBIT_E); \
- FP_PACK_RAW_E_flo.bits.frac0 = X##_f[0]; \
- FP_PACK_RAW_E_flo.bits.frac1 = X##_f[1]; \
- FP_PACK_RAW_E_flo.bits.exp = X##_e; \
- FP_PACK_RAW_E_flo.bits.sign = X##_s; \
- \
- (val) = FP_PACK_RAW_E_flo.flt; \
- } \
- while (0)
- # define FP_PACK_RAW_EP(val, X) \
- do \
- { \
- if (!FP_INHIBIT_RESULTS) \
- { \
- union _FP_UNION_E *FP_PACK_RAW_EP_flo \
- = (union _FP_UNION_E *) (val); \
- \
- if (X##_e) \
- X##_f[1] |= _FP_IMPLBIT_E; \
- else \
- X##_f[1] &= ~(_FP_IMPLBIT_E); \
- FP_PACK_RAW_EP_flo->bits.frac0 = X##_f[0]; \
- FP_PACK_RAW_EP_flo->bits.frac1 = X##_f[1]; \
- FP_PACK_RAW_EP_flo->bits.exp = X##_e; \
- FP_PACK_RAW_EP_flo->bits.sign = X##_s; \
- } \
- } \
- while (0)
- # define FP_UNPACK_E(X, val) \
- do \
- { \
- FP_UNPACK_RAW_E (X, (val)); \
- _FP_UNPACK_CANONICAL (E, 4, X); \
- } \
- while (0)
- # define FP_UNPACK_EP(X, val) \
- do \
- { \
- FP_UNPACK_RAW_EP (X, (val)); \
- _FP_UNPACK_CANONICAL (E, 4, X); \
- } \
- while (0)
- # define FP_UNPACK_SEMIRAW_E(X, val) \
- do \
- { \
- FP_UNPACK_RAW_E (X, (val)); \
- _FP_UNPACK_SEMIRAW (E, 4, X); \
- } \
- while (0)
- # define FP_UNPACK_SEMIRAW_EP(X, val) \
- do \
- { \
- FP_UNPACK_RAW_EP (X, (val)); \
- _FP_UNPACK_SEMIRAW (E, 4, X); \
- } \
- while (0)
- # define FP_PACK_E(val, X) \
- do \
- { \
- _FP_PACK_CANONICAL (E, 4, X); \
- FP_PACK_RAW_E ((val), X); \
- } \
- while (0)
- # define FP_PACK_EP(val, X) \
- do \
- { \
- _FP_PACK_CANONICAL (E, 4, X); \
- FP_PACK_RAW_EP ((val), X); \
- } \
- while (0)
- # define FP_PACK_SEMIRAW_E(val, X) \
- do \
- { \
- _FP_PACK_SEMIRAW (E, 4, X); \
- FP_PACK_RAW_E ((val), X); \
- } \
- while (0)
- # define FP_PACK_SEMIRAW_EP(val, X) \
- do \
- { \
- _FP_PACK_SEMIRAW (E, 4, X); \
- FP_PACK_RAW_EP ((val), X); \
- } \
- while (0)
- # define FP_ISSIGNAN_E(X) _FP_ISSIGNAN (E, 4, X)
- # define FP_NEG_E(R, X) _FP_NEG (E, 4, R, X)
- # define FP_ADD_E(R, X, Y) _FP_ADD (E, 4, R, X, Y)
- # define FP_SUB_E(R, X, Y) _FP_SUB (E, 4, R, X, Y)
- # define FP_MUL_E(R, X, Y) _FP_MUL (E, 4, R, X, Y)
- # define FP_DIV_E(R, X, Y) _FP_DIV (E, 4, R, X, Y)
- # define FP_SQRT_E(R, X) _FP_SQRT (E, 4, R, X)
- # define FP_FMA_E(R, X, Y, Z) _FP_FMA (E, 4, 8, R, X, Y, Z)
- /* Square root algorithms:
- We have just one right now, maybe Newton approximation
- should be added for those machines where division is fast.
- This has special _E version because standard _4 square
- root would not work (it has to start normally with the
- second word and not the first), but as we have to do it
- anyway, we optimize it by doing most of the calculations
- in two UWtype registers instead of four. */
- # define _FP_SQRT_MEAT_E(R, S, T, X, q) \
- do \
- { \
- (q) = (_FP_W_TYPE) 1 << (_FP_W_TYPE_SIZE - 1); \
- _FP_FRAC_SRL_4 (X, (_FP_WORKBITS)); \
- while (q) \
- { \
- T##_f[1] = S##_f[1] + (q); \
- if (T##_f[1] <= X##_f[1]) \
- { \
- S##_f[1] = T##_f[1] + (q); \
- X##_f[1] -= T##_f[1]; \
- R##_f[1] += (q); \
- } \
- _FP_FRAC_SLL_2 (X, 1); \
- (q) >>= 1; \
- } \
- (q) = (_FP_W_TYPE) 1 << (_FP_W_TYPE_SIZE - 1); \
- while (q) \
- { \
- T##_f[0] = S##_f[0] + (q); \
- T##_f[1] = S##_f[1]; \
- if (T##_f[1] < X##_f[1] \
- || (T##_f[1] == X##_f[1] \
- && T##_f[0] <= X##_f[0])) \
- { \
- S##_f[0] = T##_f[0] + (q); \
- S##_f[1] += (T##_f[0] > S##_f[0]); \
- _FP_FRAC_DEC_2 (X, T); \
- R##_f[0] += (q); \
- } \
- _FP_FRAC_SLL_2 (X, 1); \
- (q) >>= 1; \
- } \
- _FP_FRAC_SLL_4 (R, (_FP_WORKBITS)); \
- if (X##_f[0] | X##_f[1]) \
- { \
- if (S##_f[1] < X##_f[1] \
- || (S##_f[1] == X##_f[1] \
- && S##_f[0] < X##_f[0])) \
- R##_f[0] |= _FP_WORK_ROUND; \
- R##_f[0] |= _FP_WORK_STICKY; \
- } \
- } \
- while (0)
- # define FP_CMP_E(r, X, Y, un, ex) _FP_CMP (E, 4, (r), X, Y, (un), (ex))
- # define FP_CMP_EQ_E(r, X, Y, ex) _FP_CMP_EQ (E, 4, (r), X, Y, (ex))
- # define FP_CMP_UNORD_E(r, X, Y, ex) _FP_CMP_UNORD (E, 4, (r), X, Y, (ex))
- # define FP_TO_INT_E(r, X, rsz, rsg) _FP_TO_INT (E, 4, (r), X, (rsz), (rsg))
- # define FP_FROM_INT_E(X, r, rs, rt) _FP_FROM_INT (E, 4, X, (r), (rs), rt)
- # define _FP_FRAC_HIGH_E(X) (X##_f[2])
- # define _FP_FRAC_HIGH_RAW_E(X) (X##_f[1])
- # define _FP_FRAC_HIGH_DW_E(X) (X##_f[4])
- #else /* not _FP_W_TYPE_SIZE < 64 */
- union _FP_UNION_E
- {
- XFtype flt;
- struct _FP_STRUCT_LAYOUT
- {
- # if __BYTE_ORDER == __BIG_ENDIAN
- _FP_W_TYPE pad : (_FP_W_TYPE_SIZE - 1 - _FP_EXPBITS_E);
- unsigned sign : 1;
- unsigned exp : _FP_EXPBITS_E;
- _FP_W_TYPE frac : _FP_W_TYPE_SIZE;
- # else
- _FP_W_TYPE frac : _FP_W_TYPE_SIZE;
- unsigned exp : _FP_EXPBITS_E;
- unsigned sign : 1;
- # endif
- } bits;
- };
- # define FP_DECL_E(X) _FP_DECL (2, X)
- # define FP_UNPACK_RAW_E(X, val) \
- do \
- { \
- union _FP_UNION_E FP_UNPACK_RAW_E_flo; \
- FP_UNPACK_RAW_E_flo.flt = (val); \
- \
- X##_f0 = FP_UNPACK_RAW_E_flo.bits.frac; \
- X##_f1 = 0; \
- X##_e = FP_UNPACK_RAW_E_flo.bits.exp; \
- X##_s = FP_UNPACK_RAW_E_flo.bits.sign; \
- } \
- while (0)
- # define FP_UNPACK_RAW_EP(X, val) \
- do \
- { \
- union _FP_UNION_E *FP_UNPACK_RAW_EP_flo \
- = (union _FP_UNION_E *) (val); \
- \
- X##_f0 = FP_UNPACK_RAW_EP_flo->bits.frac; \
- X##_f1 = 0; \
- X##_e = FP_UNPACK_RAW_EP_flo->bits.exp; \
- X##_s = FP_UNPACK_RAW_EP_flo->bits.sign; \
- } \
- while (0)
- # define FP_PACK_RAW_E(val, X) \
- do \
- { \
- union _FP_UNION_E FP_PACK_RAW_E_flo; \
- \
- if (X##_e) \
- X##_f0 |= _FP_IMPLBIT_E; \
- else \
- X##_f0 &= ~(_FP_IMPLBIT_E); \
- FP_PACK_RAW_E_flo.bits.frac = X##_f0; \
- FP_PACK_RAW_E_flo.bits.exp = X##_e; \
- FP_PACK_RAW_E_flo.bits.sign = X##_s; \
- \
- (val) = FP_PACK_RAW_E_flo.flt; \
- } \
- while (0)
- # define FP_PACK_RAW_EP(fs, val, X) \
- do \
- { \
- if (!FP_INHIBIT_RESULTS) \
- { \
- union _FP_UNION_E *FP_PACK_RAW_EP_flo \
- = (union _FP_UNION_E *) (val); \
- \
- if (X##_e) \
- X##_f0 |= _FP_IMPLBIT_E; \
- else \
- X##_f0 &= ~(_FP_IMPLBIT_E); \
- FP_PACK_RAW_EP_flo->bits.frac = X##_f0; \
- FP_PACK_RAW_EP_flo->bits.exp = X##_e; \
- FP_PACK_RAW_EP_flo->bits.sign = X##_s; \
- } \
- } \
- while (0)
- # define FP_UNPACK_E(X, val) \
- do \
- { \
- FP_UNPACK_RAW_E (X, (val)); \
- _FP_UNPACK_CANONICAL (E, 2, X); \
- } \
- while (0)
- # define FP_UNPACK_EP(X, val) \
- do \
- { \
- FP_UNPACK_RAW_EP (X, (val)); \
- _FP_UNPACK_CANONICAL (E, 2, X); \
- } \
- while (0)
- # define FP_UNPACK_SEMIRAW_E(X, val) \
- do \
- { \
- FP_UNPACK_RAW_E (X, (val)); \
- _FP_UNPACK_SEMIRAW (E, 2, X); \
- } \
- while (0)
- # define FP_UNPACK_SEMIRAW_EP(X, val) \
- do \
- { \
- FP_UNPACK_RAW_EP (X, (val)); \
- _FP_UNPACK_SEMIRAW (E, 2, X); \
- } \
- while (0)
- # define FP_PACK_E(val, X) \
- do \
- { \
- _FP_PACK_CANONICAL (E, 2, X); \
- FP_PACK_RAW_E ((val), X); \
- } \
- while (0)
- # define FP_PACK_EP(val, X) \
- do \
- { \
- _FP_PACK_CANONICAL (E, 2, X); \
- FP_PACK_RAW_EP ((val), X); \
- } \
- while (0)
- # define FP_PACK_SEMIRAW_E(val, X) \
- do \
- { \
- _FP_PACK_SEMIRAW (E, 2, X); \
- FP_PACK_RAW_E ((val), X); \
- } \
- while (0)
- # define FP_PACK_SEMIRAW_EP(val, X) \
- do \
- { \
- _FP_PACK_SEMIRAW (E, 2, X); \
- FP_PACK_RAW_EP ((val), X); \
- } \
- while (0)
- # define FP_ISSIGNAN_E(X) _FP_ISSIGNAN (E, 2, X)
- # define FP_NEG_E(R, X) _FP_NEG (E, 2, R, X)
- # define FP_ADD_E(R, X, Y) _FP_ADD (E, 2, R, X, Y)
- # define FP_SUB_E(R, X, Y) _FP_SUB (E, 2, R, X, Y)
- # define FP_MUL_E(R, X, Y) _FP_MUL (E, 2, R, X, Y)
- # define FP_DIV_E(R, X, Y) _FP_DIV (E, 2, R, X, Y)
- # define FP_SQRT_E(R, X) _FP_SQRT (E, 2, R, X)
- # define FP_FMA_E(R, X, Y, Z) _FP_FMA (E, 2, 4, R, X, Y, Z)
- /* Square root algorithms:
- We have just one right now, maybe Newton approximation
- should be added for those machines where division is fast.
- We optimize it by doing most of the calculations
- in one UWtype registers instead of two, although we don't
- have to. */
- # define _FP_SQRT_MEAT_E(R, S, T, X, q) \
- do \
- { \
- (q) = (_FP_W_TYPE) 1 << (_FP_W_TYPE_SIZE - 1); \
- _FP_FRAC_SRL_2 (X, (_FP_WORKBITS)); \
- while (q) \
- { \
- T##_f0 = S##_f0 + (q); \
- if (T##_f0 <= X##_f0) \
- { \
- S##_f0 = T##_f0 + (q); \
- X##_f0 -= T##_f0; \
- R##_f0 += (q); \
- } \
- _FP_FRAC_SLL_1 (X, 1); \
- (q) >>= 1; \
- } \
- _FP_FRAC_SLL_2 (R, (_FP_WORKBITS)); \
- if (X##_f0) \
- { \
- if (S##_f0 < X##_f0) \
- R##_f0 |= _FP_WORK_ROUND; \
- R##_f0 |= _FP_WORK_STICKY; \
- } \
- } \
- while (0)
- # define FP_CMP_E(r, X, Y, un, ex) _FP_CMP (E, 2, (r), X, Y, (un), (ex))
- # define FP_CMP_EQ_E(r, X, Y, ex) _FP_CMP_EQ (E, 2, (r), X, Y, (ex))
- # define FP_CMP_UNORD_E(r, X, Y, ex) _FP_CMP_UNORD (E, 2, (r), X, Y, (ex))
- # define FP_TO_INT_E(r, X, rsz, rsg) _FP_TO_INT (E, 2, (r), X, (rsz), (rsg))
- # define FP_FROM_INT_E(X, r, rs, rt) _FP_FROM_INT (E, 2, X, (r), (rs), rt)
- # define _FP_FRAC_HIGH_E(X) (X##_f1)
- # define _FP_FRAC_HIGH_RAW_E(X) (X##_f0)
- # define _FP_FRAC_HIGH_DW_E(X) (X##_f[2])
- #endif /* not _FP_W_TYPE_SIZE < 64 */
|