Code C ++ pour tester la conjecture de Collatz plus rapidement que l’assemblage écrit à la main – pourquoi?

J’ai écrit ces deux solutions pour Project Euler Q14 , en assemblage et en C ++. Ils sont la même approche de force brute identique pour tester la conjecture de Collatz . La solution d’assemblage a été assemblée avec

nasm -felf64 p14.asm && gcc p14.o -o p14 

Le C ++ a été compilé avec

 g++ p14.cpp -o p14 

Assembly, p14.asm

 section .data fmt db "%d", 10, 0 global main extern printf section .text main: mov rcx, 1000000 xor rdi, rdi ; max i xor rsi, rsi ; i l1: dec rcx xor r10, r10 ; count mov rax, rcx l2: test rax, 1 jpe even mov rbx, 3 mul rbx inc rax jmp c1 even: mov rbx, 2 xor rdx, rdx div rbx c1: inc r10 cmp rax, 1 jne l2 cmp rdi, r10 cmovl rdi, r10 cmovl rsi, rcx cmp rcx, 2 jne l1 mov rdi, fmt xor rax, rax call printf ret 

C ++, p14.cpp

 #include  using namespace std; int sequence(long n) { int count = 1; while (n != 1) { if (n % 2 == 0) n /= 2; else n = n*3 + 1; ++count; } return count; } int main() { int max = 0, maxi; for (int i = 999999; i > 0; --i) { int s = sequence(i); if (s > max) { max = s; maxi = i; } } cout << maxi << endl; } 

Je connais les optimisations du compilateur pour améliorer la vitesse et tout, mais je ne vois pas beaucoup de façons d’optimiser davantage ma solution d’assemblage (ne pas parler de programme de manière mathématique).

Le code C ++ a un module à chaque terme et division tous les termes pairs, où l’assemblage n’est qu’une division par terme pair.

Mais l’assemblage prend en moyenne 1 seconde de plus que la solution C ++. Pourquoi est-ce? Je demande surtout de la curiosité.

Temps d’exécution

Mon système: Linux 64 bits sur Intel Celeron 2955U à 1,4 GHz (microarchitecture Haswell).

  • g++ (non optimisé): av 1272 ms

  • g++ -O3 avg 578 ms

  • original asm (div) avg 2650 ms

  • Asm (shr) avg 679 ms

  • @johnfound asm , assemblé avec nasm avg 501 ms

  • @hidefromkgb asm avg 200 ms

  • @hidefromkgb asm optimisé par @Peter Cordes avg 145 ms

  • @Veedrac C ++ moyenne 81 ms avec -O3 , 305 ms avec -O0

Si vous pensez qu’une instruction DIV de 64 bits est un bon moyen de diviser par deux, il n’est pas surprenant que la sortie asm du compilateur ait battu votre code écrit à la main, même avec -O0 (compilation rapide, aucune optimisation supplémentaire et stockage / rechargement en mémoire) après / avant chaque instruction C, un débogueur peut modifier les variables).

Reportez-vous au guide de l’assemblage d’optimisation d’Agner Fog pour apprendre à écrire des asm efficaces. Il a également des tables d’instructions et un guide de microarchie pour des détails spécifiques sur des processeurs spécifiques. Voir aussi le wiki du tag x86 pour plus de liens de perf.

Voir aussi cette question plus générale sur le fait de battre le compilateur avec asm écrit à la main: le langage d’assemblage en ligne est-il plus lent que le code C ++ natif? . TL: DR: oui si vous le faites mal (comme cette question).

Habituellement, vous pouvez très bien laisser le compilateur faire son travail, surtout si vous essayez d’écrire du C ++ capable de comstackr efficacement . Voir aussi l’ assemblage plus rapide que les langages compilés? . Une des réponses à ces diapositives ordonnées montre comment différents compilateurs C optimisent certaines fonctions très simples avec des astuces géniales.


 even: mov rbx, 2 xor rdx, rdx div rbx 

Sur Intel Haswell, div r64 est de 36 uops, avec une latence de 32-96 cycles et un débit de 1 à 21-74 cycles. (Plus les 2 uops pour configurer RBX et zéro RDX, mais une exécution dans le désordre peut les exécuter en avance). Les instructions à haut débit comme DIV sont microcodées, ce qui peut également provoquer des goulots d’étranglement. Dans ce cas, la latence est le facteur le plus pertinent car elle fait partie d’une chaîne de dépendances en boucle.

shr rax, 1 fait la même division non signée: c’est 1 uop, avec une latence de 1c , et peut courir 2 par cycle d’horloge.

A titre de comparaison, la division 32 bits est plus rapide, mais rest horrible par rapport aux décalages. idiv r32 est de 9 uops, 22-29c de latence et un par 8-11c sur Haswell.


Comme vous pouvez le constater en regardant la sortie asm de GCC ( explorateur du compilateur Godbolt ), elle n’utilise que des instructions Shifts . clang -O0 comstack naïvement comme vous le pensiez, même en utilisant deux fois IDIV 64 bits. (Lors de l’optimisation, les compilateurs utilisent les deux sorties d’IDIV lorsque la source effectue une division et un module avec les mêmes opérandes, s’ils utilisent IDIV du tout)

GCC n’a pas de mode totalement naïf; il transforme toujours par GIMPLE, ce qui signifie que certaines “optimisations” ne peuvent pas être désactivées . Cela inclut la reconnaissance de la division par constante et l’utilisation de décalages (puissance de 2) ou d’ un inverse multiplicatif à virgule fixe (non de puissance 2) pour éviter les IDIV (voir div_by_13 dans le lien ci-dessus).

