À un moment donné de mon programme, je calcule un diviseur entier d. A partir de là, d sera constant.

Plus tard dans le code, je diviserai plusieurs fois par d - en effectuant une division entière, puisque la valeur de d n'est pas une constante connue à la compilation.

Étant donné que la division entière est un processus relativement lent par rapport à d'autres types d'arithmétique entière, je voudrais l'optimiser. Existe-t-il un autre format dans lequel je pourrais stocker d, afin que le processus de division s'exécute plus rapidement? Peut-être une réciproque d'une certaine forme?

Je n'ai besoin de la valeur de d pour rien d'autre.

La valeur de d est n'importe quel entier 64 bits, mais s'intègre généralement assez bien dans 32 bits.

22
CygnusX1 27 juil. 2017 à 17:23

2 réponses

mise à jour - dans ma réponse originale, j'ai noté un algorithme mentionné dans un thread précédent pour le code généré par le compilateur pour diviser par constante. Le code d'assembly a été écrit pour correspondre à un document lié à partir de ce thread précédent. Le code généré par le compilateur implique des séquences légèrement différentes selon le diviseur.

Dans cette situation, le diviseur n'est pas connu avant l'exécution, un algorithme commun est donc souhaité. L'exemple de la réponse de geza montre un algorithme commun, qui pourrait être intégré dans le code d'assemblage avec GCC, mais Visual Studio ne prend pas en charge l'assemblage en ligne en mode 64 bits. Dans le cas de Visual Studio, il existe un compromis entre le code supplémentaire impliqué si vous utilisez des éléments intrinsèques et l'appel d'une fonction écrite en assembly. Sur mon système (Intel 3770k 3,5 GHz), j'ai essayé d'appeler une seule fonction qui fait | mul add adc shr |, et j'ai également essayé d'utiliser un pointeur pour fonctionner pour utiliser éventuellement des séquences plus courtes | mul shr | ou | shr (1) mul shr | selon le diviseur, mais cela a fourni peu ou pas de gain, selon le compilateur. La surcharge principale dans ce cas est l'appel (par rapport à | mul add adc shr |). Même avec la surcharge d'appel, la séquence | call mul add adc shr ret | en moyenne environ 4 fois plus vite que diviser sur mon système.

Notez que le code source lié à libdivide dans la réponse de geza n'a pas de routine commune qui peut gérer un diviseur == 1. La séquence commune de libdivide est multiplier, soustraire, déplacer (1), ajouter, déplacer, par rapport à l'exemple de séquence c ++ de geza de multiplier, ajouter, adc, décalage.


Ma réponse originale: l'exemple de code ci-dessous utilise l'algorithme décrit dans un fil de discussion précédent.

Pourquoi GCC utilise-t-il la multiplication par un nombre étrange dans l'implémentation de la division entière?

Ceci est un lien vers le document mentionné dans l'autre fil:

http://gmplib.org/~tege/divcnst-pldi94.pdf

L'exemple de code ci-dessous est basé sur le document pdf et est destiné à Visual Studio, utilisant ml64 (assembleur 64 bits), fonctionnant sous Windows (système d'exploitation 64 bits). Le code avec les étiquettes main ... et dcm ... est le code pour générer un pré-décalage (rspre, nombre de bits de fin de zéro dans le diviseur), un multiplicateur et un post-décalage (rspost). Le code avec les étiquettes dct ... est le code pour tester la méthode.

        includelib      msvcrtd
        includelib      oldnames

sw      equ     8                       ;size of word

        .data
arrd    dq      1                       ;array of test divisors
        dq      2
        dq      3
        dq      4
        dq      5
        dq      7
        dq      256
        dq      3*256
        dq      7*256
        dq      67116375
        dq      07fffffffffffffffh      ;max divisor
        dq      0
        .data?

        .code
        PUBLIC  main

main    PROC
        push    rbp
        push    rdi
        push    rsi
        sub     rsp,64                  ;allocate stack space
        mov     rbp,rsp
        lea     rsi,arrd                ;set ptr to array of divisors
        mov     [rbp+6*sw],rsi
        jmp     main1

main0:  mov     [rbp+0*sw],rbx          ;[rbp+0*sw] = rbx = divisor = d
        cmp     rbx,1                   ;if d <= 1, q=n or overflow
        jbe     main1
        bsf     rcx,rbx                 ;rcx = rspre
        mov     [rbp+1*sw],rcx          ;[rbp+1*sw] = rspre
        shr     rbx,cl                  ;rbx = d>>rsc
        bsr     rcx,rbx                 ;rcx = floor(log2(rbx))
        mov     rsi,1                   ;rsi = 2^floor(log2(rbx))
        shl     rsi,cl
        cmp     rsi,rbx                 ;br if power of 2
        je      dcm03
        inc     rcx                     ;rcx = ceil(log2(rcx)) = L = rspost
        shl     rsi,1                   ;rsi = 2^L
