Pourquoi GCC n'optimise-t-il pas un * a * a * a * a * a à (a * a * a) * (a * a * a)?


Answers

Lambdageek souligne correctement que parce que l'associativité ne tient pas pour les nombres à virgule flottante, l'optimisation d' a*a*a*a*a*a à (a*a*a)*(a*a*a) peut changer la valeur. C'est pourquoi il est interdit par C99 (sauf si spécifiquement autorisé par l'utilisateur, via un flag de compilation ou un pragma). En général, l'hypothèse est que le programmeur a écrit ce qu'elle a fait pour une raison, et le compilateur devrait respecter cela. Si vous voulez (a*a*a)*(a*a*a) , écrivez cela.

Cela peut être une douleur à écrire, cependant; pourquoi le compilateur ne peut-il pas faire ce que vous considérez être la bonne chose quand vous utilisez pow(a,6) ? Parce que ce serait la mauvaise chose à faire. Sur une plate-forme avec une bonne bibliothèque mathématique, pow(a,6) est significativement plus précis que a*a*a*a*a*a ou (a*a*a)*(a*a*a) . Juste pour fournir quelques données, j'ai couru une petite expérience sur mon Mac Pro, en mesurant la pire erreur en évaluant un ^ 6 pour tous les nombres flottants à simple précision entre [1,2]:

worst relative error using    powf(a, 6.f): 5.96e-08
worst relative error using (a*a*a)*(a*a*a): 2.94e-07
worst relative error using     a*a*a*a*a*a: 2.58e-07

Utiliser pow au lieu d'un arbre de multiplication réduit l'erreur liée à un facteur de 4 . Les compilateurs ne devraient pas (et ne font généralement pas) des "optimisations" qui augmentent l'erreur à moins d'être autorisés à le faire par l'utilisateur (par exemple via -ffast-math ).

Notez que GCC fournit __builtin_powi(x,n) comme alternative à pow( ) , ce qui devrait générer un arbre de multiplication en ligne. Utilisez-le si vous voulez perdre de la précision en termes de performances, mais ne voulez pas activer le calcul rapide.

Question

Je fais de l'optimisation numérique sur une application scientifique. Une chose que j'ai remarquée est que GCC optimisera le pow(a,2) appel pow(a,2) en le compilant en a*a , mais le pow(a,6) appel pow(a,6) n'est pas optimisé et appellera réellement la fonction pow la bibliothèque, ce qui ralentit énormément la performance. (En revanche, Intel C ++ Compiler , icc exécutable, va éliminer l'appel de la bibliothèque pour pow(a,6) .)

Ce que je suis curieux, c'est que quand j'ai remplacé pow(a,6) par a*a*a*a*a*a utilisant GCC 4.5.1 et les options " -O3 -lm -funroll-loops -msse4 ", il utilise 5 instructions mulsd :

movapd  %xmm14, %xmm13
mulsd   %xmm14, %xmm13
mulsd   %xmm14, %xmm13
mulsd   %xmm14, %xmm13
mulsd   %xmm14, %xmm13
mulsd   %xmm14, %xmm13

alors que si j'écris (a*a*a)*(a*a*a) , cela produira

movapd  %xmm14, %xmm13
mulsd   %xmm14, %xmm13
mulsd   %xmm14, %xmm13
mulsd   %xmm13, %xmm13

ce qui réduit le nombre d'instructions de multiplication à 3. icc a un comportement similaire.

Pourquoi les compilateurs ne reconnaissent-ils pas cette astuce d'optimisation?