gcc -Os (optimiser pour la taille) utilise IDIV pour la division non-power-of-2, malheureusement même dans les cas où le code inverse multiplicatif est seulement légèrement plus grand mais beaucoup plus rapide.


Aider le compilateur

(résumé pour ce cas: utilisez uint64_t n )

Tout d’abord, il n’est intéressant que d’examiner la sortie optimisée du compilateur. ( -O3 ). -O0 vitesse est fondamentalement sans signification.

Regardez votre sortie asm (sur Godbolt, ou voir comment supprimer “le bruit” de la sortie d’assemblage GCC / clang? ). Lorsque le compilateur ne fait pas du code optimal en premier lieu: écrire votre source C / C ++ d’une manière qui guide le compilateur à créer un meilleur code est généralement la meilleure approche . Vous devez savoir, et savoir ce qui est efficace, mais vous appliquez ces connaissances indirectement. Les compilateurs sont également une bonne source d’idées: parfois, clang fait quelque chose de cool, et vous pouvez tenir gcc à faire la même chose: voyez cette réponse et ce que j’ai fait avec la boucle non déroulée dans le code de @Veedrac ci-dessous.)

Cette approche est portable, et dans 20 ans, un futur compilateur pourra le comstackr de manière efficace sur du matériel futur (x86 ou non), en utilisant peut-être une nouvelle extension ISA ou une vectorisation automatique. Les manuscrits x86-64 manuscrits d’il y a 15 ans ne seraient généralement pas optimisés pour Skylake. Par exemple, la macro-fusion compare & branch n’existait pas à l’époque. Ce qui est optimal pour les asm fabriqués manuellement pour une microarchitecture peut ne pas être optimal pour les autres processeurs actuels et futurs. Les commentaires sur la réponse de @ johnfound discutent des différences majeures entre AMD Bulldozer et Intel Haswell, qui ont un impact important sur ce code. Mais en théorie, g++ -O3 -march=bdver3 et g++ -O3 -march=skylake fera le bon choix. (Ou -march=native .) Ou -mtune=... juste pour régler, sans utiliser les instructions que d’autres processeurs pourraient ne pas prendre en charge.

Mon sentiment est que guider le compilateur pour qu’il fonctionne bien pour un processeur actuel dont vous vous souciez ne devrait pas être un problème pour les futurs compilateurs. Ils espèrent mieux que les compilateurs actuels trouver des moyens de transformer le code et trouver un moyen de fonctionner pour les futurs CPU. Quoi qu’il en soit, les futurs x86 ne seront probablement pas terribles pour tout ce qui est bon sur le x86 actuel, et le futur compilateur évitera tout piège spécifique à l’asm tout en implémentant quelque chose comme le transfert de données de votre source C, s’il ne voit pas quelque chose de mieux.

Asm écrit à la main est une boîte noire pour l’optimiseur, de sorte que la propagation constante ne fonctionne pas lorsque l’inline fait d’une entrée une constante de compilation. D’autres optimisations sont également affectées. Lisez https://gcc.gnu.org/wiki/DontUseInlineAsm avant d’utiliser asm. (Et éviter les asm en ligne de style MSVC: les entrées / sorties doivent traverser la mémoire, ce qui ajoute une surcharge .)

Dans ce cas : votre n a un type signé et gcc utilise la séquence SAR / SHR / ADD qui donne l’arrondi correct. (IDIV et décalage arithmétique “round” différemment pour les entrées négatives, voir la saisie manuelle de SAR insn set ). (IDK si gcc a essayé et a échoué à prouver que n ne peut pas être négatif, ou quoi. Le débordement signé est un comportement indéfini, donc il aurait dû être capable de le faire.)

Vous devriez avoir utilisé uint64_t n , il ne peut donc que SHR. Et donc, il est portable sur les systèmes où la long est seulement de 32 bits (par exemple, x86-64 Windows).


BTW, la sortie asm optimisée de gcc semble plutôt bonne (en utilisant unsigned long n ) : la boucle interne qu’elle insère dans main() fait ceci:

  # from gcc5.4 -O3 plus my comments # edx= count=1 # rax= uint64_t n .L9: # do{ lea rcx, [rax+1+rax*2] # rcx = 3*n + 1 mov rdi, rax shr rdi # rdi = n>>1; test al, 1 # set flags based on n%2 (aka n&1) mov rax, rcx cmove rax, rdi # n= (n%2) ? 3*n+1 : n/2; add edx, 1 # ++count; cmp rax, 1 jne .L9 #}while(n!=1) cmp/branch to update max and maxi, and then do the next n 

La boucle interne est sans twig et le chemin critique de la chaîne de dépendances en boucle est le suivant:

  • LEA à 3 composants (3 cycles)
  • cmov (2 cycles sur Haswell, 1c sur Broadwell ou plus tard).

Total: 5 cycles par itération, goulot d’étranglement de latence . L’exécution en désordre prend en charge tout le rest en parallèle (en théorie: je n’ai pas testé avec les compteurs de perf pour voir si ça fonctionne vraiment à 5c / iter).

L’entrée FLAGS de cmov (produite par TEST) est plus rapide à produire que l’entrée RAX (à partir de LEA-> MOV), elle n’est donc pas sur le chemin critique.

