5 * Core division routines
8 * Copyright (C) 1992,1993,1994
9 * W. Metzenthen, 22 Parker St, Ormond, Vic 3163,
10 * Australia. E-mail billm@vaxc.cc.monash.edu.au
11 * All rights reserved.
13 * This copyright notice covers the redistribution and use of the
14 * FPU emulator developed by W. Metzenthen. It covers only its use
15 * in the 386BSD, FreeBSD and NetBSD operating systems. Any other
16 * use is not permitted under this copyright.
18 * Redistribution and use in source and binary forms, with or without
19 * modification, are permitted provided that the following conditions
21 * 1. Redistributions of source code must retain the above copyright
22 * notice, this list of conditions and the following disclaimer.
23 * 2. Redistributions in binary form must include information specifying
24 * that source code for the emulator is freely available and include
26 * a) an offer to provide the source code for a nominal distribution
28 * b) list at least two alternative methods whereby the source
29 * can be obtained, e.g. a publically accessible bulletin board
30 * and an anonymous ftp site from which the software can be
32 * 3. All advertising materials specifically mentioning features or use of
33 * this emulator must acknowledge that it was developed by W. Metzenthen.
34 * 4. The name of W. Metzenthen may not be used to endorse or promote
35 * products derived from this software without specific prior written
38 * THIS SOFTWARE IS PROVIDED ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES,
39 * INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY
40 * AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL
41 * W. METZENTHEN BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
42 * EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
43 * PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
44 * PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
45 * LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
46 * NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
47 * SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
50 * The purpose of this copyright, based upon the Berkeley copyright, is to
51 * ensure that the covered software remains freely available to everyone.
53 * The software (with necessary differences) is also available, but under
54 * the terms of the GNU copyleft, for the Linux operating system and for
55 * the djgpp ms-dos extender.
57 * W. Metzenthen June 1994.
64 /*---------------------------------------------------------------------------+
65 | Kernel for the division routines. |
67 | void reg_u_div(FPU_REG *a, FPU_REG *a, |
68 | FPU_REG *dest, unsigned int control_word) |
70 | Does not compute the destination exponent, but does adjust it. |
71 +---------------------------------------------------------------------------*/
73 #include <gnu/i386/fpemul/fpu_asm.h>
76 /* #define dSIGL(x) (x) */
77 /* #define dSIGH(x) 4(x) */
83 Result: accum_3:accum_2:accum_1:accum_0
84 Overflow flag: ovfl_flag
105 .globl _divide_kernel
115 movl PARAM1,%esi /* pointer to num */
116 movl PARAM2,%ebx /* pointer to denom */
117 movl PARAM3,%edi /* pointer to answer */
119 #ifdef DENORM_OPERAND
124 call _denormal_operand
133 call _denormal_operand
138 #endif DENORM_OPERAND
142 /* testl $0x80000000, SIGH(%esi) *//* Dividend */
144 testl $0x80000000, SIGH(%ebx) /* Divisor*/
148 /* Check if the divisor can be treated as having just 32 bits */
150 jnz L_Full_Division /* Can't do a quick divide */
152 /* We should be able to zip through the division here */
153 movl SIGH(%ebx),%ecx /* The divisor */
154 movl SIGH(%esi),%edx /* Dividend */
155 movl SIGL(%esi),%eax /* Dividend */
158 setaeb ovfl_flag /* Keep a record */
161 subl %ecx,%edx /* Prevent the overflow */
164 /* Divide the 64 bit number by the 32 bit denominator */
168 /* Work on the remainder of the first division */
173 /* Work on the remainder of the 64 bit division */
177 testb $255,ovfl_flag /* was the num > denom ? */
180 /* Do the shifting here */
181 /* increase the exponent */
184 /* shift the mantissa right one bit */
185 stc /* To set the ms bit */
191 jmp LRound_precision /* Do the rounding as required*/
194 /*---------------------------------------------------------------------------+
195 | Divide: Return arg1/arg2 to arg3. |
197 | This routine does not use the exponents of arg1 and arg2, but does |
198 | adjust the exponent of arg3. |
200 | The maximum returned value is (ignoring exponents) |
201 | .ffffffff ffffffff |
202 | ------------------ = 1.ffffffff fffffffe |
203 | .80000000 00000000 |
204 | and the minimum is |
205 | .80000000 00000000 |
206 | ------------------ = .80000000 00000001 (rounded) |
207 | .ffffffff ffffffff |
209 +---------------------------------------------------------------------------*/
213 /* Save extended dividend in local register*/
219 movl %eax,accum_1 /* zero the extension */
220 movl %eax,accum_0 /* zero the extension */
222 movl SIGL(%esi),%eax /* Get the current num */
225 /*----------------------------------------------------------------------*/
226 /* Initialization done */
227 /* Do the first 32 bits */
230 cmpl SIGH(%ebx),%edx /* Test for imminent overflow */
238 /* The dividend is greater or equal, would cause overflow */
239 setaeb ovfl_flag /* Keep a record */
242 sbbl SIGH(%ebx),%edx /* Prevent the overflow */
247 /* At this point, we have a dividend < divisor, with a record of
248 adjustment in ovfl_flag */
250 /* We will divide by a number which is too large */
255 /* here we need to divide by 100000000h,
256 i.e., no division at all.. */
261 divl %ecx /* Divide the numerator by the augmented
265 movl %eax,result_2 /* Put the result in the answer */
267 mull SIGH(%ebx) /* mul by the ms dw of the denom */
269 subl %eax,accum_2 /* Subtract from the num local reg */
272 movl result_2,%eax /* Get the result back */
273 mull SIGL(%ebx) /* now mul the ls dw of the denom */
275 subl %eax,accum_1 /* Subtract from the num local reg */
278 je LDo_2nd_32_bits /* Must check for non-zero result here */
284 /* need to subtract another once of the denom */
285 incl result_2 /* Correct the answer */
289 subl %eax,accum_1 /* Subtract from the num local reg */
294 jne L_bugged_1 /* Must check for non-zero result here */
297 /*----------------------------------------------------------------------*/
298 /* Half of the main problem is done, there is just a reduced numerator
300 /* Work with the second 32 bits, accum_0 not used from now on */
302 movl accum_2,%edx /* get the reduced num */
305 /* need to check for possible subsequent overflow */
308 ja LPrevent_2nd_overflow
313 LPrevent_2nd_overflow:
314 /* The numerator is greater or equal, would cause overflow */
315 /* prevent overflow */
321 incl result_2 /* Reflect the subtraction in the answer */
324 je L_bugged_2 /* Can't bump the result to 1.0 */
328 cmpl $0,%ecx /* augmented denom msw*/
329 jnz LSecond_div_not_1
331 /* %ecx == 0, we are dividing by 1.0 */
336 divl %ecx /* Divide the numerator by the denom ms dw */
339 movl %eax,result_1 /* Put the result in the answer */
341 mull SIGH(%ebx) /* mul by the ms dw of the denom */
343 subl %eax,accum_1 /* Subtract from the num local reg */
350 movl result_1,%eax /* Get the result back */
351 mull SIGL(%ebx) /* now mul the ls dw of the denom */
353 subl %eax,accum_0 /* Subtract from the num local reg */
354 sbbl %edx,accum_1 /* Subtract from the num local reg */
368 /* need to subtract another once of the denom */
371 subl %eax,accum_0 /* Subtract from the num local reg */
380 addl $1,result_1 /* Correct the answer */
384 jc L_bugged_2 /* Must check for non-zero result here */
387 /*----------------------------------------------------------------------*/
388 /* The division is essentially finished here, we just need to perform
389 tidying operations. */
390 /* deal with the 3rd 32 bits */
392 movl accum_1,%edx /* get the reduced num */
395 /* need to check for possible subsequent overflow */
396 cmpl SIGH(%ebx),%edx /* denom*/
398 ja LPrevent_3rd_overflow
400 cmpl SIGL(%ebx),%eax /* denom */
403 LPrevent_3rd_overflow:
404 /* prevent overflow */
410 addl $1,result_1 /* Reflect the subtraction in the answer */
415 /* This is a tricky spot, there is an overflow of the answer */
416 movb $255,ovfl_flag /* Overflow -> 1.000 */
419 /* Prepare for rounding.
420 // To test for rounding, we just need to compare 2*accum with the
426 jz LRound_ovfl /* The accumulator contains zero.*/
432 jc LRound_large /* No need to compare, denom smaller */
438 movl $0x70000000,%eax /* Denom was larger */
444 movl $0x80000000,%eax /* Remainder was exactly 1/2 denom */
448 movl $0xff000000,%eax /* Denom was smaller */
451 /* We are now ready to deal with rounding, but first we must get
452 the bits properly aligned */
453 testb $255,ovfl_flag /* was the num > denom ? */
458 /* shift the mantissa right one bit */
459 stc /* Will set the ms bit */
464 /* Round the result as required */
466 decl EXP(%edi) /* binary point between 1st & 2nd bits */
475 /* The logic is wrong if we got here */
477 pushl EX_INTERNAL|0x202
483 pushl EX_INTERNAL|0x203
489 pushl EX_INTERNAL|0x204