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