De même, le MOV-> SHR qui produit l’entrée RDI de CMOV est hors du chemin critique, car il est également plus rapide que le LEA. MOV sur IvyBridge et les versions ultérieures ne présentent aucune latence (gérée au moment du changement de registre). (Cela prend toujours un uop, et un slot dans le pipeline, donc ce n’est pas gratuit, juste une latence nulle). Le MOV supplémentaire dans la chaîne dep de LEA fait partie du goulot d’étranglement sur les autres processeurs.

Le fichier cmp / jne ne fait pas non plus partie du chemin critique: il ne porte pas de boucle, car les dépendances de contrôle sont traitées avec une prédiction de twig + une exécution spéculative, contrairement aux dépendances de données sur le chemin critique.


Battre le compilateur

GCC a fait du bon travail ici. Il pourrait enregistrer un octet de code en utilisant inc edx au lieu de add edx, 1 , car personne ne se soucie de P4 et de ses fausses dépendances pour les instructions de modification de drapeau partiel.

Il pourrait également enregistrer toutes les instructions MOV, et le TEST: SHR définit CF = le bit décalé, nous pouvons donc utiliser cmovc au lieu de test / cmovz .

  ### Hand-optimized version of what gcc does .L9: #do{ lea rcx, [rax+1+rax*2] # rcx = 3*n + 1 shr rax, 1 # n>>=1; CF = n&1 = n%2 cmovc rax, rcx # n= (n&1) ? 3*n+1 : n/2; inc edx # ++count; cmp rax, 1 jne .L9 #}while(n!=1) 

Voir la réponse de @ johnfound pour un autre truc astucieux: supprimez le CMP en effectuant un twigment sur le résultat du flag SHR et en l’utilisant pour CMOV: zéro uniquement si n était 1 (ou 0) pour commencer. (Fait amusant: SHR avec nombre! = 1 sur Nehalem ou antérieur provoque un blocage si vous lisez les résultats de l’indicateur . C’est comme ça qu’ils l’ont rendu unique. L’encodage spécial shift-by-1 est correct, cependant.)

Eviter MOV n’aide pas du tout à la latence sur Haswell ( Le MOV de x86 peut-il vraiment être “libre”? Pourquoi ne puis-je pas le reproduire du tout? ). Cela aide de manière significative sur les processeurs comme Intel pré-IvB et AMD Bulldozer-family, où le MOV n’est pas une latence nulle. Les instructions MOV gaspillées par le compilateur affectent le chemin critique. Les complexes BD-LEA et CMOV sont tous deux à faible latence (respectivement 2c et 1c), ce qui en fait une fraction plus importante de la latence. En outre, les goulots d’étranglement de débit deviennent un problème, car il ne comporte que deux canaux ALU entiers. Voir la réponse de @ johnfound , où il a des résultats de timing d’un processeur AMD.

Même sur Haswell, cette version peut aider un peu en évitant certains retards occasionnels où un uop non critique vole un port d’exécution de l’un sur le chemin critique, retardant ainsi l’exécution d’un cycle. (Ceci s’appelle un conflit de ressources). Il enregistre également un registre, ce qui peut aider lors de la réalisation de plusieurs valeurs n en parallèle dans une boucle entrelacée (voir ci-dessous).

La latence de LEA dépend du mode d’adressage , sur les processeurs de la famille Intel SnB. 3c pour 3 composants ( [base+idx+const] , qui prend deux ajouts distincts), mais seulement 1c avec 2 composants ou moins (un add). Certains processeurs (comme Core2) effectuent même un LEA à 3 composants en un seul cycle, mais pas SnB-family. Pire encore, la famille Intel SnB standardise les latences de sorte qu’il n’y a pas de raccourcis 2c , sinon le LEA à 3 composants ne serait que 2c comme Bulldozer. (Le LEA à 3 composants est également plus lent sur AMD, mais pas autant).

Donc lea rcx, [rax + rax*2] / inc rcx est seulement une latence de 2 c, plus rapide que lea rcx, [rax + rax*2 + 1] , sur les processeurs de la famille Intel SnB comme Haswell. Pause sur BD, et pire sur Core2. Il en coûte un uop supplémentaire, ce qui ne vaut normalement pas la peine de sauver la latence 1c, mais la latence est le principal goulot d’étranglement et Haswell dispose d’un pipeline suffisamment large pour gérer le débit supplémentaire.

Ni gcc, ni icc, ni clang (sur godbolt) n’utilisaient la sortie CF de SHR, toujours en utilisant un ET ou un TEST . Compilateurs stupides. : P Ce sont des machines complexes, mais un être humain intelligent peut souvent les battre avec des problèmes à petite échelle. (Bien entendu, des milliers de fois plus de temps pour y penser!) Les compilateurs n’utilisent pas d’algorithmes exhaustifs pour rechercher toutes les manières possibles de faire les choses, car cela prendrait trop de temps lors de l’optimisation ils font de leur mieux. Ils ne modélisent pas non plus le pipeline dans la microarchitecture cible, du moins pas de la même manière que l’ IACA ou d’autres outils d’parsing statique; ils utilisent simplement des heuristiques.)


Le déroulement simple de la boucle ne vous aidera pas ; cette boucle gêne la latence d’une chaîne de dépendances transscope en boucle, et non pas le temps / débit de la boucle. Cela signifie qu’il ferait bien avec l’hyperthreading (ou tout autre type de SMT), puisque le CPU a beaucoup de temps pour intercaler des instructions à partir de deux threads. Cela impliquerait de paralléliser la boucle dans main , mais c’est correct parce que chaque thread peut juste vérifier une plage de n valeurs et produire une paire d’entiers en conséquence.

