Est-il possible de rouler une version nettement plus rapide de sqrt

Question

Dans une application je profilage, je trouve que dans certains scénarios cette fonction est en mesure de prendre plus de 10% du temps total d'exécution.

J'ai vu la discussion au cours des années de mises en œuvre plus rapide sqrt utilisant la ruse sournoise à virgule flottante, mais je ne sais pas si ces choses sont obsolètes sur les processeurs modernes.

MSVC ++ 2008 compilateur est utilisé, pour référence ... mais je suppose que sqrt ne va pas ajouter de frais généraux bien.

Voir aussi ici pour discussion similaire sur cette est une méthode largement utilisée, mais est-il réellement beaucoup plus rapide? Combien de cycles est de toute façon RACINE ces jours-ci?

Answer 1

Oui, il est possible même sans supercherie:

1) précision de sacrifice pour la vitesse:. L'algorithme de sqrt est itérative, re-mettre en œuvre avec moins d'itérations

2) tables de recherche:. Soit juste pour le point de départ de l'itération, ou combiné avec interpolation pour vous aider à tout le chemin

3) la mise en cache: êtes-vous sqrting toujours le même ensemble limité de valeurs? le cas échéant, la mise en cache peut bien fonctionner. J'ai trouvé cela utile dans les applications graphiques où la même chose est calculée pour beaucoup de formes de la même taille, de sorte que les résultats peuvent être utilement mises en cache.

Answer 2

Il y a une grande table de comparaison ici: http://assemblyrequired.crashworks.org/timing-square-root/

Longue histoire courte, les ssqrts de SSE2 est d'environ 2x plus rapide que FPU fsqrt, et une approximation + itération est d'environ 4x plus rapide que celle (8x globale).

En outre, si vous essayez de prendre une sqrt simple précision, assurez-vous que est en fait ce que vous obtenez. Je l'ai entendu parler d'au moins un compilateur qui transformerait l'argument flottant à un double, appelez sqrt double précision, puis reconvertir à flotter.

Answer 3

Vous êtes très probablement d'obtenir des améliorations plus rapides en changeant votre algorithmes qu'en changeant leur implémentations : Essayez d'appel sqrt() moins au lieu de faire des appels plus rapidement. (Et si vous pensez que cela est impossible - les améliorations pour sqrt() que vous mentionnez ne sont que: l'amélioration de la algorithme utilisé pour calculer une racine carrée.)

Comme il est utilisé très souvent, il est probable que votre mise en œuvre de la bibliothèque standard de sqrt() est presque optimale pour le cas général. Sauf si vous avez un domaine restreint (par exemple, si vous avez besoin moins de précision) où l'algorithme peut prendre des raccourcis, il est très peu probable que quelqu'un arrive avec une implémentation qui est plus rapide.

Notez que, puisque cette fonction utilise 10% de votre temps d'exécution, même si vous parvenez à trouver une mise en œuvre qui ne prend 75% du temps de std::sqrt(), ce encore apportera seulement votre temps d'exécution par 2,5% . Pour la plupart des applications utilisateurs ne remarquent même pas cela, sauf si elles utilisent une montre pour mesurer.

Answer 4

Quelle précision besoin de votre sqrt-vous être? Vous pouvez obtenir une approximation raisonnable très rapidement: voir l'excellent Quake3 inverse fonction racine carrée d'inspiration (notez que le code sous GPL, vous voudrez peut-être pas intégrer directement).

Answer 5

Je ne sais pas si vous avez fixé, mais je l'ai lu à ce sujet avant, et il semble que le plus rapide chose à faire est de remplacer la fonction sqrt avec une version de montage en ligne;

vous pouvez voir une description d'une charge d'alternatives

Voici un moyen rapide avec une table de seulement 8 Ko. L'erreur est d'environ 0,5% du résultat. Vous pouvez facilement agrandir la table, ce qui réduit l'erreur. Runs environ 5 fois plus rapide que la sqrt régulière ()

// LUT for fast sqrt of floats. Table will be consist of 2 parts, half for sqrt(X) and half for sqrt(2X).
const int nBitsForSQRTprecision = 11;                       // Use only 11 most sagnificant bits from the 23 of float. We can use 15 bits instead. It will produce less error but take more place in a memory. 
const int nUnusedBits   = 23 - nBitsForSQRTprecision;       // Amount of bits we will disregard
const int tableSize     = (1 << (nBitsForSQRTprecision+1)); // 2^nBits*2 because we have 2 halves of the table.
static short sqrtTab[tableSize]; 
static unsigned char is_sqrttab_initialized = FALSE;        // Once initialized will be true

