c++ - recorrer - División entera repetida por un valor constante de tiempo de ejecución
recorrer listas en scheme (3)
El libro "Hacker''s delight" tiene "Capítulo 10: División de enteros por constante" que abarca 74 páginas. Puede encontrar todos los ejemplos de código de forma gratuita en este directorio: http://www.hackersdelight.org/hdcode.htm En su caso, las Figs. 10-1., 10-2 y 10-3 son lo que quieres.
El problema de dividir por una constante d es equivalente a mutiplying por c = 1 / d . Estos algoritmos calculan esa constante para usted. Una vez que tienes c , calculas el dividendo como tal:
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;
}
En algún punto de mi programa calculo un divisor entero d . Desde ese punto en adelante, d será constante.
Más adelante en el código dividiré entre eso d varias veces, realizando una división entera, ya que el valor de d no es una constante conocida en tiempo de compilación.
Dado que la división de enteros es un proceso relativamente lento en comparación con otro tipo de aritmética de enteros, me gustaría optimizarlo. ¿Hay algún formato alternativo en el que pueda almacenar d , para que el proceso de división se realice más rápido? Tal vez un recíproco de alguna forma?
No necesito el valor de d para nada más.
El valor de d es cualquier entero de 64 bits, pero por lo general encaja bastante bien en 32 bits.
Hay una biblioteca para esto: libdivide :
libdivide es una biblioteca de código abierto para optimizar la división de enteros
libdivide le permite reemplazar costosas divisiones de enteros con multiplicaciones y cambios de bits comparativamente baratos. Los compiladores suelen hacer esto, pero solo cuando el divisor se conoce en el momento de la compilación. libdivide le permite aprovecharla en tiempo de ejecución. El resultado es que la división de enteros puede ser más rápida, mucho más rápida. Además, libdivide le permite dividir un vector SSE2 por una constante de tiempo de ejecución, lo cual es especialmente bueno porque SSE2 no tiene instrucciones de división de enteros.
libdivide es gratuito y de código abierto con una licencia permisiva. El nombre "libdivide" es un poco chistoso, ya que no hay una biblioteca propiamente dicha: el código está empaquetado completamente como un único archivo de encabezado, con una API de C y C ++.
Puedes leer sobre el algoritmo detrás de él en este blog ; por ejemplo, esta entry .
Básicamente, el algoritmo subyacente es el mismo que los compiladores utilizan para optimizar la división por una constante, excepto que permite que estas optimizaciones de reducción de la fuerza se realicen en tiempo de ejecución.
Nota: puedes crear una versión aún más rápida de libdivide. La idea es que para cada divisor, siempre puedes crear un triplete ( mul / add / shift ), por lo que esta expresión da el resultado: ( num * mul + add ) >> shift (multiplicar es una multiplicación amplia aquí). ¡Curiosamente, este método podría superar la versión del compilador para una división constante para varias microbenchmarks!
Aquí está mi implementación (esto no se puede compilar de forma inmediata, pero se puede ver el algoritmo general):
struct Divider_u32 {
u32 mul;
u32 add;
s32 shift;
void set(u32 divider);
};
void Divider_u32::set(u32 divider) {
s32 l = indexOfMostSignificantBit(divider);
if (divider&(divider-1)) {
u64 m = static_cast<u64>(1)<<(l+32);
mul = static_cast<u32>(m/divider);
u32 rem = static_cast<u32>(m)-mul*divider;
u32 e = divider-rem;
if (e<static_cast<u32>(1)<<l) {
mul++;
add = 0;
} else {
add = mul;
}
shift = l;
} else {
if (divider==1) {
mul = 0xffffffff;
add = 0xffffffff;
shift = 0;
} else {
mul = 0x80000000;
add = 0;
shift = l-1;
}
}
}
u32 operator/(u32 v, const Divider_u32 &div) {
u32 t = static_cast<u32>((static_cast<u64>(v)*div.mul+div.add)>>32)>>div.shift;
return t;
}
actualización : en mi respuesta original, observé un algoritmo mencionado en un hilo anterior para el código generado por el compilador para dividir por constante. El código de ensamblaje se escribió para coincidir con un documento vinculado a ese hilo anterior. El código generado por el compilador involucra secuencias ligeramente diferentes dependiendo del divisor.
En esta situación, el divisor no se conoce hasta el tiempo de ejecución, por lo que se desea un algoritmo común. El ejemplo en la respuesta de geza muestra un algoritmo común, que podría estar incluido en el código de ensamblaje con GCC, pero Visual Studio no admite el ensamblaje en línea en el modo de 64 bits. En el caso de Visual Studio, hay un intercambio entre el código adicional involucrado si se usa intrínseco, en lugar de llamar a una función escrita en ensamblador. En mi sistema (Intel 3770k 3.5ghz) intenté llamar a una sola función que hace | mul add adc shr |, y también intenté usar un puntero para funcionar para usar opcionalmente secuencias más cortas | mul shr | o | shr (1) mul shr | dependiendo del divisor, pero esto proporcionó poca o ninguna ganancia, dependiendo del compilador. La sobrecarga principal en este caso es la llamada (versus | mul add adc shr |). Incluso con la sobrecarga de llamadas, la secuencia | call mul add adc shr ret | un promedio de aproximadamente 4 veces más rápido que la división en mi sistema.
Tenga en cuenta que el enlace al código fuente para libdivide en la respuesta de geza no tiene una rutina común que pueda manejar un divisor == 1. La secuencia común de libdivide es multiplicar, restar, shift (1), add, shift, versus la secuencia c ++ de geza de ejemplo de multiplicar, agregar, adc, shift.
Mi respuesta original: el código de ejemplo siguiente utiliza el algoritmo descrito en un hilo anterior.
¿Por qué GCC usa la multiplicación por un número extraño al implementar la división de enteros?
Este es un enlace al documento mencionado en el otro hilo:
http://gmplib.org/~tege/divcnst-pldi94.pdf
El código de ejemplo a continuación se basa en el documento pdf y está destinado a Visual Studio, utilizando ml64 (ensamblador de 64 bits), ejecutándose en Windows (sistema operativo de 64 bits). El código con etiquetas main ... y dcm ... es el código para generar un preshift (rspre, número de bits de cero finales en el divisor), multiplicador y posthift (rspost). El código con etiquetas dct ... es el código para probar el método.
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