L’entrelacement à la main dans un seul fil peut aussi être viable . Peut-être calculer la séquence pour une paire de nombres en parallèle, puisque chacun ne prend que quelques registres, et ils peuvent tous mettre à jour le même max / maxi . Cela crée plus de parallélisme au niveau de l’instruction .

L’astuce consiste à décider d’attendre que toutes les n valeurs aient atteint 1 avant d’obtenir une autre paire de valeurs n de départ, ou de sortir et d’obtenir un nouveau sharepoint départ pour une seule qui a atteint la condition finale, sans toucher les registres du autre séquence. Il est probablement préférable de faire en sorte que chaque chaîne travaille sur des données utiles, sinon vous devrez incrémenter conditionnellement son compteur.


Vous pourriez peut-être même faire cela avec des éléments SSE compactés pour incrémenter conditionnellement le compteur pour les éléments vectoriels où n n’a pas encore atteint 1 . Et puis, pour masquer la latence encore plus longue d’une implémentation à incrément conditionnel SIMD, vous devez garder plus de vecteurs de valeurs n dans l’air. Peut-être que la valeur avec un vecteur 256b (4x uint64_t ).

Je pense que la meilleure stratégie pour détecter un 1 “collant” consiste à masquer le vecteur de tout-en-un que vous ajoutez pour incrémenter le compteur. Donc, après avoir vu un 1 dans un élément, le vecteur-incrément aura un zéro et + = 0 est un no-op.

Idée non testée pour la vectorisation manuelle

 # starting with YMM0 = [ n_d, n_c, n_b, n_a ] (64-bit elements) # ymm4 = _mm256_set1_epi64x(1): increment vector # ymm5 = all-zeros: count vector .inner_loop: vpaddq ymm1, ymm0, xmm0 vpaddq ymm1, ymm1, xmm0 vpaddq ymm1, ymm1, set1_epi64(1) # ymm1= 3*n + 1. Maybe could do this more efficiently? vprllq ymm3, ymm0, 63 # shift bit 1 to the sign bit vpsrlq ymm0, ymm0, 1 # n /= 2 # There may be a better way to do this blend, avoiding the bypass delay for an FP blend between integer insns, not sure. Probably worth it vpblendvpd ymm0, ymm0, ymm1, ymm3 # variable blend controlled by the sign bit of each 64-bit element. I might have the source operands backwards, I always have to look this up. # ymm0 = updated n in each element. vpcmpeqq ymm1, ymm0, set1_epi64(1) vpandn ymm4, ymm1, ymm4 # zero out elements of ymm4 where the compare was true vpaddq ymm5, ymm5, ymm4 # count++ in elements where n has never been == 1 vptest ymm4, ymm4 jnz .inner_loop # Fall through when all the n values have reached 1 at some point, and our increment vector is all-zero vextracti128 ymm0, ymm5, 1 vpmaxq .... crap this doesn't exist # Actually just delay doing a horizontal max until the very very end. But you need some way to record max and maxi. 

Vous pouvez et devez implémenter cela avec des éléments insortingnsèques, au lieu de l’Asm écrit à la main.


Amélioration algorithmique / implémentation:

En plus de mettre en œuvre la même logique avec un système d’exploitation plus efficace, recherchez des moyens de simplifier la logique ou d’éviter les tâches inutiles. p.ex. mémoriser pour détecter les terminaisons communes aux séquences. Ou mieux encore, regardez 8 bits à la fois (réponse de gnasher)

@EOF souligne que tzcnt (ou bsf ) pourrait être utilisé pour faire plusieurs itérations n/=2 en une seule étape. C’est probablement mieux que la vectorisation SIMD, car aucune instruction SSE ou AVX ne peut le faire. Il est néanmoins compatible avec la réalisation de plusieurs scalaires en parallèle dans différents registres entiers.

Ainsi, la boucle pourrait ressembler à ceci:

 goto loop_entry; // C++ structured like the asm, for illustration only do { n = n*3 + 1; loop_entry: shift = _tzcnt_u64(n); n >>= shift; count += shift; } while(n != 1); 

Cela peut faire beaucoup moins d’itérations, mais les décalages à comptage variable sont lents sur les processeurs Intel SnB-family sans BMI2. 3 uops, 2c de latence. (Ils ont une dépendance d’entrée sur les FLAGS car count = 0 signifie que les flags ne sont pas modifiés. Ils gèrent cela comme une dépendance de données et prennent plusieurs UOP car un UOP ne peut avoir que 2 entrées (pre-HSW / BDW)). C’est le genre auquel se réfèrent les gens qui se plaignent de la conception délirante de CISC de x86. Cela rend les processeurs x86 plus lents qu’ils ne l’auraient été si l’ISA était conçue à partir de zéro, même de la même manière. (c.-à-d. cela fait partie de la “taxe x86” qui coûte de la vitesse / puissance.) SHRX / SHLX / SARX (BMI2) sont une grande victoire (latence 1 uop / 1c).

Il place également tzcnt (3c sur Haswell et versions ultérieures) sur le chemin critique, ce qui allonge considérablement la latence totale de la chaîne de dépendances en boucle. Cela supprime tout besoin de CMOV, ou prépare un registre contenant n>>1 , cependant. La réponse de @Veedrac surmonte tout cela en reportant le tzcnt / shift pour des itérations multiples, ce qui est très efficace (voir ci-dessous).