Les fonctions de bibliothèque telles que "pow" sont généralement soigneusement conçues pour produire l'erreur minimale possible (dans le cas générique). On obtient généralement des fonctions d'approximation avec des splines (selon le commentaire de Pascal, l'implémentation la plus courante semble être l'utilisation de l' algorithme de Remez )

fondamentalement l'opération suivante:

pow(x,y);

a une erreur inhérente d'approximativement la même ampleur que l'erreur dans une multiplication ou division unique .

Alors que l'opération suivante:

float a=someValue;
float b=a*a*a*a*a*a;

a une erreur inhérente qui est supérieure à plus de 5 fois l'erreur d'une seule multiplication ou division (parce que vous combinez 5 multiplications).

Le compilateur devrait être très attentif au type d'optimisation qu'il est en train de faire:

  1. Si l'optimisation de pow(a,6) en a*a*a*a*a*a peut améliorer les performances, mais réduire drastiquement la précision des nombres à virgule flottante.
  2. si on optimise a*a*a*a*a*a pour pow(a,6) cela peut en fait réduire la précision car "a" était une valeur spéciale qui permet une multiplication sans erreur (une puissance de 2 ou un petit nombre entier)
  3. si l'optimisation de pow(a,6) à (a*a*a)*(a*a*a) ou (a*a)*(a*a)*(a*a) il peut y avoir une perte de précision par rapport à la fonction de pow .

En général, vous savez que pour des valeurs flottantes arbitraires, "pow" a une meilleure précision que n'importe quelle fonction que vous pourriez éventuellement écrire, mais dans certains cas particuliers, plusieurs multiplications peuvent avoir une meilleure précision et performance, c'est au développeur de choisir ce qui convient le mieux. finalement commentant le code afin que personne d'autre "optimise" ce code.

La seule chose qui a du sens (opinion personnelle, et apparemment un choix dans GCC sans optimisation ou compilateur particulier) à optimiser devrait remplacer "pow (a, 2)" par "a * a". Ce serait la seule chose sensée qu'un vendeur de compilateur devrait faire.




gcc peut réellement faire cette optimisation, même pour les nombres à virgule flottante. Par exemple,

double foo(double a) {
  return a*a*a*a*a*a;
}

devient

foo(double):
    mulsd   %xmm0, %xmm0
    movapd  %xmm0, %xmm1
    mulsd   %xmm0, %xmm1
    mulsd   %xmm1, %xmm0
    ret

avec -O -funsafe-math-optimizations . Cette réorganisation viole IEEE-754, cependant, il nécessite donc le drapeau.

Les entiers signés, comme Peter Cordes l'a souligné dans un commentaire, peuvent effectuer cette optimisation sans optimisations -funsafe-math-optimizations puisqu'ils sont exactement là où il n'y a pas de débordement et s'il y a débordement, vous obtenez un comportement indéfini. Donc vous obtenez

foo(long):
    movq    %rdi, %rax
    imulq   %rdi, %rax
    imulq   %rdi, %rax
    imulq   %rax, %rax
    ret

avec juste -O . Pour les entiers non-signés, c'est encore plus simple car ils fonctionnent avec des puissances de 2 et peuvent donc être réordonnés librement même en cas de débordement.




Aucune affiche n'a encore mentionné la contraction des expressions flottantes (norme ISO C, 6.5p8 et 7.12.2). Si le pragma FP_CONTRACT est défini sur ON , le compilateur peut considérer une expression comme a*a*a*a*a*a comme une opération unique, comme si elle était évaluée exactement avec un seul arrondi. Par exemple, un compilateur peut le remplacer par une fonction de puissance interne à la fois plus rapide et plus précise. Ceci est particulièrement intéressant car le comportement est partiellement contrôlé par le programmeur directement dans le code source, alors que les options du compilateur fournies par l'utilisateur final peuvent parfois être utilisées incorrectement.

L'état par défaut du pragma FP_CONTRACT est défini par l'implémentation, de sorte qu'un compilateur est autorisé à effectuer de telles optimisations par défaut. Ainsi, un code portable qui doit suivre strictement les règles IEEE 754 devrait explicitement le mettre sur OFF .

Si un compilateur ne supporte pas ce pragma, il doit être prudent en évitant une telle optimisation, dans le cas où le développeur a choisi de le mettre sur OFF .

GCC ne supporte pas ce pragma, mais avec les options par défaut, il suppose qu'il est sur ON ; donc pour les cibles avec un FMA matériel, si l'on veut empêcher la transformation a*b+c en fma (a, b, c), il faut prévoir une option telle que -ffp-contract=off (pour paramétrer explicitement le pragma à OFF ) ou -std=c99 (pour indiquer à GCC de se conformer à une version standard de C, ici C99, donc suivez le paragraphe ci-dessus). Dans le passé, cette dernière option n'empêchait pas la transformation, ce qui signifiait que GCC ne se conformait pas sur ce point: https://gcc.gnu.org/bugzilla/show_bug.cgi?id=37845




Fortran (conçu pour le calcul scientifique) a un opérateur de puissance intégré, et pour autant que je sache, les compilateurs Fortran optimiseront généralement l'augmentation des puissances entières d'une manière similaire à ce que vous décrivez. C / C ++ n'a malheureusement pas d'opérateur de puissance, seulement la fonction de bibliothèque pow() . Cela n'empêche pas les compilateurs intelligents de traiter pow spécialement et de le calculer plus rapidement pour des cas particuliers, mais il semble qu'ils le font moins souvent ...

Il y a quelques années, j'essayais de rendre plus pratique le calcul optimal des puissances entières et j'ai trouvé ce qui suit. C'est C ++, et non C, et cela dépend toujours du fait que le compilateur soit un peu plus intelligent sur la façon d'optimiser / inline les choses. Quoi qu'il en soit, j'espère que vous trouverez cela utile en pratique:

template<unsigned N> struct power_impl;

template<unsigned N> struct power_impl {
    template<typename T>
    static T calc(const T &x) {
        if (N%2 == 0)
            return power_impl<N/2>::calc(x*x);
        else if (N%3 == 0)
            return power_impl<N/3>::calc(x*x*x);
        return power_impl<N-1>::calc(x)*x;
    }
};

template<> struct power_impl<0> {
    template<typename T>
    static T calc(const T &) { return 1; }
};

template<unsigned N, typename T>
inline T power(const T &x) {
    return power_impl<N>::calc(x);
}

Clarification pour les curieux: ceci ne trouve pas la façon optimale de calculer les puissances, mais trouver la solution optimale est un problème NP-complet et cela ne vaut pas mieux pour les petits pouvoirs (contrairement à pow ), il n'y a pas de raison de agiter avec le détail.

Ensuite, utilisez-le simplement comme power<6>(a) .

Cela permet de taper facilement des puissances (pas besoin d'épeler 6 pars avec parens), et vous permet d'avoir ce type d'optimisation sans -ffast-math au cas où vous auriez quelque chose de dépendant de la précision comme la sommation compensée (un exemple où la commande des opérations est essentiel).

Vous pouvez probablement aussi oublier que c'est du C ++ et juste l'utiliser dans le programme C (s'il compile avec un compilateur C ++).

J'espère que cela peut être utile.

MODIFIER:

C'est ce que je reçois de mon compilateur:

Pour a*a*a*a*a*a ,

    movapd  %xmm1, %xmm0
    mulsd   %xmm1, %xmm0
    mulsd   %xmm1, %xmm0
    mulsd   %xmm1, %xmm0
    mulsd   %xmm1, %xmm0
    mulsd   %xmm1, %xmm0

Pour (a*a*a)*(a*a*a) ,

    movapd  %xmm1, %xmm0
    mulsd   %xmm1, %xmm0
    mulsd   %xmm1, %xmm0
    mulsd   %xmm0, %xmm0

Pour la power<6>(a) ,

    mulsd   %xmm0, %xmm0
    movapd  %xmm0, %xmm1
    mulsd   %xmm0, %xmm1
    mulsd   %xmm0, %xmm1



Parce qu'un nombre à virgule flottante de 32 bits - tel que 1.024 - n'est pas 1.024. Dans un ordinateur, 1.024 est un intervalle: de (1.024-e) à (1.024 + e), où "e" représente une erreur. Certaines personnes ne réalisent pas cela et croient également que * dans un * a représente la multiplication de nombres de précision arbitraires sans qu'il y ait d'erreurs attachées à ces nombres. La raison pour laquelle certaines personnes ne réalisent pas cela est peut-être les calculs mathématiques qu'ils ont pratiqués dans les écoles primaires: travailler uniquement avec des nombres idéaux sans erreurs, et croire qu'il est acceptable de simplement ignorer "e" tout en effectuant la multiplication. Ils ne voient pas le "e" implicite dans "float a = 1.2", "a * a * a" et les codes C similaires.

Si la majorité des programmeurs reconnaissent (et peuvent exécuter) l'idée que C expression a * a * a * a * a * a ne fonctionne pas réellement avec des nombres idéaux, le compilateur GCC serait alors LIBRE pour optimiser "a * a * a * a * a * a "dire" t = (a * a); t * t * t "qui nécessite un plus petit nombre de multiplications. Mais malheureusement, le compilateur GCC ne sait pas si le programmeur qui écrit le code pense que "a" est un nombre avec ou sans erreur. Et GCC ne fera que ce à quoi ressemble le code source - parce que c'est ce que GCC voit à l'œil nu.

... une fois que vous savez quel genre de programmeur vous êtes, vous pouvez utiliser le commutateur "-fast-math" pour dire à GCC que "Hey, GCC, je sais ce que je fais!". Cela permettra à GCC de convertir un * a * a * a * a * a en un texte différent - il a l'air différent d'un * a * a * a * a * a - mais calcule toujours un nombre dans l'intervalle d'erreur de a * a * a * a * a * a. C'est OK, puisque vous savez déjà que vous travaillez avec des intervalles, pas des nombres idéaux.




Links