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