1 .file "reg_u_div.S"
2 /*---------------------------------------------------------------------------+ 3 | reg_u_div.S | 4 | | 5 | Core division routines | 6 | | 7 | Copyright (C) 1992,1993 | 8 | W. Metzenthen, 22 Parker St, Ormond, Vic 3163, | 9 | Australia. E-mail apm233m@vaxc.cc.monash.edu.au | 10 | | 11 | | 12 +---------------------------------------------------------------------------*/ 13
14 /*---------------------------------------------------------------------------+ 15 | Kernel for the division routines. | 16 | | 17 | void reg_u_div(FPU_REG *a, FPU_REG *a, | 18 | FPU_REG *dest, unsigned int control_word) | 19 | | 20 | Does not compute the destination exponent, but does adjust it. | 21 +---------------------------------------------------------------------------*/ 22
23 #include "exception.h"
24 #include "fpu_asm.h"
25 #include "control_w.h"
26
27
28 // #define dSIGL(x) (x) 29 // #define dSIGH(x) 4(x) 30
31
32 .data
33 /* 34 Local storage: 35 Result: accum_3:accum_2:accum_1:accum_0 36 Overflow flag: ovfl_flag 37 */ 38 .align 2,0
39 accum_3:
40 .long 0
41 accum_2:
42 .long 0
43 accum_1:
44 .long 0
45 accum_0:
46 .long 0
47 result_1:
48 .long 0
49 result_2:
50 .long 0
51 ovfl_flag:
52 .byte 0
53
54
55 .text
56 .align 2,144
57
58 .globl _reg_u_div
59
60 .globl _divide_kernel
61
62 _reg_u_div:
63 pushl %ebp
64 movl %esp,%ebp
65
66 pushl %esi
67 pushl %edi
68 pushl %ebx
69
70 movl PARAM1,%esi /* pointer to num */ 71 movl PARAM2,%ebx /* pointer to denom */ 72 movl PARAM3,%edi /* pointer to answer */ 73
74 #ifdefDENORM_OPERAND 75 movl EXP(%esi),%eax
76 cmpl EXP_UNDER,%eax
77 jg xOp1_not_denorm
78
79 call _denormal_operand
80 orl %eax,%eax
81 jnz FPU_Arith_exit
82
83 xOp1_not_denorm:
84 movl EXP(%ebx),%eax
85 cmpl EXP_UNDER,%eax
86 jg xOp2_not_denorm
87
88 call _denormal_operand
89 orl %eax,%eax
90 jnz FPU_Arith_exit
91
92 xOp2_not_denorm:
93 #endif DENORM_OPERAND
94
95 _divide_kernel:
96 #ifdefPARANOID 97 // testl $0x80000000, SIGH(%esi) // Dividend 98 // je L_bugged 99 testl $0x80000000, SIGH(%ebx) // Divisor 100 je L_bugged
101 #endif PARANOID
102
103 /* Check if the divisor can be treated as having just 32 bits */ 104 cmpl $0,SIGL(%ebx)
105 jnz L_Full_Division /* Can't do a quick divide */ 106
107 /* We should be able to zip through the division here */ 108 movl SIGH(%ebx),%ecx /* The divisor */ 109 movl SIGH(%esi),%edx /* Dividend */ 110 movl SIGL(%esi),%eax /* Dividend */ 111
112 cmpl %ecx,%edx
113 setaeb ovfl_flag /* Keep a record */ 114 jb L_no_adjust
115
116 subl %ecx,%edx /* Prevent the overflow */ 117
118 L_no_adjust:
119 /* Divide the 64 bit number by the 32 bit denominator */ 120 divl %ecx
121 movl %eax,result_2
122
123 /* Work on the remainder of the first division */ 124 xorl %eax,%eax
125 divl %ecx
126 movl %eax,result_1
127
128 /* Work on the remainder of the 64 bit division */ 129 xorl %eax,%eax
130 divl %ecx
131
132 testb $255,ovfl_flag /* was the num > denom ? */ 133 je L_no_overflow
134
135 /* Do the shifting here */ 136 /* increase the exponent */ 137 incl EXP(%edi)
138
139 /* shift the mantissa right one bit */ 140 stc /* To set the ms bit */ 141 rcrl result_2
142 rcrl result_1
143 rcrl %eax
144
145 L_no_overflow:
146 jmp LRound_precision // Do the rounding as required 147
148
149 /*---------------------------------------------------------------------------+ 150 | Divide: Return arg1/arg2 to arg3. | 151 | | 152 | This routine does not use the exponents of arg1 and arg2, but does | 153 | adjust the exponent of arg3. | 154 | | 155 | The maximum returned value is (ignoring exponents) | 156 | .ffffffff ffffffff | 157 | ------------------ = 1.ffffffff fffffffe | 158 | .80000000 00000000 | 159 | and the minimum is | 160 | .80000000 00000000 | 161 | ------------------ = .80000000 00000001 (rounded) | 162 | .ffffffff ffffffff | 163 | | 164 +---------------------------------------------------------------------------*/ 165
166
167 L_Full_Division:
168 // Save extended dividend in local register 169 movl SIGL(%esi),%eax
170 movl %eax,accum_2
171 movl SIGH(%esi),%eax
172 movl %eax,accum_3
173 xorl %eax,%eax
174 movl %eax,accum_1 /* zero the extension */ 175 movl %eax,accum_0 /* zero the extension */ 176
177 movl SIGL(%esi),%eax /* Get the current num */ 178 movl SIGH(%esi),%edx
179
180 /*----------------------------------------------------------------------*/ 181 /* Initialization done */ 182 /* Do the first 32 bits */ 183
184 movb $0,ovfl_flag
185 cmpl SIGH(%ebx),%edx /* Test for imminent overflow */ 186 jb LLess_than_1
187 ja LGreater_than_1
188
189 cmpl SIGL(%ebx),%eax
190 jb LLess_than_1
191
192 LGreater_than_1:
193 /* The dividend is greater or equal, would cause overflow */ 194 setaeb ovfl_flag /* Keep a record */ 195
196 subl SIGL(%ebx),%eax
197 sbbl SIGH(%ebx),%edx /* Prevent the overflow */ 198 movl %eax,accum_2
199 movl %edx,accum_3
200
201 LLess_than_1:
202 /* At this point, we have a dividend < divisor, with a record of 203 adjustment in ovfl_flag */ 204
205 /* We will divide by a number which is too large */ 206 movl SIGH(%ebx),%ecx
207 addl $1,%ecx
208 jnc LFirst_div_not_1
209
210 /* here we need to divide by 100000000h, 211 i.e., no division at all.. */ 212 mov %edx,%eax
213 jmp LFirst_div_done
214
215 LFirst_div_not_1:
216 divl %ecx /* Divide the numerator by the augmented 217 denom ms dw */ 218
219 LFirst_div_done:
220 movl %eax,result_2 /* Put the result in the answer */ 221
222 mull SIGH(%ebx) /* mul by the ms dw of the denom */ 223
224 subl %eax,accum_2 /* Subtract from the num local reg */ 225 sbbl %edx,accum_3
226
227 movl result_2,%eax /* Get the result back */ 228 mull SIGL(%ebx) /* now mul the ls dw of the denom */ 229
230 subl %eax,accum_1 /* Subtract from the num local reg */ 231 sbbl %edx,accum_2
232 sbbl $0,accum_3
233 je LDo_2nd_32_bits /* Must check for non-zero result here */ 234
235 #ifdefPARANOID 236 jb L_bugged_1
237 #endif PARANOID
238
239 /* need to subtract another once of the denom */ 240 incl result_2 /* Correct the answer */ 241
242 movl SIGL(%ebx),%eax
243 movl SIGH(%ebx),%edx
244 subl %eax,accum_1 /* Subtract from the num local reg */ 245 sbbl %edx,accum_2
246
247 #ifdefPARANOID 248 sbbl $0,accum_3
249 jne L_bugged_1 /* Must check for non-zero result here */ 250 #endif PARANOID
251
252 /*----------------------------------------------------------------------*/ 253 /* Half of the main problem is done, there is just a reduced numerator 254 to handle now */ 255 /* Work with the second 32 bits, accum_0 not used from now on */ 256 LDo_2nd_32_bits:
257 movl accum_2,%edx /* get the reduced num */ 258 movl accum_1,%eax
259
260 /* need to check for possible subsequent overflow */ 261 cmpl SIGH(%ebx),%edx
262 jb LDo_2nd_div
263 ja LPrevent_2nd_overflow
264
265 cmpl SIGL(%ebx),%eax
266 jb LDo_2nd_div
267
268 LPrevent_2nd_overflow:
269 /* The numerator is greater or equal, would cause overflow */ 270 /* prevent overflow */ 271 subl SIGL(%ebx),%eax
272 sbbl SIGH(%ebx),%edx
273 movl %edx,accum_2
274 movl %eax,accum_1
275
276 incl result_2 /* Reflect the subtraction in the answer */ 277
278 #ifdefPARANOID 279 je L_bugged_2 /* Can't bump the result to 1.0 */ 280 #endif PARANOID
281
282 LDo_2nd_div:
283 cmpl $0,%ecx // augmented denom msw 284 jnz LSecond_div_not_1
285
286 /* %ecx == 0, we are dividing by 1.0 */ 287 mov %edx,%eax
288 jmp LSecond_div_done
289
290 LSecond_div_not_1:
291 divl %ecx /* Divide the numerator by the denom ms dw */ 292
293 LSecond_div_done:
294 movl %eax,result_1 /* Put the result in the answer */ 295
296 mull SIGH(%ebx) /* mul by the ms dw of the denom */ 297
298 subl %eax,accum_1 /* Subtract from the num local reg */ 299 sbbl %edx,accum_2
300
301 #ifdefPARANOID 302 jc L_bugged_2
303 #endif PARANOID
304
305 movl result_1,%eax /* Get the result back */ 306 mull SIGL(%ebx) /* now mul the ls dw of the denom */ 307
308 subl %eax,accum_0 /* Subtract from the num local reg */ 309 sbbl %edx,accum_1 /* Subtract from the num local reg */ 310 sbbl $0,accum_2
311
312 #ifdefPARANOID 313 jc L_bugged_2
314 #endif PARANOID
315
316 jz LDo_3rd_32_bits
317
318 #ifdefPARANOID 319 cmpl $1,accum_2
320 jne L_bugged_2
321 #endif PARANOID
322
323 /* need to subtract another once of the denom */ 324 movl SIGL(%ebx),%eax
325 movl SIGH(%ebx),%edx
326 subl %eax,accum_0 /* Subtract from the num local reg */ 327 sbbl %edx,accum_1
328 sbbl $0,accum_2
329
330 #ifdefPARANOID 331 jc L_bugged_2
332 jne L_bugged_2
333 #endif PARANOID
334
335 addl $1,result_1 /* Correct the answer */ 336 adcl $0,result_2
337
338 #ifdefPARANOID 339 jc L_bugged_2 /* Must check for non-zero result here */ 340 #endif PARANOID
341
342 /*----------------------------------------------------------------------*/ 343 /* The division is essentially finished here, we just need to perform 344 tidying operations. */ 345 /* deal with the 3rd 32 bits */ 346 LDo_3rd_32_bits:
347 movl accum_1,%edx /* get the reduced num */ 348 movl accum_0,%eax
349
350 /* need to check for possible subsequent overflow */ 351 cmpl SIGH(%ebx),%edx // denom 352 jb LRound_prep
353 ja LPrevent_3rd_overflow
354
355 cmpl SIGL(%ebx),%eax // denom 356 jb LRound_prep
357
358 LPrevent_3rd_overflow:
359 /* prevent overflow */ 360 subl SIGL(%ebx),%eax
361 sbbl SIGH(%ebx),%edx
362 movl %edx,accum_1
363 movl %eax,accum_0
364
365 addl $1,result_1 /* Reflect the subtraction in the answer */ 366 adcl $0,result_2
367 jne LRound_prep
368 jnc LRound_prep
369
370 /* This is a tricky spot, there is an overflow of the answer */ 371 movb $255,ovfl_flag /* Overflow -> 1.000 */ 372
373 LRound_prep:
374 // Prepare for rounding. 375 // To test for rounding, we just need to compare 2*accum with the 376 // denom. 377 movl accum_0,%ecx
378 movl accum_1,%edx
379 movl %ecx,%eax
380 orl %edx,%eax
381 jz LRound_ovfl // The accumulator contains zero. 382
383 // Multiply by 2 384 clc
385 rcll $1,%ecx
386 rcll $1,%edx
387 jc LRound_large // No need to compare, denom smaller 388
389 subl SIGL(%ebx),%ecx
390 sbbl SIGH(%ebx),%edx
391 jnc LRound_not_small
392
393 movl $0x70000000,%eax // Denom was larger 394 jmp LRound_ovfl
395
396 LRound_not_small:
397 jnz LRound_large
398
399 movl $0x80000000,%eax // Remainder was exactly 1/2 denom 400 jmp LRound_ovfl
401
402 LRound_large:
403 movl $0xff000000,%eax // Denom was smaller 404
405 LRound_ovfl:
406 /* We are now ready to deal with rounding, but first we must get 407 the bits properly aligned */ 408 testb $255,ovfl_flag /* was the num > denom ? */ 409 je LRound_precision
410
411 incl EXP(%edi)
412
413 /* shift the mantissa right one bit */ 414 stc // Will set the ms bit 415 rcrl result_2
416 rcrl result_1
417 rcrl %eax
418
419 // Round the result as required 420 LRound_precision:
421 decl EXP(%edi) /* binary point between 1st & 2nd bits */ 422
423 movl %eax,%edx
424 movl result_1,%ebx
425 movl result_2,%eax
426 jmp FPU_round
427
428
429 #ifdefPARANOID 430 /* The logic is wrong if we got here */ 431 L_bugged:
432 pushl EX_INTERNAL|0x202
433 call EXCEPTION
434 pop %ebx
435 jmp L_exit
436
437 L_bugged_1:
438 pushl EX_INTERNAL|0x203
439 call EXCEPTION
440 pop %ebx
441 jmp L_exit
442
443 L_bugged_2:
444 pushl EX_INTERNAL|0x204
445 call EXCEPTION
446 pop %ebx
447 jmp L_exit
448
449 L_exit:
450 popl %ebx
451 popl %edi
452 popl %esi
453
454 leave
455 ret
456 #endif PARANOID