// Table of precalculated sqrt() for future fast calculation. Approximates the exact with an error of about 0.5%
// Note: To access the bits of a float in C quickly we must misuse pointers.
// More info in: http://en.wikipedia.org/wiki/Single_precision
void build_fsqrt_table(void){
    unsigned short i;
    float f;
    UINT32 *fi = (UINT32*)&f;

    if (is_sqrttab_initialized)
        return;

    const int halfTableSize = (tableSize>>1);
    for (i=0; i < halfTableSize; i++){
         *fi = 0;
         *fi = (i << nUnusedBits) | (127 << 23); // Build a float with the bit pattern i as mantissa, and an exponent of 0, stored as 127

         // Take the square root then strip the first 'nBitsForSQRTprecision' bits of the mantissa into the table
         f = sqrtf(f);
         sqrtTab[i] = (short)((*fi & 0x7fffff) >> nUnusedBits);

         // Repeat the process, this time with an exponent of 1, stored as 128
         *fi = 0;
         *fi = (i << nUnusedBits) | (128 << 23);
         f = sqrtf(f);
         sqrtTab[i+halfTableSize] = (short)((*fi & 0x7fffff) >> nUnusedBits);
    }
    is_sqrttab_initialized = TRUE;
}

// Calculation of a square root. Divide the exponent of float by 2 and sqrt() its mantissa using the precalculated table.
float fast_float_sqrt(float n){
    if (n <= 0.f) 
        return 0.f;                           // On 0 or negative return 0.
    UINT32 *num = (UINT32*)&n;
    short e;                                  // Exponent
    e = (*num >> 23) - 127;                   // In 'float' the exponent is stored with 127 added.
    *num &= 0x7fffff;                         // leave only the mantissa 

    // If the exponent is odd so we have to look it up in the second half of the lookup table, so we set the high bit.
    const int halfTableSize = (tableSize>>1);
    const int secondHalphTableIdBit = halfTableSize << nUnusedBits;
    if (e & 0x01) 
        *num |= secondHalphTableIdBit;  
    e >>= 1;                                  // Divide the exponent by two (note that in C the shift operators are sign preserving for signed operands

    // Do the table lookup, based on the quaternary mantissa, then reconstruct the result back into a float
    *num = ((sqrtTab[*num >> nUnusedBits]) << nUnusedBits) | ((e + 127) << 23);
    return n;
}

Answer 6

Answer 7

Answer 8

Voici un moyen rapide avec une table de seulement 8 Ko. L'erreur est d'environ 0,5% du résultat. Vous pouvez facilement agrandir la table, ce qui réduit l'erreur. Runs environ 5 fois plus rapide que la sqrt régulière ()

// LUT for fast sqrt of floats. Table will be consist of 2 parts, half for sqrt(X) and half for sqrt(2X).
const int nBitsForSQRTprecision = 11;                       // Use only 11 most sagnificant bits from the 23 of float. We can use 15 bits instead. It will produce less error but take more place in a memory. 
const int nUnusedBits   = 23 - nBitsForSQRTprecision;       // Amount of bits we will disregard
const int tableSize     = (1 << (nBitsForSQRTprecision+1)); // 2^nBits*2 because we have 2 halves of the table.
static short sqrtTab[tableSize]; 
static unsigned char is_sqrttab_initialized = FALSE;        // Once initialized will be true

// Table of precalculated sqrt() for future fast calculation. Approximates the exact with an error of about 0.5%
// Note: To access the bits of a float in C quickly we must misuse pointers.
// More info in: http://en.wikipedia.org/wiki/Single_precision
void build_fsqrt_table(void){
    unsigned short i;
    float f;
    UINT32 *fi = (UINT32*)&f;

    if (is_sqrttab_initialized)
        return;

    const int halfTableSize = (tableSize>>1);
    for (i=0; i < halfTableSize; i++){
         *fi = 0;
         *fi = (i << nUnusedBits) | (127 << 23); // Build a float with the bit pattern i as mantissa, and an exponent of 0, stored as 127

         // Take the square root then strip the first 'nBitsForSQRTprecision' bits of the mantissa into the table
         f = sqrtf(f);
         sqrtTab[i] = (short)((*fi & 0x7fffff) >> nUnusedBits);

         // Repeat the process, this time with an exponent of 1, stored as 128
         *fi = 0;
         *fi = (i << nUnusedBits) | (128 << 23);
         f = sqrtf(f);
         sqrtTab[i+halfTableSize] = (short)((*fi & 0x7fffff) >> nUnusedBits);
    }
    is_sqrttab_initialized = TRUE;
}

// Calculation of a square root. Divide the exponent of float by 2 and sqrt() its mantissa using the precalculated table.
float fast_float_sqrt(float n){
    if (n <= 0.f) 
        return 0.f;                           // On 0 or negative return 0.
    UINT32 *num = (UINT32*)&n;
    short e;                                  // Exponent
    e = (*num >> 23) - 127;                   // In 'float' the exponent is stored with 127 added.
    *num &= 0x7fffff;                         // leave only the mantissa 

    // If the exponent is odd so we have to look it up in the second half of the lookup table, so we set the high bit.
    const int halfTableSize = (tableSize>>1);
    const int secondHalphTableIdBit = halfTableSize << nUnusedBits;
    if (e & 0x01) 
        *num |= secondHalphTableIdBit;  
    e >>= 1;                                  // Divide the exponent by two (note that in C the shift operators are sign preserving for signed operands

    // Do the table lookup, based on the quaternary mantissa, then reconstruct the result back into a float
    *num = ((sqrtTab[*num >> nUnusedBits]) << nUnusedBits) | ((e + 127) << 23);
    return n;
}

Answer 9

Answer 10

Answer 11