1 // SPDX-License-Identifier: GPL-2.0-or-later
2 /*
3 * Linux/PA-RISC Project (http://www.parisc-linux.org/)
4 *
5 * Floating-point emulation code
6 * Copyright (C) 2001 Hewlett-Packard (Paul Bame) <bame@debian.org>
7 */
8 /*
9 * BEGIN_DESC
10 *
11 * File:
12 * @(#) pa/spmath/dfdiv.c $Revision: 1.1 $
13 *
14 * Purpose:
15 * Double Precision Floating-point Divide
16 *
17 * External Interfaces:
18 * dbl_fdiv(srcptr1,srcptr2,dstptr,status)
19 *
20 * Internal Interfaces:
21 *
22 * Theory:
23 * <<please update with a overview of the operation of this file>>
24 *
25 * END_DESC
26 */
27
28
29 #include "float.h"
30 #include "dbl_float.h"
31
32 /*
33 * Double Precision Floating-point Divide
34 */
35
36 int
dbl_fdiv(dbl_floating_point * srcptr1,dbl_floating_point * srcptr2,dbl_floating_point * dstptr,unsigned int * status)37 dbl_fdiv (dbl_floating_point * srcptr1, dbl_floating_point * srcptr2,
38 dbl_floating_point * dstptr, unsigned int *status)
39 {
40 register unsigned int opnd1p1, opnd1p2, opnd2p1, opnd2p2;
41 register unsigned int opnd3p1, opnd3p2, resultp1, resultp2;
42 register int dest_exponent, count;
43 register boolean inexact = FALSE, guardbit = FALSE, stickybit = FALSE;
44 boolean is_tiny;
45
46 Dbl_copyfromptr(srcptr1,opnd1p1,opnd1p2);
47 Dbl_copyfromptr(srcptr2,opnd2p1,opnd2p2);
48 /*
49 * set sign bit of result
50 */
51 if (Dbl_sign(opnd1p1) ^ Dbl_sign(opnd2p1))
52 Dbl_setnegativezerop1(resultp1);
53 else Dbl_setzerop1(resultp1);
54 /*
55 * check first operand for NaN's or infinity
56 */
57 if (Dbl_isinfinity_exponent(opnd1p1)) {
58 if (Dbl_iszero_mantissa(opnd1p1,opnd1p2)) {
59 if (Dbl_isnotnan(opnd2p1,opnd2p2)) {
60 if (Dbl_isinfinity(opnd2p1,opnd2p2)) {
61 /*
62 * invalid since both operands
63 * are infinity
64 */
65 if (Is_invalidtrap_enabled())
66 return(INVALIDEXCEPTION);
67 Set_invalidflag();
68 Dbl_makequietnan(resultp1,resultp2);
69 Dbl_copytoptr(resultp1,resultp2,dstptr);
70 return(NOEXCEPTION);
71 }
72 /*
73 * return infinity
74 */
75 Dbl_setinfinity_exponentmantissa(resultp1,resultp2);
76 Dbl_copytoptr(resultp1,resultp2,dstptr);
77 return(NOEXCEPTION);
78 }
79 }
80 else {
81 /*
82 * is NaN; signaling or quiet?
83 */
84 if (Dbl_isone_signaling(opnd1p1)) {
85 /* trap if INVALIDTRAP enabled */
86 if (Is_invalidtrap_enabled())
87 return(INVALIDEXCEPTION);
88 /* make NaN quiet */
89 Set_invalidflag();
90 Dbl_set_quiet(opnd1p1);
91 }
92 /*
93 * is second operand a signaling NaN?
94 */
95 else if (Dbl_is_signalingnan(opnd2p1)) {
96 /* trap if INVALIDTRAP enabled */
97 if (Is_invalidtrap_enabled())
98 return(INVALIDEXCEPTION);
99 /* make NaN quiet */
100 Set_invalidflag();
101 Dbl_set_quiet(opnd2p1);
102 Dbl_copytoptr(opnd2p1,opnd2p2,dstptr);
103 return(NOEXCEPTION);
104 }
105 /*
106 * return quiet NaN
107 */
108 Dbl_copytoptr(opnd1p1,opnd1p2,dstptr);
109 return(NOEXCEPTION);
110 }
111 }
112 /*
113 * check second operand for NaN's or infinity
114 */
115 if (Dbl_isinfinity_exponent(opnd2p1)) {
116 if (Dbl_iszero_mantissa(opnd2p1,opnd2p2)) {
117 /*
118 * return zero
119 */
120 Dbl_setzero_exponentmantissa(resultp1,resultp2);
121 Dbl_copytoptr(resultp1,resultp2,dstptr);
122 return(NOEXCEPTION);
123 }
124 /*
125 * is NaN; signaling or quiet?
126 */
127 if (Dbl_isone_signaling(opnd2p1)) {
128 /* trap if INVALIDTRAP enabled */
129 if (Is_invalidtrap_enabled()) return(INVALIDEXCEPTION);
130 /* make NaN quiet */
131 Set_invalidflag();
132 Dbl_set_quiet(opnd2p1);
133 }
134 /*
135 * return quiet NaN
136 */
137 Dbl_copytoptr(opnd2p1,opnd2p2,dstptr);
138 return(NOEXCEPTION);
139 }
140 /*
141 * check for division by zero
142 */
143 if (Dbl_iszero_exponentmantissa(opnd2p1,opnd2p2)) {
144 if (Dbl_iszero_exponentmantissa(opnd1p1,opnd1p2)) {
145 /* invalid since both operands are zero */
146 if (Is_invalidtrap_enabled()) return(INVALIDEXCEPTION);
147 Set_invalidflag();
148 Dbl_makequietnan(resultp1,resultp2);
149 Dbl_copytoptr(resultp1,resultp2,dstptr);
150 return(NOEXCEPTION);
151 }
152 if (Is_divisionbyzerotrap_enabled())
153 return(DIVISIONBYZEROEXCEPTION);
154 Set_divisionbyzeroflag();
155 Dbl_setinfinity_exponentmantissa(resultp1,resultp2);
156 Dbl_copytoptr(resultp1,resultp2,dstptr);
157 return(NOEXCEPTION);
158 }
159 /*
160 * Generate exponent
161 */
162 dest_exponent = Dbl_exponent(opnd1p1) - Dbl_exponent(opnd2p1) + DBL_BIAS;
163
164 /*
165 * Generate mantissa
166 */
167 if (Dbl_isnotzero_exponent(opnd1p1)) {
168 /* set hidden bit */
169 Dbl_clear_signexponent_set_hidden(opnd1p1);
170 }
171 else {
172 /* check for zero */
173 if (Dbl_iszero_mantissa(opnd1p1,opnd1p2)) {
174 Dbl_setzero_exponentmantissa(resultp1,resultp2);
175 Dbl_copytoptr(resultp1,resultp2,dstptr);
176 return(NOEXCEPTION);
177 }
178 /* is denormalized, want to normalize */
179 Dbl_clear_signexponent(opnd1p1);
180 Dbl_leftshiftby1(opnd1p1,opnd1p2);
181 Dbl_normalize(opnd1p1,opnd1p2,dest_exponent);
182 }
183 /* opnd2 needs to have hidden bit set with msb in hidden bit */
184 if (Dbl_isnotzero_exponent(opnd2p1)) {
185 Dbl_clear_signexponent_set_hidden(opnd2p1);
186 }
187 else {
188 /* is denormalized; want to normalize */
189 Dbl_clear_signexponent(opnd2p1);
190 Dbl_leftshiftby1(opnd2p1,opnd2p2);
191 while (Dbl_iszero_hiddenhigh7mantissa(opnd2p1)) {
192 dest_exponent+=8;
193 Dbl_leftshiftby8(opnd2p1,opnd2p2);
194 }
195 if (Dbl_iszero_hiddenhigh3mantissa(opnd2p1)) {
196 dest_exponent+=4;
197 Dbl_leftshiftby4(opnd2p1,opnd2p2);
198 }
199 while (Dbl_iszero_hidden(opnd2p1)) {
200 dest_exponent++;
201 Dbl_leftshiftby1(opnd2p1,opnd2p2);
202 }
203 }
204
205 /* Divide the source mantissas */
206
207 /*
208 * A non-restoring divide algorithm is used.
209 */
210 Twoword_subtract(opnd1p1,opnd1p2,opnd2p1,opnd2p2);
211 Dbl_setzero(opnd3p1,opnd3p2);
212 for (count=1; count <= DBL_P && (opnd1p1 || opnd1p2); count++) {
213 Dbl_leftshiftby1(opnd1p1,opnd1p2);
214 Dbl_leftshiftby1(opnd3p1,opnd3p2);
215 if (Dbl_iszero_sign(opnd1p1)) {
216 Dbl_setone_lowmantissap2(opnd3p2);
217 Twoword_subtract(opnd1p1,opnd1p2,opnd2p1,opnd2p2);
218 }
219 else {
220 Twoword_add(opnd1p1, opnd1p2, opnd2p1, opnd2p2);
221 }
222 }
223 if (count <= DBL_P) {
224 Dbl_leftshiftby1(opnd3p1,opnd3p2);
225 Dbl_setone_lowmantissap2(opnd3p2);
226 Dbl_leftshift(opnd3p1,opnd3p2,(DBL_P-count));
227 if (Dbl_iszero_hidden(opnd3p1)) {
228 Dbl_leftshiftby1(opnd3p1,opnd3p2);
229 dest_exponent--;
230 }
231 }
232 else {
233 if (Dbl_iszero_hidden(opnd3p1)) {
234 /* need to get one more bit of result */
235 Dbl_leftshiftby1(opnd1p1,opnd1p2);
236 Dbl_leftshiftby1(opnd3p1,opnd3p2);
237 if (Dbl_iszero_sign(opnd1p1)) {
238 Dbl_setone_lowmantissap2(opnd3p2);
239 Twoword_subtract(opnd1p1,opnd1p2,opnd2p1,opnd2p2);
240 }
241 else {
242 Twoword_add(opnd1p1,opnd1p2,opnd2p1,opnd2p2);
243 }
244 dest_exponent--;
245 }
246 if (Dbl_iszero_sign(opnd1p1)) guardbit = TRUE;
247 stickybit = Dbl_allp1(opnd1p1) || Dbl_allp2(opnd1p2);
248 }
249 inexact = guardbit | stickybit;
250
251 /*
252 * round result
253 */
254 if (inexact && (dest_exponent > 0 || Is_underflowtrap_enabled())) {
255 Dbl_clear_signexponent(opnd3p1);
256 switch (Rounding_mode()) {
257 case ROUNDPLUS:
258 if (Dbl_iszero_sign(resultp1))
259 Dbl_increment(opnd3p1,opnd3p2);
260 break;
261 case ROUNDMINUS:
262 if (Dbl_isone_sign(resultp1))
263 Dbl_increment(opnd3p1,opnd3p2);
264 break;
265 case ROUNDNEAREST:
266 if (guardbit && (stickybit ||
267 Dbl_isone_lowmantissap2(opnd3p2))) {
268 Dbl_increment(opnd3p1,opnd3p2);
269 }
270 }
271 if (Dbl_isone_hidden(opnd3p1)) dest_exponent++;
272 }
273 Dbl_set_mantissa(resultp1,resultp2,opnd3p1,opnd3p2);
274
275 /*
276 * Test for overflow
277 */
278 if (dest_exponent >= DBL_INFINITY_EXPONENT) {
279 /* trap if OVERFLOWTRAP enabled */
280 if (Is_overflowtrap_enabled()) {
281 /*
282 * Adjust bias of result
283 */
284 Dbl_setwrapped_exponent(resultp1,dest_exponent,ovfl);
285 Dbl_copytoptr(resultp1,resultp2,dstptr);
286 if (inexact)
287 if (Is_inexacttrap_enabled())
288 return(OVERFLOWEXCEPTION | INEXACTEXCEPTION);
289 else Set_inexactflag();
290 return(OVERFLOWEXCEPTION);
291 }
292 Set_overflowflag();
293 /* set result to infinity or largest number */
294 Dbl_setoverflow(resultp1,resultp2);
295 inexact = TRUE;
296 }
297 /*
298 * Test for underflow
299 */
300 else if (dest_exponent <= 0) {
301 /* trap if UNDERFLOWTRAP enabled */
302 if (Is_underflowtrap_enabled()) {
303 /*
304 * Adjust bias of result
305 */
306 Dbl_setwrapped_exponent(resultp1,dest_exponent,unfl);
307 Dbl_copytoptr(resultp1,resultp2,dstptr);
308 if (inexact)
309 if (Is_inexacttrap_enabled())
310 return(UNDERFLOWEXCEPTION | INEXACTEXCEPTION);
311 else Set_inexactflag();
312 return(UNDERFLOWEXCEPTION);
313 }
314
315 /* Determine if should set underflow flag */
316 is_tiny = TRUE;
317 if (dest_exponent == 0 && inexact) {
318 switch (Rounding_mode()) {
319 case ROUNDPLUS:
320 if (Dbl_iszero_sign(resultp1)) {
321 Dbl_increment(opnd3p1,opnd3p2);
322 if (Dbl_isone_hiddenoverflow(opnd3p1))
323 is_tiny = FALSE;
324 Dbl_decrement(opnd3p1,opnd3p2);
325 }
326 break;
327 case ROUNDMINUS:
328 if (Dbl_isone_sign(resultp1)) {
329 Dbl_increment(opnd3p1,opnd3p2);
330 if (Dbl_isone_hiddenoverflow(opnd3p1))
331 is_tiny = FALSE;
332 Dbl_decrement(opnd3p1,opnd3p2);
333 }
334 break;
335 case ROUNDNEAREST:
336 if (guardbit && (stickybit ||
337 Dbl_isone_lowmantissap2(opnd3p2))) {
338 Dbl_increment(opnd3p1,opnd3p2);
339 if (Dbl_isone_hiddenoverflow(opnd3p1))
340 is_tiny = FALSE;
341 Dbl_decrement(opnd3p1,opnd3p2);
342 }
343 break;
344 }
345 }
346
347 /*
348 * denormalize result or set to signed zero
349 */
350 stickybit = inexact;
351 Dbl_denormalize(opnd3p1,opnd3p2,dest_exponent,guardbit,
352 stickybit,inexact);
353
354 /* return rounded number */
355 if (inexact) {
356 switch (Rounding_mode()) {
357 case ROUNDPLUS:
358 if (Dbl_iszero_sign(resultp1)) {
359 Dbl_increment(opnd3p1,opnd3p2);
360 }
361 break;
362 case ROUNDMINUS:
363 if (Dbl_isone_sign(resultp1)) {
364 Dbl_increment(opnd3p1,opnd3p2);
365 }
366 break;
367 case ROUNDNEAREST:
368 if (guardbit && (stickybit ||
369 Dbl_isone_lowmantissap2(opnd3p2))) {
370 Dbl_increment(opnd3p1,opnd3p2);
371 }
372 break;
373 }
374 if (is_tiny) Set_underflowflag();
375 }
376 Dbl_set_exponentmantissa(resultp1,resultp2,opnd3p1,opnd3p2);
377 }
378 else Dbl_set_exponent(resultp1,dest_exponent);
379 Dbl_copytoptr(resultp1,resultp2,dstptr);
380
381 /* check for inexact */
382 if (inexact) {
383 if (Is_inexacttrap_enabled()) return(INEXACTEXCEPTION);
384 else Set_inexactflag();
385 }
386 return(NOEXCEPTION);
387 }
388