root/kernel/FPU-emu/reg_u_div.S

/* [previous][next][first][last][top][bottom][index][help] */
   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 #ifdef PARANOID
  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 #ifdef PARANOID
 135         jc      xL_bugged2
 136 #endif PARANOID
 137 
 138 xL_no_round:
 139         jmp     xL_done
 140 
 141 
 142 #ifdef PARANOID
 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 #ifdef PARANOID
 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 #ifdef PARANOID
 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 #ifdef PARANOID
 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 #ifdef PARANOID
 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 #ifdef PARANOID
 322         jc      L_bugged
 323 #endif PARANOID
 324 
 325         jz      L35             /* Just deal with rounding now */
 326 
 327 #ifdef PARANOID
 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 #ifdef PARANOID
 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 #ifdef PARANOID
 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 #ifdef PARANOID
 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 #ifdef PARANOID
 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

/* [previous][next][first][last][top][bottom][index][help] */