Nous pouvons utiliser en toute sécurité BSF ou TZCNT de manière interchangeable, car n ne peut jamais être nul à ce point. Le code machine de TZCNT décode en tant que BSF sur les CPU qui ne supportent pas BMI1. (Les préfixes sans signification sont ignorés, de sorte que REP BSF s’exécute en tant que BSF).

TZCNT fonctionne beaucoup mieux que BSF sur les processeurs AMD qui le supportent, il peut donc être utile d’utiliser REP BSF , même si vous ne vous souciez pas de définir ZF si l’entrée est nulle plutôt que la sortie. Certains compilateurs le font quand vous utilisez __builtin_ctzll même avec -mno-bmi .

Ils fonctionnent de la même manière sur les processeurs Intel, alors sauvegardez simplement l’octet si c’est tout ce qui compte. TZCNT sur Intel (pré-Skylake) a toujours une fausse dépendance sur l’opérande de sortie supposé être en écriture seule, tout comme BSF, pour supporter le comportement non documenté que BSF avec input = 0 laisse sa destination non modifiée. Vous devez donc contourner ce problème à moins d’optimiser uniquement pour Skylake, donc il n’y a rien à gagner avec l’octet REP supplémentaire. (Intel va souvent au-delà de ce qu’exige le manuel ISA x86, pour éviter de casser du code largement utilisé qui dépend de quelque chose qu’il ne devrait pas ou qui est interdit rétroactivement. Par exemple, Windows 9x ne présuppose aucune lecture anticipée des entrées TLB . lorsque le code a été écrit, avant que Intel ne mette à jour les règles de gestion du TLB .)

Quoi qu’il en soit, LZCNT / TZCNT sur Haswell ont le même faux dep que POPCNT: voir cette Q & R. C’est pourquoi dans la sortie asm de gcc pour le code de @ Veedrac, vous voyez qu’il brise la chaîne dep avec xor-zeroing sur le registre qu’il va utiliser comme destination de TZCNT, quand il n’utilise pas dst = src. Puisque TZCNT / LZCNT / POPCNT ne laisse jamais sa destination indéfinie ou non modifiée, cette fausse dépendance sur la sortie des processeurs Intel est purement une limitation / un bogue de performance. Vraisemblablement, il vaut la peine que certains transistors / puissances les fassent se comporter comme d’autres uops qui vont à la même unité d’exécution. Le seul avantage logiciel visible est l’interaction avec une autre limitation microarchitecturale: ils peuvent micro-fusionner un opérande de mémoire avec un mode d’adressage indexé sur Haswell, mais sur Skylake où Intel a supprimé la fausse dépendance pour LZCNT / TZCNT, ils “se décollent”. modes d’adressage indexés tandis que POPCNT peut toujours fusionner n’importe quel mode addr.


Améliorations des idées / codes provenant d’autres réponses:

La réponse de @ hidefromkgb a une belle observation: vous êtes assuré de pouvoir effectuer un changement de position correct après un 3n + 1. Vous pouvez calculer cela plus efficacement encore que de simplement omettre les contrôles entre les étapes. L’implémentation asm dans cette réponse est cassée, mais (dépend de OF, qui n’est pas défini après SHRD avec un compte> 1), et lente: ROR rdi,2 est plus rapide que SHRD rdi,rdi,2 et utilise deux instructions CMOV sur le chemin critique est plus lent qu’un test supplémentaire qui peut être exécuté en parallèle.

J’ai mis en ordre / amélioré C (qui guide le compilateur pour produire de meilleures données), et testé + travail plus rapide asm (dans les commentaires en dessous du C) sur Godbolt: voir le lien dans la réponse de @ hidefromkgb . (Cette réponse a atteint la limite des 30 carats à partir des grandes URL Godbolt, mais les raccourcis peuvent pourrir et sont trop longs pour goo.gl de toute façon.)

A également amélioré l’impression en sortie pour convertir en une chaîne et en faire une write() au lieu d’écrire un caractère à la fois. Cela minimise l’impact sur le chronométrage du programme entier avec perf stat ./collatz (pour enregistrer les compteurs de performance), et j’ai désembrouillé certaines des asm non critiques.


@ Le code de Veedrac

J’ai eu un très petit accélérateur de changement de vitesse autant que nous soaps besoin de faire, et vérifier pour continuer la boucle. De 7.5s pour la limite = 1e8 à 7.275s, sur Core2Duo (Merom), avec un facteur de déroulement de 16.

code + commentaires sur Godbolt . N’utilisez pas cette version avec clang; il fait quelque chose de stupide avec la boucle de report. Utiliser un compteur tmp k puis l’append pour count plus tard change ce que fait clang, mais cela blesse légèrement gcc.

Voir la discussion dans les commentaires: le code de Veedrac est excellent sur les CPU avec BMI1 (c.-à-d. Pas Celeron / Pentium)

Affirmer que le compilateur C ++ peut produire un code plus optimal qu’un programmeur compétent en langage assembleur est une très mauvaise erreur. Et surtout dans ce cas. L’humain peut toujours rendre le code meilleur que le compilateur, et cette situation particulière illustre bien cette affirmation.

La différence de temps que vous voyez est que le code d’assemblage dans la question est très loin d’être optimal dans les boucles internes.

(Le code ci-dessous est 32 bits, mais peut être facilement converti en 64 bits)

Par exemple, la fonction de séquence ne peut être optimisée que pour 5 instructions:

  .seq: inc esi ; counter lea edx, [3*eax+1] ; edx = 3*n+1 shr eax, 1 ; eax = n/2 cmovc eax, edx ; if CF eax = edx jnz .seq ; jmp if n<>1 

