arch/x86/math-emu/reg_u_div.S

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