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