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