op-1.h 9.2 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303
  1. /* Software floating-point emulation.
  2. Basic one-word fraction declaration and manipulation.
  3. Copyright (C) 1997,1998,1999 Free Software Foundation, Inc.
  4. This file is part of the GNU C Library.
  5. Contributed by Richard Henderson (rth@cygnus.com),
  6. Jakub Jelinek (jj@ultra.linux.cz),
  7. David S. Miller (davem@redhat.com) and
  8. Peter Maydell (pmaydell@chiark.greenend.org.uk).
  9. The GNU C Library is free software; you can redistribute it and/or
  10. modify it under the terms of the GNU Library General Public License as
  11. published by the Free Software Foundation; either version 2 of the
  12. License, or (at your option) any later version.
  13. The GNU C Library is distributed in the hope that it will be useful,
  14. but WITHOUT ANY WARRANTY; without even the implied warranty of
  15. MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
  16. Library General Public License for more details.
  17. You should have received a copy of the GNU Library General Public
  18. License along with the GNU C Library; see the file COPYING.LIB. If
  19. not, write to the Free Software Foundation, Inc.,
  20. 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. */
  21. #ifndef __MATH_EMU_OP_1_H__
  22. #define __MATH_EMU_OP_1_H__
  23. #define _FP_FRAC_DECL_1(X) _FP_W_TYPE X##_f=0
  24. #define _FP_FRAC_COPY_1(D,S) (D##_f = S##_f)
  25. #define _FP_FRAC_SET_1(X,I) (X##_f = I)
  26. #define _FP_FRAC_HIGH_1(X) (X##_f)
  27. #define _FP_FRAC_LOW_1(X) (X##_f)
  28. #define _FP_FRAC_WORD_1(X,w) (X##_f)
  29. #define _FP_FRAC_ADDI_1(X,I) (X##_f += I)
  30. #define _FP_FRAC_SLL_1(X,N) \
  31. do { \
  32. if (__builtin_constant_p(N) && (N) == 1) \
  33. X##_f += X##_f; \
  34. else \
  35. X##_f <<= (N); \
  36. } while (0)
  37. #define _FP_FRAC_SRL_1(X,N) (X##_f >>= N)
  38. /* Right shift with sticky-lsb. */
  39. #define _FP_FRAC_SRS_1(X,N,sz) __FP_FRAC_SRS_1(X##_f, N, sz)
  40. #define __FP_FRAC_SRS_1(X,N,sz) \
  41. (X = (X >> (N) | (__builtin_constant_p(N) && (N) == 1 \
  42. ? X & 1 : (X << (_FP_W_TYPE_SIZE - (N))) != 0)))
  43. #define _FP_FRAC_ADD_1(R,X,Y) (R##_f = X##_f + Y##_f)
  44. #define _FP_FRAC_SUB_1(R,X,Y) (R##_f = X##_f - Y##_f)
  45. #define _FP_FRAC_DEC_1(X,Y) (X##_f -= Y##_f)
  46. #define _FP_FRAC_CLZ_1(z, X) __FP_CLZ(z, X##_f)
  47. /* Predicates */
  48. #define _FP_FRAC_NEGP_1(X) ((_FP_WS_TYPE)X##_f < 0)
  49. #define _FP_FRAC_ZEROP_1(X) (X##_f == 0)
  50. #define _FP_FRAC_OVERP_1(fs,X) (X##_f & _FP_OVERFLOW_##fs)
  51. #define _FP_FRAC_CLEAR_OVERP_1(fs,X) (X##_f &= ~_FP_OVERFLOW_##fs)
  52. #define _FP_FRAC_EQ_1(X, Y) (X##_f == Y##_f)
  53. #define _FP_FRAC_GE_1(X, Y) (X##_f >= Y##_f)
  54. #define _FP_FRAC_GT_1(X, Y) (X##_f > Y##_f)
  55. #define _FP_ZEROFRAC_1 0
  56. #define _FP_MINFRAC_1 1
  57. #define _FP_MAXFRAC_1 (~(_FP_WS_TYPE)0)
  58. /*
  59. * Unpack the raw bits of a native fp value. Do not classify or
  60. * normalize the data.
  61. */
  62. #define _FP_UNPACK_RAW_1(fs, X, val) \
  63. do { \
  64. union _FP_UNION_##fs _flo; _flo.flt = (val); \
  65. \
  66. X##_f = _flo.bits.frac; \
  67. X##_e = _flo.bits.exp; \
  68. X##_s = _flo.bits.sign; \
  69. } while (0)
  70. #define _FP_UNPACK_RAW_1_P(fs, X, val) \
  71. do { \
  72. union _FP_UNION_##fs *_flo = \
  73. (union _FP_UNION_##fs *)(val); \
  74. \
  75. X##_f = _flo->bits.frac; \
  76. X##_e = _flo->bits.exp; \
  77. X##_s = _flo->bits.sign; \
  78. } while (0)
  79. /*
  80. * Repack the raw bits of a native fp value.
  81. */
  82. #define _FP_PACK_RAW_1(fs, val, X) \
  83. do { \
  84. union _FP_UNION_##fs _flo; \
  85. \
  86. _flo.bits.frac = X##_f; \
  87. _flo.bits.exp = X##_e; \
  88. _flo.bits.sign = X##_s; \
  89. \
  90. (val) = _flo.flt; \
  91. } while (0)
  92. #define _FP_PACK_RAW_1_P(fs, val, X) \
  93. do { \
  94. union _FP_UNION_##fs *_flo = \
  95. (union _FP_UNION_##fs *)(val); \
  96. \
  97. _flo->bits.frac = X##_f; \
  98. _flo->bits.exp = X##_e; \
  99. _flo->bits.sign = X##_s; \
  100. } while (0)
  101. /*
  102. * Multiplication algorithms:
  103. */
  104. /* Basic. Assuming the host word size is >= 2*FRACBITS, we can do the
  105. multiplication immediately. */
  106. #define _FP_MUL_MEAT_1_imm(wfracbits, R, X, Y) \
  107. do { \
  108. R##_f = X##_f * Y##_f; \
  109. /* Normalize since we know where the msb of the multiplicands \
  110. were (bit B), we know that the msb of the of the product is \
  111. at either 2B or 2B-1. */ \
  112. _FP_FRAC_SRS_1(R, wfracbits-1, 2*wfracbits); \
  113. } while (0)
  114. /* Given a 1W * 1W => 2W primitive, do the extended multiplication. */
  115. #define _FP_MUL_MEAT_1_wide(wfracbits, R, X, Y, doit) \
  116. do { \
  117. _FP_W_TYPE _Z_f0, _Z_f1; \
  118. doit(_Z_f1, _Z_f0, X##_f, Y##_f); \
  119. /* Normalize since we know where the msb of the multiplicands \
  120. were (bit B), we know that the msb of the of the product is \
  121. at either 2B or 2B-1. */ \
  122. _FP_FRAC_SRS_2(_Z, wfracbits-1, 2*wfracbits); \
  123. R##_f = _Z_f0; \
  124. } while (0)
  125. /* Finally, a simple widening multiply algorithm. What fun! */
  126. #define _FP_MUL_MEAT_1_hard(wfracbits, R, X, Y) \
  127. do { \
  128. _FP_W_TYPE _xh, _xl, _yh, _yl, _z_f0, _z_f1, _a_f0, _a_f1; \
  129. \
  130. /* split the words in half */ \
  131. _xh = X##_f >> (_FP_W_TYPE_SIZE/2); \
  132. _xl = X##_f & (((_FP_W_TYPE)1 << (_FP_W_TYPE_SIZE/2)) - 1); \
  133. _yh = Y##_f >> (_FP_W_TYPE_SIZE/2); \
  134. _yl = Y##_f & (((_FP_W_TYPE)1 << (_FP_W_TYPE_SIZE/2)) - 1); \
  135. \
  136. /* multiply the pieces */ \
  137. _z_f0 = _xl * _yl; \
  138. _a_f0 = _xh * _yl; \
  139. _a_f1 = _xl * _yh; \
  140. _z_f1 = _xh * _yh; \
  141. \
  142. /* reassemble into two full words */ \
  143. if ((_a_f0 += _a_f1) < _a_f1) \
  144. _z_f1 += (_FP_W_TYPE)1 << (_FP_W_TYPE_SIZE/2); \
  145. _a_f1 = _a_f0 >> (_FP_W_TYPE_SIZE/2); \
  146. _a_f0 = _a_f0 << (_FP_W_TYPE_SIZE/2); \
  147. _FP_FRAC_ADD_2(_z, _z, _a); \
  148. \
  149. /* normalize */ \
  150. _FP_FRAC_SRS_2(_z, wfracbits - 1, 2*wfracbits); \
  151. R##_f = _z_f0; \
  152. } while (0)
  153. /*
  154. * Division algorithms:
  155. */
  156. /* Basic. Assuming the host word size is >= 2*FRACBITS, we can do the
  157. division immediately. Give this macro either _FP_DIV_HELP_imm for
  158. C primitives or _FP_DIV_HELP_ldiv for the ISO function. Which you
  159. choose will depend on what the compiler does with divrem4. */
  160. #define _FP_DIV_MEAT_1_imm(fs, R, X, Y, doit) \
  161. do { \
  162. _FP_W_TYPE _q, _r; \
  163. X##_f <<= (X##_f < Y##_f \
  164. ? R##_e--, _FP_WFRACBITS_##fs \
  165. : _FP_WFRACBITS_##fs - 1); \
  166. doit(_q, _r, X##_f, Y##_f); \
  167. R##_f = _q | (_r != 0); \
  168. } while (0)
  169. /* GCC's longlong.h defines a 2W / 1W => (1W,1W) primitive udiv_qrnnd
  170. that may be useful in this situation. This first is for a primitive
  171. that requires normalization, the second for one that does not. Look
  172. for UDIV_NEEDS_NORMALIZATION to tell which your machine needs. */
  173. #define _FP_DIV_MEAT_1_udiv_norm(fs, R, X, Y) \
  174. do { \
  175. _FP_W_TYPE _nh, _nl, _q, _r, _y; \
  176. \
  177. /* Normalize Y -- i.e. make the most significant bit set. */ \
  178. _y = Y##_f << _FP_WFRACXBITS_##fs; \
  179. \
  180. /* Shift X op correspondingly high, that is, up one full word. */ \
  181. if (X##_f < Y##_f) \
  182. { \
  183. R##_e--; \
  184. _nl = 0; \
  185. _nh = X##_f; \
  186. } \
  187. else \
  188. { \
  189. _nl = X##_f << (_FP_W_TYPE_SIZE - 1); \
  190. _nh = X##_f >> 1; \
  191. } \
  192. \
  193. udiv_qrnnd(_q, _r, _nh, _nl, _y); \
  194. R##_f = _q | (_r != 0); \
  195. } while (0)
  196. #define _FP_DIV_MEAT_1_udiv(fs, R, X, Y) \
  197. do { \
  198. _FP_W_TYPE _nh, _nl, _q, _r; \
  199. if (X##_f < Y##_f) \
  200. { \
  201. R##_e--; \
  202. _nl = X##_f << _FP_WFRACBITS_##fs; \
  203. _nh = X##_f >> _FP_WFRACXBITS_##fs; \
  204. } \
  205. else \
  206. { \
  207. _nl = X##_f << (_FP_WFRACBITS_##fs - 1); \
  208. _nh = X##_f >> (_FP_WFRACXBITS_##fs + 1); \
  209. } \
  210. udiv_qrnnd(_q, _r, _nh, _nl, Y##_f); \
  211. R##_f = _q | (_r != 0); \
  212. } while (0)
  213. /*
  214. * Square root algorithms:
  215. * We have just one right now, maybe Newton approximation
  216. * should be added for those machines where division is fast.
  217. */
  218. #define _FP_SQRT_MEAT_1(R, S, T, X, q) \
  219. do { \
  220. while (q != _FP_WORK_ROUND) \
  221. { \
  222. T##_f = S##_f + q; \
  223. if (T##_f <= X##_f) \
  224. { \
  225. S##_f = T##_f + q; \
  226. X##_f -= T##_f; \
  227. R##_f += q; \
  228. } \
  229. _FP_FRAC_SLL_1(X, 1); \
  230. q >>= 1; \
  231. } \
  232. if (X##_f) \
  233. { \
  234. if (S##_f < X##_f) \
  235. R##_f |= _FP_WORK_ROUND; \
  236. R##_f |= _FP_WORK_STICKY; \
  237. } \
  238. } while (0)
  239. /*
  240. * Assembly/disassembly for converting to/from integral types.
  241. * No shifting or overflow handled here.
  242. */
  243. #define _FP_FRAC_ASSEMBLE_1(r, X, rsize) (r = X##_f)
  244. #define _FP_FRAC_DISASSEMBLE_1(X, r, rsize) (X##_f = r)
  245. /*
  246. * Convert FP values between word sizes
  247. */
  248. #define _FP_FRAC_CONV_1_1(dfs, sfs, D, S) \
  249. do { \
  250. D##_f = S##_f; \
  251. if (_FP_WFRACBITS_##sfs > _FP_WFRACBITS_##dfs) \
  252. { \
  253. if (S##_c != FP_CLS_NAN) \
  254. _FP_FRAC_SRS_1(D, (_FP_WFRACBITS_##sfs-_FP_WFRACBITS_##dfs), \
  255. _FP_WFRACBITS_##sfs); \
  256. else \
  257. _FP_FRAC_SRL_1(D, (_FP_WFRACBITS_##sfs-_FP_WFRACBITS_##dfs)); \
  258. } \
  259. else \
  260. D##_f <<= _FP_WFRACBITS_##dfs - _FP_WFRACBITS_##sfs; \
  261. } while (0)
  262. #endif /* __MATH_EMU_OP_1_H__ */