;       jz      main1                   ;d > 08000000000000000h, use compare
        add     rcx,[rbp+1*sw]          ;rcx = L+rspre
        cmp     rcx,8*sw                ;if d > 08000000000000000h, use compare
        jae     main1
        mov     rax,1                   ;[rbp+3*sw] = 2^(L+rspre)
        shl     rax,cl
        mov     [rbp+3*sw],rax
        sub     rcx,[rbp+1*sw]          ;rcx = L
        xor     rdx,rdx
        mov     rax,rsi                 ;hi N bits of 2^(N+L)
        div     rbx                     ;rax == 1
        xor     rax,rax                 ;lo N bits of 2^(N+L)
        div     rbx
        mov     rdi,rax                 ;rdi = mlo % 2^N
        xor     rdx,rdx
        mov     rax,rsi                 ;hi N bits of 2^(N+L) + 2^(L+rspre)
        div     rbx                     ;rax == 1
        mov     rax,[rbp+3*sw]          ;lo N bits of 2^(N+L) + 2^(L+rspre)
        div     rbx                     ;rax = mhi % 2^N
        mov     rdx,rdi                 ;rdx = mlo % 2^N
        mov     rbx,8*sw                ;rbx = e = # bits in word
dcm00:  mov     rsi,rdx                 ;rsi = mlo/2
        shr     rsi,1
        mov     rdi,rax                 ;rdi = mhi/2
        shr     rdi,1
        cmp     rsi,rdi                 ;break if mlo >= mhi
        jae     short dcm01
        mov     rdx,rsi                 ;rdx = mlo/2
        mov     rax,rdi                 ;rax = mhi/2
        dec     rbx                     ;e -= 1
        loop    dcm00                   ;loop if --shpost != 0
dcm01:  mov     [rbp+2*sw],rcx          ;[rbp+2*sw] = shpost
        cmp     rbx,8*sw                ;br if N+1 bit multiplier
        je      short dcm02
        xor     rdx,rdx
        mov     rdi,1                   ;rax = m = mhi + 2^e
        mov     rcx,rbx
        shl     rdi,cl
        or      rax,rdi
        jmp     short dct00

dcm02:  mov     rdx,1                   ;rdx = 2^N
        dec     qword ptr [rbp+2*sw]    ;dec rspost
        jmp     short dct00

dcm03:  mov     rcx,[rbp+1*sw]          ;rcx = rsc
        jmp     short dct10

;       test    rbx = n, rdx = m bit N, rax = m%(2^N)
;               [rbp+1*sw] = rspre, [rbp+2*sw] = rspost

dct00:  mov     rdi,rdx                 ;rdi:rsi = m
        mov     rsi,rax
        mov     rbx,0fffffffff0000000h  ;[rbp+5*sw] = rbx = n
dct01:  mov     [rbp+5*sw],rbx
        mov     rdx,rdi                 ;rdx:rax = m
        mov     rax,rsi
        mov     rcx,[rbp+1*sw]          ;rbx = n>>rspre
        shr     rbx,cl
        or      rdx,rdx                 ;br if 65 bit m
        jnz     short dct02
        mul     rbx                     ;rdx = (n*m)>>N
        jmp     short dct03

dct02:  mul     rbx
        sub     rbx,rdx
        shr     rbx,1
        add     rdx,rbx
dct03:  mov     rcx,[rbp+2*sw]          ;rcx = rspost
        shr     rdx,cl                  ;rdx = q = quotient
        mov     [rbp+4*sw],rdx          ;[rbp+4*sw] = q
        xor     rdx,rdx                 ;rdx:rax = n
        mov     rax,[rbp+5*sw]
        mov     rbx,[rbp+0*sw]          ;rbx = d
        div     rbx                     ;rax = n/d
        mov     rdx,[rbp+4*sw]          ;br if ok
        cmp     rax,rdx                 ;br if ok
        je      short dct04
        nop                             ;debug check
dct04:  mov     rbx,[rbp+5*sw]
        inc     rbx
        jnz     short dct01
        jmp     short main1

;       test    rbx = n, rcx = rsc
dct10:  mov     rbx,0fffffffff0000000h  ;rbx = n
dct11:  mov     rsi,rbx                 ;rsi = n
        shr     rsi,cl                  ;rsi = n>>rsc
        xor     edx,edx
        mov     rax,rbx
        mov     rdi,[rbp+0*sw]          ;rdi = d
        div     rdi
        cmp     rax,rsi                 ;br if ok
        je      short dct12
        nop
dct12:  inc     rbx
        jnz     short dct11

main1:  mov     rsi,[rbp+6*sw]          ;rsi ptr to divisor
        mov     rbx,[rsi]               ;rbx = divisor = d
        add     rsi,1*sw                ;advance ptr
        mov     [rbp+6*sw],rsi
        or      rbx,rbx
        jnz     main0                   ;br if not end table

        add     rsp,64                  ;restore regs
        pop     rsi
        pop     rdi
        pop     rbp
        xor     rax,rax
        ret     0

main    ENDP
        END
2
rcgldr 4 août 2017 à 23:42

Le livre "Le plaisir du pirate" a "Chapitre 10: Division entière par constante" s'étendant sur 74 pages. Vous pouvez trouver tous les exemples de code gratuitement dans ce répertoire: http://www.hackersdelight.org/hdcode.htm Dans votre cas, les fig. 10-1., 10-2 et 10-3 sont ce que vous voulez.

Le problème de la division par une constante d équivaut à une multiplication par c = 1 / d . Ces algorithmes calculent une telle constante pour vous. Une fois que vous avez c , vous calculez le dividende en tant que tel:

int divideByMyConstant(int dividend){
  int c = MAGIC; // Given by the algorithm

  // since 1/d < 1, c is actually (1<<k)/d to fit nicely ina 32 bit int
  int k = MAGIC_SHIFT; //Also given by the algorithm

  long long tmp = (long long)dividend * c; // use 64 bit to hold all the precision...

  tmp >>= k; // Manual floating point number =)

  return (int)tmp;
}
6
ZeroUltimax 31 juil. 2017 à 13:55