Efficiente funzione del piano intero in C++




performance x86-64 (4)

Il casting su int è notoriamente lento.

Forse hai vissuto sotto una roccia dal x86-64, o per il resto mancato che questo non sia stato vero per un po 'su x86. :)

SSE / SSE2 hanno un'istruzione da convertire con troncamento (anziché la modalità di arrotondamento predefinita). L'ISA supporta questa operazione in modo efficiente proprio perché la conversione con la semantica C non è rara nelle codebase reali. Il codice x86-64 utilizza registri SSE / SSE2 XMM per la matematica FP scalare, non x87, a causa di questo e di altre cose che lo rendono più efficiente. Anche il moderno codice a 32 bit utilizza i registri XMM per la matematica scalare.

Durante la compilazione per x87 (senza SSE3 fisttp ), i compilatori dovevano cambiare la modalità di arrotondamento x87 al troncamento, FP archiviare in memoria, quindi modificare nuovamente la modalità di arrotondamento. (E poi ricaricare il numero intero dalla memoria, tipicamente da un locale in pila, se si fa qualcosa in più con esso.) X87 è stato terribile per questo.

Sì, è stato terribilmente lento, ad esempio nel 2006, quando è stato scritto il link nella risposta di @ Kirjain, se avevi ancora una CPU a 32 bit o se stavi usando una CPU x86-64 per eseguire codice a 32 bit.

La conversione con una modalità di arrotondamento diversa da quella troncata o predefinita (quella più vicina) non è direttamente supportata, e fino a SSE4.1 roundps / roundpd tua migliore scommessa è stata quella dei numeri magici come nel collegamento del 2006 dalla risposta di @ Kirjain.

Alcuni trucchi interessanti, ma solo per il double -> 32 bit interi. Improbabile che valga la pena di espandersi per double se si dispone di float .

O più in generale, è sufficiente aggiungere un numero di grande magnitudine per attivare l'arrotondamento, quindi sottrarlo di nuovo per tornare all'intervallo originale. Questo può funzionare per il float senza espandersi per double , ma non sono sicuro di quanto sia facile far funzionare il floor .

Ad ogni modo, la soluzione più ovvia qui è _mm256_floor_ps() e _mm256_cvtps_epi32 ( vroundps e vcvtps2dq ). Una versione non AVX di questo può funzionare con SSE4.1.

Non sono sicuro di poter fare ancora meglio; Se avessi un enorme array da elaborare (e non riuscissi a inserire questo lavoro con altri lavori), potresti impostare la modalità di arrotondamento MXCSR su "verso -Inf" (piano) e utilizzare semplicemente vcvtps2dq (che usa la modalità di arrotondamento corrente). ). Quindi ripristinalo. Ma probabilmente è meglio bloccare la conversione in cache o eseguirla al volo mentre si generano i dati, presumibilmente da altri calcoli FP che richiedono la modalità di arrotondamento FP impostata su Arresto predefinito.

roundps / pd / ss / sd è 2 uops su CPU Intel, ma solo 1 uop (per 128-bit lane) su AMD Ryzen. cvtps2dq è anche 1 uop. la conversione doppia -> int include anche un shuffle. La conversione scalare FP-> int (che copia in un registro intero) di solito costa anche un extra-uop.

Quindi c'è spazio per la possibilità che i trucchi del numero magico siano una vittoria in alcuni casi; forse vale la pena di investigare se _mm256_floor_ps() + cvt fa parte di un collo di bottiglia critico (o più probabilmente se hai doppio e vuoi int32).

@ Cássio Renan int foo = floorf(f) verrà automaticamente auto-compilato se compilato con gcc -O3 -fno-trapping-math (o -ffast-math ), con -march= qualcosa che ha SSE4.1 o AVX. https://godbolt.org/z/ae_KPv

Questo è forse utile se lo stai usando con altri codici scalari che non sono vettorizzati manualmente. Soprattutto se speri che il compilatore auto-vettori tutto il resto.

Voglio definire una funzione efficiente del piano intero, ovvero una conversione da float o double che esegue il troncamento verso il meno infinito.

Possiamo supporre che i valori siano tali che non si verifichi un overflow di numeri interi. Finora ho alcune opzioni

  • casting per int; ciò richiede una gestione speciale dei valori negativi, in quanto il cast viene troncato verso zero;

    I= int(F); if (I < 0 && I != F) I--;
  • fondere il risultato di floor in int;

    int(floor(F));
  • lanciare a int con un grande spostamento per ottenere i positivi (questo può restituire risultati errati per valori elevati);

    int(F + double(0x7fffffff)) - 0x7fffffff;

Il casting su int è notoriamente lento. Così sono se le prove. Non ho cronometrato la funzione del piano, ma ho visto post che sostengono che è anche lento.

Riesci a pensare a migliori alternative in termini di velocità, precisione o intervallo consentito? Non ha bisogno di essere portatile. Gli obiettivi sono recenti architetture x86 / x64.


