Rappelons que les différentes mémoires sont des éléments essentiels des ordinateurs. On distingue en général deux sortes de mémoires
Toutes les données (nombres, caractères, noms, adresses, etc.) manipulées doivent être stockées dans l’ordinateur. A cet effet, toutes ces mémoires sont divisées en petites cases appelées octets ; ces octets peuvent stocker une suite de 8 chiffres égaux à 0 ou à 1 (des bits), soit 256 éléments différents. Ces octets sont numérotés dans chaque mémoire de 0 jusqu’à la capacité de la mémoire (moins un). L’ordinateur se charge d’effectuer une traduction entre les données entrées et un ensemble d’octets, c’est l’opération de codage. Typiquement
A partir de là, on voit que toute donnée (la liste des élèves d’une classe et leurs notes au cours de l’année, l’annuaire du téléphone) a une taille que l’on peut exprimer en octets. Cette taille dépend de la méthode de codage utilisée (par exemple MS-DOS code les caractères sur un octet alors que Windows les code sur deux octets) mais obéit à quelques règles simples
Nous allons dans la suite de cette section étudier le rapport entre la taille des données et le temps d’exécution d’un traitement concernant ces données. Il nous suffit donc d’étudier le rapport entre le nombre N d’objets élémentaires à traiter et le temps d’exécution.
Rappelons que si l’on dispose de deux fonctions f et g définies sur N et à valeurs
réelles positives, on dit que f(n) = O(g(n)) s’il existe une constante K ≥ 0 telle que
∀n N, f(n) ≤ Kg(n).
Considérons une donnée constituée d’un certain nombre d’ objets et un certain traitement concernant ces données. Soit T(N) le temps d’exécution de ce traitement lorsqu’il concerne N objets. Ce temps d’exécution dépend évidemment de la machine utilisée et du codage des données, mais il possède certains invariants : en particulier, si l’on utilise deux machines de performances différentes, les temps T1(N) et T2(N) sont en gros proportionnels (le facteur de proportionnalité mesurant justement la différence de rapidité des deux machines).
Ceci va nous permettre de différencier certains types d’algorithmes, du point de vue de leur complexité, c’est à dire de la rapidité de croissance de T(N) quand N devient grand. On dit qu’un algorithme de traitement de données est à complexité
Ce temps d’exécution dépend bien entendu du nombre d’opérations élémentaires effectuées sur ces données et de la plus ou moins grande rapidité de ces opérations élémentaires. Parmi ces opérations élémentaires on peut citer
Dans l’évaluation de ce temps d’exécution, on sera amené à préciser quelles sont les opérations élémentaires à prendre en considération : si l’on effectue à peu près autant d’additions que de multiplications, on peut négliger les additions qui sont beaucoup plus rapides et ne compter que les multiplications, les échanges en mémoire interne peuvent souvent être négligés vis à vis d’échanges entre la mémoire interne et la mémoire externe, etc.
Prenons une machine (très rapide) capable d’accomplir 108 opérations par seconde, et imaginons des algorithmes qui effectuent un traitement donné dans un temps T(N). Nous donnons ci dessous un tableau des temps de traitement que l’on peut espérer pour N éléments suivant la complexité de l’algorithme (où s,m,h,j,a désignent respectivement la seconde, la minute, l’heure, le jour, l’année et u l’âge estimé de l’univers, et où les nombres totalement inimaginables n’ont pas été affichés)
log 2N | N | N log N | N2 | N4 | 2N | 4N | N! | |
N = 5 | 2.10-8s | 5.10-8s | 8.10-8s | 3.10-7s | 6.10-6s | 3.10-7s | 1.10-5s | 1.10-6s |
N = 10 | 2.10-8s | 1.10-7s | 2.10-7s | 1.10-6s | 1.10-4s | 1.10-5s | 1.10-2s | 4.10-2s |
N = 15 | 3.10-8s | 2.10-7s | 4.10-7s | 2.10-6s | 5.10-4s | 3.10-4s | 1.101s | 2h |
N = 20 | 3.10-8s | 2.10-7s | 6.10-7s | 4.10-6s | 2.10-3s | 1.10-2s | 2h | 6.102a |
N = 30 | 3.10-8s | 3.10-7s | 1.10-6s | 9.10-6s | 8.10-3s | 1.101s | 3.102a | 6.106u |
N = 50 | 4.10-8s | 5.10-7s | 2.10-6s | 3.10-5s | 6.10-2s | 1.102j | 2.104u | 6.1038u |
N = 100 | 5.10-8s | 1.10-6s | 5.10-6s | 1.10-4s | 1s | 2.104u | 4.1034u | |
N = 500 | 6.10-8s | 5.10-6s | 3.10-5s | 3.10-3s | 1.101m | |||
N = 1000 | 7.10-8s | 1.10-5s | 7.10-5s | 1.10-2s | 2.100h | |||
N = 5000 | 9.10-8s | 5.10-5s | 4.10-4s | 3.10-1s | 7.101j | |||
N = 10000 | 9.10-8s | 1.10-4s | 9.10-4s | 1s | 3a | |||
N = 50000 | 1.10-7s | 5.10-4s | 5.10-3s | 3.101s | 2.103a | |||
On constate que les algorithmes hyperexponentiels deviennent impraticables dès que N dépasse une dizaine, que les algorithmes exponentiels dépassent les limites acceptables dès que N atteint la trentaine. Ces deux types d’algorithmes sont donc incapables de traiter raisonnablement même de petits volumes de données.
On constate également, que pour rechercher un nom parmi 50000 personnes, entre le temps d’exécution d’une recherche linéaire en O(N) et celui d’une recherche dichotomique en O(log 2N), il existe un facteur 1000. De même, pour trier 10000 nombres réels (ce qui est courant en synthèse d’images avec une grille 100 × 100), on gagne un facteur 1000 en utilisant un tri en O(N log 2N) comme le tri fusion, par rapport à un tri en O(N2) comme le tri par insertion : sur la même machine, l’un des tris pourra prendre 4 secondes alors que l’autre tri demandera une heure !
Remarque Il convient de nuancer un peu la dernière argumentation. Les algorithmes performants sont souvent complexes à écrire et cette complexité peut se ressentir sur les temps d’exécution pour de petits volumes de données. Autrement dit, pour un volume de données de l’ordre de quelques unités, il arrive souvent que les algorithmes bêtes soient plus rapides que les algorithmes intelligents. Entre un algorithme bête en N2∕2 et un algorithme intelligent en 3N log 2N, la différence ne se fera sentir en faveur du second qu’à partir de N = 30. Il peut parfois être nécessaire d’adapter l’algorithme utilisé au volume des données à traiter. N’oublions pas également que les algorithmes sophistiqués utilisent souvent plus de mémoire, ce qui peut être un inconvénient même si celle-ci n’est plus une denrée rare ; il faut savoir établir un compromis entre rapidité et occupation mémoire.
La méthode de programmation généralement désignée par diviser pour régner procède d’une philosophie dont nous allons tenter de donner une idée. Supposons que nous ayions à traiter un problème algorithmique concernant un objet de taille N (cet objet pouvant être lui même un ensemble d’objets). On peut procéder en trois étapes
Cette méthode est bien connue des adeptes de Machiavel : pour détruire un ensemble d’ennemis, les diviser en sous ensembles que l’on écrase un à un, la fusion des résultats ne posant alors guère de difficultés.
Beaucoup d’algorithmes performants se fondent sur ces principes, très souvent avec p = 2. C’est le seul cas que nous examinerons par la suite, mais le lecteur n’aura pas de mal à adapter raisonnements et méthodes à des cas plus complexes. Dans ce cas la méthode devient donc
Remarque En classe de Mathématiques Supérieures, on pourra sauter les démonstrations de cette section et se contenter de retenir les énoncés des deux théorèmes ci-dessous. Ces résultats seront repris dans un cadre un peu plus général en classe de Mathématiques Spéciales.
Appelons T(N) le temps nécessaire pour traiter un objet de taille N, P(N) le temps nécessaire pour partitionner l’objet de taille N en 2 sous-objets de taille approximative N∕2 d’une manière adéquate et F(N) le temps nécessaire pour fusionner deux solutions de taille N∕2 en une solution de taille N. On voit que la fonction T(N) va vérifier
Supposons en particulier que N est de la forme 2k et posons t(k) = T(2k). On aura alors
Nous nous intéresserons tout particulièrement au cas où P et F sont des fonctions polynomiales de N. Quitte à les majorer, on peut alors supposer que P(N) = αNp et F(N) = βNp, auquel cas on obtient une suite t(k) vérifiant
avec γ = α + β. Posons xk = 2-kt(k). La suite xk vérifie la relation de récurrence
Si p = 1 (c’est à dire si le temps nécessaire à la partition et à la fusion est proportionnel à la taille N des données ce qui est un cas courant), on obtient xk = x1 + γ(k - 1) ~ γk, donc t(k) = 2kxk ~ γk2k. Autrement dit T(N) ~ γN log 2N lorsque N est de la forme 2k. Si N n’est pas une puissance de 2, choisissons k de telle sorte que 2k-1 < N < 2k. On a alors T(2k-1) ≤ T(N) ≤ T(2k) où k est égal à la partie entière de log 2N équivalente à log 2N et 2k = 2 ⋅ 2k-1 < 2N, ce qui nous donnera une majoration du type T(N) ≤ δN log 2N.
Si p > 1, on a
avec δ = γ. On a alors t(k) = 2kxk ~ δ2kp. On obtient donc T(N) ~ δNp
lorsque N est de la forme 2k. Lorsque N n’est pas de la forme 2k, on choisit k de telle sorte que
2k-1 < N < 2k. On a alors
puisque 2k = 22k-1 < 2N.
Théorème 5.2.1 On suppose que les temps de partition et de fusion sont des O(Np) pour un problème de taille N. Alors le temps de calcul dans une méthode ”diviser pour régner” est
En fait les méthodes diviser pour régner sont parfois un peu plus complexes et ne permettent pas tout à fait de traiter les deux sous problèmes de taille N∕2 de manière indépendante. Comme nous le verrons par la suite, il peut être nécessaire en fait de traiter q sous problèmes de taille N∕2 ce qui conduit, pour le temps de calcul T0(N) à des relations de récurrence du type T0(N) = qT0(N∕2) + P(N) + F(N). Supposons de nouveau que P(N) + F(N) ≤ γNp. Alors T0(N) est clairement majoré par la solution de la récurrence T(N) = γT(N∕2) + γNp. Posons alors t(k) = T(2k). Nous avons donc
Soit ω = log 2q et xk = 2-ωktk. On a alors
soit encore
On en déduit que xk = x0 + γ ∑ i=1k2i(p-ω).
Si p > ω, soit 2p > q, ceci nous donne
soit encore
Par contre, si p < ω, soit 2p < q, la suite xk admet une limite. On a donc xk = O(1), ce qui montre que
Enfin, si p = ω, soit 2p = q, on a xk = x0 + γ(k - 1) = O(k). On en déduit que
Les mêmes méthodes par encadrement que dans le cas q = 2 montre que ces estimations restent valables lorsque N n’est pas de la forme 2k. On a donc
Théorème 5.2.2 On suppose que dans une méthode ”diviser pour régner”, la méthode conduise à résoudre q problèmes de taille N∕2 avec des processus de partition et de fusion en O(Np), si bien que le temps de calcul T(N) vérifie une relation
Alors, en posant ω = log 2q
Nous allons évaluer le temps de calcul T(N) de la puissance xN d’un entier x en fonction de l’exposant N. Nous avons déjà vu deux méthodes pour calculer cette puissance
Dans les deux méthodes, le temps de calcul est proportionnel à N. Nous allons donc utiliser
une méthode diviser pour régner pour diminuer le temps de calcul. Pour cela nous remarquons
que xN = (x2)p si N = 2p est pair et xN = x × (x2)p si N = 2p + 1 est impair. Nous
obtenons ainsi T(N) = T(N∕2) + α si N est pair et T(N) = T() + 2α si N est
impair, où α est le temps nécessaire à une multiplication. Dans tous les cas, on a
donc T(N) = T(N∕2) + O(1), ce qui d’après le théorème précédent nous donne
T(N) = O(log 2N).
La fonction puissance ainsi définie est naturellement récursive et on obtient une fonction Caml
|
Une version non récursive est un peu plus difficile à construire. Définissons une suite xn par
x0 = x et xn+1 = (xn)2. On a facilement par récurrence ∀n ≥ 0, xn = x2n
. Décomposons alors
n en base 2 sous la forme n = m020 + m121 + + mk2k avec mi
{0,1}. On a alors
xn = x0m0x1m1…xkmk. Or on sait calculer les mi : il suffit de définir les ni et mi par récurrence
par
où l’on désigne par p div q le quotient de la division entière de p par q. On stocke alors dans trois références d’une part ni, d’autre part xi et enfin x0m0x1m1…ximi. Ceci conduit à la fonction suivante :
|
Proposition 5.2.1 La fonction puiss2it calcule xn en un temps O(log 2N)
Démonstration: L’assertion sur le temps est évidente puisque p étant à itération divisé par 2 et étant initialisé à n au départ, le corps de la boucle ne peut être effectué qu’au plus log 2N fois. Ceci montre en même temps que la boucle s’achève bien.
Il nous suffit donc d’exhiber un invariant de la boucle. Or on constate que l’on a l’invariant suivant à la fin de la i-ième itération
y = xi, p = ni et res = x0m0x1m1…xi-1mi-1
C’est vrai après la 0-ième itération, c’est à dire à l’initialisation de la boucle puisque y = x = x0, p = n = n0 et res = 1. De plus, si cette condition est vérifiée après la i-ième itération, on a alors p mod 2 = ni mod2 = mi et en distinguant suivant que cet entier vaut 1 ou 0, on obtient après l’affectation de res,
Après les affectations de y et p, on a y = xi2 = xi+1 et p = ni mod2 = ni+1, donc la propriété est encore vérifiée après la (i + 1)-ième itération.
A la sortie de la boucle, on a alors ni = 0, donc i = k + 1 et alors res = x0m0x1m1…xkmk = xn, ce que l’on voulait montrer.
Les méthodes (récursive ou itérative) qui permettent de calculer la puissance n-ième d’un entier s’adaptent sans aucune difficulté au calcul de la puissance n-ième de tout objet pour lequel cela peut avoir un sens (nombre réel, polynôme, matrice, etc). Revenons par exemple au calcul des nombres de Fibonacci Fn. On vérifie immédiatement que
d’où bien entendu
Le calcul de la puissance n-ième de la matrice pouvant être effectué en un temps O(log 2N), il en est donc de même du calcul de Fn. De la méthode récursive brutale d’ordre 2 à la méthode récursive dichotomique (diviser pour régner), nous sommes donc passés d’un temps de calcul en O(2n) à un temps de calcul en O(log 2n) ; de performances intermédiaires, la méthode itérative ou la méthode récursive d’ordre 1 ont un temps de calcul en O(n).
Pour simplifier, nous travaillerons avec des polynômes à coefficients entiers. Si
P(X) = ∑
i=0naiXi Z[X], nous stockerons les coefficients dans un tableau de taille n + 1, la
valeur de ai étant stockée dans l’élément de numéro i du tableau. Pour des raisons de
commodité, nous définirons une fonction qui renvoie le degré d’un polynôme (égal à la
longueur du tableau diminuée de 1) et une autre fonction qui initialise un polynôme
P(X) = ∑
i=0naiXi avec tous les ai égaux à 0. Pour cela nous utiliserons les deux fonctions
prédéfinies de Caml : la fonction vect_length qui renvoie la longueur d’un tableau et la
fonction make_vect qui renvoie un tableau de longueur donné avec tous les éléments initialisés
à une valeur donnée :
|
La somme de deux polynômes est facile à définir à l’aide de la fonction suivante :
|
Note de programmation On a ici utilisé une petite astuce de programmation en rendant la fonction récursive. Ceci nous permet, dans le cas où le degré de p est strictement plus petit que le degré de q, d’intervertir facilement les paramètres en se contentant de rappeler la fonction add_poly avec les paramètres q et p, grâce à la commutativité de l’addition. Ceci améliore sensiblement la lisibilité du programme sans ajouter de complexité supplémentaire. En fait ici, la récursivité n’est que de façade.
La structure même des deux boucles indexées successives montre la proposition suivante :
Proposition 5.2.2 La fonction add_poly calcule la somme de deux polynômes p et q en un temps O(max(degp,degq)).
En ce qui concerne le produit de deux polynômes, dans une première méthode, nous utiliserons le schéma de calcul suivant
qui est le plus facile à programmer. Il suffit donc de faire parcourir à i l’intervalle [0,m], puis de faire parcourir à j l’intervalle [0,n] en ajoutant le produit aibj au terme de degré i + j du polynôme produit. Ce polynôme produit aura été préalablement initialisé à 0.
|
L’imbrication des deux boucles indexées montre que le temps de calcul du produit de deux polynômes de degrés respectifs m et n est proportionnel à (m + 1)(n + 1). Pour deux polynômes de degrés inférieurs ou égaux à N, le temps de calcul par cette méthode est donc un O(N2) (d’où la qualification d’algorithme quadratique).
Prenons maintenant deux polynômes P(X) et Q(X) de degrés strictement inférieurs à N = 2k = 2M et posons P(X) = P1(X) + XMP2(X), Q(X) = Q1(X) + XMQ2(X) (partition). On a donc
Il semble qu’il faille calculer 4 produits de polynômes de degré strictement inférieurs à M = N∕2 pour calculer le produit P(X)Q(X). Comme le temps nécessaire à l’addition de polynômes de degré N est un O(N), on trouve un temps de calcul T(N) = 4T(N∕2) + O(N). Avec les notations du paragraphe précédent, on a donc p = 1 et q = 4 ; on a 2p < q et donc T(N) = O(Nlog 24) = O(N2). Nous n’avons donc rien gagné par rapport à la méthode brutale.
En fait un peu d’astuce va nous éviter de calculer un des produits. Il suffit de remarquer que
autrement dit, si nous posons R1 = P1Q1, R2 = (P1 + P2)(Q1 + Q2), R3 = Q1Q2, nous avons (fusion)
Nous n’avons plus que trois produits de polynômes de degrés strictement inférieurs à M = N∕2 à effectuer pour calculer le produit P(X)Q(X). Le temps de calcul nécessaire pour effectuer la partition P(X) = P1(X) + XMP2(X), Q(X) = Q1(X) + XMQ2(X) et celui nécessaire pour opérer la fusion
étant clairement un O(N), le temps de calcul nécessaire pour calculer par cette méthode le produit de deux polynômes de degrés strictement inférieurs à N = 2k vérifiera
ce qui nous donne grâce au théorème du paragraphe précédent
avec log 23 ~ 1.58, ce qui est un gain non négligeable.
La méthode de calcul ci-dessus, encore appelée méthode de Karatsuba, est clairement récursive, la récursivité portant sur l’entier N = 2k. Pour k = 0, on a N = 1, auquel cas les deux polynômes sont constants, leur produit est alors évident ce qui initialise la récursion.
Les polynômes de degré strictement inférieur à N = 2k seront stockés dans des tableaux de longueur N. Nous utiliserons quelques fonctions utilitaires :
|
La fonction fusion reçoit en paramètre trois polynômes r1, r2 et r3 de degrés strictement inférieurs à N et l’entier N = 2M, et renvoie le polynôme R1(X) + R2(X)XM + R3(X)X2M grâce à un calcul évident.
|
A partir de là, il suffit de suivre la formulation mathématique pour obtenir la procédure de calcul du produit de deux polynômes :
|
Un petit programme nous a permis de calculer le nombre de multiplications de nombres entiers effectuées au cours du calcul de (1 + X)2n par la méthode quadratique et par la méthode diviser pour régner (sachant que la multiplication est de loin l’opération la plus lente parmi toutes celles qui sont effectuées sur les coefficients)
n | méthode quadratique | diviser pour régner |
3 | 38 | 39 |
4 | 119 | 120 |
5 | 408 | 363 |
6 | 1497 | 1092 |
7 | 5722 | 3279 |
Le gain est peu sensible (voire inexistant) pour des polynômes de petit degré et ne commence à devenir sensible que pour n = 6 (résultat de degré 26 = 64). Pour un résultat de degré 128, on voit que la deuxième fonction est environ deux fois plus rapide que la fonction initiale, ce qui n’est pas négligeable1.
Note de programmation Notre deuxième fonction n’est évidemment pas universelle. Elle suppose a priori que les polynômes sont tous deux entrés dans des tableaux de taille N = 2k. Si l’on veut en faire une fonction d’usage général, pour multiplier deux polynômes p et q de degrés respectifs m et n, il faut commencer par déterminer le plus petit entier k tel que m < 2k et n < 2k, ensuite recopier nos polynômes dans deux tableaux de taille 2k, appliquer la fonction mult_poly et enfin recopier le tableau obtenu (de taille 2k+1) dans un polynôme de degré m+n en supprimant les 2k+1 -(m+n+1) derniers coefficients qui sont nuls. Le lecteur vérifiera que ceci ne change pas la complexité de la fonction dont le temps d’exécution est donc un O(max(m,n)log 23).
Soit b un entier strictement positif appelé base de numération. On sait que tout nombre entier
n peut s’écrire de manière unique sous la forme n = akbk + ak-1bk-1 + … + a0 avec
0 ≤ ai ≤ b - 1. Introduisons alors le polynôme Pn(X) = akXk + ak-1Xk-1 + … + a0.
On a n = Pn(b). On en déduit que si m et n sont deux entiers, on peut écrire
mn = Pm(b)Pn(b) = Q(b) où Q(X) est le produit Pm(X)Pn(X). Le polynôme à coefficients
entiers Q(X) s’écrit sous la forme Q(X) = αlXl + … + α0, seulement il ne vérifie
plus la condition 0 ≤ αi ≤ b - 1. Pour calculer les chiffres du produit mn dans la
base b, il nous faut donc renormaliser le polynôme Q(X). Pour cela, il nous suffit de
remarquer que si l’on prend i [0,k] et qi
Z, le polynôme Q(X) et le polynôme
Qi(X) = αlXl + … + (αi+1 + qi)Xi+1 + (αi - bqi)Xi + … + α0, on a Qi(b) = Q(b) ; en effet
L’usage d’une base de numération permet de dépasser les limitations du langage de programmation en termes de taille des entiers manipulés. Il suffit en effet de remplacer l’entier n = akbk + ak-1bk-1 + … + a0 par un tableau noté à la Caml [|a0;…;ak|], qui coïncide justement avec la notation utilisée pour le polynôme Pn(X) = akXk + ak-1Xk-1 + … + a0 vérifiant n = Pn(b). Il nous suffit donc de définir une opération de renormalisation d’un polynôme relativement à la base b pour adapter les procédures de calcul sur les polynômes aux calculs sur les grands nombres entiers stockés sous formes de tableaux en base b.
Note de programmation La procédure de renormalisation d’un polynôme P par
rapport à la base b suit l’algorithme décrit ci-dessus en ce qui concerne tous les coefficients,
sauf celui de plus haut degré. Il n’est en effet pas possible en Caml de faire croître le tableau
au fur et à mesure des besoins pour tenir compte de l’augmentation du degré. Mais il suffit
de remarquer qu’une fois le polynôme P(x) écrit sous la forme P(X) = βlXl + … + β0 avec
∀i [0,l - 1], 0 ≤ βi ≤ b - 1, si l’écriture de βl en base b est βl = γλbλ + … + γ0, alors le
polynôme normalisé de P(X) n’est autre que
Nous utiliserons donc une fonction table_decimale qui renvoie dans un tableau les chiffres de l’écriture en base b d’un entier n. Pour cela nous stockons au fur et à mesure les chiffres dans une liste qui croît dynamiquement suivant l’algorithme évident par divisions euclidiennes successives : le coefficient de b0 est le reste de la division euclidienne de n par b, puis on remplace n par son quotient par b et on recommence jusqu’à obtenir n = 0. Une fois cette liste construite, on la retourne à l’aide de la fonction prédéfinie rev et on la transforme en tableau à l’aide de la fonction prédéfinie vect_of_list.
|
Il nous suffit maintenant de renormaliser le polynôme suivant la base b par l’algorithme suivant :
|
Il est alors facile de définir l’addition et la multiplication de deux nombres entiers écrits en base b dans deux tableaux :
|
A titre d’exemples, nous allons calculer 100 ! en base 10 ; nous utiliserons une petite fonction d’affichage qui se contente de juxtaposer dans le bon ordre (c’est-à-dire l’ordre inverse) les éléments du tableau
|
On constate que les polynômes manipulés par cette méthode peuvent être de très grands degrés (environ 160 pour celui décrivant l’écriture de 100! en base 10) et que la méthode diviser pour régner prend tout son intérêt pour de tels polynômes.
La recherche d’éléments dans un ensemble est une opération essentielle en algorithmique. Dans le cas où l’ensemble n’est muni d’aucune structure, la seule possibilité est une recherche exhaustive : parcourir tout les éléments de l’ensemble jusqu’à trouver l’élément choisi. Clairement, si l’ensemble a N éléments, le temps de recherche est un O(N). Le placage d’une structure sur l’ensemble permet souvent de réduire de manière drastique les temps de recherche.
Nous allons montrer que si l’ensemble est un tableau ordonné, le temps de recherche peut
être réduit en un O(log 2N) par une méthode de type diviser pour régner. Considérons en effet
un ensemble X = {a0,…,aN-1} ordonné de telle sorte que a0 ≤ a1 ≤… ≤ aN-1 et soit b un
élément de même type que les ai. Nous voulons savoir si b appartient à X et éventuellement
calculer un indice i tel que ai = b. Donnons nous pour cela m {0,…,N - 1}. Si b ≤ am, alors
b
X si et seulement si b
{a0,…,am} ; si par contre b > am, alors b
X si et seulement si
b
{am+1,…,aN-1}.
En choisissant un m de l’ordre de N∕2, on voit que le temps de recherche T(N) d’un élément b dans un ensemble à N éléments vérifiera une récurrence du type T(N) = T(N∕2) + O(1). Il s’agit ici d’une méthode diviser pour régner où la partition s’effectue en temps constant (la recherche du milieu m et la comparaison à am) et où la fusion est triviale (s’effectue en temps nul). On sait d’après un paragraphe précédent que dans ce cas T(N) = O(log 2N).
La recherche dans un tableau trié peut donc se faire de manière récursive (nous avons
choisi de renvoyer l’indice i de l’élément ai tel que ai = b si b X, et - 1 si b
X) ; pour cela
nous définissons une procédure récursive cherche qui effectue le gros du travail : elle reçoit en
paramètre les indices i et j qui sont les bornes du sous-tableau ai,ai+1,…,aj dans lequel doit
s’effectuer la recherche.
La procédure de recherche dichotomique se contente de recevoir les paramètres b et X, de créer l’environnement nécessaire à le recherche et d’appeler la procédure cherche avec comme paramètres les limites du tableau initial.
|
Proposition 5.3.1 La fonction de recherche binaire récursive effectue bien la recherche d’un élément b dans un tableau trié de N éléments en un temps O(log 2N).
Démonstration: l’assertion sur le temps de calcul a déjà été démontrée. Il
nous suffit donc de démontrer que la fonction cherche renvoie bien le résultat de la
recherche de b dans le sous-tableau ai ≤… ≤ aj. Nous allons le faire par récurrence
sur j -i. C’est clair, si j -i = 0. Si j -i ≥ 1, m est la partie entière de si
bien que i ≤ m ≤
< j.
En conclusion, l’appel cherche 0 (N-1) renvoie bien le résultat de la recherche de b dans {a0,…,aN-1} donc dans X.
Note de programmation La vérification de la relation i ≤ m ≤ < j est une
étape essentielle de la démonstration. Il serait en effet catastrophique d’avoir un croisement
des indices c’est à dire qu’à un moment donné on puisse appeler la fonction cherche i j
avec i > j ; c’est une des erreurs fréquentes dans les méthodes de dichotomie.
Il est clair que la récursivité ci-dessus est terminale (la fonction se termine par l’appel récursif). Ceci nous permet de construire facilement une version itérative de la recherche en introduisant
|
Proposition 5.3.2 La fonction de recherche binaire itérative effectue bien la recherche d’un élément b dans un tableau trié de N éléments en un temps O(log 2N).
Démonstration: l’assertion sur le temps de calcul a déjà été
démontrée. Il nous faut donc démontrer que la boucle se termine et que la
fonction renvoie bien le résultat de la recherche. Pour cela nous allons exhiber
un invariant de la boucle booléenne. Posons f(i,j) = -1 si b{ai,…,aj} et
f(i,j) = k si k est le plus petit indice tel que b = ak avec i ≤ k ≤ j. Soit
in et jn les valeurs des références i et j à l’entrée dans la n-ième itération.
Montrons alors que la condition in < jn et la valeur de f(in,jn) est un
invariant de la boucle booléenne à l’entrée dans l’itération. A la première
exécution du corps de la boucle, on a i1 = 0, j1 = N - 1 et puisque
ce corps est exécuté c’est que 0 < N - 1 ; quant à f(i1,j1) il vaut tout
simplement f(0,N - 1) (c’est à dire la valeur recherchée). Au début de la
n-ième itération, on sait que in < jn (puisque le corps de la boucle est
exécuté) ; on montre comme dans la démonstration de la fonction récursive
que in ≤ m ≤
< jn. On a donc in ≤ m < m + 1 ≤ jn. Si b ≤ am,
on a alors f(in,jn) = f(in,m) = f(in+1,jn+1) et à l’itération suivante, soit
in+1 = jn+1 auquel cas on sortira de la boucle, soit in+1 < jn+1. Si par
contre, b > am, on a alors f(in,jn) = f(m + 1,jn) = f(in+1,jn+1) et à
l’itération suivante, soit in+1 = jn+1 auquel cas on sortira de la boucle, soit
in+1 < jn+1.
La stricte décroissance de la suite d’entiers naturels jn - in montre que la boucle s’achève au bout d’un nombre fini d’itérations. A ce moment on a in = jn et f(in,jn) = f(i0,j0) = f(0,N - 1) puisque cette valeur est un invariant de la boucle. Or f(in,jn) = f(in,in) vaut in si b = ain et -1 sinon. On a donc bien calculé f(0,N - 1), qui était l’entier recherché.
Remarque Le fait que le tableau soit trié est bien évidemment essentiel. Mais il est aussi essentiel que la structure soit un tableau, ce qui permet d’accéder à l’élément médian am en un temps constant. Si au lieu d’un tableau on travaillait avec une liste, cette recherche ne pourrait se faire qu’en un temps O(N) (puisqu’on ne peut accéder aux éléments que séquentiellement) et la recherche dichotomique perdrait tout son intérêt. C’est toute la différence entre une structure à accès direct (où l’on peut accéder directement à tous les éléments) et une structure à accès séquentiel. Remarquons également qu’il est stupide de commencer par trier un tableau avant d’effectuer une recherche dans ce tableau puisque le tri lui-même demande un temps en O(N log 2N) au minimum ; dans ce cas une banale recherche exhaustive en O(N) est plus performante. Le tri, puis la recherche dichotomique dans un tableau ne peuvent se justifier que si l’on effectue de nombreuses recherches dans ce tableau. Si l’on effectue p recherches dans un tableau à N éléments,
En gros, on peut dire que dès que p devient de l’ordre de log 2N, il peut être intéressant d’effectuer un tri préalable.
Supposons qu’un enseignant se trouve devant deux paquets de copies, chacun des paquets étant trié par ordre croissant de notes. Il cherche à réunir les deux paquets en un seul trié par ordre croissant. La solution est évidente : il prend à chaque étape la copie ayant la plus petite note parmi les sommets des deux paquets jusqu’à avoir épuisé l’un des deux.
Autrement dit la fusion de deux tableaux triés par ordre croissant peut se faire grâce à la fonction suivante, où il faut simplement prendre garde à l’épuisement possible d’un des deux tableaux :
|
Dans le cas où l’on dispose d’une place supplémentaire à la fin des deux tableaux a et b, on peut y disposer une sentinelle dont on est certain qu’elle est supérieure à tous les éléments des deux tableaux. Ceci permet alors de se passer des deux tests de débordement puisque lorsque i = m, on est certain que ai = am ≥ bj et lorsque j = n, on est certain que bj = bn ≥ ai. Cependant l’usage d’une telle sentinelle n’est pas évident dans le cas général : ceci oblige à trouver un élément plus grand que tous les éléments des deux tableaux et ensuite à réserver une place supplémentaire dans chacun d’eux pour y stocker la sentinelle.
Dans tous les cas, on constate que le temps d’exécution de la fusion de deux tableaux triés de longueurs m et n est un O(m + n) puisque le corps de la boucle indexée est exécuté m + n fois.
Une fois l’idée de la fusion de deux tableaux triés acquise, la procédure de tri par fusion d’un tableau ne pose plus de difficulté. Etant donné un tableau de taille N, on le partage en deux tableaux de taille (approximative) N∕2 que l’on trie récursivement et que l’on fusionne ensuite. Bien entendu les appels récursifs se terminent lorsque le tableau est de taille 1 auquel cas il est tout trié.
La partition du tableau se fait en temps constant et la fusion des deux tableaux en un
temps O( +
) = O(N) si bien que le temps T(N) d’exécution du tri par fusion sur un
tableau de N éléments vérifie T(N) = 2T(N∕2) + O(N) + O(1) = 2T(N∕2) + O(N). On en
déduit que T(N) = O(N log 2N).
Il nous faut utiliser une procédure de fusion de deux sous-tableaux consécutifs d’un même tableau, plutôt que de fusion de deux tableaux indépendants. La procédure de fusion des deux tableaux triés {al,al+1,…,am} et {am+1,am+2,…,au} peut s’écrire de la même manière que ci-dessus :
|
En fait, plutôt que de réinitialiser à chaque fusion un nouveau tableau b, il est plus simple de définir un tableau auxiliaire ce qui conduit simplement à supprimer le b=make_vect (u-l+1) 0 à condition de disposer dans l’environnement de la fonction de fusion d’un tableau auxiliaire b. Ceci conduit à une procédure de tri d’un tableau par fusion sous la forme :
|
Note de programmation L’auteur dans une première version avait fait une erreur difficilement décelable en remplaçant l’instruction let b=make_vect N 0 par let b=a ; il oubliait que dans une telle liaison,dans la mesure où un tableau est une structure modifiable en place, le tableau b devenait physiquement égal au tableau a ; toute modification du tableau b (et il y en a) se répercutait sur le tableau a avec toutes les conséquences que l’on peut imaginer. C’est un point à ne pas oublier dans la programmation en Caml.
Proposition 5.3.3 Le tri par fusion trie un tableau de N éléments en un temps O(N log 2N).
Démonstration: Nous allons commencer par montrer que la procédure de fusion fusionne bien deux sous-tableaux triés al ≤ al+1 ≤ … ≤ am et am+1 ≤ am+2 ≤ … ≤ au en un tableau b0 ≤ b1 ≤ … ≤ bl-u qui est ensuite recopié en place dans le tableau a.
L’achèvement de la procédure de fusion est garanti par le fait qu’il s’agit d’une boucle indexée. Il nous suffit donc d’exhiber un invariant de la boucle. Pour cela nous allons utiliser un invariant quadruple à l’entrée dans le corps de la boucle
et la famille b0,…,bk-1 est à une permutation près la famille
C’est clair avant la première itération puisqu’alors i = l, j = m + 1 et k = 0. Supposons que ce soit vrai avant la n-ième itération. Les trois premières conditions sont clairement conservées par le corps de la boucle : les deux premières puisque i n’est pas modifié si i = m + 1 et j n’est pas modifié si j = u + 1, la troisième puisque lors de l’exécution du corps de la boucle, soit i, soit j est augmenté de 1 et ce de manière exclusive, et qu’ensuite k est augmenté de 1 pour l’exécution suivante du corps de la boucle, ce qui conserve l’égalité i + j = l + m + 1 + k pour la (n + 1)-ième itération. Si maintenant ai ≤ bj alors le corps de la boucle affecte ai à bk et augmente i de 1, si bien que l’on a encore b0 ≤… ≤ bk-1 ≤ bk = ai-1 ≤ min(ai,bj) et la quatrième condition est bien vérifiée. De manière similaire, la quatrième condition est bien conservée si bj ≤ ai. Avant l’itération d’indice k = u - l, on a donc
et la famille b0,…,bu-l-1 est à une permutation près la famille
Comme i + j = u + m + 1 avec i ≤ m + 1 et j ≤ u + 1, on a soit i = m + 1 et j = u, soit j = u + 1 et i = m. Dans le premier cas, le corps de la boucle se contente d’affecter au à bk et d’augmenter i, donc après la sortie de la boucle on a b0 ≤… ≤ bk = au et cette famille est à une permutation près la famille
ce qui est ce que l’on souhaitait. Le raisonnement est similaire dans le second cas.
Une fois que l’on sait que la procédure de fusion effectue bien son travail, la démonstration de la procédure tri est une banale récurrence sur u-l. Si u = l la procédure ne fait rien et le tableau est tout trié, donc elle le trie bien. Si u-l ≥ 1, la procédure trie deux sous tableaux de taille strictement inférieure, ce qu’elle fait bien par hypothèse de récurrence, puis fusionne les deux tableaux, ce qui est correct.
L’assertion sur le temps d’exécution a déjà été justifiée. Ceci achève la démonstration.
Il procède des mêmes idées que le tri fusion, en partant d’une autre méthode de partition. Considérons un tableau a0,…,aN-1 de N éléments. Effectuons alors une permutation du tableau pour obtenir un nouveau tableau b0,…,bN-1 tel que a0 = bm avec
Il est clair que si maintenant on trie les sous-tableaux b0,…,bm-1 et bm+1,…,bN-1, le tableau obtenu sera trié.
Contrairement au tri fusion dans lequel la partition est triviale et la fusion plus délicate, dans le tri rapide, la fusion est inexistante (les sous-tableaux étant triés en place) et la partition cruciale. De plus, contrairement au tri fusion, le tri rapide n’a pas besoin de tableau auxiliaire pour effectuer des recopies, il utilise donc moins de mémoire que le tri fusion. Par contre, on verra que dans certains cas très particuliers, il n’a de rapide que le nom !
Pour effectuer la partition d’un sous-tableau al,…,au nous allons effectuer une boucle
indexée par i [l + 1,u] avec un deuxième indice m de manière à maintenir l’invariant
suivant
On initialise la boucle avec i = m = l. Lorsqu’on augmente i de 1, deux cas sont à considérer
Après exécution de la totalité de la boucle, on aura donc la situation
Il suffit alors d’échanger al et am pour avoir la situation
qui est exactement ce que l’on souhaitait.
L’implémentation en Caml de cette procédure de partition peut se faire sous la forme suivante qui partitionne le tableau et retourne l’entier m où est situé le pivot :
|
Il suffit ensuite d’insérer cette procédure de partition à l’intérieur d’une procédure de tri récursive en tenant compte de ce qu’un sous-tableau qui a au plus un élément est tout trié, et d’appeler cette procédure récursive pour le tableau entier.
|
Proposition 5.3.4 Le tri rapide trie un tableau de N éléments en un temps O(N log 2N) dans le cas moyen et O(N2) dans le cas le plus défavorable (tableau déjà trié).
Démonstration: La construction même de la procédure de partition à
l’aide d’un invariant de boucle montre qu’elle effectue bien son travail, c’est à
dire que lors du retour partition l u, l’indice m a été modifié de telle sorte
que l ≤ m ≤ u, ∀i [l,m-1], ai ≤ am et ∀i
[m+1,u], ai > am. Montrons
maintenant par récurrence sur u - l que la procédure tri l u trie bien le
sous-tableau de l à u. C’est clair si u-l ≤ 0, puisque la procédure ne fait rien
et que le sous-tableau est tout trié. Si u - l ≥ 1, alors la procédure tri l u
commence par appeler la procédure de partition. A la sortie de la procédure,
on a donc l ≤ m ≤ u, ∀i
[l,m - 1], ai ≤ am et ∀i
[m + 1,u], ai > am ;
ensuite la procédure s’appelle récursivement sur les deux sous-tableaux de
tailles strictement inférieures al,…,am-1 et am+1,…,au ; par hypothèse de
récurrence la procédure trie bien ces sous-tableaux et à la fin de la procédure
on a donc
ce qui achève la démonstration par récurrence.
Le temps d’exécution de la procédure tri l u dépend essentiellement de la taille des sous-tableaux al,…,am-1 et am+1,…,au. Si on note t(l,u) le temps d’exécution de la procédure tri l u, comme le temps d’exécution de la procédure de partition est clairement proportionnel à u - l, on aura t(l,u) = t(l,m - 1) + t(m + 1,l) + α(l - u). Le cas le plus défavorable est le cas m = l ou m = u c’est à dire où l’élément al est le plus petit ou le plus grand élément du sous-tableau. Dans ce cas, l’un des sous-tableaux est vide, et l’autre est de taille u - l. Si ceci se produit tout au long du tri, c’est à dire si le tableau initial est déjà trié par ordre croissant ou décroissant, la formule de récurrence concernant le temps T(N) de tri du tableau sera T(N) = T(N - 1) + O(N) et on sait alors que T(N) = 0(N2).
Si par contre, le tableau n’est en aucune façon ordonné et donc
parfaitement aléatoire, on peut espérer2
que m est de l’ordre de , auquel cas les deux sous-tableaux sont de
taille moitié. La formule de récurrence devient alors T(N) = 2T(N∕2)+O(N)
et on obtient alors T(N) = O(N log 2N).