Le code entier ressemble à:

 include "%lib%/freshlib.inc" @BinaryType console, compact options.DebugMode = 1 include "%lib%/freshlib.asm" start: InitializeAll mov ecx, 999999 xor edi, edi ; max xor ebx, ebx ; max i .main_loop: xor esi, esi mov eax, ecx .seq: inc esi ; counter lea edx, [3*eax+1] ; edx = 3*n+1 shr eax, 1 ; eax = n/2 cmovc eax, edx ; if CF eax = edx jnz .seq ; jmp if n<>1 cmp edi, esi cmovb edi, esi cmovb ebx, ecx dec ecx jnz .main_loop OutputValue "Max sequence: ", edi, 10, -1 OutputValue "Max index: ", ebx, 10, -1 FinalizeAll stdcall TerminateAll, 0 

Pour comstackr ce code, FreshLib est nécessaire.

Dans mes tests (processeur AMD A4-1200 à 1 GHz), le code ci-dessus est environ quatre fois plus rapide que le code C ++ de la question (lorsqu’il est compilé avec -O0 : 430 ms vs 1900 ms) et plus de deux fois plus rapide (430 ms vs 830 ms) lorsque le code C ++ est compilé avec -O3 .

La sortie des deux programmes est la même: séquence max = 525 sur i = 837799.

Pour plus de performance: une simple modification consiste à observer qu’après n = 3n + 1, n sera pair, vous pouvez donc diviser par 2 immédiatement. Et n ne sera pas 1, vous n’avez donc pas besoin de le tester. Donc, vous pourriez enregistrer quelques déclarations si et écrire:

 while (n % 2 == 0) n /= 2; if (n > 1) for (;;) { n = (3*n + 1) / 2; if (n % 2 == 0) { do n /= 2; while (n % 2 == 0); if (n == 1) break; } } 

Voici un gros gain: Si vous regardez les 8 bits les plus bas de n, toutes les étapes jusqu’à ce que vous les divisiez par 2 huit fois sont complètement déterminées par ces huit bits. Par exemple, si les huit derniers bits sont 0x01, c’est en binary que votre numéro est ???? 0000 0001 alors les prochaines étapes sont:

 3n+1 -> ???? 0000 0100 / 2 -> ???? ?000 0010 / 2 -> ???? ??00 0001 3n+1 -> ???? ??00 0100 / 2 -> ???? ???0 0010 / 2 -> ???? ???? 0001 3n+1 -> ???? ???? 0100 / 2 -> ???? ???? ?010 / 2 -> ???? ???? ??01 3n+1 -> ???? ???? ??00 / 2 -> ???? ???? ???0 / 2 -> ???? ???? ???? 

Donc, toutes ces étapes peuvent être prédites, et 256k + 1 est remplacé par 81k + 1. Quelque chose de similaire se produira pour toutes les combinaisons. Donc, vous pouvez faire une boucle avec une grande déclaration de changement:

 k = n / 256; m = n % 256; switch (m) { case 0: n = 1 * k + 0; break; case 1: n = 81 * k + 1; break; case 2: n = 81 * k + 1; break; ... case 155: n = 729 * k + 425; break; ... } 

Lancez la boucle jusqu’à ce que n ≤ 128, car à ce point, n pourrait devenir 1 avec moins de huit divisions par 2, et faire huit pas ou plus à la fois vous ferait rater le point où vous atteignez 1 pour la première fois. Continuez ensuite la boucle “normale” ou préparez un tableau indiquant le nombre d’étapes supplémentaires à atteindre 1.

PS Je soupçonne fortement que la suggestion de Peter Cordes le rendrait encore plus rapide. Il n’y aura aucune twig conditionnelle, sauf une, et celle-ci sera prédite correctement, sauf lorsque la boucle se termine. Donc, le code serait quelque chose comme

 static const unsigned int multipliers [256] = { ... } static const unsigned int adders [256] = { ... } while (n > 128) { size_t lastBits = n % 256; n = (n >> 8) * multipliers [lastBits] + adders [lastBits]; } 

En pratique, vous mesureriez si le traitement des 9, 10, 11, 12 derniers bits de n à la fois serait plus rapide. Pour chaque bit, le nombre d’entrées dans la table doublerait et j’excelle un ralentissement lorsque les tables ne rentrent plus dans le cache L1.

PPS. Si vous avez besoin du nombre d’opérations: Dans chaque itération, nous faisons exactement huit divisions par deux et un nombre variable d’opérations (3n + 1). Une méthode évidente pour compter les opérations serait donc un autre tableau. Mais nous pouvons réellement calculer le nombre d’étapes (en fonction du nombre d’itérations de la boucle).

Nous pourrions redéfinir légèrement le problème: Remplacer n par (3n + 1) / 2 si impair, et remplacer n par n / 2 si pair. Ensuite, chaque itération fera exactement 8 étapes, mais vous pourriez considérer cette sortingche 🙂 Supposons donc qu’il y ait eu r opérations n <- 3n + 1 et s opérations n <- n / 2. Le résultat sera tout à fait exactement n '= n * 3 ^ r / 2 ^ car n <- 3n + 1 signifie n <- 3n * (1 + 1 / 3n). En prenant le logarithme, nous trouvons r = (s + log2 (n '/ n)) / log2 (3).

Si nous faisons la boucle jusqu’à n ≤ 1 000 000 et que nous avons une table précalculée, combien d’itérations sont nécessaires à partir de n’importe quel sharepoint départ n ≤ 1 000 000, puis calculons comme ci-dessus, arrondi à l’entier le plus proche.

On a rather unrelated note: more performance hacks!

  • [the first «conjecture» has been finally debunked by @ShreevatsaR; removed]

  • When traversing the sequence, we can only get 3 possible cases in the 2-neighborhood of the current element N (shown first):

    1. [even] [odd]
    2. [odd] [even]
    3. [even] [even]

    To leap past these 2 elements means to compute (N >> 1) + N + 1 , ((N << 1) + N + 1) >> 1 and N >> 2 , respectively.

    Let`s prove that for both cases (1) and (2) it is possible to use the first formula, (N >> 1) + N + 1 .

    Case (1) is obvious. Case (2) implies (N & 1) == 1 , so if we assume (without loss of generality) that N is 2-bit long and its bits are ba from most- to least-significant, then a = 1 , and the following holds:

     (N << 1) + N + 1: (N >> 1) + N + 1: b10 b1 b1 b + 1 + 1 ---- --- bBb0 bBb 

    where B = !b . Right-shifting the first result gives us exactly what we want.

    QED: (N & 1) == 1 ⇒ (N >> 1) + N + 1 == ((N << 1) + N + 1) >> 1 .

    As proven, we can traverse the sequence 2 elements at a time, using a single ternary operation. Another 2× time reduction.

The resulting algorithm looks like this:

 uint64_t sequence(uint64_t size, uint64_t *path) { uint64_t n, i, c, maxi = 0, maxc = 0; for (n = i = (size - 1) | 1; i > 2; n = i -= 2) { c = 2; while ((n = ((n & 3)? (n >> 1) + n + 1 : (n >> 2))) > 2) c += 2; if (n == 2) c++; if (c > maxc) { maxi = i; maxc = c; } } *path = maxc; return maxi; } int main() { uint64_t maxi, maxc; maxi = sequence(1000000, &maxc); printf("%llu, %llu\n", maxi, maxc); return 0; } 

Here we compare n > 2 because the process may stop at 2 instead of 1 if the total length of the sequence is odd.

[EDIT:]

Let`s translate this into assembly!

 MOV RCX, 1000000; DEC RCX; AND RCX, -2; XOR RAX, RAX; MOV RBX, RAX; @main: XOR RSI, RSI; LEA RDI, [RCX + 1]; @loop: ADD RSI, 2; LEA RDX, [RDI + RDI*2 + 2]; SHR RDX, 1; SHRD RDI, RDI, 2; ror rdi,2 would do the same thing CMOVL RDI, RDX; Note that SHRD leaves OF = undefined with count>1, and this doesn't work on all CPUs. CMOVS RDI, RDX; CMP RDI, 2; JA @loop; LEA RDX, [RSI + 1]; CMOVE RSI, RDX; CMP RAX, RSI; CMOVB RAX, RSI; CMOVB RBX, RCX; SUB RCX, 2; JA @main; MOV RDI, RCX; ADD RCX, 10; PUSH RDI; PUSH RCX; @itoa: XOR RDX, RDX; DIV RCX; ADD RDX, '0'; PUSH RDX; TEST RAX, RAX; JNE @itoa; PUSH RCX; LEA RAX, [RBX + 1]; TEST RBX, RBX; MOV RBX, RDI; JNE @itoa; POP RCX; INC RDI; MOV RDX, RDI; @outp: MOV RSI, RSP; MOV RAX, RDI; SYSCALL; POP RAX; TEST RAX, RAX; JNE @outp; LEA RAX, [RDI + 59]; DEC RDI; SYSCALL; 

Use these commands to comstack:

 nasm -f elf64 file.asm ld -o file file.o 

See the C and an improved/bugfixed version of the asm by Peter Cordes on Godbolt . (editor’s note: Sorry for putting my stuff in your answer, but my answer hit the 30k char limit from Godbolt links + text!)

C++ programs are translated to assembly programs during the generation of machine code from the source code. It would be virtually wrong to say assembly is slower than C++. Moreover, the binary code generated differs from comstackr to comstackr. So a smart C++ comstackr may produce binary code more optimal and efficient than a dumb assembler’s code.

However I believe your profiling methodology has certain flaws. The following are general guidelines for profiling:

  1. Make sure your system is in its normal/idle state. Stop all running processes (applications) that you started or that use CPU intensively (or poll over the network).
  2. Your datasize must be greater in size.
  3. Your test must run for something more than 5-10 seconds.
  4. Do not rely on just one sample. Perform your test N times. Collect results and calculate the mean or median of the result.

From comments:

But, this code never stops (because of integer overflow) !?! Yves Daoust

For many numbers it will not overflow.

If it will overflow – for one of those unlucky initial seeds, the overflown number will very likely converge toward 1 without another overflow.

Still this poses interesting question, is there some overflow-cyclic seed number?

Any simple final converging series starts with power of two value (obvious enough?).

2^64 will overflow to zero, which is undefined infinite loop according to algorithm (ends only with 1), but the most optimal solution in answer will finish due to shr rax producing ZF=1.

Can we produce 2^64? If the starting number is 0x5555555555555555 , it’s odd number, next number is then 3n+1, which is 0xFFFFFFFFFFFFFFFF + 1 = 0 . Theoretically in undefined state of algorithm, but the optimized answer of johnfound will recover by exiting on ZF=1. The cmp rax,1 of Peter Cordes will end in infinite loop (QED variant 1, “cheapo” through undefined 0 number).

