Note: This website is archived. For up-to-date information about D projects and development, please visit wiki.dlang.org.

Changeset 884

Show
Ignore:
Timestamp:
12/15/08 05:56:03 (16 years ago)
Author:
walter
Message:

improved long division speed

Files:

Legend:

Unmodified
Added
Removed
Modified
Copied
Moved

  • branches/phobos-1.x/phobos/internal/llmath.d

    r295 r884  
    11// llmath.d 
    2 // Copyright (C) 1993-2003 by Digital Mars, www.digitalmars.com 
     2// Copyright (C) 1993-2008 by Digital Mars, http://www.digitalmars.com 
    33// All Rights Reserved 
    44// Written by Walter Bright 
    55 
    66// Compiler runtime support for 64 bit longs 
    77 
    88extern (C): 
    99 
    1010 
    1111/*************************************** 
    1212 * Unsigned long divide. 
    1313 * Input: 
    1414 *  [EDX,EAX],[ECX,EBX] 
    1515 * Output: 
    1616 *  [EDX,EAX] = [EDX,EAX] / [ECX,EBX] 
    1717 *  [ECX,EBX] = [EDX,EAX] % [ECX,EBX] 
    18  *  ESI,EDI destroyed 
    1918 */ 
    2019 
    2120void __ULDIV__() 
    2221{ 
    2322    asm 
    2423    { 
    2524    naked           ; 
    2625    test    ECX,ECX     ; 
    2726    jz  uldiv       ; 
    2827 
     28    // if ECX > EDX, then quotient is 0 and remainder is [EDX,EAX] 
     29    cmp ECX,EDX     ; 
     30    ja  quo0        ; 
     31 
     32    test    ECX,ECX     ; 
     33    js  Lleft       ; 
     34 
     35    /* We have n>d, and know that n/d will fit in 32 bits. 
     36         * d will be left justified if we shift it left s bits. 
     37     * [d1,d0] <<= s 
     38     * [n2,n1,n0] = [n1,n0] << s 
     39     * 
     40     * Use one divide, by this reasoning: 
     41     * ([n2,n1]<<32 + n0)/(d1<<32 + d0) 
     42     * becomes: 
     43     * ([n2,n1]<<32)/(d1<<32 + d0) + n0/(d1<<32 + d0) 
     44     * The second divide is always 0. 
     45     * Ignore the d0 in the first divide, which will yield a quotient 
     46     * that might be too high by 1 (because d1 is left justified). 
     47     * We can tell if it's too big if: 
     48     *  q*[d1,d0] > [n2,n1,n0] 
     49     * which is: 
     50     *  q*[d1,d0] > [[q*[d1,0]+q%[d1,0],n1,n0] 
     51     * If we subtract q*[d1,0] from both sides, we get: 
     52     *  q*d0 > [[n2,n1]%d1,n0] 
     53     * So if it is too big by one, reduce q by one to q'=q-one. 
     54     * Compute remainder as: 
     55     *  r = ([n1,n0] - q'*[d1,d0]) >> s 
     56     * Again, we can subtract q*[d1,0]: 
     57     *  r = ([n1,n0] - q*[d1,0] - (q'*[d1,d0] - q*[d1,0])) >> s 
     58     *  r = ([[n2,n1]%d1,n0] + (q*[d1,0] - (q - one)*[d1,d0])) >> s 
     59     *  r = ([[n2,n1]%d1,n0] + (q*[d1,0] - [d1 *(q-one),d0*(1-q)])) >> s 
     60     *  r = ([[n2,n1]%d1,n0] + [d1 *one,d0*(one-q)])) >> s 
     61     */ 
     62 
    2963    push    EBP     ; 
    30  
    31     // left justify [ECX,EBX] and leave count of shifts + 1 in EBP 
    32  
    33     mov EBP,1       ;   // at least 1 shift 
    34     test    ECX,ECX     ;   // left justified? 
    35     js  L1      ;   // yes 
    36     jnz L2      ; 
    37     add EBP,8       ; 
    38     mov CH,CL       ; 
    39     mov CL,BH       ; 
    40     mov BH,BL       ; 
    41     xor BL,BL       ;   // [ECX,EBX] <<= 8 
    42     test    ECX,ECX     ; 
    43     js  L1      ; 
    44     even            ; 
    45 L2: inc EBP     ;   // another shift 
    46     shl EBX,1       ; 
    47     rcl ECX,1       ;   // [ECX,EBX] <<= 1 
    48     jno L2      ;   // not left justified yet 
    49  
    50 L1: mov ESI,ECX     ; 
    51     mov EDI,EBX     ;   // [ESI,EDI] = divisor 
    52  
    53     mov ECX,EDX     ; 
    54     mov EBX,EAX     ;   // [ECX,EBX] = [EDX,EAX] 
    55     xor EAX,EAX     ; 
    56     cdq         ;   // [EDX,EAX] = 0 
    57     even            ; 
    58 L4: cmp ESI,ECX     ;   // is [ECX,EBX] > [ESI,EDI]? 
    59     ja  L3      ;   // yes 
    60     jb  L5      ;   // definitely less than 
    61     cmp EDI,EBX     ;   // check low order word 
    62     ja  L3      ; 
    63 L5: sub EBX,EDI     ; 
    64     sbb ECX,ESI     ;   // [ECX,EBX] -= [ESI,EDI] 
    65     stc         ;   // rotate in a 1 
    66 L3: rcl EAX,1       ; 
    67     rcl EDX,1       ;   // [EDX,EAX] = ([EDX,EAX] << 1) + C 
    68     shr ESI,1       ; 
    69     rcr EDI,1       ;   // [ESI,EDI] >>= 1 
    70     dec EBP     ;   // control count 
    71     jne L4      ; 
    72     pop EBP     ; 
    73     ret         ; 
    74  
    75 div0:   mov EAX,-1      ; 
    76     cwd         ;   // quotient is -1 
    77 //  xor ECX,ECX     ; 
    78 //  mov EBX,ECX     ;   // remainder is 0 (ECX and EBX already 0) 
     64    push    ESI     ; 
     65    push    EDI     ; 
     66 
     67    mov ESI,EDX     ; 
     68    mov EDI,EAX     ; 
     69    mov EBP,ECX     ; 
     70 
     71    bsr EAX,ECX     ;   // EAX is now 30..0 
     72    xor EAX,0x1F    ;   // EAX is now 1..31 
     73    mov CH,AL       ; 
     74    neg EAX     ; 
     75    add EAX,32      ; 
     76    mov CL,AL       ; 
     77 
     78    mov EAX,EBX     ; 
     79    shr EAX,CL      ; 
     80    xchg    CH,CL       ; 
     81    shl EBP,CL      ; 
     82    or  EBP,EAX     ; 
     83    shl EBX,CL      ; 
     84 
     85    mov EDX,ESI     ; 
     86    xchg    CH,CL       ; 
     87    shr EDX,CL      ; 
     88 
     89    mov EAX,EDI     ; 
     90    shr EAX,CL      ; 
     91    xchg    CH,CL       ; 
     92    shl EDI,CL      ; 
     93    shl ESI,CL      ; 
     94    or  EAX,ESI     ; 
     95 
     96    div EBP     ; 
     97    push    EBP     ; 
     98    mov EBP,EAX     ; 
     99    mov ESI,EDX     ; 
     100 
     101    mul EBX     ; 
     102    cmp EDX,ESI     ; 
     103    ja  L1      ; 
     104    jb  L2      ; 
     105    cmp EAX,EDI     ; 
     106    jbe L2      ; 
     107L1: dec EBP     ; 
     108    sub EAX,EBX     ; 
     109    sbb EDX,0[ESP]  ; 
     110L2: 
     111    add ESP,4       ; 
     112    sub EDI,EAX     ; 
     113    sbb ESI,EDX     ; 
     114    mov EAX,ESI     ; 
     115    xchg    CH,CL       ; 
     116    shl EAX,CL      ; 
     117    xchg    CH,CL       ; 
     118    shr EDI,CL      ; 
     119    or  EDI,EAX     ; 
     120    shr ESI,CL      ; 
     121    mov EBX,EDI     ; 
     122    mov ECX,ESI     ; 
     123    mov EAX,EBP     ; 
     124    xor EDX,EDX     ; 
     125 
     126    pop EDI     ; 
     127    pop ESI     ; 
    79128    pop EBP     ; 
    80129    ret         ; 
    81130 
    82131uldiv:  test    EDX,EDX     ; 
    83132    jnz D3      ; 
    84133    // Both high words are 0, we can use the DIV instruction 
    85134    div EBX     ; 
    86135    mov EBX,EDX     ; 
    87136    mov EDX,ECX     ;   // EDX = ECX = 0 
    88137    ret         ; 
    …… 
    93142    mov EAX,EDX     ; 
    94143    xor EDX,EDX     ; 
    95144    div EBX     ; 
    96145    xchg    ECX,EAX     ; 
    97146    div EBX     ; 
    98147    // ECX,EAX = result 
    99148    // EDX = remainder 
    100149    mov EBX,EDX     ; 
    101150    mov EDX,ECX     ; 
    102151    xor ECX,ECX     ; 
     152    ret         ; 
     153 
     154quo0:   // Quotient is 0 
     155    // Remainder is [EDX,EAX] 
     156    mov EBX,EAX     ; 
     157    mov ECX,EDX     ; 
     158    xor EAX,EAX     ; 
     159    xor EDX,EDX     ; 
     160    ret         ; 
     161 
     162Lleft:  // The quotient is 0 or 1 and EDX >= ECX 
     163    cmp EDX,ECX     ; 
     164    ja  quo1        ;   // [EDX,EAX] > [ECX,EBX] 
     165    // EDX == ECX 
     166    cmp EAX,EBX     ; 
     167    jb  quo0        ; 
     168 
     169quo1:   // Quotient is 1 
     170    // Remainder is [EDX,EAX] - [ECX,EBX] 
     171    sub EAX,EBX     ; 
     172    sbb EDX,ECX     ; 
     173    mov EBX,EAX     ; 
     174    mov ECX,EDX     ; 
     175    mov EAX,1       ; 
     176    xor EDX,EDX     ; 
    103177    ret         ; 
    104178    } 
    105179} 
    106180 
    107181 
    108182/*************************************** 
    109183 * Signed long divide. 
    110184 * Input: 
    111185 *  [EDX,EAX],[ECX,EBX] 
    112186 * Output: