De l’algorithme au code

Un algorithme est une idée : « regarder au milieu, jeter la moitié qui ne peut pas contenir la cible, recommencer ». Un programme, ce sont des milliers de décisions précises que l’idée laisse ouvertes. Où l’intervalle commence-t-il et où finit-il ? Et si la liste est vide ? Quel type contient un indice ? Que renvoie « pas trouvé » ? La plupart des bogues se logent dans ces décisions, pas dans l’idée.

Ce chapitre, le dernier du niveau des algorithmes, emmène un algorithme, la recherche dichotomique, décrite dans le chapitre recherche et tri, depuis l’idée jusqu’au code machine qui s’exécute sur un vrai CPU. En chemin, il montre les outils qui rendent la traduction fiable : un pseudo-code précis, un invariant explicite, des tests qui essaient tous les cas limites, et un coup d’œil à ce que le compilateur a réellement produit.

De l’idée au pseudo-code

La première étape consiste à écrire l’algorithme assez précisément pour que chaque décision soit prise, sans s’engager dans un langage. Un bon pseudo-code de la recherche dichotomique rend l’intervalle explicite :

search(a[0..n), x):            -- a is sorted; find x, or report that it's absent
    lo ← 0, hi ← n              -- the part of a still to be searched is a[lo..hi)
    while lo < hi:
        mid ← a point with lo ≤ mid < hi, near the middle
        if a[mid] < x:  lo ← mid + 1     -- x can only be to the right of mid
        if a[mid] > x:  hi ← mid         -- x can only be to the left of mid
        if a[mid] = x:  return mid
    return "absent"

La notation a[lo..hi) désigne un intervalle semi-ouvert : il inclut lo et exclut hi. C’est un choix, et c’est ce choix qui rend tout le reste simple. Le tableau entier est a[0..n) ; l’intervalle vide, c’est lo = hi ; le nombre d’éléments restants est hi − lo. Mélanger intervalles semi-ouverts et fermés dans la même boucle est la source la plus courante de bogues de recherche dichotomique.

Du pseudo-code au C

Traduire en C force les décisions restantes :

• « absent » doit être une valeur. Renvoyer −1 fonctionne parce qu’un indice valide n’est jamais négatif.
• Les indices ont besoin d’un type. int convient pour des tableaux de moins de 2³¹ éléments ; au-delà, les programmeurs C utilisent size_t, et le type va se révéler important.
• « près du milieu » a besoin d’une formule. (lo + hi) / 2 paraît évident. Elle contient un bogue.

int bsearch_int(const int *a, int n, int x) {
    int lo = 0, hi = n;              /* invariant: if x is in a, it is in a[lo..hi) */
    while (lo < hi) {
        int mid = lo + (hi - lo) / 2;
        if (a[mid] < x) lo = mid + 1;
        else if (a[mid] > x) hi = mid;
        else return mid;
    }
    return -1;
}

Choisir la structure de données

La recherche dichotomique ne fonctionne que sur des données gardées triées dans un tableau : elle a besoin de sauter à n’importe quelle position en O(1), ce qu’une liste chaînée ne sait pas faire. L’algorithme dicte donc la structure, et la structure a ses propres coûts : insérer dans un tableau trié oblige à décaler tout ce qui suit le nouvel élément, en O(n).

C’est pourquoi le choix dépend de ce que le programme fait d’autre avec les données :

Choisir la structure fait partie du choix de l’algorithme. La meilleure recherche ne sert à rien si garder les données triées coûte plus qu’elle ne rapporte.

Les invariants : pourquoi la boucle est juste

Comment savoir que la boucle est correcte ? Pas en essayant quelques exemples, mais en énonçant ce qui reste vrai à chaque tour (l’invariant de boucle) et en vérifiant trois choses :

1. Il est vrai au départ. Avant la boucle, lo = 0 et hi = n : l’intervalle est le tableau entier, donc si x est dans a, il est dans a[lo..hi). Vrai.
2. Chaque tour le préserve. Si a[mid] < x, alors, comme le tableau est trié, tous les éléments jusqu’à mid inclus sont inférieurs à x, donc x ne peut être que dans a[mid+1..hi), d’où lo = mid + 1. Si a[mid] > x, alors x ne peut être que dans a[lo..mid), d’où hi = mid, et non mid − 1, puisque hi est exclu.
3. À la sortie de la boucle, il donne la réponse. La boucle se termine avec lo = hi : l’intervalle est vide, donc d’après l’invariant x n’est pas dans a, et renvoyer −1 est juste.

Il faut aussi que la boucle se termine : hi − lo est un entier qui diminue strictement à chaque tour (mid est dans [lo, hi), donc lo = mid + 1 comme hi = mid réduisent l’intervalle), et il ne peut pas devenir négatif. L’invariant nous dit que le programme est juste s’il s’arrête ; l’intervalle qui rétrécit nous dit qu’il s’arrête.

Chaque ligne du code découle de l’invariant. C’est le vrai usage des invariants : non pas une documentation écrite après coup, mais l’outil qui vous dit, pendant que vous écrivez, s’il faut mid ou mid − 1, < ou <=.

Tester les cas limites

Le raisonnement attrape la plupart des bogues ; les tests attrapent le reste, et les bogues de la recherche dichotomique se cachent aux limites : un tableau vide, un tableau d’un élément, une cible avant le premier élément, après le dernier, entre deux éléments, égale au premier ou au dernier. Une petite fonction de recherche a un petit espace d’entrées, donc on peut le tester en entier jusqu’à une certaine taille, et comparer avec une version trop simple pour être fausse : ici, la recherche linéaire.

La démo ci-dessous essaie toutes les tailles de tableau de 0 à 16, remplis de 1, 3, 5…, et toutes les cibles de 0 à 2n + 1, de sorte que chaque élément et chaque trou soient cherchés. Chaque case est un test, comparé à la recherche linéaire :

La version correcte réussit les 306 tests, avec 964 sondages en tout. Cassez-la maintenant exprès, comme un vrai programmeur le ferait par accident :

• hi = mid - 1 confond intervalles fermés et semi-ouverts : 56 échecs sur 306. La recherche saute des éléments situés juste à côté d’un milieu écarté. Cliquez sur une case rouge pour voir lequel elle a sauté.
• lo <= hi laisse la boucle tourner une fois de plus quand la plage est vide. Quand lo = hi et a[mid] > x, l’affectation hi = mid ne change rien, et la boucle tourne indéfiniment : 136 des tests ne se terminent jamais. Et chaque recherche d’une cible plus grande que le dernier élément finit par lire a[n], un élément après la fin du tableau.

Tester exhaustivement les petits cas est l’une des techniques les plus efficaces qui soient : c’est peu coûteux, cela trouve exactement les cas limites qu’on oublie, et les bogues d’algorithmes comme celui-ci se montrent presque toujours sur de petites entrées.

Le débordement que personne n’avait testé

Voici ce que donnent les deux formules du milieu avec des int de 32 bits quand lo vaut 1,5 milliard et hi 2 milliards :

(lo + hi) / 2      = -397483648
lo + (hi - lo) / 2 = 1750000000

La somme lo + hi, 3,5 milliards, ne tient pas dans un int de 32 bits, dont le maximum est 2 147 483 647. Elle reboucle sur un nombre négatif (le chapitre sur les types de données explique pourquoi), et le « milieu » devient −397 483 648. La formule sûre ne calcule jamais rien de plus grand que hi.

C’est le bogue qui est resté environ neuf ans dans Arrays.binarySearch de Java, jusqu’à ce que Joshua Bloch le signale en 2006, après qu’il eut frappé un vrai programme. Personne ne l’avait attrapé parce que le tester demande un tableau de plus d’un milliard d’éléments, et que personne ne teste avec ça. Nous, si. Un programme C peut réserver 8 Go d’espace d’adressage pour 2 milliards d’int avec mmap ; les pages ne sont jamais touchées, sauf la trentaine que lit la recherche dichotomique, donc elles ne coûtent presque aucune mémoire réelle. Y chercher une valeur plus grande que tous les éléments pousse lo vers la fin, et la somme déborde. Sur ce Mac :

Aucun des deux résultats qui « marchent » n’est un sursis. La version non optimisée lit de la mémoire qui n’appartient pas au tableau, et que cela plante ou non dépend de l’endroit où la randomisation de l’espace d’adressage a placé les choses lors de cette exécution. Quant à la version optimisée : en C, le débordement d’un entier signé est un comportement indéfini, donc le compilateur a le droit de supposer qu’il n’arrive jamais, et clang, en déduisant que lo + hi n’est donc jamais négatif, a calculé le milieu avec un décalage non signé, qui se trouve donner la bonne réponse ici. Un autre compilateur, une autre version ou un autre code autour peuvent faire un autre choix. Le sanitizer de comportement indéfini (-fsanitize=undefined) est l’outil qui signale le bogue dans tous les cas. Une raison de plus de lancer les tests avec les sanitizers activés.

Ce qu’en fait le compilateur

Voici la version correcte compilée par clang avec -O2 pour x86-64 :

bsearch_int:
    test    esi, esi              ; n <= 0? then not found
    jle     .notfound
    xor     ecx, ecx              ; lo = 0
    jmp     .loop
.right:
    inc     eax                   ; mid + 1
    mov     ecx, eax              ; lo = mid + 1
    mov     eax, esi              ; hi unchanged
.left:
    mov     esi, eax              ; hi = mid (or unchanged)
    cmp     ecx, eax
    jge     .notfound             ; while (lo < hi)
.loop:
    mov     eax, esi
    sub     eax, ecx              ; hi - lo
    shr     eax                   ; / 2, as an unsigned shift
    add     eax, ecx              ; + lo = mid
    cmp     dword ptr [rdi + 4*rax], edx    ; a[mid] vs x
    jl      .right
    jg      .left
    ret                           ; equal: return mid
.notfound:
    mov     eax, -1
    ret

(Commentaires ajoutés ; étiquettes renommées.) La division par 2 est devenue un simple shr. Pour un int signé, / 2 demande normalement des instructions supplémentaires pour arrondir vers zéro les nombres négatifs, mais le compilateur a prouvé que hi − lo n’est jamais négatif dans la boucle, le même fait que notre invariant. La comparaison a[mid] < x est un seul cmp avec un opérande mémoire à index mis à l’échelle, comme dans le chapitre sur les modes d’adressage, et le if à trois branches se réduit à deux sauts conditionnels sur les mêmes indicateurs. Pour l’ARM 64 bits, clang calcule tout le milieu en deux instructions, la seconde décalant et additionnant d’un coup :

    sub     w8, w1, w9            // hi - lo
    add     w8, w9, w8, lsr #1    // lo + (hi - lo) / 2
    ldr     w10, [x0, w8, uxtw #2]        // a[mid]

Comment il s’exécute

Sur la vraie machine, le grand O de l’algorithme se retrouve presque exactement. En cherchant 2 000 valeurs aléatoires (dont la moitié absentes) dans un tableau trié d’un million d’int :

Environ 2 000 fois plus rapide, comme le prévoit le chapitre sur la complexité : 20 sondages au lieu de centaines de milliers de comparaisons. Mais 20 sondages en une centaine de nanosecondes, cela fait environ 5 ns par sondage, bien plus que la poignée d’instructions que chacun exécute. Les sondages sautent d’un bout à l’autre d’un tableau de 4 Mo, donc les premiers ratent en général les caches, et la direction de chaque pas est aléatoire, donc les sauts conditionnels sont mal prédits environ une fois sur deux. Les implémentations rapides remplacent les sauts par des déplacements conditionnels et rangent le tableau dans un ordre favorable au cache : le même algorithme, adapté à la machine qui est en dessous.

Le relais vers les langages de programmation

C’est le chemin que prend tout programme : un algorithme, rendu précis par un pseudo-code et des invariants, traduit dans un langage dont les détails (types entiers, débordements, comportement indéfini) peuvent casser une idée juste, puis compilé en instructions dont la vitesse dépend des caches et des prédicteurs de branchement.

Le niveau suivant, en descendant, est le langage de programmation lui-même. On décrit souvent les langages de haut niveau comme traduits par des compilateurs, parfois interprétés, et Java comme compilé en bytecode ensuite interprété. C’est incomplet depuis longtemps : la machine virtuelle Java compile en code machine le bytecode souvent exécuté, pendant que le programme tourne. Le chapitre langages compilés, interprétés et compilés à la volée commence là, avec la même boucle exécutée de trois façons.

À retenir

• La plupart des bogues se logent dans les décisions qu’un algorithme laisse ouvertes : bornes des intervalles, entrées vides, types, valeurs de retour.
• Écrivez un pseudo-code aux intervalles explicites ; utilisez des intervalles semi-ouverts [lo, hi) de façon cohérente.
• L’invariant de boucle, « si x est dans a, il est dans a[lo..hi) », dicte chaque ligne ; un hi − lo qui diminue strictement prouve la terminaison.
• Testez exhaustivement les petites entrées contre une version trivialement correcte : 306 cas ont attrapé hi = mid − 1 (56 échecs) et lo <= hi (ne se termine jamais).
• (lo + hi) / 2 déborde au-delà de 2³¹ ; sur un vrai tableau de 2 milliards d’éléments, il a planté ou lu hors du tableau en -O0, a marché par chance en -O2, et a été signalé par -fsanitize=undefined. Utilisez lo + (hi − lo) / 2.
• Le compilateur transforme / 2 en un seul décalage parce qu’il peut prouver que la valeur est positive ou nulle ; sur le vrai CPU, la recherche dichotomique a été 2 000 fois plus rapide que la recherche linéaire, et son coût restant, ce sont les défauts de cache et les erreurs de prédiction de branchement.