How about some more complex number, which will create cycle without 0 ? Frankly, I’m not sure, my Math theory is too hazy to get any serious idea, how to deal with it in serious way. But intuitively I would say the series will converge to 1 for every number : 0 < number, as the 3n+1 formula will slowly turn every non-2 prime factor of original number (or intermediate) into some power of 2, sooner or later. So we don't need to worry about infinite loop for original series, only overflow can hamper us.

So I just put few numbers into sheet and took a look on 8 bit truncated numbers.

There are three values overflowing to 0 : 227 , 170 and 85 ( 85 going directly to 0 , other two progressing toward 85 ).

But there’s no value creating cyclic overflow seed.

Funnily enough I did a check, which is the first number to suffer from 8 bit truncation, and already 27 is affected! It does reach value 9232 in proper non-truncated series (first truncated value is 322 in 12th step), and the maximum value reached for any of the 2-255 input numbers in non-truncated way is 13120 (for the 255 itself), maximum number of steps to converge to 1 is about 128 (+-2, not sure if “1” is to count, etc…).

Interestingly enough (for me) the number 9232 is maximum for many other source numbers, what’s so special about it? :-O 9232 = 0x2410 … hmmm.. no idea.

Unfortunately I can’t get any deep grasp of this series, why does it converge and what are the implications of truncating them to k bits, but with cmp number,1 terminating condition it’s certainly possible to put the algorithm into infinite loop with particular input value ending as 0 after truncation.

But the value 27 overflowing for 8 bit case is sort of alerting, this looks like if you count the number of steps to reach value 1 , you will get wrong result for majority of numbers from the total k-bit set of integers. For the 8 bit integers the 146 numbers out of 256 have affected series by truncation (some of them may still hit the correct number of steps by accident maybe, I’m too lazy to check).

You did not post the code generated by the comstackr, so there’ some guesswork here, but even without having seen it, one can say that this:

 test rax, 1 jpe even 

… has a 50% chance of mispredicting the branch, and that will come expensive.

The comstackr almost certainly does both computations (which costs neglegibly more since the div/mod is quite long latency, so the multiply-add is “free”) and follows up with a CMOV. Which, of course, has a zero percent chance of being mispredicted.

Even without looking at assembly, the most obvious reason is that /= 2 is probably optimized as >>=1 and many processors have a very quick shift operation. But even if a processor doesn’t have a shift operation, the integer division is faster than floating point division.

Edit: your milage may vary on the “integer division is faster than floating point division” statement above. The comments below reveal that the modern processors have prioritized optimizing fp division over integer division. So if someone were looking for the most likely reason for the speedup which this thread’s question asks about, then comstackr optimizing /=2 as >>=1 would be the best 1st place to look.


On an unrelated note , if n is odd, the expression n*3+1 will always be even. So there is no need to check. You can change that branch to

 { n = (n*3+1) >> 1; count += 2; } 

So the whole statement would then be

 if (n & 1) { n = (n*3 + 1) >> 1; count += 2; } else { n >>= 1; ++count; } 

As a generic answer, not specifically directed at this task: In many cases, you can significantly speed up any program by making improvements at a high level. Like calculating data once instead of multiple times, avoiding unnecessary work completely, using caches in the best way, and so on. These things are much easier to do in a high level language.

Writing assembler code, it is possible to improve on what an optimising comstackr does, but it is hard work. And once it’s done, your code is much harder to modify, so it is much more difficult to add algorithmic improvements. Sometimes the processor has functionality that you cannot use from a high level language, inline assembly is often useful in these cases and still lets you use a high level language.

In the Euler problems, most of the time you succeed by building something, finding why it is slow, building something better, finding why it is slow, and so on and so on. That is very, very hard using assembler. A better algorithm at half the possible speed will usually beat a worse algorithm at full speed, and getting the full speed in assembler isn’t sortingvial.

For the Collatz problem, you can get a significant boost in performance by caching the “tails”. This is a time/memory trade-off. See: memoization ( https://en.wikipedia.org/wiki/Memoization ). You could also look into dynamic programming solutions for other time/memory trade-offs.

Example python implementation:

 import sys inner_loop = 0 def collatz_sequence(N, cache): global inner_loop l = [ ] stop = False n = N tails = [ ] while not stop: inner_loop += 1 tmp = n l.append(n) if n <= 1: stop = True elif n in cache: stop = True elif n % 2: n = 3*n + 1 else: n = n // 2 tails.append((tmp, len(l))) for key, offset in tails: if not key in cache: cache[key] = l[offset:] return l def gen_sequence(l, cache): for elem in l: yield elem if elem in cache: yield from gen_sequence(cache[elem], cache) raise StopIteration if __name__ == "__main__": le_cache = {} for n in range(1, 4711, 5): l = collatz_sequence(n, le_cache) print("{}: {}".format(n, len(list(gen_sequence(l, le_cache))))) print("inner_loop = {}".format(inner_loop)) 

The simple answer:

  • doing a MOV RBX, 3 and MUL RBX is expensive; just ADD RBX, RBX twice

  • ADD 1 is probably faster than INC here

  • MOV 2 and DIV is very expensive; just shift right

  • 64-bit code is usually noticeably slower than 32-bit code and the alignment issues are more complicated; with small programs like this you have to pack them so you are doing parallel computation to have any chance of being faster than 32-bit code

If you generate the assembly listing for your C++ program, you can see how it differs from your assembly.