Dai un'occhiata ai numeri magici . L'algoritmo proposto sulla pagina web dovrebbe essere molto più efficiente della semplice trasmissione. Non l'ho mai usato da solo, ma questo è il confronto delle prestazioni che offrono sul sito (xs_ToInt e xs_CRoundToInt sono le funzioni proposte):

Performing 10000000 times:
simple cast           2819 ms i.e. i = (long)f;
xs_ToInt              1242 ms i.e. i = xs_ToInt(f); //numerically same as above
bit-twiddle(full)     1093 ms i.e. i = BitConvertToInt(f); //rounding from Fluid
fistp                  676 ms i.e. i = FISTToInt(f); //Herf, et al x86 Assembly rounding 
bit-twiddle(limited)   623 ms i.e. i = FloatTo23Bits(f); //Herf, rounding only in the range (0...1]  
xs_CRoundToInt         609 ms i.e. i = xs_CRoundToInt(f); //rounding with "magic" numbers

Inoltre, xs_ToInt è apparentemente modificato in modo che le prestazioni migliorino:

Performing 10000000 times:
simple cast convert   3186 ms i.e. fi = (f*65536);
fistp convert         3031 ms i.e. fi = FISTToInt(f*65536);
xs_ToFix               622 ms i.e. fi = xs_Fix<16>::ToFix(f);

Breve spiegazione di come funziona il metodo dei "numeri magici":

"Fondamentalmente, per aggiungere due numeri in virgola mobile, il processore" allinea "i punti decimali dei numeri in modo che possa facilmente aggiungere i bit. Ciò avviene" normalizzando "i numeri in modo tale che i bit più significativi vengano conservati , cioè il numero più piccolo "normalizza" per corrispondere a quello più grande, quindi il principio della conversione del "numero magico" che usa xs_CRoundToInt () è questo: aggiungiamo un numero in virgola mobile abbastanza grande (un numero così grande che ci sono cifre significative solo fino al punto decimale, e nessuno dopo di esso) a quello che stai convertendo in modo tale che: (a) il numero viene normalizzato dal processore al suo intero equivalente e (b) l'aggiunta dei due non cancella l'integrale bit significativi nel numero che stavi cercando di convertire (cioè XX00 + 00YY = XXYY). "

La citazione è presa dalla stessa pagina web.


Perché non usare solo questo:

#include <cmath>

auto floor_(float const x) noexcept
{
  int const t(x);

  return t - (t > x);
}

Se lo stai facendo in batch, il compilatore può autorizzarlo automaticamente, se sai cosa stai facendo. Ad esempio, ecco una piccola implementazione che autorizza la conversione di float in numeri interi, su GCC:

#include <cmath>

// Compile with -O3 and -march=native to see autovectorization
__attribute__((optimize("-fno-trapping-math")))
void testFunction(float* input, int* output, int length) {
  // Assume the input and output are aligned on a 32-bit boundary.
  // Of course, you have  to ensure this when calling testFunction, or else
  // you will have problems.
  input = static_cast<float*>(__builtin_assume_aligned(input, 32));
  output = static_cast<int*>(__builtin_assume_aligned(output, 32));

  // Also assume the length is a multiple of 32.
  if (length & 31) __builtin_unreachable();

  // Do the conversion
  for (int i = 0; i < length; ++i) {
    output[i] = floor(input[i]);
  }
}

Questo è l'assembly generato per x86-64 (con istruzioni AVX512):

testFunction(float*, int*, int):
        test    edx, edx
        jle     .L5
        lea     ecx, [rdx-1]
        xor     eax, eax
.L3:
        # you can see here that the conversion was vectorized
        # to a vrndscaleps (that will round the float appropriately)
        # and a vcvttps2dq (thal will perform the conversion)
        vrndscaleps     ymm0, YMMWORD PTR [rdi+rax], 1
        vcvttps2dq      ymm0, ymm0
        vmovdqa64       YMMWORD PTR [rsi+rax], ymm0
        add     rax, 32
        cmp     rax, rdx
        jne     .L3
        vzeroupper
.L5:
        ret

Se il tuo obiettivo non supporta AVX512, continuerà a essere auto-ottimizzato usando le istruzioni SSE4.1, supponendo che tu abbia quelle. Questo è l'output con -O3 -msse4.1 :

testFunction(float*, int*, int):
        test    edx, edx
        jle     .L1
        shr     edx, 2
        xor     eax, eax
        sal     rdx, 4
.L3:
        roundps xmm0, XMMWORD PTR [rdi+rax], 1
        cvttps2dq       xmm0, xmm0
        movaps  XMMWORD PTR [rsi+rax], xmm0
        add     rax, 16
        cmp     rax, rdx
        jne     .L3
.L1:
        ret

Guardalo vivere su Godbolt







floor