ALGORITMIQUE

FASCICULE DE COURS

EUPI

A NNÉE U NIVERSITAIRE

2018-2019
L2 SPI
 Table des matières

1 Introduction à l’algorithmique (Ref. [5]) . . . . . . . . . . . . . . . . . . . . . . . . . 7

1.1 Présentation 7
1.2 Exemple de problème 7
1.2.1 Première Solution . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8
1.2.2 Deuxième solution . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9
1.2.3 Les approximations . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10
1.2.4 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10

2 Les briques de base des algorithmes en Matlab (Ref[1], [2]) . . . . 13

2.1 Introduction à Matlab 13
2.2 Prise en main 13
2.2.1 La fenêtre de commande . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 13
2.2.2 Le répertoire de travail . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 14
2.2.3 L’éditeur . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 14
2.2.4 Le workspace . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 14
2.3 Les variables 15
2.3.1 Les types numériques . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 16
2.3.2 Le type logique . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 17
2.3.3 Les tableaux et matrices . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 17
2.3.4 Les chaînes de caractères . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 19
2.3.5 Les structures . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 19
2.4 Les opérateurs 20
2.4.1 L’addition et la soustraction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 20
2.4.2 Multiplication, division et puissance . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 20
2.4.3 Opérations relationnelles . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 21
2.4.4 Opérations logiques . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 22

3
2.4.5 Quelques autres fonctions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 23
2.5 Les fichiers m 24
2.6 Les fonctions 24
2.7 Les structures conditionnelles 26
2.7.1 Instruction "If" (Si) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 26
2.7.2 Instruction "Switch" . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 26
2.8 Les boucles 26
2.8.1 Instruction "For" (Pour) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 26
2.8.2 Instruction "While" (Tant que) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 27
2.9 Exemple 27

3 Architecture des ordinateurs (Ref. [3], [4]) . . . . . . . . . . . . . . . . . . . . . . 33

3.1 Introduction 33
3.2 Le modèle "Little Man Computer" 33
3.3 Architecture interne des processeurs 35
3.3.1 L’unité de commande . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 35
3.3.2 L’unité de traitement . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 36
3.3.3 Exemple . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 36
3.4 Les pipelines 37
3.5 La mémoire cache 39
3.6 Conclusion 42

4 Algorithmes et performances . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 43
4.1 Introduction (Ref. [6]) 43
4.2 La notation en O 44
4.3 Analyse meilleur cas, pire cas, cas moyen 45
4.3.1 Pire cas . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 45
4.3.2 Meilleur cas . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 46
4.3.3 Cas moyen . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 46
4.3.4 Famille de performance . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 46
4.4 Exemples d’algorithmes 47
4.4.1 Divination (Ref. [5] . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 47
4.5 Anagrammes (Ref. [6] 49
4.5.1 Premier algorithme . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 49
4.5.2 Deuxième solution . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 50
4.5.3 Troisième algorithme . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 51
4.6 Conclusion 52

5 La recherche et le tri (Ref. [6]) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 53

5.1 Introduction 53
5.2 Quelques algorithmes de recherche 53
5.2.1 La recherche séquentielle . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 53
5.2.2 La recherche séquentielle triée . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 54
5.2.3 La recherche dichotomique . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 55
5.2.4 Les tables de hachage . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 57

4
5.3 Quelques algorithmes de tri 62
5.3.1 Le tri a bulle . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 62
5.3.2 Le tri par sélection . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 64
5.3.3 Le tri par insertion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 65
5.3.4 Le tri shell . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 66
5.3.5 Le tri fusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 67
5.3.6 Le tri rapide (Quick-sort) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 69
5.4 Conclusion 71

6 Introduction au calcul numérique (Ref. [7]) . . . . . . . . . . . . . . . . . . . . . 73

6.1 Introduction 73
6.2 Erreur, stabilité et precision 73
6.3 Exemple : circuit RLC 75
6.3.1 Calcul de l’intégrale . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 78
6.3.2 Calcul de la dérivée . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 79
6.3.3 Retour au circuit RLC . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 79
6.3.4 Calcul du sinus . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 82
6.3.5 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 84

5
 Introduction à l’algorithmique (Ref. [5])

1.1 Présentation
L’informatique a fortement évolué ces dernières années avec l’augmentation en puissance des
micro-processeurs tandis les réseaux de communications ont pris une place prépondérante dans la
société. Les langages de programmation et les plate-formes de développement disponibles ont suivi
ces mutations en s’adaptant et en se multipliant. C’est ce cadre complexe et dynamique que les
informaticiens doivent appréhender.
Cet environnent difficile tend à cacher les principes de base régissant la conception de pro-
gramme informatique. Nous allons en aborder quelques un dans ce cours d’introduction à l’algo-
rithmique.
L’étude des algorithmes est largement antérieure aux premiers ordinateurs. Les premiers dont
on a retrouvé la trace datent des Babyloniens au IIIème millénaire avant JC pour la résolution d’équa-
tions simples. Le mot algorithme lui même vient du nom du mathématicien Al-Khwârizmî (IXème
siècle) auteur d’un ouvrage sur la résolution de systèmes d’équations linéaires et quadratiques.
Pour ce chapitre, nous allons étudier un problème de géométrie algorithmique. C’est une branche
qui a été longtemps poussée par l’infographie mais elle est désormais utilisée dans de nombreux
domaines comme la robotique, la vision ou encore la génération de maillage.

1.2 Exemple de problème

Nous allons étudier un des problèmes connus de la géométrie algorithmique : il s’agit de
déterminer l’enveloppe convexe (la forme convexe minimale) entourant un ensemble de points
P sur un plan à deux dimensions. Un exemple avec 15 points est donné sur la Fig. 1.1. Nous
nous attacherons, dans cette partie du cours, à écrire un algorithme de résolution pour n’importe
quelle série de n points Pi constituant P. Il devra trouver l’ensemble Li de points de P composant
l’enveloppe convexe. Ce ensemble devra être trié afin que toute les séries de points (Li−1 , Li , Li+1 )
forment un angle ouvert.

F IGURE 1.1 – 15 points aléatoirement placés dans un plan.

La solution du problème de la Fig. 1.1 est visuellement évidente. Nous pouvons déterminer les
points Li et tracer les segments de l’enveloppe (Fig. 1.2) directement.

7
 Néanmoins, cette méthode ne constitue pas un algorithme. Pour cela, nous devons écrire
une succession d’instructions qui permettent de trouver efficacement l’enveloppe convexe pour
n’importe quelle série de points. Ce problème n’est pas simple, et il ne se rattache pas de manière
évidente aux algorithmes de tri et de recherche que nous verrons plus tard.

F IGURE 1.2 – Solution de l’exemple

La résolution d’un problème implique de le comprendre puis d’utiliser ses propriétés. Ainsi et
pour notre exemple, nous devons commencer par définir ce qui caractérise les points composant
l’enveloppe convexe.

1.2.1 Première Solution

Pour ce premier algorithme, nous allons utiliser la propriété suivante : les points Li de l’ensemble
L composant l’enveloppe convexe ne font partie d’aucun triangle Pi , Pj , Pk de points de P. Vous
pouvez le vérifier sur l’exemple de la Fig. 1.2 en traçant tous les triangles possibles. A partir de
cette idée, nous pouvons écrire la suite d’instructions composant l’algorithme. Soit, en langage
naturel :

Pour tous les points p0 de P faire:

Pour tous les points p1 de {P-p0} faire
Pour tous les points p2 de {P-p0-p1} faire
Pour tous les points p3 de {P-p0-p1-p2} faire
Si p3 est dans le triangle (p0,p1,p2) alors
Marquer p3 comme ne faisant pas partie de l'enveloppe
Créer un ensemble L de points avec les points non marqués

Nous pouvons remarquer que cette étape nécessite 4 boucles imbriqués sur l’ensemble des n points
de P. On dira que la complexité de cet algorithme croît en O(n4 ), nous en reparlerons plus tard. A
la fin de cette étape, nous avons l’ensemble des éléments de L, néanmoins il ne sont pas triés (un
segment Li , Li+1 ne fait pas forcément partie de l’enveloppe convexe.
Ainsi, il nous reste a construire l’enveloppe à partir des points composants L. Pour cela, nous
allons partir d’un point, celui le plus à gauche et balayer les points par angle croissant (Fig. 1.3). Le
résultat de cette opération nous permet de former l’enveloppe convexe :

Choisir comme point L_0 de L le point le plus à gauche

Pour tous les points L_i de {L-L_0}:
Calculer l'angle (L_0,L_i) avec la verticale
Trier les L_i par angle

8
 F IGURE 1.3 – Tri des points Li

1.2.2 Deuxième solution

L’algorithme précédent n’est pas très performant. En effet, nous avions remarqué que le temps
de calcul nécessaire pour obtenir l’enveloppe convexe croît en O(n4 ). Nous allons construire un
autre algorithme basé sur une autre propriété de l’enveloppe convexe : le plus bas de P en fait
toujours partie (P8 dans l’exemple de la Fig. 1.2).
La première solution finissait par le tri des éléments de L par angle. Nous allons désormais
commencer par cette étape à partir du point le plus bas. Le résultat pour l’exemple est donné sur la
Fig. 1.4. Nous pouvons immédiatement observer que les deux points avec les angles les plus grands
sont également sur l’enveloppe (P0 et P6 ).
Pour résumer la première étape de l’algorithme est :
Définir L_0 comme le point de P le plus bas
Ajouter L_0 à L
Pour tous les points pi de {P-L_0}:
Calculer l'angle ai de (L_0,pi) avec la verticale
Trier l'ensemble P par angle
Définir L_1 et L_n à partir des plus grands angles

F IGURE 1.4 – Angles par rapport au point le plus bas

A la fin de cette étape, nous disposons d’un morceau (Ln , L0 , L1 ) de l’enveloppe. L’enveloppe
totale est obtenue en parcourant par angle croissant les points triés Pi de P. Le point est ajouté si
l’angle qu’il forme est ouvert avec l’enveloppe partielle, sinon il est retiré. Soit :

9
 Initialisation: L_j est L_1
Pour tous les points pi de {P-L_0-L_1}:
Calculer l'angle (L_(j-1),L_j,pi)
Si l'angle est fermé alors
Ajouter pi à L et L_j devient pi
sinon
Enlever L_j de L et L_j devient L_(j-1)
Recommencer tant que (L_(j-1),L_j,pi) est ouvert

Pour ce nouvel algorithme, nous n’avons que deux boucles imbriquées en plus du tri. Il sera
donc, dans le pire cas, en O(n2 ).

1.2.3 Les approximations

Les deux algorithmes précédents calculent toujours la solution exacte au problème proposé.
Néanmoins, dans de nombreux cas, une solution approchée mais plus rapide à obtenir est suffisante.
Reprenons l’exemple de la Fig. 1.2. Si nous découpons le plan en 3 bandes, nous pouvons
observer qu’une approximation de l’enveloppe convexe est obtenue en prenant le point le plus haut
et le plus bas de chaque bande (Fig. 1.5). Ce qui nous permet d’écrire l’algorithme suivant :

Découper le plan en M bandes

Pour tous les points pi de {P}:
Déterminer le maximum et le minimum de chaque bande
Déterminer le point le plus à droite et le plus à gauche
Ajouter les points trouvés à L
Pour tous les points L_i de {L}:
Enlever les points L_i qui forment un angle (Li-1,Li,Li+1) ouvert

F IGURE 1.5 – Approximation de l’enveloppe.

Nous n’avons ici qu’une seule boucle sur les points. La performance de cet algorithme sera
donc en O(n).

1.2.4 Conclusion
A partir d’un problème géométrique, nous avons pu voir trois algorithmes très différents
permettant sa résolution. Il en existe de nombreux autres beaucoup plus performants dans la
littérature scientifique avec des degrés de complexités variés.
Incidemment, nous n’avons pas encore parlé de langage ni d’implémentation. En effet, les
algorithmes que nous avons défini en font abstraction. Pour résumé, nous avons cherché à :

10
— Comprendre le problème posé et l’objectif à atteindre ;
— Trouver des approches de résolution ;
— Rechercher un algorithme respectant les contraintes de complexité (temps de calcul, mémoire)
et de précision (approximations) ;
— Puis il faudra implémenter la solution trouvée dans un langage adapté.

11
 Les briques de base des algorithmes en Matlab (Ref[1], [2])

2.1 Introduction à Matlab

Matlab est à la fois un langage de programmation et un environnement de développement. Ce
logiciel dispose d’une interface utilisateur que nous verrons dans la première partie de ce chapitre
avant d’aborder la syntaxe de son langage.
Le nom Matlab vient de l’anglais MATrix LABoratory. En effet, il a tout d’abord été crée
pour pouvoir traiter facilement les problèmes de calcul matriciel. Matlab est désormais largement
répandu et est utilisé dans de nombreux domaines comme :
— Le calcul numérique ;
— Les statistiques, les probabilités ;
— Le traitement du signal ;
— Le traitement d’image ;
— L’automatisme...
Le langage Matlab est facilement accessible au programmeur débutant. C’est un langage
interprété (c’est à dire qu’il n’y a pas d’étape de compilation) et faiblement typé (le type des
variables n’est pas explicitement donné). Les performances brutes d’un programme écrit en Matlab
sont globalement inférieures à d’autres langages compilés comme le C ou le Fortran. De plus,
même si la notion de type ou la gestion de la mémoire est caché à l’utilisateur, elles sont présentes
et il est utile d’en tenir compte.

2.2 Prise en main

Après le démarrage Matlab, son IDE apparaît. La Fig. 2.1 présente celle de la version R2017a.
IDE signifie Integrated Development Environment ou encore Environnement de Développement
Intégré en français. Elle se compose de plusieurs parties et nous allons décrire certaines.

2.2.1 La fenêtre de commande

Cette fenêtre permet d’accéder à l’interpréteur de commande de Matlab. Il permet d’exécuter
directement une fonction ou une opération. Pour cela, elle doit être écrite sur la ligne comportant
des chevrons (»). par exemple, on peut réaliser des opérations arithmétiques simples. Le résultat
sera donné directement sur la ligne suivante (sans chevrons) :
1 >> 30+2
2 ans =
3 32

La ligne de commande permet également d’utiliser directement des fonctions :

1 >> l o g ( 1 5 )
2 ans =
3 2.7081

Le nombre de chiffre affiché est indépendant de la précision du calcul. La précision utilisée par
Matlab est bien plus importante. Nous pouvons ainsi, par exemple, demander un affichage étendu
du résultat :

13
 F IGURE 2.1 – L’IDE de Matlab

1 >> f o r m a t l o n g
2 >> l o g ( 1 5 )
3 ans =
4 2.708050201102210

Les point-virgules à la fin d’une ligne de commande bloquent l’affichage du résultat.

2.2.2 Le répertoire de travail

L’IDE permet également de gérer les fichiers comportant des programmes ou des sauvegardes de
l’environnement. Cette partie de l’IDE comporte un explorateur de fichier et permet de sélectionner
le répertoire de travail actif.

2.2.3 L’éditeur
Matlab dispose d’un éditeur de texte intégré. Il permet d’écrire des programmes ou des fonctions
qui pourront être utilisé. Il permet de visualiser les mots clés du langages (comme par exemple
"for").
Cet éditeur permet également d’exécuter le programme écrit, de l’exécuter pas à pas ou d’insérer
des points d’arrêt dans le code. Ces deux dernières options sont utiles pour trouver les erreurs
de programmation en figeant le programme dans un état particulier. En observant le contenu des
variables, il est alors possible de comprendre ce qui cause le dysfonctionnement observé.

2.2.4 Le workspace
Cette partie de l’IDE affiche toutes les variables en mémoire de Matlab. Elle permet d’en
visualiser les noms, leur contenu et leur type.
Une variable est créé lors d’une affectation avec l’opérateur "=". Par exemple :

14
 1 >> a = 1

Nous pouvons effacer une variable avec la commande :

1 >> c l e a r a

Ou plus simplement effacer toutes les variables :

1 >> c l e a r

2.3 Les variables

Dès que le programme devient plus complexe qu’une opération simple, il devient nécessaire
de passer par des éléments intermédiaires qui retiennent la valeur d’une opération. Ce sont les
variables.
Une variable a toujours un type, c’est-à-dire qu’elle peut être une valeur entière, une valeur
décimale, un nombre complexe, une chaîne de caractère, etc. C’est toujours le cas même si cette
déclaration n’est pas explicite dans les langages faiblement typés (dont Matlab fait partie).
Avec Matlab, une variable est créé ou réaffectée avec l’opérateur "=". Par exemple, le code
suivant créera une variable de type double puis réaffectera cette variable dans un type chaîne de
caractère :

1 >> a = 1 ;
2 >> a = ’ t o t o ’ ;

Nous n’avons jamais indiqué que nous souhaitions que a soit du type double ou chaîne de caractère,
car c’est un langage faiblement typé. Matlab s’est chargé de toute cette partie à notre place. Il
indique le type d’une variable avec la commande "whos" :

1 >> a = 1 ;
2 >> whos ( ’ a ’ )
3 Name Size Bytes Class Attributes
4 a 1 x1 8 double

Le choix du nom d’une variable est laissé à l’utilisateur. Néanmoins, il y a quelques règles à
respecter pour celui-ci :
— Il peut avoir n’importe quelle taille ;
— Il peut contenir des lettres minuscules, majuscules et des chiffres ;
— Il ne doit pas contenir de caractères spéciaux (comme de la ponctuation) ;
— Il doit commencer par une lettre.
Il est de bon ton de choisir des noms parlants qui ne soient ni trop longs ni trop courts. On
peut utiliser des majuscules pour matérialiser le nom une variable composée de plusieurs mots
(masseDeLaFarine par exemple).
Nous allons maintenant voir plus en détails quelques types de variables gérés par Matlab.
L’ensemble des types disponibles peut être obtenu avec l’aide intégrée en utilisant la commande :

1 >> doc ’ D a t a Types ’

15
2.3.1 Les types numériques
Lorsqu’une variable numérique est créé sous Matlab, elle est automatiquement du type double.
C’est un format de donnée qui occupe 8 octets 1 (64 bits) et permet de représenter les nombres en
virgule flottante. Il comporte une fraction f composée de 52 bits b, un exposant e sur 11 bits et un
signe s sur 1 bit (Fig. 2.2). La valeur réelle mémorisée dans la variable v est alors obtenue grâce à
la formule suivante :
!
52
v = (−1)s 1 + ∑ b52−i 2−i × 2e−1023 (2.1)
i=1

F IGURE 2.2 – Le format double

Ce format permet donc de représenter tous les nombres positifs (et négatifs) compris entre
environ 2.22 · 10−308 et 1.8 · 10308 . Des codes spéciaux permettent d’encoder le +0, −0, −∞, +∞
et NaN (Not a Number, lors d’une division par 0 par exemple).
Matlab est capable de gérer de nombreux autres formats numériques : précision simple (single),
et les entiers signés (int) et non signés (uint) sur 8, 16, 32 et 64 bits. Ces format de données n’offrent
évidemment pas la versatilité du format double. Par exemple, le type uint8 permet uniquement de
gérer les entiers entre 0 et 255. Néanmoins, leur utilisation permet de gagner un espace mémoire
conséquent (1 byte pour un uint8 au lieu de 8 bytes pour un double) et de diminuer le temps de
calcul en facilitant les opérations pour le processeur. Nous pouvons néanmoins remarquer que la
plupart des fonctions intégrées en Matlab sont écrites pour des entrées de type double, ce qui peut
induire un temps d’exécution supérieur pour les autres formats de données.
Le changement de type d’une variable fait appel au nom du type correspondant. Par exemple,
pour changer une variable a de double précision à single précision puis en int16 (entier sur 16 bits),
nous pouvons écrire le code suivant :
1 >> a = 1 0 ;
2 >> whos ( ’ a ’ )
3 Name Size Bytes Class Attributes
4 a 1 x1 8 double
5 >> a = s i n g l e ( a ) ;
6 >> whos ( ’ a ’ )
7 Name Size Bytes Class Attributes
8 a 1 x1 4 single
9 >> a = i n t 1 6 ( a ) ;
10 >> whos ( ’ a ’ )
11 Name Size Bytes Class Attributes
12 a 1 x1 2 int16

Nous noterons également que Matlab connait un type complexe. Celui-ci est créé par l’utilisation
de la forme a + ib ou du mot clé "complex". Par exemple :
1 >> a = 10 i + 1 ;
1. bytes en anglais

16
 2 >> whos ( ’ a ’ )
3 Name Size Bytes Class Attributes
4 a 1 x1 16 double complex
5 >> a = complex ( 1 0 , 1 ) ;
6 >> whos ( ’ a ’ )
7 Name Size Bytes Class Attributes
8 a 1 x1 16 double complex

Nous pouvons résumer les différents types numériques dans le tableau suivant :
Type Bytes Valeur min Valeur max Pas minimal
double 8 ≈ −1.8 · 10308 ≈ 1.8 · 10308 ≈ 2.22 · 10−308
single 4 ≈ −3.4 · 1038 ≈ 3.4 · 1038 ≈ 1.2 · 10−38
int8 1 -128 127 1
int16 2 -32768 32767 1
int32 4 -2147483648 2147483647 1
int64 8 -9223372036854775808 9223372036854775807 1
uint8 1 0 255 1
uint16 2 0 65535 1
uint32 4 0 4294967295 1
uint64 8 0 18446744073709551615 1
Les formats flottants (single et double) induisent des erreurs inhérentes à ces représentations. Nous
pouvons voir, par exemple, dans l’exemple suivant que le nombre 3.98 est mémorisé avec une
erreur de l’ordre de 10−6 dans le format single :
1 >> f o r m a t l o n g
2 >> a = s i n g l e ( 3 . 8 9 )
3 a =
4 single
5 3.8900001

2.3.2 Le type logique

Le type de donnée logical est particulièrement simple car il n’a que deux valeur permises : 0 ou
1 (ou true et false). Il est utilisé, par exemple, lors de d’opération de comparaison. Par exemple :
1 >> a = f a l s e
2 ans =
3 logical
4 0

2.3.3 Les tableaux et matrices

Matlab est avant tout un langage qui a été développé pour traiter et manipuler des tableaux
et des matrices. La plupart des fonctions intégrées sont spécifiquement créés pour pouvoir être
appliquées à ce type de donnée. Cette notion est au cœur de Matlab et nous pouvons même observer,
dans les exemples précédents, que les variables ont étés définies comme des tableaux 1 × 1.
Notons qu’une matrice est un cas particulier de la notion de tableau. Seul cette première a un
sens mathématique. Néanmoins, comme Matlab ne fait pas la différence, nous parlerons uniquement
des tableaux.
Pour créer un tableau, nous pouvons utiliser l’opérateur d’affectation
 "=". Les colonnes sont
séparées par un espace " ". Par exemple, pour créer le tableau t = 1 2 3 , nous pouvons écrire :

17
 1 >> t = [ 1 2 3 ]
2 t =
3 1 2 3

Les éléments sur

 différentes
 lignes sont séparées par un point virgule (" ;"). Par exemple, pour créer
1 2 3
le tableau t = , nous pouvons écrire :
4 5 6
1 >> t = [ 1 2 3 ; 4 5 6 ]
2 t =
3 1 2 3
4 4 5 6

Matlab permet également de construire

 une matrice en concaténant des matrices existantes. Par
A B
exemple pour obtenir, M = , nous pouvons écrire :
C D
1 >> A = [ 1 1 1 ; 1 1 1 ] ;
2 >> B = [ 2 ; 2 ] ;
3 >> C = [ 3 3 3 ] ;
4 >> D = [ 4 ] ;
5 >> M = [A B ; C D]
6 M =
7 1 1 1 2
8 1 1 1 2
9 3 3 3 4

Matlab facilite grandement la manipulation de tableaux avec ses fonctions de manipulation. On

peut aisément accéder ou modifier un élément en indiquant son indice entre parenthèse. On notera
que le premier élement d’un tableau a l’indice 1. Par exemple :
1 >> M = [ 1 2 3 ] ;
2 >> M( 2 )
3 ans =
4 2
5 >> M = [ 4 ; 5 ; 6 ] ;
6 >> M( 2 )
7 ans =
8 5
9 >> M = [ 1 2 3 ; 4 5 6 ] ;
10 >> M( 2 , 1 )
11 ans =
12 4
13 >> M( 2 , 1 ) = 7 ;
14 >> M
15 M =
16 1 2 3
17 7 5 6

Le " :" permet de spécifier un ensemble d’indices que l’on souhaite afficher ou modifier. Il peut
s’écrire sous la forme indiceMin:indiceMax ou indiceMin:pas:indiceMax. Par exemple :
1 >> M( 2 , 1 : 2 )

18
 2 ans =
3 7 5
4 >> M = [ 1 2 3 4 5 6 7 8 ] ;
5 >> M( 1 : 7 )
6 ans =
7 1 2 3 4 5 6 7
8 >> M( 1 : 2 : 7 )
9 ans =
10 1 3 5 7
11 >> M( 1 : 2 : 7 ) = 1 0 ;
12 >> M
13 M =
14 10 2 10 4 10 6 10 8

Pour terminer, Matlab possède de nombreuses fonctions qui permettent de générer rapidement un
tableau. Nous nous limiterons ici à ces quelques exemples :
— zeros(N,M) : Crée un tableau de 0 de taille N × M ;
— ones(N,M) : Crée un tableau de 1 de taille N × M ;
— eye(N) : Crée un tableau de 0 avec des 1 sur la diagonale de taille N × N ;
— rand(N,M) : Crée un tableau de taille N × M rempli de nombres aléatoires compris entre 0
et 1 ;
— L’opérateur "N :P :M" : Crée un tableau de taille (N − M)/P rempli des points entre N et M
avec un pas P.
Comme vous le verrez en étudiant d’autres langages, ces fonctions de créations et de manipula-
tion de tableaux sont particulièrement puissantes et sont une des raisons du succès de Matlab.

2.3.4 Les chaînes de caractères

En Matlab, les chaînes de caractères sont simplement des tableaux de caractères. Elles se
déclarent grâce à des simples quotes. Les fonctions utilisables sur des tableaux le sont donc
également sur des chaînes de caractère. Par exemple :
1 >> c = ’ A m i c a l e m e n t ’ ;
2 >> c ( 3 )
3 ans =
4 i
5 >> c ( 3 : 7 )
6 ans =
7 icale

2.3.5 Les structures

Le dernier type donnée que nous étudierons ici est le type structure. C’est un type composite
qui permet d’agréger plusieurs types de données dans le même conteneur. Elles se définissent avec
le mot clé "struct". En voici quelques exemples :
1 >> s = s t r u c t ;
2 >> s . nom = ’ A l b e r t ’ ;
3 >> s . a g e = 3 4 ;
4 >> s
5 s =
6 s t r u c t with f i e l d s :

19
 7 nom : ’ A l b e r t ’
8 a g e : 34
9 >> b = s t r u c t ( ’nom ’ , { ’ A l b e r t ’ , ’ Rene ’ } , ’ a g e ’ , { 1 2 , 3 4 } )
10 b =
11 1 x2 s t r u c t a r r a y w i t h f i e l d s :
12 nom
13 age
14 >> b ( 1 )
15 ans =
16 s t r u c t with f i e l d s :
17 nom : ’ A l b e r t ’
18 a g e : 12

2.4 Les opérateurs

Nous avons vu que Matlab traite de la même façon les scalaires et les données sous forme
tableaux. Nous allons étudier dans cette partie quelques unes des opérations mathématiques dispo-
nibles.

2.4.1 L’addition et la soustraction

A partir de deux tableaux A et B de taille identique d’éléments Ai j et Bi j , l’opérateur addition
("+") et soustraction ("-") crée un nouveau tableau C de même taille et d’éléments Ci j = Ai j + Bi j
(resp. Ci j = Ai j − Bi j ). Par exemple :
1 >> A = [ 1 2 ] ;
2 >> B = [ 3 4 ] ;
3 >> C = A+B
4 C =
5 4 6

2.4.2 Multiplication, division et puissance

Dans le cas de la multiplication, division et puissance, on dispose d’opérateurs différent suivant
que nous souhaitons la faire terme à terme ou au sens des matrices. Dans le cas d’opérations terme à
terme, les opérateurs s’écrivent .∗ (multiplication), ./ (division) et .ˆ (puissance). Nous obtiendrons
ainsi, par exemple :
1 >> A = [ 1 2 ] ;
2 >> B = [ 1 4 ] ;
3 >> C = A. ∗ B
4 C =
5 1 8
6 >> C = B . / A
7 C =
8 1 2
9 >> C = A. ^ 2
10 C =
11 1 4

Pour résumer, pour deux tableaux de même taille A et B, lorsque :

20
 — C = A. ∗ B alors Ci j = Ai j Bi j
— C = A./B alors Ci j = Ai j /Bi j
— C = A.ˆn alors Ci j = Anij
B
— C = A.ˆB alors Ci j = Ai ji j
Ces opérateurs existent également au sens des matrices. Il s’écrivent alors ∗ (multiplication), /
(division) et ˆ (puissance). En résumé, lorsque :
— C = A ∗ B alors C = ∑k Ai j Bi j
— C = A/B alors C = AB−1
— C = A\B alors C = A−1 B
— C = Aˆn alors C = An

2.4.3 Opérations relationnelles

Les opérations relationnelles permettent de faire des comparaisons logiques entre deux éléments.
Ces éléments peuvent être des valeurs numériques ou des tableaux. Les opérateurs relationnels
renvoient une valeur de type booléen, c’est à dire la valeur true (1) lorsque la condition est
vérifiée et false (0) dans le cas contraires. Les principaux opérateurs disponibles sont, pour deux
valeurs/tableaux de même taille A et B :
— L’égalité, A == B ;
— La différence A ∼= B ;
— La supériorité A > B ;
— La supériorité ou l’égalité A >= B ;
— L’infériorité A < B ;
— L’infériorité ou l’égalité A <= B.
Prenons quelques exemples :
1 >> 9 == 9
2 ans =
3 logical
4 1
5 >> 9 > 12
6 ans =
7 logical
8 0
9 >> [ 1 2 3 4 ] <= [ 2 1 3 8 ]
10 ans =
11 1 x4 l o g i c a l a r r a y
12 1 0 1 1
13 >> 2 == [ 1 2 3 ]
14 ans =
15 1 x3 l o g i c a l a r r a y
16 0 1 0

Notons que, puisque les valeurs flottantes ne sont que des approximations, la plus simple des
opérations de comparaison peut devenir suspecte [6]. Reprenons l’exemple introductif :
1 >> f o r m a t l o n g
2 >> a = s i n g l e ( 3 . 8 9 )
3 a =
4 single
5 3.8900001

21
 Ainsi, si nous comparons la variable a à 3.89, le résultat sera faux (la précision de la comparaison
est double) :
1 >> a == 3 . 8 9
2 ans =
3 logical
4 0

Nous pouvons également observer le même phénomène pour le résultat de fonctions mathéma-
tiques :
1 >> s i n ( p i ) == 0
2 ans =
3 logical
4 0

Une solution pour surmonter cette difficulté est de définir un δ petit et de vérifier la condition
suivante |a − 3.89| < δ . La solution n’est pas toujours aussi simple. Par exemple considérons
trois points (a,b) (c,d) (e,f), nous pouvons déterminer leur position en fonction du signe de s =
(c − a) · ( f − b) − (d − b) · (e − a) :
— Ils sont colinéaires si s = 0 ;
— Ils tournent vers la gauche si s < 0 ;
— Ils tournent vers la droite si s > 0.
Pour cet exemple, nous prendrons trois points sur la droite y = 5x. Observons ce que donne le
résultat de ce calcul en Matlab en précision single et double :
1 >> a = 1 / 3 ; b = 5 / 3 ; c = 3 3 ; d = 1 6 5 ; e = 1 9 ; f = 9 5 ;
2 >> ( c−a ) ∗ ( f−b ) −(d−b ) ∗ ( e−a )
3 ans =
4 −4.5475 e −13
5 >> a= s i n g l e ( 1 / 3 ) ; b= s i n g l e ( 5 / 3 ) ; c= s i n g l e ( 3 3 ) ; d= s i n g l e ( 1 6 5 ) ; e=
single (19) ; f= single (95) ;
6 >> ( c−a ) ∗ ( f−b ) −(d−b ) ∗ ( e−a )
7 ans =
8 single
9 4 . 8 8 2 8 e −04

Nous avons donc une conclusion différente (et fausse) dans les deux cas. Ces exemples montrent
qu’il faut faire attention lors des calculs en virgule flottante.

2.4.4 Opérations logiques

Les valeurs du type logical ont leur propre jeu d’opérateur. Ils permettent de combiner des
conditions logiques. Nous disposons des opérateurs suivants pour deux valeurs logiques ou deux
tableaux logiques de même taille A et B :
— L’opérateur Et qui s’écrit A&B ou and(A, B) ;
— L’opérateur Ou qui s’écrit A|B ou or(A, B) ;
— L’opérateur Ou exclusif qui s’écrit xor(A, B) ;
— L’opérateur Non qui s’écrit ∼ A ou not(A).
Par exemple :
1 >> [ 1 0 1 ] | [ 1 1 0 ]
2 ans =

22
 3 1 x3 l o g i c a l a r r a y
4 1 1 1
5 >> 0 & 1
6 ans =
7 logical
8 0
9 >> [ 1 0 1 ] & 0
10 ans =
11 1 x3 l o g i c a l a r r a y
12 0 0 0

Les opérateurs Et et Ou disposent en sus de versions court-circuit qui s’écrivent respectivement &&
et ||. Lorsque on utilise ces opérateurs, la première condition est d’abord évaluée puis la seconde si
la première condition le permet. Nous pouvons ainsi, par exemple, vérifier qu’une opérande est non
nulle avant une division par 0. Par exemple :
1 >> s = ( a ~= 0 ) && ( b / a > 1 )

2.4.5 Quelques autres fonctions

Nous aurons besoin de quelques fonctions supplémentaires dans ce cours. En voici une liste
non exhaustive :
— La transposition, elle inverses les ligne et les colonnes d’une matrice. Son opérateur est "’".
Par exemple B = A0 ;
— La taille d’un tableau. Elle est obtenue grâce à la fonction size ;
— La longueur d’un tableau suivant sa plus grande dimension est obtenue grâce à la fonction
length ;
— Le redimensionnement grâce à la fonction reshape. Cette opération permet de modifier un
tableau de taille N × M en un tableau de taille N 0 × M 0 si NM = N 0 M 0 ;
— L’affichage. La fonction disp permet d’afficher une chaîne de caractère.
Par exemple :
1 >> A = [ 1 2 3 4 5 6 ]
2 A =
3 1 2 3 4 5 6
4 >> B = r e s h a p e (A, 2 , 3 )
5 B =
6 1 3 5
7 2 4 6
8 >> B’
9 ans =
10 1 2
11 3 4
12 5 6
13 >> s i z e ( B )
14 ans =
15 2 3
16 >> l e n g t h ( B )
17 ans =
18 3

23
2.5 Les fichiers m
Depuis le début de ce chapitre, nous avons toujours utilisé la fenêtre de commande pour créer les
variables et les modifier. C’est une approche pratique pour démarrer et tester des idées rapidement
mais elle n’est pas adaptée à la réalisation de programmes.
Pour cela Matlab dispose de fichiers de commandes. Ce sont des fichiers textes portant l’exten-
sion .m. Ils comportent des commandes Matlab identiques à celles qui peuvent être utilisées dans la
fenêtre de commande. Les commentaires sont indiqués grâce au caractère %. Prenons l’exemple
suivant, enregistré dans le fichier Prog1.m :
1 %% P r o g 1 .m
2 % Exemple de f i c h i e r m
3

4 A = [ 1 2 3 ] ; % d e f i n i t i o n de A
5 B = A − 1; % d e f i n i t i o n de B
6 disp ( ’ Resultat ’ )
7 B == 2

Les commentaires en début du fichier sont affichés lors de de l’appel de la fonction help :
1 >> h e l p p r o g 1
2 P r o g 1 .m
3 Exemple de f i c h i e r m

Le programme est appelé en utilisant son nom, sans l’extension .m, dans l’interface de commande :
1 >> p r o g 1
2 Resultat
3 ans =
4 1 x3 l o g i c a l a r r a y
5 0 0 1

Par défaut, le programme est cherché dans le répertoire actuel, celui défini dans la fenêtre current
folder. Nous pouvons ajouter d’autres répertoires de recherche avec la fonction addpath. La façon
d’indiquer ces répertoires dépend du système d’exploitation. Par exemple sous Microsoft Windows :
1 >> a d d p a t h ( ’C : \ Documents and S e t t i n g s \ A l b e r t \ a l g o \ m a t l a b \ ’ )

2.6 Les fonctions

Les fonctions en Matlab sont des fichiers commandes particuliers. Ce sont des programmes
qui permettent de réaliser des opérations répétitives et d’éviter de recourir au copier-coller. Elles
permettent ainsi d’avoir une plus grand ré-utilisabilité des codes développés.
Les informations fournies à la fonction sont appelées les arguments ou les paramètres de la
fonction. Les informations retournées par la fonctions sont appelées le résultat.
En Matlab, une fonction utilise :
— Le mot clé function ;
— Un nom de fonction (qui doit être identique au nom du fichier .m) ;
— Des paramètres d’entrées i1 , i2 , ..., in (ces paramètres doivent être présents lors de l’appel de
la fonction), ils ne peuvent pas être modifiées par la fonction ;
— Un ensemble de résultats o1 , o2 , ..., on (seul le premier résultat est obligatoirement retourné) ;
— La fonction se termine par le mot clé end.
Par exemple, nous pouvons définir une fonction testfun comme suit :

24
 1 f u n c t i o n [ o1 , o2 ] = t e s t f u n ( i 1 , i 2 )
2 %TESTFUN F o n c t i o n de t e s t
3 % Exemple p o u r l e c o u r s
4 tmp = i 1 + i 2 ;
5

6 o1 = tmp ;
7 o2 = struct ;
8 o2 . somme = tmp ;
9 o2 . d i f f e r e n c e = i 1 − i 2 ;
10 end

Notons que la variable tmp utilisée dans l’exemple n’est visible que dans le cadre de la fonction
testfun. C’est une variable locale, elle est détruite à la fin de l’appel de la fonction. Elle n’existe que
dans ce contexte particulier. Si une autre partie du programme utilise une variable de même nom,
elles sont traitées comme des variables différentes.
Voyons quelques exemples d’appel à cette fonction :

1 >> [ a , b ] = t e s t f u n ( 1 , 2 )
2 a =
3 3
4 b =
5 s t r u c t with f i e l d s :
6 somme : 3
7 d i f f e r e n c e : −1
8 >> a = t e s t f u n ( 5 , 3 )
9 a =
10 8
11 >> tmp
12 U n d e f i n e d f u n c t i o n o r v a r i a b l e ’ tmp ’ .

De la même manière qu’il existe des variables locales à une fonction, il existe des sous-fonctions
qui ne sont visibles que dans le contexte de la fonction mère. Ce sont des fonctions décrites dans le
même fichier mais avec un nom différent. Par exemple :

1 f u n c t i o n [ o1 , o2 ] = t e s t f u n 2 ( i 1 , i 2 )
2 %TESTFUN2 F o n c t i o n de t e s t
3 % Exemple p o u r l e c o u r s
4 o1 = mysomme ( i 1 , i 2 ) ;
5 o2 = struct ;
6 o2 . somme = mysomme ( i 1 , i 2 ) ;
7 o2 . d i f f e r e n c e = i 1 − i 2 ;
8 end
9

10 f u n c t i o n s = mysomme ( i 1 , i 2 )
11 s = i1+i2 ;
12 end

25
 2.7 Les structures conditionnelles
Nous avons vu, dans la partie précédente, des fonctions comportant une suite d’assignements.
Néanmoins, nous souhaitons pouvoir contrôler le comportement d’un programme en fonction des
valeurs de variables ou de paramètres. Pour cela, nous avons besoin de structures spécifiques. Elles
se présentent sous la forme de blocs se terminant par le mot clé end.

2.7.1 Instruction "If" (Si)

L’instruction Si utilise le mot clé if. Il doit être suivi d’un type logique (une variable ou une
condition). Le bloc d’instructions sera uniquement exécuté si la condition vaut 1 ou true. Lorsque
la condition vaut 0 (ou false), d’autre blocs peuvent être définis avec le mot clé elseif (suivi d’une
autre condition logique) ou du mot clé else si toutes les conditions sont fausses.
Par exemple :
1 i f a >= 5
2 d i s p ( ’ Bloc 1 ’ ) ;
3 elseif b < 3
4 d i s p ( ’ Bloc 2 ’ ) ;
5 else
6 d i s p ( ’ Bloc 3 ’ ) ;
7 end

Nous pouvons remarquer que le bloc répondant à la condition b < 3 ne pourra être exécuté que si
a < 5.

2.7.2 Instruction "Switch"

La structure switch est également une structure conditionnelle, de même nature que le if. La
différence principale est que, dans ce cas, c’est la valeur d’une variable qui est testée. Le mot clé
case permet de définir le bloc à exécuter et le mot clé otherwise le comportement par défaut.
Par exemple :
1 switch a
2 case {1 ,3 ,5}
3 d i s p ( ’ Bloc 1 ’ ) ;
4 case 7
5 d i s p ( ’ Bloc 2 ’ ) ;
6 otherwise
7 d i s p ( ’ Bloc 3 ’ ) ;
8 end

2.8 Les boucles

Les boucles permettent d’itérer sur certaines instructions du programmes. Nous allons en voir
deux, l’instruction for et l’instruction while.

2.8.1 Instruction "For" (Pour)

L’instruction for permet d’exécuter un bloc d’instruction un nombre fini de fois, défini en début
de structure. Elle se compose du mot clé for suivi d’une liste de valeurs. Elle est souvent utilisée
sous la forme :

26
 for variable = début:pas:fin
bloc d'instructions
end
Nous pouvons noter qu’en Matlab la liste de valeur n’est pas forcément une liste de valeurs entières
mais peut être de n’importe quel type. Par exemple :
1 f o r i =1:10
2 disp ( i )
3 end
4

5 for i =0.2:0.3:1.4
6 disp ( i )
7 end
8

9 for i = { ’ toto ’ , ’ tata ’}

10 disp ( i )
11 end

2.8.2 Instruction "While" (Tant que)

La boucle while permet d’exécuter un bloc d’instructions tant qu’une condition est vraie, sans
connaître a priori le nombre de fois où elle sera exécutée. Elle se compose du mot clé while suivi
d’une condition logique. Par exemple :
1 while t e s t <0.2
2 test = test ∗0.1;
3 end

2.9 Exemple
Nous avons désormais tous les outils pour implémenter un algorithme. Prenons un des algo-
rithmes de détermination de l’enveloppe convexe du premier chapitre :

Découper le plan en M bandes

Pour tous les points pi de {P}:
Déterminer le maximum et le minimum de chaque bande
Déterminer le point le plus à droite et le plus à gauche
Ajouter les points trouvés à L
Pour tous les points L_i de {L}:
Enlever les points L_i qui forment un angle (Li-1,Li,Li+1) ouvert
Nous supposerons que les points pi de P sont présents sous la forme d’un tableau de structure
avec deux entrées doubles x et y. L’enveloppe des points li sera retournée sous la même forme.
Les coordonnées des bandes sont des paramètres de la fonction. Finalement, pour simplifier, nous
supposerons qu’il y a des points dans chaque bande. Nous pouvons alors écrire l’implémentation
suivante dans le fichier convexHullApprox.m :
1 f u n c t i o n [ hull , h u l l t , pointsMin , pointsMax ] = convexHullApprox (
points , bandes )
2 %% F o n c t i o n de d e t e r m i n a t i o n de l ’ e n v e l o p p e c o n v e x e
3 % U t i l i s e des approximations

27
 4 % E n t r e e : p o i n t s : T a b l e a u de s t r u c t u r e s de p o i n t s ( x , y )
5 % b a n d e s : Bandes p o u r l ’ a p p r o x i m a t i o n ( t a b l e a u de
v a l e u r s de x )
6 % S o r t i e : h u l l : T a b l e a u s t r u c t u r e s de p o i n t s ( x , y ) de l ’
e n v e l l o p e convexe
7

9 %% I n i t i a l i s a t i o n
10 % Points
11 pointDroite = points (1) ;
12 pointGauche = p o i n t s (1) ;
13 pointHaut = points (1) ;
14 pointBas = points (1) ;
15

16 % D e t e r m i n a t i o n de l ’ e n v e l o p p e max
17 for unPoint = points
18 i f unPoint . x > pointDroite . x
19 pointDroite = unPoint ;
20 end
21 i f unPoint . x < pointGauche . x
22 pointGauche = unPoint ;
23 end
24 i f unPoint . y > pointHaut . y
25 pointHaut = unPoint ;
26 end
27 i f unPoint . y < pointBas . y
28 pointBas = unPoint ;
29 end
30 end
31

32 % I n i t i a l i s a t i o n de l ’ e n v e l o p p e
33 pointsMax = s t r u c t ( ’x ’ ,{} , ’y ’ ,{}) ;
34 pointsMin = s t r u c t ( ’x ’ ,{} , ’y ’ ,{}) ;
35 nBandes = l e n g t h ( bandes ) +1;
36 f o r i = 1 : nBandes
37 pointsMax ( i ) = pointBas ;
38 pointsMin ( i ) = pointHaut ;
39 end
40

41 % C a l c u l d e s p o i n t s min e t max de c h a q u e b a n d e
42 for unPoint = points
43 i f unPoint . x < bandes ( 1 )
44 bandeC = 1 ;
45 e l s e i f u n P o i n t . x >= b a n d e s ( nBandes −1)
46 bandeC = nBandes ;
47 else
48 f o r i = 1 : nBandes −2
49 i f u n P o i n t . x >= b a n d e s ( i ) && u n P o i n t . x < b a n d e s ( i + 1 )
50 bandeC = i + 1 ;

28
51 end
52 end
53 end
54 i f p o i n t s M a x ( bandeC ) . y < u n P o i n t . y
55 p o i n t s M a x ( bandeC ) = u n P o i n t ;
56 end
57 i f p o i n t s M i n ( bandeC ) . y > u n P o i n t . y
58 p o i n t s M i n ( bandeC ) = u n P o i n t ;
59 end
60 end
61

62 % C o n s t r u c t i o n de l ’ e n v e l o p p e
63 i f ( ptComp ( p o i n t s M i n ( 1 ) , p o i n t s M a x ( 1 ) ) )
64 p o i n t s M i n = p o i n t s M i n ( 2 : nBandes −1) ;
65 end
66 i f ( ptComp ( p o i n t s M i n ( l e n g t h ( p o i n t s M i n ) ) , p o i n t s M a x ( nBandes −1) ) )
67 p o i n t s M i n = p o i n t s M i n ( 1 : nBandes −2) ;
68 end
69 i f ( ptComp ( p o i n t G a u c h e , p o i n t s M a x ( 1 ) ) | | ptComp ( p o i n t G a u c h e ,
pointsMin (1) ) )
70 h u l l t = [ pointsMax ] ;
71 else
72 h u l l t = [ pointGauche , pointsMax ] ;
73 end
74 i f ( ptComp ( p o i n t D r o i t e , p o i n t s M a x ( nBandes −1) ) | | . . .
75 ptComp ( p o i n t D r o i t e , p o i n t s M i n ( l e n g t h ( p o i n t s M i n ) ) ) )
76 h u l l t = [ h ul l t , f l i p l r ( pointsMin ) ] ;
77 else
78 h u l l t = [ h ul l t , pointDroite , f l i p l r ( pointsMin ) ] ;
79 end
80

81 % On ne g a r d e que l e s p o i n t s f o r m a n t l ’ e n v e l o p p e c o n v e x e
82 hull = hullt (1) ;
83 ip = 1;
84 ih = 1;
85 w h i l e ( i h <= l e n g t h ( h u l l t ) )
86 j = ih +1;
87 k = j +1;
88 if j > length ( hullt )
89 j = j −l e n g t h ( h u l l t ) ;
90 end
91

92 if k > length ( hullt )

93 k = k−l e n g t h ( h u l l t ) ;
94 end
95 v1vectv2 = vectprod ( h u l l t ( ip ) , h u l l t ( j ) , h u l l t ( k ) ) ;
96 i f v 1 v e c t v 2 > 0 % Le p o i n t e s t s u r l ’ e n v e l o p p e c o n v e x e
97 ip = j ;
98 % On v e r i f i e que l e p o i n t a j o u t e p r e c e d e m m e n t e s t b i e n

29
99 % s u r l ’ e n v e l o p p e convexe , s i n o n on l ’ e n l e v e
100 lh = length ( hull ) ;
101 i f l h >=2
102 v 1 v e c t v 2 = v e c t p r o d ( h u l l ( l h −1) , h u l l ( l h ) , h u l l t ( j ) ) ;
103 w h i l e v 1 v e c t v 2 < 0 && l h >= 2
104 h u l l = h u l l ( 1 : l h −1) ;
105 lh = lh − 1;
106 end
107 end
108 hull = [ hull hullt ( j ) ];
109 end
110 ih = ih + 1;
111 end
112

113 end
114

115 %% F o n c t i o n de c o m p a r a i s o n de deux p o i n t s
116 f u n c t i o n r e s = ptComp ( p t 1 , p t 2 )
117 r e s = ( ( p t 1 . x == p t 2 . x ) && ( p t 1 . y == p t 2 . y ) ) ;
118 end
119

120 %% F o n c t i o n p e r m e t t a n t de d e t e r m i n e r s i l ’ a n g l e e s t o u v e r t ou
ferme
121 f u n c t i o n r e s = v e c t p r o d ( pt1 , pt2 , pt3 )
122 v1 . x= ( p t 1 . x−p t 2 . x ) ;
123 v1 . y= ( p t 1 . y−p t 2 . y ) ;
124 v2 . x= ( p t 3 . x−p t 2 . x ) ;
125 v2 . y= ( p t 3 . y−p t 2 . y ) ;
126 r e s = v1 . x∗ v2 . y−v1 . y ∗ v2 . x ;
127 end

A la vue de cette implémentation, nous pouvons voir que l’algorithme précédent n’était que grossiè-
rement détaillé car de nombreux cas particulier n’étaient pas traités. Néanmoins, ce programme a
toujours la performance de l’algorithme.
Pour tester ce programme, nous allons générer 15 points aléatoirement répartis dans [0,10]
dans quatre bandes [0,2[, [2,5[, [5,8[ et [8,10[. Les résultats des différentes étapes seront affichés
dans des figures. Le code ci-dessous réalise ces opérations. Ce dernier code utilise de nombreuses
fonctions que nous n’avons pas vu dans ce chapitre pour l’affichage de courbes.

1 vx = r a n d ( 1 , 1 5 ) ∗ 1 0 ;
2 v x l = m a t 2 c e l l ( vx , 1 , o n e s ( 1 , numel ( vx ) ) ) ;
3 vy = r a n d ( 1 , 1 5 ) ∗ 1 0 ;
4 v y l = m a t 2 c e l l ( vy , 1 , o n e s ( 1 , numel ( vy ) ) ) ;
5

6 mespoints = s t r u c t ( ’ x ’ , vxl , ’ y ’ , vyl ) ;

7 mesbandes = [2 5 8 ] ;
8

9 [ h u l l , h u l l t , ptMin , ptMax ] = c o n v e x H u l l A p p r o x ( m e s p o i n t s , m e s b a n d e s )
;
10

30
11 figure ;
12 a x i s ( [ 0 10 0 1 0 ] ) ;
13 h o l d on ;
14 f o r i = 1 : l e n g t h ( ptMin )−1
15 l i n e ( [ ptMin ( i ) . x , ptMin ( i + 1 ) . x ] , [ ptMin ( i ) . y , ptMin ( i + 1 ) . y ] , ’
color ’ , ’ r ’ )
16 end
17 f o r i = 1 : l e n g t h ( ptMax )−1
18 l i n e ( [ ptMax ( i ) . x , ptMax ( i + 1 ) . x ] , [ ptMax ( i ) . y , ptMax ( i + 1 ) . y ] , ’
color ’ , ’g ’ )
19 end
20 s c a t t e r ( vx , vy )
21 f o r i = 1 : l e n g t h ( mesbandes )
22 l i n e ( [ mesbandes ( i ) , mesbandes ( i ) ] , [ 0 , 1 0 ] )
23 end
24 t i t l e ( ’ P o i n t s minimum e t maximum de c h a q u e b a n d e ’ ) ;
25

26

27 figure ;
28 a x i s ( [ 0 10 0 1 0 ] ) ;
29 h o l d on ;
30 f o r i = 1 : l e n g t h ( h u l l t )−1
31 l i n e ( [ h u l l t ( i ) . x , h u l l t ( i +1) . x ] , [ h u l l t ( i ) . y , h u l l t ( i +1) . y ] )
32 end
33 line ([ hullt (1) . x , hullt ( length ( hullt ) ) . x ] ,[ hullt (1) . y , hullt (
length ( hullt ) ) . y ])
34 s c a t t e r ( vx , vy )
35 f o r i = 1 : l e n g t h ( mesbandes )
36 l i n e ( [ mesbandes ( i ) , mesbandes ( i ) ] , [ 0 , 1 0 ] )
37 end
38 t i t l e ( ’ Enveloppe avec t o u s l e s p o i n t s ’ ) ;
39

40 figure ;
41 a x i s ( [ 0 10 0 1 0 ] ) ;
42 h o l d on ;
43 f o r i = 1 : l e n g t h ( h u l l )−1
44 l i n e ( [ h u l l ( i ) . x , h u l l ( i +1) . x ] , [ h u l l ( i ) . y , h u l l ( i +1) . y ] )
45 end
46 line ([ hull (1) . x , hull ( length ( hull ) ) . x ] ,[ hull (1) . y , hull ( length (
hull ) ) . y ])
47 s c a t t e r ( vx , vy )
48 f o r i = 1 : l e n g t h ( mesbandes )
49 l i n e ( [ mesbandes ( i ) , mesbandes ( i ) ] , [ 0 , 1 0 ] )
50 end
51 t i t l e ( ’ Enveloppe convexe ’ ) ;

La Fig. 2.3 montre les différentes étapes d’exécution de l’algorithme. Tout d’abord, il cherche
les points maximum et minimum de chaque bande. Les points sont regroupés en une seule enveloppe.
Finalement, l’enveloppe convexe est construite en ne gardant que les angles ouverts.

31
F IGURE 2.3 – Approximation de l’enveloppe.

32
 Architecture des ordinateurs (Ref. [3], [4])

3.1 Introduction
Les ordinateurs ne comprennent pas directement le code Matlab. Les fichiers .m comportent une
liste d’instructions que nous aimerions voir réalisée par la machine. Pour pouvoir être exécutés, les
programmes doivent être sous la forme d’un code machine qui pourra être interprété par l’élément
central : le micro-processeur.
Le code machine n’est pas unique et dépend du l’achitecture du système sur lequel il va être
utilisé (c’est pourquoi, par exemple, les programmes exécutables sur un Mac ne le sont pas sur un
PC). Ce code n’est pas très lisible comme nous pouvons le voir sur cette exemple réalisé pour un
ordinateur du type "Little Man Computer" :
1 901
2 334
3 901
4 134
5 902
6 000

Il peut être mit sous une forme plus lisible, appelé langage assembleur. Par exemple, le code
précédent peut se mettre sous la forme :
1 INP
2 STA 34
3 INP
4 ADD 34
5 OUT
6 HLT

Il est probable que la lecture de celui-ci ne vous aide pas beaucoup. Nous avons désormais à
disposition des langages évolués qui nous permettent, après une étape de compilation, de générer un
code assembleur. Matlab utilise un langage interprété, c’est à dire que le programme est transformés
en une suite d’instruction qui est ensuite digérée par une machine logicielle afin de finir en langage
machine. Ces langages nous permettent d’obtenir une forme facilement lisible pour nous. Par
exemple, le code précédent peut se mettre sous la forme :
1 a = i n p u t ( ’ ’ ) ; % Demande un p r e m i e r nombre
2 d i s p ( a+ i n p u t ( ’ ’ ) ) % Demande un s e c o n d nombre e t a f f i c h e l a somme

Dans cette partie du cours, nous allons étudier quelques principes du fonctionnement d’un
micro-processeur afin de mieux appréhender la manière dont les programmes fonctionnent.

3.2 Le modèle "Little Man Computer"

Nous allons commencer ce chapitre par étudier le modèle du Little Man Computer (LMC).
C’est un modèle simplifié utilisé pour l’apprentissage de l’architecture des ordinateurs.
Un petit homme ("Little man") est enfermé dans une pièce. Cette pièce comporte (Fig. 3.1) :

33
 — 100 boîtes numérotées de 0 à 99 représentant la mémoire de l’ordinateur, chacune pouvant
recevoir un nombre décimal ;
— Deux boîtes nommées inbox et outbox permettant au petit homme de recevoir ou d’envoyer à
l’extérieur des nombres décimaux ;
— Son bureau est un modèle du processeur, il comporte :
— Une boîte appelée accumulateur ;
— Une calculatrice qui écrit son résultat dans l’accumulateur ;
— Une boîte appelée compteur programme indiquant le numéro de la boîte de la mémoire
contenant l’instruction que doit effectuer le petit homme ;
— Une boîte appelée registre d’instruction qui comporte l’instruction que réalise le petit
homme.

F IGURE 3.1 – L’environnement de travail du petit homme, chargé avec le programme introduit au
début de ce chapitre.

Au début le compteur programme indique la boite mémoire d’adresse 0. Dès que le petit homme a
copié son instruction dans le registre d’instruction, ce compteur est incrémenté automatiquement
de 1.

Les instructions sont sous forme d’un code décimal à 3 chiffres. C’est son langage machine. Les
instructions peuvent être écrites avec un assembleur mais alors elles doivent êtres compilées pour
pouvoir être comprises par le petit homme. Le tableau suivant résume les instructions disponibles,
le code machine et le code assembleur correspondant :

34
 Code numérique Code assembleur Description
Ajoute le contenu de la boîte mémoire xx à la valeur
1xx ADD contenue dans l’accumulateur, le résultat est placé
dans l’accumulateur
Soustrait le contenu de la boîte mémoire xx à la
2xx SUB valeur contenue dans l’accumulateur, le résultat est
placé dans l’accumulateur
Enregistre le contenu de l’accumulateur dans la boîte
3xx STA
mémoire xx
5xx LDA Charge le contenu de la boite xx dans l’accumulateur
Fixe le compteur programme à la valeur xx, la boite
6xx BRA
mémoire xx contiendra donc la prochaine instruction
Fixe le compteur programme à la valeur xx si le
7xx BRZ
contenu de l’accumulateur est 0
Fixe le compteur programme à la valeur xx si le
8xx BRP
contenu de l’accumulateur est positif
Transfère le contenu de la boite inbox dans l’accu-
901 INP
mulateur
Transfère le contenu de l’accumulateur dans la boite
902 OUT
d’envoi (outbox)
000 HLT ou COB Fin du programme (pause café - COB)

3.3 Architecture interne des processeurs

Nous avons vu le modèle du LMC, mais dans la réalité un processeur est d’abord un morceau
de silicium sur lequel sont gravés plusieurs millions de transistors. Il fonctionne en logique binaire,
c’est-à-dire avec deux états qui peuvent être 0 ou 1 en binaire et qui sont représentés par 0 V ou par
0.8 V en interne.
Nous pouvons les décomposer les micro-processeurs en deux grandes parties : l’unité de
commande et l’unité de traitement. Ces deux unités dialoguent avec l’extérieur en utilisant le bus
de données (pour transmettre des données) et le bus d’adresses (pour savoir où sont les données).

F IGURE 3.2 – Unité de traitement et unité de commande

3.3.1 L’unité de commande

L’unité de commande s’occupe de séquencer le déroulement des opérations. Elle effectue une
recherche en mémoire de l’instruction, le décodage de l’instruction et finalement, elle pilote son
exécution.

35
 Pour cela, elle utilise un compteur de programme (PC) qui indique en permanence l’adresse de
la prochaine instruction à exécuter. Il peut écrire sur le bus d’adresse. L’instruction correspondante
est ensuite chargée dans le registre 1 d’instruction via le bus de données. Finalement, l’instruction
est traitée après son décodage par une unité appelée le séquenceur.

3.3.2 L’unité de traitement

L’unité de traitement regroupe les circuits qui sont nécessaire à l’exécution des instructions.
Pour cela, elle dispose d’accumulateurs qui sont des registres de travail permettant de mémoriser
des opérandes ou le résultat d’une opération arithmétique. Ces opérandes peuvent être utilisées
dans l’unité arithmétique et logique (UAL) qui permet d’effectuer les opérations logiques (ou, et...)
et arithmétiques (addition, soustraction, etc). L’opération réalisée est choisie par le séquenceur de
la partie commande.

3.3.3 Exemple
A partir des informations données dans les deux sections précédentes, nous pouvons détailler la
Fig. 3.2 pour obtenir la Fig 3.3. Reprenons l’exemple donné au début du chapitre :

F IGURE 3.3 – Unité de traitement et unité de commande

1. Instruction "901"
Le compteur programme (PC) contient 0. Il pointe donc sur le contenu de la boite mémoire
d’adresse 0 ("901"). "901" est chargé dans le registre d’instruction. Le compteur programme
pointe désormais sur l’instruction suivante ("334"). Le décodeur convertit le code "901" en
une demande de lecture d’un nombre qui sera sur le bus de donnée après sa saisie. Le nombre
est chargé dans l’accumulateur.
2. Instruction "334"
Le registre d’instruction est chargé avec l’instruction "334". PC pointe désormais sur l’ins-
truction "901". L’instruction "334" comporte deux parties, le code 3 indique un chargement et
34 est l’adresse de la boite mémoire dans lequel le contenu de l’accumulateur sera mémorisé.
3. Instruction "901"
Le registre d’instruction est chargé avec l’instruction "901". PC pointe maintenant vers

1. Un registre est une toute petite mémoire

36
 l’instruction "134". Un nouveau nombre est demandé à l’utilisateur et mémorisé dans l’accu-
mulateur.
4. Instruction "134"
Le registre d’instruction est chargé avec l’instruction "134". PC pointe vers l’instruction
"902". L’instruction "134" est décodée, le code 1 indique une opération d’addition entre
l’accumulateur et le contenu de la mémoire à l’adresse 34. Le décodeur sélectionne l’opération
addition dans l’UAL et demande à chercher le contenu de la mémoire 34.
5. Instruction "902"
Le contenu de l’accumulateur est affiché à l’utilisateur.
6. Instruction "000"
Fin du programme.
Cette architecture avec un seul accès à une mémoire contenant les données et le programme
est appelé architecture de Von Neuman. C’est l’architecture adoptée par les processeurs d’usage
général.

3.4 Les pipelines

Dans notre exemple, nous avons vu le fonctionnement basique d’une opération de calcul :
1. Charger les instructions depuis la mémoire ;
2. Charger les opérandes depuis la mémoire ;
3. Effectuer les calculs ;
4. Mémoriser les résultats.
Pendant qu’une de ces étapes est effectuée, le reste du processeur est inoccupé. L’architecture
pipeline s’attaque à ce problème.
Ce concept est directement issu du travail à la chaîne dans les lignes de montages. Imaginons
une voiture qui doit subir trois opérations qui prennent un temps T avant d’être prête (Fig. 3.4).
Une voiture nécessitera 3T pour passer dans toutes les étapes de la chaînes. Ainsi, si l’on fait passer
les voitures une à une, il faudra 12T pour quatre voitures. Dans le cas d’une architecture pipeline,
une nouvelle voiture est traitée dès qu’un poste est libre. Les voitures nécessiterons toujours 3T
pour être prête mais il ne faudra plus que 6T pour quatre voitures.

F IGURE 3.4 – Chaîne sans architecture pipeline (a) et avec architecture pipeline (b)

Dans un ordinateur, l’exécution d’une instruction est découpée en petits morceaux appelés
étages du pipeline. Plusieurs instructions se chevauchent donc en même temps. Comme nous
l’avons vu pour les voiture, le temps d’exécution d’une instruction reste le même mais le débit est
augmenté. L’architecture pipeline est invisible pour le programmeur.
Par exemple, un Pentium 4 Prescott a 31 étages de pipeline.
L’utilisation d’un pipeline n’est pas sans poser de nouveaux problèmes. Ils peuvent être de trois
types :

37
 — Aléa de structure : L’implémentation empêche une certaine combinaison d’instructions
(ressources matérielles demandées par plusieurs étages) ;
— Aléa de données : Le résultat d’une opération dépend des opérations précédentes, par
exemple :
1 a = b + c;
2 c = a ∗5;

— Aléa de contrôle : L’exécution d’une instruction du type "if" ou "while" nécessite d’évaluer
la condition pour savoir quelle branche charger dans le pipeline. Par exemple :
1 a = b + c;
2 i f a == 5
3 c = a ∗ 5;
4 else
5 c = 4;
6 end

Pour résoudre les aléas de données, le compilateur réordonne les opérations ou insère des blancs
dans le pipeline.
Les aléa de contrôles sont plus compliqués à prendre en compte. Lorsque l’étage de décodeur
d’instruction ne peut pas calculer l’adresse du saut, il doit la prédire (par exemple en supposant que
la condition d’une boucle "while" est plus souvent vraie que fausse). Si le choix n’est pas correct,
alors il faudra vider le pipeline, ce qui introduira une latence supplémentaire.
Les aléa de contrôle ont incidemment une influence sur les performances d’un code Matlab.
Prenons le code suivant :
1 tic
2 j = 0;
3 f o r i = 1:100000000
4 if i > 1
5 j = 0;
6 else
7 j = 0;
8 end
9 end
10 toc

Matlab nous indique que ce programme à un temps d’exécution de 0.39 secondes sur la machine sur
laquelle il a été utilisée. Un second code, avec exactement la même fonctionnalité mais la condition
inverse met quand à lui 0.26 secondes.
1 tic
2 j = 0;
3 f o r i = 1:100000000
4 i f i <= 1
5 j = 0;
6 else
7 j = 0;
8 end
9 end
10 toc

38
 Matlab est un langage interprété, et non pas compilé, aussi il est difficile d’interpréter le résultat.
Le même code n’aurait pas fonctionné dans un langage tel que le C, il aurait été supprimé par le
compilateur car il ne fait rien.
Nous pouvons néanmoins voir que la condition if dans la boucle while est particulièrement
difficile à évaluer pour la prédiction d’embranchement.

3.5 La mémoire cache

La mémoire principale d’un ordinateur doit avoir une taille importante tout en gardant un temps
d’accès court. Ce sont deux contraintes fortement antagonistes. L’accès à la mémoire est un goulot
d’étranglement, en effet, nous avons vu que le micro-processeur doit y accéder pour connaître ses
instructions ainsi que pour récupérer les opérandes nécessaire à leur exécution. Au fur et à mesure
de l’évolution des mémoires et des processeurs l’écart de performance entre les accès à la mémoire
et la vitesse d’exécution du processeur n’a fait que s’accentuer (Fig. 3.5).

F IGURE 3.5 – Évolution comparées de performances des accès à la mémoire et des performances
des processeurs.

Afin d’améliorer ses performances, une mémoire beaucoup plus petite et rapide que la mémoire
principale est intégrée directement au cœur processeur, c’est la mémoire cache. Ainsi, si dans les
années 1990, la plupart des micro-processeur n’avaient pas de cache, en 2001, il en ont pratiquement
tous au moins 2. Si nous reprenons le modèle du LMC, nous pourrions voir la mémoire cache
comme quelques boîtes placées à proximité du bureau pour que le petit homme n’ait ni à marcher
ni à chercher longtemps (Fig. 3.6).
La mémoire cache utilise deux grands principes :
— La localité temporelle : une donnée récemment utilisée à tendance à être réutilisée ensuite.
Par exemple, c’est le cas de la variable a dans le code suivant :
1 a = 3;
2 b = a + c;

— La localité spatiale : les données dans la même zone mémoire qu’une donnée utilisée ont
tendance à être utilisées ensuite. C’est le cas des tableaux, comme dans l’exemple suivant :

39
 F IGURE 3.6 – L’environnement de travail du petit homme avec une mémoire cache.

1 f o r i = 1:100
2 a[ i ] = i ;
3 end

La mémoire cache utilise des blocs. Ces blocs sont des morceaux de mémoire de taille comprise
entre 4 et 128 octets. La taille de ces blocs doit être identique dans la mémoire principale et la
mémoire cache. Lorsque l’unité de calcul demande une donnée dans la mémoire, elle fournit le
numéro du bloc et la position de la donnée dans le bloc. Cette donnée est tout d’abord cherchée
dans le cache. Si elle n’est pas présente (défaut de cache), cette donnée est rapatriée de la mémoire
centrale et le bloc mémoire contenant cette donnée est recopié dans la mémoire cache.
La mémoire cache est gérée en interne par le micro-processeur. Il s’occupe de l’accès à la
mémoire, de gérer la consistance des données entre le cache et la mémoire centrale, et de choisir
les données qui peuvent être remplacées. Néanmoins, il est à la charge du programmeur de prendre
en considération la mémoire cache pour optimiser son code. Par exemple prenons le code Matlab
suivant :
1 tic
2 a = zeros (10000 ,10000) ;
3 f o r i = 1:10000
4 f o r j = 1:10000
5 a ( i , j ) = rand ( ) ;
6 end
7 end
8 toc
9

10

11 tic
12 b = zeros (10000 ,10000) ;
13 f o r j = 1:10000

40
14 f o r i = 1:10000
15 b ( i , j ) = rand ( ) ;
16 end
17 end
18 toc

Un tableau à deux dimensions N, M est mémorisé sous la forme d’un tableau à une dimension de
taille N × M en mémoire. Ainsi, le code utilisant b(i, j) et respectant la localité spatiale a mis 1.9
secondes à s’exécuter tandis que le code utilisant a(i, j) a mis 2.2 secondes. Cette différence peut
ne pas vous sembler énorme. En effet, les mauvaises performances de la gestion des boucles dans
Matlab cache une partie du problème. Si nous faisons la même expérience avec un code en C :
1 # i n c l u d e < t i m e . h>
2 # i n c l u d e < s t d i o . h>
3

4 i n t main ( )
5 {
6 double a [1000∗1000];
7 int i , j ;
8

9 clock_t s t a r t = clock () , d i f f ;
10

11 f o r ( i = 0 ; i < 1 0 0 0 ; i ++)
12 {
13 f o r ( j = 0 ; j < 1 0 0 0 ; j ++)
14 {
15 a [ i +1000∗ j ] = i ;
16 }
17 }
18 d i f f = clock () − s t a r t ;
19 }

Et :
1 # i n c l u d e < t i m e . h>
2 # i n c l u d e < s t d i o . h>
3

4 i n t main ( )
5 {
6 double b [1000∗1000];
7 int i , j ;
8

9 clock_t s t a r t = clock () , d i f f ;
10

11 f o r ( j = 0 ; j < 1 0 0 0 ; j ++)
12 {
13 f o r ( i = 0 ; i < 1 0 0 0 ; i ++)
14 {
15 b [ i +1000∗ j ] = i ;
16 }
17 }

41
 18 d i f f = clock () − s t a r t ;
19 }

Le premier code nécessite 6 ms sur mon ordinateur, et le second 3 ms. La prise en compte de la
présence du cache a induit un gain d’un facteur 2 dans les performances - juste en inversant deux
indices.

3.6 Conclusion
Nous avons vu, dans cette partie, un exemple simplifié d’ordinateur (Little Man Computer).
Nous l’avons utilisé pour comprendre les principaux éléments permettant d’exécuter un programme.
A partir de ce modèle, nous avons abordé les notions de pipeline et de mémoire cache. Ces éléments
sont gérés par le micro-processeur, néanmoins, les connaître est important pour améliorer les
performances de ses programmes.
Nota Bene : Le modèle LMC utilisait une boite particulière appelée "Compteur programme".
Son contenu était l’adresse d’une des 100 boîtes dans laquelle le petit homme devait lire ses instruc-
tions. Ces variables qui contiennent l’adresse d’un contenu et non pas le contenu lui même sont
appelés pointeur. C’est une notion très importante pour beaucoup de langages de programmation.

42
 Algorithmes et performances

4.1 Introduction (Ref. [6])

Un des paramètres fondamentaux pour sélectionner et choisir un algorithme est la vitesse avec
laquelle il est capable de résoudre le problème posé. Nous allons étudier, au cours de ce chapitre,
quelques méthodes qui nous permettront de les comparer.
Commençons par étudier quelques exemples d’algorithmes permettant de calculer la somme des
n premiers entiers. Le premier algorithme, que nous appellerons algo_1, utilise un accumulateur qui
est initialisé à 0. L’algorithme itère ensuite n fois, en ajoutant à chaque fois la valeur de l’itérateur à
l’accumulateur. Nous pouvons l’écrire ainsi :

Initialiser la somme à 0
Pour i variant de 1 à n faire
Ajouter i à somme

Un algorithme est une succession d’instruction, c’est une méthode permettant de calculer la réponse
au problème posé à partir d’un vecteur d’entrée. A partir de celui-ci, nous devons construire
une représentation dans un langage de programmation. Il existe de nombreuses implémentations
possibles d’un même algorithme, dépendant du programmeur et du langage utilisé. Dans ce cours,
nous utilisons Matlab. Une première implémentation possible de l’algo_1 est, par exemple (impl_1) :

somme = 0;
for i = 1:n
somme = somme + i;
end

Nous pouvons réaliser la même fonction, dans le même langage, en utilisant des noms de variables
différents et des étapes supplémentaires. Le programme suivant réalise exactement la même fonction
(impl_2) :

lisa = 0;
for michel = 1:n
muriel = michel;
lisa = lisa + muriel;
end

Nous pouvons également avoir une représentation de l’algorithme dans un langage différent.
L’exemple suivant utilise le C (impl_3) :

somme = 0;
for (i = 1; i < n ; i++)
{
somme += i;
}

La caractérisation de ses implémentation dépend de nos critères. Ainsi l’imp_1 est une meilleure
implémentation que impl_2 car elle est plus lisible. L’impl_3 sera sans doute, sur la même machine,

43
 beaucoup plus rapide que les solutions utilisant Matlab. Néanmoins dans cette partie, ces critères
sont basés sur l’implémentation de l’algorithme et non pas sur l’algorithme lui même.
Pour cela, nous allons plutôt nous baser sur les ressources utilisées. Nous voulons pouvoir dire
que cet algorithme est plus performant parce qu’il utilise moins de ressources que tel autre. De cette
perspective, les implémentations précédentes sont similaires car elles se basent le même algorithme
pour résoudre le problème.
Les ressources disponibles sur un ordinateur sont la mémoire et la puissance de calcul. Pour
commencer, nous nous focaliserons sur le temps d’exécution nécessaire à un algorithme. C’est à
dire que nous allons enregistrer et analyser le temps nécessaire à un programme pour calculer le
résultat. Cette opération est réalisée en Matlab avec les fonctions tic et toc.
Pour comparer notre premier algorithme à quelque chose, nous devons en définir un deuxième.
Celui-ci (algo_2) utilise la formule suivante :
n
n(n + 1)
∑i = 2
(4.1)
i=0

L’algorithme résultant n’est pas itératif et calcule directement le résultat. Une implémentation en
Matlab pourrait être (impl_4) :

somme = n * (n + 1)/2;
Voyons comment cet algorithme se comporte. Observons comment le temps d’exécution évolue en
fonction du la valeur de n :
n impl_1 impl_4
103 5 µs 2 µs
104 35 µs 2 µs
105 344 µs 2 µs
106 3.5 ms 2 µs
107 34.8 ms 2 µs
Nous observons que le temps d’exécution texec du premier algorithme semble croître linéairement
avec la valeur de n tandis que le second algorithme garde un temps d’exécution constant. Néanmoins,
si nous utilisons un autre ordinateur ou un autre langage pour exécuter ce programme, les résultats
trouvés seront différents.
Nous devons trouver une manière de caractériser les algorithmes qui ne doit dépendre ni de la
machine utilisée, ni du langage utilisé.

4.2 La notation en O
Afin d’avoir une mesure de l’efficacité d’un algorithme, nous devons définir une quantité
mesurable qui dépende explicitement de la taille du problème étudié. Elle dépend de ce que nous
voulons comparer. Cela peut être le nombre d’itération fait dans les boucles, le nombre d’affectations
ou encore le temps de calcul. C’est cette quantité mesurable (1 seconde, 1 itération...) qui sera la
base de la mesure de la performance.
Sur notre exemple, nous allons quantifier le nombre d’affectation que l’algorithme nécessite.
Ainsi, à chaque fois que la valeur d’une variable est modifiée, le nombre d’affectation est incrémen-
tée de 1. Dans ce cas, l’algorithme nécessite un pas pour l’initialisation (somme = 0) suivi de n pas
pour l’affectation de la variable i et finalement n pas pour réaliser la somme (somme = somme + i).
Soit T la fonction de n représentant l’évolution du nombre de pas. Nous pouvons écrire :

T (n) = 1 + 2 · n (4.2)

44
 Le paramètre n est appelé taille du problème. En utilisant la fonction T (n), nous pouvons déduire
que le nombre de pas, et donc le temps d’exécution sera plus important pour n = 108 que pour
n = 100 et que la progression sera linéaire.
Le nombre exact d’opérations ou de pas est de moins en moins important à mesure que la valeur
de n est conséquente. Ainsi, pour comparer différents algorithmes, nous utiliserons seulement un
ordre de grandeur, c’est à dire uniquement le terme de T (n) qui croît le plus rapidement. Nous
utiliserons pour cela la notation en O( f (n)) où f (n) est une fonction simple comme, par exemple,
f (n) = n ou f (n) = n2 . Ainsi, pour reprendre l’exemple précédent, nous pouvons dire que le
nombre d’itération croît en O(n).
Cette notation ne doit pas nous faire oublier que :
— La valeur des constantes est importante (c’est pour cela que nous changeons d’ordinateur
pour des versions plus puissantes) ;
— Suivant les problèmes, la valeur de n n’est pas toujours grande.
Notons que, lorsque l’algorithme est trop complexe ou que son implémentation est inconnue,
on ne peut pas déterminer simplement la fonction T (n) à partir du nombre d’itération ou d’af-
fectation. Pour cela, nous pouvons utiliser un programme (qui peut être commercial) qui testera
une implémentation à partir d’un ensemble de données générées aléatoirement. Une régression est
ensuite effectué sur les temps d’exécution afin de déterminer la fonction qui correspond le mieux
aux résultats obtenus. Dans notre exemple, nous avons obtenu :
T impl_1 (n) = 3.4 · 10−11 × n + 5 · 10−6 (4.3)
impl_4 −6
T (n) ≈ 2 · 10 (4.4)
Ce qui nous permet de retrouver que l’algorithme 1 est en O(n) et le 2 en O(1).

4.3 Analyse meilleur cas, pire cas, cas moyen

Bien que cela ne soit pas visible dans l’exemple précédent, les performances d’un algorithme
dépendent non seulement de la taille du problème mais également des données.
Pour voir cet effet, examinons un algorithme de recherche séquentiel. A partir d’une série de
donnée D de taille n, nous souhaitons trouver la position du premier élément égal à une valeur Ds .
L’algorithme pourrait être :
i = 0
Tant que D[i] est différent de Ds et que i < n faire
i = i + 1
Retourner la valeur de i

4.3.1 Pire cas

Nous pouvons observer dans l’algorithme de recherche séquentiel que le nombre d’opérations
effectuées dépend directement de la position de la donnée recherchée dans la série de donnée D. Il
y aura une seule itération dans le cas où la donnée recherchée Ds est en première position (i = 0) et
n itération si elle est en dernière position (i = n − 1), c’est le pire cas.
Ainsi, le pire cas correspond aux jeux de données qui donnent un temps d’exécution maximum
pour l’algorithme étudié. Il permet de borner et d’estimer son comportement dans n’importe quelle
situation. Cette analyse permet de comprendre et d’analyser les raisons qui empêchent l’algorithme
de se comporter de manière normale. Dans l’exemple de la recherche séquentielle, nous pouvons
observer que les performances de l’algorithme sont liées à la position de la donnée recherchée.
Ainsi, nous pouvons essayer d’améliorer cet algorithme en traitant ces cas.
De manière plus formelle, si Dn est un ensemble de données di de taille n et t une fonction
caractérisant le travail effectué par l’algorithme, alors le pire cas est la valeur maximale de t(di )

45
 pour tous les di de Dn . La fonction Twc (n) de croissance de t en fonction de n définit la performance
dans les pires cas de l’algorithme.
Ainsi, si nous considérons le nombre d’affectations dans le pire cas pour l’algorithme de
recherche séquentielle, nous obtenons :

Twc (n) = n (4.5)

4.3.2 Meilleur cas

De la même façon que nous avons défini le pire cas à partir du jeu de données donnant un temps
d’exécution maximum, le meilleur cas sera obtenu lorsque le temps d’exécution est minimal. Sur
l’algorithme en exemple, il sera obtenu pour toutes les données dont la donnée recherchée est en
première position.

4.3.3 Cas moyen

Considérons un ensemble Dn de données di de taille n et t une fonction caractérisant le travail
effectué par l’algorithme. Nous pouvons associer à chaque donnée di , une probabilité P(di ) qu’elle
soit présentée en entrée de l’algorithme telle que :

∑ Ṗ(di ) = 1 (4.6)
di ∈Dn

On définit alors la performance de l’algorithme dans le cas moyen par la fonction T (n) :

T (n) = ∑ t(di ) · P(di ) (4.7)

di ∈Dn

Pour une recherche séquentielle sur un ensemble de données aléatoires de taille n, la valeur
recherchée sera en moyenne trouvée à la moitié du domaine de recherche. Ainsi, si nous considérons
le nombre d’affectation comme mesure de base :
1
T (n) = 1 + n (4.8)
2

4.3.4 Famille de performance

Nous pouvons évaluer et comparer la performance de différents algorithmes en étudiant l’évolu-
tion de la fonction T (n) avec la taille n du problème. Pour rappel, cette évaluation est liée au temps
de calcul. Nous pourrions également considérer la place mémoire nécessaire à l’algorithme pour
effectuer son traitement avec une approche similaire.
Nous pouvons classiquement classifier les algorithmes en grandes familles suivant leur perfor-
mance (Fig. 4.1) :
— Constante : O(1)
— Logarithmique : O(log(n))
— Linéaire : O(n)
— Log-linéaire : O(nlog(n))
— Quadratique : O(n2 )
— Exponentielle : O(2n )
Lorsque un algorithme est composé de plusieurs tâches indépendantes, sa performance est
liée aux traitements qui nécessitent le plus de ressources. Par exemple, si un algorithme peut se
décomposer en deux tâches, la première en O(n) et la seconde en O(2n ), la performance globale de
l’algorithme sera en O(2n ).

46
 F IGURE 4.1 – Familles de performance

4.4 Exemples d’algorithmes

4.4.1 Divination (Ref. [5]
Commençons par un jeu connu. Un premier joueur pense très fort à un nombre d entier compris
entre 1 et n. Un second joueur doit deviner ce nombre en un minimum de proposition prop, le
premier joueur lui indique si prop est inférieur ou supérieur à d. En tant que second joueur, nous
savons que proposer aléatoirement des valeurs est assez peu efficace. Nous allons donc commencer
par proposer la valeur médiane du domaine puis continuer à chaque fois en proposant la valeur
médiane du domaine restant. Nous pouvons résumer cette méthode avec l’algorithme suivant :

N_maximum = n
N_minimum = 1
prop = E((N_maximum+N_minimum)/2)
Tant que prop n'est pas égal à d faire
Si prop < d alors
N_minimum = prop + 1
Sinon
N_maximum = prop - 1
prop = E((N_maximum+N_minimum)/2)

Lorsque n est suffisamment petit, nous pouvons étudier tous les cas possibles en fonction de
la valeur à deviner d. Étudions, dans le tableau suivant, le comportement de l’algorithme pour n = 8 :

47
 Valeur à Valeur inférée par
deviner (d) l’algorithme (valeur)
1 4 2 1
2 4 2
3 4 2 3
4 4
5 4 6 5
6 4 6
7 4 6 7
8 4 6 7 8
Nous pouvons observer que cet algorithme diminue le domaine de recherche par deux à chaque
itération. Ainsi, dans le pire cas, l’algorithme nécessitera T (n) = 1 + E(log2 (n)) itérations pour
trouver la solution (où E représente la fonction partie entière).
Ainsi, nous retrouvons qu’il faut au maximum T (8) = 4 itérations pour n = 8. Nous pouvons
également voir comment l’algorithme se comporte pour des valeurs de n plus grand. Par exemple
dans le cas où n = 20, nous trouvons T (20) = 1 + E(4.32) = 5 itérations et pour n = 106 , seulement
T (106 ) = 20 itérations dans le pire cas. Nous pouvons observer que les algorithmes avec une per-
formance logarithmique sont extrêmement performant. Dans l’exemple présenté, une augmentation
de n de 5 ordre de grandeurs n’a impliqué que 15 itérations supplémentaires.
Pour finir, une implémentation possible de cet algorithme en Matlab pourrait être :
1 f u n c t i o n i = d i v i n a t i o n (D, N)
2 Nmin = 1 ;
3 Nmax = N ;
4 p r o p = f l o o r ( ( Nmax+Nmin ) / 2 ) ;
5 i = 1;
6 w h i l e D ~= p r o p
7 i = i +1;
8 if valeur > D
9 Nmax = prop −1;
10 else
11 Nmin = p r o p + 1 ;
12 end
13 p r o p = f l o o r ( ( Nmax+Nmin ) / 2 ) ;
14 end
15 end

Nous utilisons ici également Matlab pour calculer le nombre d’itérations dans le pire cas (en
prenant D = Nmax) et dans le cas moyen en utilisant un jeu de 1000 essais aléatoires. Le résultat
est présenté Fig. 4.2.
A partir de ces données brutes, nous allons chercher la meilleure fonction s’écrivant sous la
forme T (n) = A + B · log2 (n) représentant ces points. Cette opération peut se faire également sous
Matlab à partir des données brutes resMoy calculées au points Nvar :
1 ft = f i t t y p e ( ’ a + b∗ l o g ( x ) / l o g ( 2 ) ’ , ’ c o e f f i c i e n t s ’ , { ’ a ’ , ’ b ’ } )
2 m y f i t = f i t ( Nvar ’ , resMoy ’ , f t )

Nous obtenons les résultats suivants dans le pire cas (Twc (n)) et dans le meilleur cas (T (n)) :
Twc (n) = 0.82 + 0.95 · log2 (n) (4.9)
T (n) = −0.41 + 0.94 · log2 (n) (4.10)

48
 F IGURE 4.2 – Tracé de la fonction T (n) dans le pire et le cas moyen

Nous pouvons observer que l’on retrouve la performance en O(log(n)).

Nous pouvons reproduire la même expérience en nous basant sur le temps de calcul nécessaire
à l’algorithme pour trouver la solution. Le résultat est présenté sur la Fig. 4.3. Le temps d’exécution
d’un appel à la fonction est trop court pour pouvoir être mesuré avec suffisamment de précision.
C’est pour cela que chaque point représenté sur la courbe est la somme du temps d’exécution
obtenu pour 10000 appels.
La Fig. 4.3 a bien le comportement attendu. Les pics et variations présents sont des artefacts
liés à l’occupation et la charge de l’ordinateur sur lequel l’expérience a été menée.

4.5 Anagrammes (Ref. [6]

Dans ce second exemple, nous allons étudier les performances de différents algorithmes
permettant de détecter les anagrammes. Une chaîne de caractère est l’anagramme d’une autre si elles
sont un simple réarrangement de caractères. Par exemple "aube" et "beau" sont des anagrammes.
Nous supposerons que les anagrammes que nous recherchons sont uniquement composés des 26
lettres de l’alphabet et que les deux chaînes de caractères font la même taille n.
L’objectif est d’écrire une fonction qui renvoie une valeur booléenne (0 ou 1) indiquant si les
deux chaînes sont des anagrammes.

4.5.1 Premier algorithme

Dans un anagramme, chaque caractère d’une chaîne s1 doit avoir son jumeau dans la seconde
chaîne s2. Nous allons chercher un algorithme qui utilise directement cette propriété. Pour cela, il
itèrera sur chaque caractère de s1 puis sur ceux de s2 afin de trouver son homologue. Le caractère
de s2 trouvé sera remplacé par le caractère "-" afin de ne pas être compté deux fois. Si l’algorithme
trouve un caractère dans s2 pour chaque caractère de s1 alors les deux chaînes sont des anagrammes.
Cet algorithme peut s’écrire de la façon suivante en Matlab :
1 f u n c t i o n r e s = a n a 1 ( s1 , s 2 )
2

3 c2 = s2 ;
4 res = true ;
5 pos1 = 0;
6

7 w h i l e p o s 1 < l e n g t h ( s 1 ) && r e s == t r u e

49
 F IGURE 4.3 – Tracé de la fonction T (n) dans le cas moyen pour 1000 appels à la fonction divination

8 pos1 = pos1 + 1 ;
9 pos2 = 0 ;
10 res = false ;
11 w h i l e p o s 2 < l e n g t h ( s 2 ) && r e s == f a l s e
12 pos2 = pos2 + 1 ;
13 i f s 1 ( p o s 1 ) == c2 ( p o s 2 )
14 res = true ;
15 c2 ( p o s 2 ) = ’− ’ ;
16 end
17 end
18 end

Pour analyser cet algorithme, nous devons remarquer que chacun des n caractères de s1 génèrera
une itération sur les n caractères de s2. Chacune des n positions de s2 sera visitée une fois pour
chaque caractère de s1 . Si l’on compte les itérations, nous aurons donc :

T (n) = n2 (4.11)

Ainsi, la performance de cet algorithme est en O(n2 ).

4.5.2 Deuxième solution

Nous pouvons également remarquer que, même si s1 et s2 sont différentes, elles sont des
anagrammes seulement si elles ont les mêmes caractères. Ainsi, si nous commençons par trier par
ordre alphabétique s1 et s2, nous devrions trouver exactement la même chaîne dans les deux cas.
Pour cela, nous allons utiliser la fonction sort de Matlab :
1 f u n c t i o n r e s = a n a 2 ( s1 , s 2 )
2

50
 3 c1 = s o r t ( s1 ) ;
4 c2 = s o r t ( s2 ) ;
5

6 res = true ;
7 pos = 1;
8

9 w h i l e p o s < l e n g t h ( s 1 ) && r e s == t r u e
10 i f c1 ( p o s ) == c2 ( p o s )
11 pos = pos + 1 ;
12 else
13 res = false ;
14 end
15 end

Il serait tentant de dire que cet algorithme est en O(n). En effet, il n’y a qu’une boucle sur les n
caractères des deux chaînes après l’étape de tri. Néanmoins, les deux appels à la fonction sort
ne sont pas gratuits. Comme nous le verrons dans un chapitre ultérieur ces algorithmes ont une
performance variant entre O(n log n) et O(n2 ). C’est donc cette étape qui domine et qui donnera la
performance de cet algorithme.

4.5.3 Troisième algorithme

Pour cette dernière solution, nous allons nous baser sur le fait que deux anagrammes ont le
même nombres de a, de b, de c, etc. Ainsi, pour vérifier que deux chaînes de caractères sont des
anagrammes, nous allons commencer par compter le nombre d’occurrence de chaque caractère.
Comme il y a 26 caractères possibles, nous allons utiliser un tableau de 26 compteurs. Chaque fois
qu’un caractère particulier sera rencontré, la cellule correspondante sera incrémenté. Si a la fin les
deux tableaux sont identiques alors les deux chaînes sont des anagrammes. Cet algorithme peut
être implémenté de la façon suivante :
1 f u n c t i o n r e s = a n a 3 ( s1 , s 2 )
2

3 c1 = zeros (26 ,1) ;

4 c2 = zeros (26 ,1) ;
5

6 f o r i = 1 : l e n g t h ( s1 )
7 pos = u i n t 1 6 ( s1 ( i ) ) − u i n t 1 6 ( ’ a ’ ) + 1 ;
8 c1 ( p o s ) = c1 ( p o s ) + 1 ;
9 end
10

11 f o r i = 1 : l e n g t h ( s2 )
12 pos = u i n t 1 6 ( s2 ( i ) ) − u i n t 1 6 ( ’ a ’ ) + 1 ;
13 c2 ( p o s ) = c2 ( p o s ) + 1 ;
14 end
15

16 res = true ;
17 i = 1;
18 w h i l e i < 26 && r e s == t r u e
19 i f c1 ( i ) == c2 ( i )
20 i = i + 1;
21 else

51
 22 res = false ;
23 end
24 end

Cette solution utilise de nombreuses boucles, mais contrairement à la première solution, elle ne sont
pas imbriquées. Les deux premières boucles font chacune n itération et la troisième au maximum
26. Ce qui nous donne finalement :

T (n) = 2 · n + 26 (4.12)

Nous avons donc finalement une solution linéaire pour la résolution de ce problème.
Avant de finir ces exemples, nous n’avons pas parlé de l’espace mémoire. Nous pouvons
remarquer, que bien que la troisième solution est la meilleure performance, elle nécessite un espace
mémoire supplémentaire pour les deux tableaux. En d’autre terme, cet algorithme a utilisé de la
mémoire pour gagner du temps. Dans ce cas, la taille supplémentaire (2 tableau des 26 entiers) est
négligeable. Il n’en aurait pas forcément été vrai si l’alphabet faisait 220 caractères.
C’est un compromis courant. Il faudra souvent déterminer le meilleur algorithme pour résoudre
un problème particulier en fonction des contraintes qui sont posées.

4.6 Conclusion
Dans ce chapitre nous avons étudié une méthode de comparaison de la performance des
algorithmes. Nous avons vu que nous pouvons l’exprimer avec la notation en O( f (n)) puis nous
l’avons appliqué à quelques exemples. Nous avons également vu que :
— Les constantes sont importantes ;
— L’occupation mémoire peut également être importante.

52
 La recherche et le tri (Ref. [6])

5.1 Introduction
Nous allons maintenant aborder deux des problèmes les plus fréquemment rencontrés : la
recherche et le tri. Cette étude nous permettra de mettre en œuvre les techniques d’analyse des
performances d’algorithme que nous avons vu au chapitre précédent.
Nous verrons, tout d’abord, les algorithmes de recherche. Cette catégorie de programme permet
de déterminer si un élément est présent ou non dans une liste. Ils renvoient un booléen suivant le
résultat de la recherche (true ou false) et, si nécessaire, la position de l’élément recherché. Matlab
dispose directement de cet algorithme. Ainsi, nous pouvons utiliser la fonction ismember, comme
dans l’exemple suivant :
1 >> x = [ 1 2 3 4 5 ] ;
2 >> ismember ( 1 , x )
3 ans =
4 logical
5 1
6 >> ismember ( 7 , x )
7 ans =
8 logical
9 0

Nous allons nous intéresser aux algorithmes sous-jacent à cette catégorie de fonctions.
Puis, dans une seconde partie, nous nous intéresseront aux algorithmes de tri. Ils traitent une
listent quelconque d’éléments afin de les ordonner par ordre croissant. Matlab dispose de la fonction
sort qui est une implémentation de l’algorithme quick-sort que nous verrons à la fin de ce chapitre.
N.B. : Comme nous nous intéressons plus aux algorithmes qu’à leur implémentation, les code
proposés dans cette section ne sont pas les plus performant en Matlab. En effet, ils n’utilisent pas
d’opérateurs sur les tableaux mais utilisent des boucles for ou while.

5.2 Quelques algorithmes de recherche

5.2.1 La recherche séquentielle
Algorithme
Nous devons tout d’abord poser le problème afin de pouvoir l’étudier. Dans cette partie, nous
cherchons une donnée (un entier, une chaîne de caractère, une adresse...) dans un ensemble de
données du même type. Le plus souvent, ces données sont accessibles sous forme linéaire. Par
exemple, en Matlab, nous pouvons accéder au différentes données d’un tableau T en utilisant
un index i et en écrivant simplement T (i). Nous allons utiliser cette propriété dans ce premier
algorithme.
L’algorithme de recherche séquentielle est basé sur un parcours séquentiel de la liste d’éléments
(par exemple dans le tableau T ), en commençant par le premier. Si la valeur correspondant à l’indice
actuel i est identique à la valeur recherché (valeur = T (i)), alors la recherche est arrêtée et la valeur
true sera retournée. Dans le cas contraire, si la valeur recherchée n’a pas été trouvée, la fonction
devra retourner la valeur false. Un exemple est présenté sur la Fig. 5.1.

53
 F IGURE 5.1 – Exemple de recherche séquentielle sur un tableau d’entier

L’implémentation de cet algorithme en Matlab se basera sur une fonction qui prendra en entrée
un tableau des valeurs tableau, l’élément recherché valeur et retournera en sortie un nombre booléen
resultat. La variable résultat est initialisée à false au début du programme et assignée à true en cas
de comparaison positive (valeur = tableau(i)).
1 function r e s u l t a t = rechseq ( tableau , valeur )
2 resultat = false ;
3 i = 1;
4 w h i l e ( i < l e n g t h ( t a b l e a u ) && r e s u l t a t == f a l s e )
5 i f v a l e u r ~= t a b l e a u ( i )
6 i = i +1;
7 else
8 resultat = true ;
9 end
10 end

Résultats
Le calcul de la performance des algorithmes de recherche nécessite de définir la quantité que
nous souhaitons mesurer. La taille du problème sera la taille de la liste d’éléments dans lequel
la recherche sera effectuée. Pour les algorithmes de recherche, nous utiliserons la nombre de
comparaisons comme grandeur mesurable. Nous supposerons, de plus, que les éléments du tableau
sont aléatoirement répartis.
Lorsque l’élément recherché n’est pas présent, le programme de recherche séquentielle vérifia
toute la liste. Il devra effectuer n comparaisons. Quelque soit le tableau, tous les cas sont équivalents
(meilleur, pire ou cas moyen), et sa performance est en O(n).
Lorsque l’élément recherché est présent, nous pouvons séparer les trois cas :
— Meilleur cas : L’élément recherché est toujours en première position, la performance est en
O(1) ;
— Pire cas : L’élément recherché est toujours en dernière position, la performance est en O(n) ;
— Cas moyen : L’élément recherché est quelque part dans le tableau, nous le trouverons en
moyenne au bout de n/2 comparaison. La performance est en O(n).

5.2.2 La recherche séquentielle triée

Algorithme
Afin d’améliorer l’algorithme précédent, nous allons préalablement trier le tableau tableau.
Nous aurons alors, pour tout i, tableau(i) < tableau(i + 1). Nous pouvons utiliser cette information
supplémentaire dans l’algorithme de recherche séquentielle : si la valeur comparé est strictement
inférieure à la valeur recherchée alors elle n’est pas dans le tableau (Fig. 5.2).
Une implémentation possible de la recherche séquentielle triée est la suivante :

54
 F IGURE 5.2 – Exemple de recherche séquentielle sur un tableau d’entier triés

1 function r e s u l t a t = rechseq ( tableauTrie , valeur )

2 resultat = false ;
3 i = 1;
4 arret = false ;
5 w h i l e ( i < l e n g t h ( t a b l e a u T r i e ) && a r r e t == f a l s e )
6 i f v a l e u r == t a b l e a u T r i e ( i )
7 arret = true ;
8 resultat = true ;
9 e l s ei f ( valeur > tableauTrie ( i ) )
10 arret = true ;
11 else
12 i = i +1;
13 end
14 end

Cette implémentation diffère assez peu de la précédente. Nous y avons ajouté une condition d’arrêt
(arret == false).
Résultats
La performance du code dans le cas où un élément est dans la liste n’est pas modifiée. Néan-
moins, nous pouvons observer désormais trois différents cas lorsqu’il n’est pas présent :
— Meilleur cas : L’élément recherché est inférieur au premier élément de la liste, la performance
est en O(1).
— Pire cas : L’élément recherché est supérieur au dernier élément de la liste, la performance est
en O(n).
— Cas moyen : L’élément recherché est compris entre deux éléments de la liste, la condition
d’arrêt sera vraie en moyenne au bout de n/2 essais. La performance est en O(n).
Cette méthode permet donc d’améliorer les résultats pour lorsque l’élément n’est pas dans le
tableau en dégradant un peu les résultats dans le cas contraire.

5.2.3 La recherche dichotomique

Algorithme
La recherche dichotomique utilise différemment le tableau trié. La recherche séquentielle (triée
ou non) part du premier terme et examine ensuite chacun des n termes suivants, un à un. Dans le
cas de la recherche par dichotomie, nous allons débuter par le milieu du tableau (position n/2).
Comme le tableau est trié deux cas sont alors possibles :
— Le nombre recherché est supérieur à la valeur comparée, la valeur recherchée doit alors être
dans la partie supérieure du tableau ;
— Le nombre recherché est inférieur à la valeur comparée, la valeur recherchée doit alors être
dans la partie inférieure du tableau ;
Après cette étape, nous cherchons la valeur à la moitié du domaine de recherche restant (position
3n/4 ou n/4). Ces étapes sont recommencées autant de fois que nécessaires jusqu’à confirmer ou
non la présence de la valeur recherchée. Nous obtenons donc le code suivant :

55
 1 function r e s u l t a t = rechdic ( tableau , valeur )
2 resultat = false ;
3 fin = length ( tableau ) ;
4 debut = 1;
5

6 w h i l e ( d e b u t <= f i n && r e s u l t a t == f a l s e )
7 milieu = f l o o r ( ( debut+ f i n ) / 2 ) ;
8 i f v a l e u r == t a b l e a u ( m i l i e u )
9 resultat = true ;
10 else
11 i f valeur < tableau ( milieu )
12 fin = milieu − 1;
13 else
14 debut = milieu + 1;
15 end
16 end
17 end

Nous pouvons observer le fonctionnement du programme sur la Fig. 5.3. Dans le premier cas
présenté, le nombre 70 est recherché. Il est comparé à la valeur placée au milieu du tableau (41)
d’indice 6 (debut = 1, fin = 11, milieu = 6). Comme 70 est supérieur à 54, la recherche continue sur
la partie supérieure du tableau (debut = 7, fin = 11, milieu = 9). La valeur recherchée est comparée
à 70 ; Le programme prend fin avec la variable resultat qui vaut true.
Dans le second cas, la recherche se passe de la même façon jusqu’à ce que la condition
debut==fin soit vraie sans que la valeur ai été trouvée dans le tableau. Dans ce cas, la variable
resultat garde sa valeur false.

F IGURE 5.3 – Exemple de recherche dichotomique sur un tableau d’entier triés

Résultats
Pour déterminer la performance de cet algorithme, nous pouvons tout d’abord remarquer que
chaque étape réduit le domaine de recherche de moitié. Ainsi, si à la première itération le taille du
problème est n, elle n’est plus que de n/2 à la seconde itération, puis de n/4 etc... Dans le pire des
cas, le domaine de recherche sera réduit à un seul élément (que ce soit l’élément recherché ou non).
Pour arriver à cet état, nous avons besoin de i itération avec E(n/2i ) = 1, soit i = log2 (n). Nous
en déduisons que la performance du tri dichotomique est en log2 (n). Nous retrouvons le résultat
obtenu pour l’algorithme de divination du chapitre précédent.
Nous avions terminé le chapitre étudiant performance des algorithme en rappelant que les
constantes avaient leur importance. C’est toujours le cas ici. L’algorithme de recherche dichoto-
mique a une performance bien meilleure que la recherche séquentielle. Néanmoins, elle nécessite

56
 une liste trié et cette opération n’est pas neutre. Ainsi, si nous ne devons effectuer que quelques
recherches ou que l’opération de tri est trop coûteuse, la recherche séquentielle donnera les meilleur
résultats.

5.2.4 Les tables de hachage

Principe
L’algorithme de recherche dichotomique nous a permit d’avoir une performance logarithmique.
Pour cela, nous avons eu besoin d’une information supplémentaire sur la tableau d’éléments : il
devait être trié. Nous allons, dans cette section, aller encore plus loin pour obtenir une performance
en O(1). C’est à dire que le résultat devra être trouvé en une seule comparaison, même si nous n’y
arriverons malheureusement pas complètement. Pour cela, nous allons avoir besoin d’une nouvelle
notion : le hachage.
Une table de hachage est un tableau dans lequel les données sont mémorisées de façon à être
retrouvées facilement. Chaque position dans la table est appelé une alvéole (slot en anglais) et peut
contenir un élément. Il y aura donc une alvéole 0, une alvéole 1, etc... Tout les éléments de la liste
sur laquelle les recherches vont être effectuées doivent être placés dans la table de hachage. Pour
connaître leur position, nous devons définir une fonction de hachage. Cette fonction prend en entrée
un élément et retourne un entier correspondant à l’alvéole qui lui correspond.
Prenons par exemple une table de hachage de taille m = 11, et une liste de valeur [26, 13, 60, 87, 9].
La fonction de hachage retenue est le reste de la division par 11 (le modulo - (mod) en Matlab). Le
modulo est souvent présent dans les fonctions de hachage car son résultat est toujours inférieur à la
taille de la table de hachage. Nous obtenons les valeurs suivantes pour les alvéoles :
valeur mod(valeur, 11)
26 4
13 2
60 5
87 10
9 9
A partir de ces résultats, nous pouvons compléter la table de hachage. Les valeurs obtenues nous
donnent le numéro de l’alvéole correspondante. La table de hachage est donc :
Alvéole 0 1 2 3 4 5 6 7 8 9 10
Valeurs 13 26 60 9 87
Nous pouvons remarquer que toutes les positions dans la table ne sont pas occupées. Par exemple,
ici, seules 5 alvéoles sur 11 sont utilisées. Le facteur de compression λ est définit comme la taille
du tableau initial n sur le nombre d’alvéoles k :
n
λ= (5.1)
k
Désormais, lorsque nous voulons vérifier qu’un élément est dans la table, nous devons calculer la
valeur de la fonction de hachage correspondante. Nous pouvons ensuite directement comparer si il
est présent dans la table de hachage. Nous disposons d’une méthode permettant de rechercher la
présence d’un élément avec une performance en O(1).
Néanmoins, vous avez déjà dû remarquer que cette technique ne fonctionne que lorsque chaque
élément de la table d’origine a une valeur de hachage différente. Dans le cas contraire, c’est à dire
que deux éléments ont la même position dans la table de hachage, nous avons un problème qu’il
nous faudra résoudre. Ce phénomène porte le nom de collision. Par exemple, si nous souhaitons
ajouter la valeur 43 à la table précédente. Ce nombre à la même valeur de hachage (10) que 87 et
devrait donc également se trouver dans l’alvéole 11. Nous ne pouvons donc pas le placer.

57
Les fonctions de hachage
Fonctions de hachage parfaites
Le choix de la fonction de hachage est crucial pour les performances de l’algorithme de
recherche par hachage. Nous avons déjà vu que son choix doit permettre de limiter le problème des
collisions. Une fonction de hachage qui associe une alvéole unique à chaque élément de la tableau
initial est appelé une fonction de hachage parfaite. Si nous connaissons préalablement le tableau,
et que la liste d’éléments ne changera jamais, nous pouvons construire une telle fonction. Dans le
cas général, ce n’est pas le cas et nous devrons gérer les collisions.
Nous pourrions penser qu’il suffit d’avoir une table de hachage assez grande pour pouvoir un
fonction de hachage parfaite. Par exemple, si la valeur de tous les éléments est comprise entre 0 à
99, il suffit d’utiliser une table de hachage avec 100 alvéoles. Néanmoins, cette méthode n’est pas
applicable dans le cas général car elle conduit à des tables de hachages trop grandes (par exemples
vos numéros de cartes d’étudiants sont sur 16 chiffres).
Nous ne pouvons donc pas ignorer les collisions. Heureusement, de nombreuses méthodes ont
été étudiées les gérer.

Fonctions de hachage sur des entiers

Pour résumer, Nous voulons une fonction de hachage qui minimise le nombre de collision, soit
facilement calculable et distribue uniformément les éléments dans la table de hachage. Nous avons
déjà vu le modulo. Nous allons en étudier quelques autres exemples.
La méthode par repliement consiste à diviser les éléments en morceaux de taille identique pour
ensuite les additionner. Une fonction modulo finale permet de garder le résultat inférieur à la taille
de la table de hachage. Imaginons, par exemple, que nos éléments soient des numéros de téléphone
comme 04 73 40 70 09. Nous pouvons les séparer en groupe de 2 chiffres (04,73,40,70,09) et
les additionner 04 + 73 + 40 + 70 + 09 = 196. Si nous gardons une table de hachage de même
dimension que dans notre exemple (11), l’alvéole correspondant au numéro de téléphone sera donc
mod(196, 11) = 9. Une légère variation de cette méthode inverse un nombre sur deux afin d’avoir
éventuellement une meilleure répartition dans la table de hachage (ce qui donnerait dans notre
exemple 04 + 37 + 40 + 07 + 09 = 97).
Une deuxième méthode est appelée le milieu du carré. Ici, nous commençons par porter le
nombre au carré avant de prendre la partie du milieu du résultat. Par exemple, pour le nombre 67,
nous obtenons 672 = 4489 et nous pouvons garder les deux chiffres du milieu (48) pour enfin fini
par le modulo mod(48, 11) = 4.

Fonctions de hachage sur des caractères

Dans tous les exemples que nous avons vu pour l’instant, nous avons utilisé des entiers positifs
en entrée de la fonction de hachage. Or ce n’est pas souvent le cas, comme par exemple, lorsque
l’on fait un recherche sur des noms.
Afin de pouvoir traiter des chaînes de caractères, nous allons revenir au cas des entiers positifs.
Par exemple, nous pouvons utiliser la position de chaque caractère dans la table ascii et utiliser
la méthode du repliement. Par exemple, pour le mot "maison", nous obtenons dans une table de
hachage de taille 11 :

Mot m a i s o n
Valeurs ascii 109 97 105 115 111 110
Valeur repliée 532
Modulo 11 4

Ce qui nous donne le programme Matlab suivant (A vérifier) :

58
1 f u n c t i o n resHach = hach ( chaine , t a i l l e T a b l e a u )
2 somme = 0
3 for i = 1: length ( chaine )
4 somme = somme + c h a i n e ( i )
5 end
6 r e s H a c h = mod ( somme , 1 1 )

Vous pouvez imaginer de nombreuses autres fonctions de hachages pour tous types d’éléments.
Il faut néanmoins se rappeler qu’elles doivent rester suffisamment simple pour le coût associé ne
rende pas son calcul plus long qu’une simple recherche dichotomique.
La résolution des collisions
Nous avons vu, qu’à part pour le cas des fonctions de hachage parfaites, il nous faut être capable
de gérer le problème des collisions. C’est-à-dire que nous devons trouver une solution lorsque deux
éléments ont la même valeur de hachage. C’est la résolution des collisions.

Adressage ouvert
La méthode la plus simple est de placer l’élément dans une autre alvéole qui n’est pas encore
utilisée. Cette solution s’appelle l’adressage ouvert. Nous pouvons, par exemple, rechercher si les
alvéoles suivantes sont libres. C’est la technique du sondage linéaire. Reprenons l’exemple du
début avec une fonction de hachage modulo, une table de hachage de taille m = 11, et un tableau
de valeur [26, 13, 60, 87, 9].

F IGURE 5.4 – Exemple de construction d’une table de hachage

La Fig. 5.4 montre quelques exemples d’insertions. Ainsi, si nous souhaitons ajouter la valeur
49, elle doit être dans l’alvéole 5 qui est déjà occupée. Elle est donc placée dans l’alvéole 6. De
même, la valeur 53 devrait être placé dans l’alvéole 9. Comme l’alvéole 9 et 10 sont occupée,
elle est placée dans l’alvéole 0. Finalement, la valeur 11 qui devrait se trouver dans l’alvéole 0 se
retrouve dans l’alvéole 1. Nous obtenons donc la table de hachage suivante :
Alvéole 0 1 2 3 4 5 6 7 8 9 10
Valeurs 53 11 13 26 60 49 9 87
Nous ne pouvons plus désormais savoir si un élément est, ou n’est pas, dans la table de hachage
avec une seule comparaison. Il nous faut utiliser la même méthode pet effectuer une recherche
linéaire tant que l’élément n’est pas trouvé ou que l’on tombe sur une alvéole vide. Par exemple,
si nous cherchons la valeur 20 (alvéole attendue 9), nous ne pouvons pas savoir qu’elle n’est pas
présente uniquement en comparant avec la valeur 9. nous devons rechercher dans les alvéoles
10,0,1,2 pour finalement être sûr de son absence avec l’alvéole vide 3.

59
 Un des problèmes du sondage linéaire est que les éléments ont tendance à former des grumeaux
dans les tables de hachage au fur et à mesure que les collisions s’accumulent. C’est à dire que les
données ont tendance à s’accumuler dans certaines zones de la table comme dans de la Fig. 5.5
avec 4 éléments qui ont comme valeur de hachage 4.

F IGURE 5.5 – Exemple de grumeau

Ces grumeaux peuvent décaler les données de leur alvéoles naturelle, comme par exemple lors
de l’insertion du nombre 11 dans la Fig. 5.4.
Afin de simplifier l’écriture, nous pouvons écrire le sondage linéaire sous une forme plus
générale en considérant une fonction de re-hachage rehash qui donne la nouvelle position npos à
partir de la position déjà occupée pos :

npos = rehash(pos) (5.2)

Dans le cas présent (sondage linéaire) et pour une table de hachage de taille tailleHach, cette
fonction peut s’écrire :

rehash : pos 7→ (pos + 1) modulo tailleHach (5.3)

La manière la plus simple d’étendre cette fonction afin d’éviter la formation de grumeaux est
d’utiliser un autre pas que 1, par exemple 4 :

rehash : pos 7→ (pos + 5) modulo tailleHach (5.4)

Il faut bien évidemment que la table de hachage ne soit pas un multiple de 5 pour pouvoir parcourir
tout ses emplacements afin de trouver une place libre. C’est pour cette raison que nous avions
utilisé une table de taille 11 (c’est un nombre premier !).
Une variation du sondage linéaire est le sondage quadratique, le pas n’est pas linéaire et croît à
chaque essai (1, 3, 5, puis 7...). Ce qui veux dire que l’algorithme cherchera une position de libre à
l’adresse pos + 1, puis pos + 4, puis pos + 9, etc. Cette méthode ajoute à chaque fois la valeur du
prochain carré.

Listes chaînées
Une méthode alternative pour gérer les collisions est de pouvoir placer plusieurs éléments dans
chaque alvéole de la table de hachage. Cette opération utilise une liste chaînée, c’est à dire une suite
de taille variable d’éléments liées entre eux. Sur la Fig. 5.6, l’exemple de la Fig. 5.4 est modifiée en
utilisant cette méthode. Ainsi, lorsqu’un nouvel élément est ajouté, sa valeur est toujours placé dans
l’alvéole correspondante. Elle est ajoutée dans la liste d’éléments. Il n’y a donc pas de problème de
de collisions.
Néanmoins, pour rechercher un élément après avoir obtenu l’alvéole correspondante avec la
fonction de hachage, faire une nouvelle recherche dans les éléments présent. La taille du problème
résultant doit être beaucoup plus petite que la taille du tableau initial pour que cette méthode soit
avantageuse.

60
 F IGURE 5.6 – Exemple de construction d’une table de hachage avec une liste chaînée

Résultats
Dans l’introduction de cette méthode, nous avions espéré trouver un algorithme de recherche
avec une performance en O(1). Nous avons vu que cette limite n’est atteinte que pour des fonctions
de hachage parfaites. Dans le pire cas, la fonction de hachage retourne la même alvéole pour toutes
les données d’entrée. La recherche devient une recherche séquentielle avec une performance en
O(n).
Dans le cas général, le phénomène des collisions dégrade les performances de l’algorithme. Les
calculs étant complexe, nous allons admettre les résultats. Ils utilisent le facteur de compression λ
définit dans l’équation 5.1. On montre que la performance de la méthode avec chaînage est :
T (n) = O (1 + λ ) (5.5)
On montre également que la performance avec un adressage ouvert et un sondage linéaire est :
 
1
T (n) = O (5.6)
1−λ

F IGURE 5.7 – Tracé de T (n) pour une table de hachage avec sondage linéaire et avec liste chaînée.

La dépendance en n est bien entendu dans le λ . Ces deux fonctions sont tracées sur la Fig. 5.7.
Sur ces courbes, nous pouvons avoir l’impression que la liste chaînée est meilleure que le sondage
linéaire. Il nous faut ici rappeler que les constantes sont importantes et que la gestion d’une liste
chaînée est beaucoup plus complexe qu’un simple tableau. Le sondage linéaire simple sera donc,
dans de nombreuses applications, plus rapide.

61
 5.3 Quelques algorithmes de tri
Les algorithmes de tri ont été un domaine de recherche très intense en informatique donnant
naissance à de nombreuses approches. L’opération de tri permet de placer des éléments suivant un
ordre défini. Par exemple, nous pouvons souhaiter trier des entiers par ordre croissant, ou encore
des cartons par poids ou par taille. Ces algorithmes sont souvent utilisés pour préparer les données
pour d’autres algorithmes (comme la recherche dichotomique).
De la même façon que pour les algorithmes de recherche, le choix d’un algorithme dépendra
de la taille de tableau à trier. S’il est petit le temps de calcul nécessaire à la mise en place de
procédure complexes ne vaudra pas la peine. Dans le cas contraires, nous voudrons avoir le plus de
raffinements possibles pour diminuer le temps de calcul. Nous allons étudier et comparer ici un
petit nombre d’algorithmes de tri de complexité croissante.
Avant de commencer, nous devons définir une quantité qui nous servira de mesure pour
comparer les différents algorithmes. Pour cela, nous allons nous baser sur deux quantités différente.
La première est le nombre de comparaisons, que nous avons déjà utilisé pour les algorithmes de
recherche. La deuxième est le nombre d’échanges entre deux éléments nécessaire à l’opération de
tri. En effet, c’est une opération coûteuse qui justifiera certain développements.

5.3.1 Le tri a bulle

Nous allons commencer par un algorithme simple : le tri à bulle. Cet algorithme parcourt la liste
d’éléments à trier en comparant à chaque fois les éléments adjacent. Si nécessaire, il les inverse
suivant leur valeur. A chaque passage, la valeur la plus haute restante est placée à sa place. Ainsi,
les valeurs remontent comme des "bulles". La liste est parcourue autant de fois que nécessaire.
Nous pouvons déjà observer que le parcours d’une liste de n éléments nécessite n − 1 compa-
raisons. Au bout de la première itération, le premier élément trié sera en place à la fin de la liste.
Ainsi, après n − 1 parcours, la liste est triée.
La Fig. 5.8 montre un exemple d’un passage de l’algorithme du tri à bulle pour une liste de 6
entiers (et donc 5 comparaisons). Nous pouvons observer qu’une fois que la valeur 91 est traitée,
elle est déplacée à chaque nouvelle comparaison jusqu’à sa place finale en bout de tableau.

F IGURE 5.8 – Exemple d’une passe du tri à bulle sur un tableau d’entier.

Nous obtenons le programme suivant :

1 function l i s t e = triBulle ( tableau )
2 n i t e r a t i o n = l e n g t h ( t a b l e a u ) −1;
3 f o r n p a s s = n i t e r a t i o n : −1:1
4 for i = 1: npass
5 i f t a b l e a u ( i ) > t a b l e a u ( i +1)

62
 6 temp = tableau ( i ) ;
7 tableau ( i ) = t a b l e a u ( i +1) ;
8 t a b l e a u ( i + 1 ) = temp ;
9 end
10 end
11 end
12 end

Cet algorithme nécessite d’échanger la valeur de deux cases du tableau tableau. Pour cela, en géné-
ral, cette opération d’échange utilise une variable intermédiaire (variable tmp dans le programme).
Nous pouvons également noter que nous avons mit le même nom de variable en entrée et en
sortie de la fonction (tableau). En effet, la valeur d’une entrée n’est pas modifiable en Matlab (c’est
ce que l’on appelle le passage par valeur). En mettant le même nom de variable en sortie, Matlab
autorise les modifications. De plus, depuis matlabR2010, la variable tableau n’est plus recopié en
mémoire lors de cette opération. C’est ce que l’on appelle un passage par référence.

Résultats
L’algorithme du tri à bulle a besoin de parcourir le liste n − 1 fois pour une liste de taille n
avant qu’elle ne soit triée. Nous pouvons en déduire que la performance T (n) est définie comme la
somme des n − 1 premiers entier :

n2 − n
T (n) = (5.7)
2

Ainsi, la performance de cet algorithme est donc en O(n2 ). En introduction, nous avons vu que
le nombre d’échange est également un paramètre important pour comparer les algorithmes de tri.
Dans le meilleur de cas, si la liste est déjà triée, il n’y en aura aucune. En moyenne, nous aurons un
échange dans la moitié des cas.
Le tri à bulle est souvent considéré comme une des pires méthode de tri. En effet, les éléments
de la liste sont échangés pas à pas jusqu’à leur position finale. Ces opérations sont très coûteuses.
Néanmoins, le tri à bulle a un avantage (en plus de sa simplicité) sur les autres algorithmes de tri :
il peut savoir qu’une liste est triée avant la fin normale de toutes les itérations. Si au court d’un
passage, aucun échange n’est fait alors la liste est trié.
Cette observation permet de modifier légèrement l’algorithme proposé afin d’en améliorer les
performances :
1 function tableau = triBulle_court ( tableau )
2 n i t e r a t i o n = l e n g t h ( t a b l e a u ) −1;
3 npass = niteration ;
4 echange = true ;
5 w h i l e n p a s s > 0 && e c h a n g e == t r u e
6 echange = f a l s e ;
7 for i = 1: npass
8 i f t a b l e a u ( i ) > t a b l e a u ( i +1)
9 temp = tableau ( i ) ;
10 tableau ( i ) = t a b l e a u ( i +1) ;
11 t a b l e a u ( i + 1 ) = temp ;
12 echange = t r u e ;
13 end
14 end
15 npass = npass − 1;

63
 16 end
17 end

5.3.2 Le tri par sélection

Contrairement au tri à bulle , le tri par sélection ne va faire qu’un échange par passage. Pour
cela, l’algorithme va uniquement rechercher la plus grande valeur non triée dans la liste d’élément,
puis il va le placer directement à sa position finale. Ainsi, après le premier passage, l’élément le
plus grand sera en place, après le second passage, le deuxième élément le plus grand sera en place,
etc.. Il nous faudra donc également n − 1 passage pour trier la liste.
La Fig. 5.10 montre un exemple de tri par sélection. Les 6 éléments de la liste sont triés en 5
passage de l’algorithme.

F IGURE 5.9 – Exemple de tri par sélection sur une liste d’entier.

L’implémentation Matlab correspondante pourrait être :

1 function tableau = triSelection ( tableau )
2 n i t e r a t i o n = l e n g t h ( t a b l e a u ) −1;
3 f o r p o s i t i o n = n i t e r a t i o n +1: −1:1
4 maxPosition = 1;
5 for i = 1: position
6 i f tableau ( i ) > l i s t e ( maxPosition )
7 maxPosition = i ;
8 end
9 end
10 temp = t a b l e a u ( p o s i t i o n ) ;
11 tableau ( p o s i t i o n ) = tableau ( maxPosition ) ;
12 t a b l e a u ( m a x P o s i t i o n ) = temp ;
13 end
14 end

La variable maxPosition mémorise l’emplacement du plus grand élément non trié. Une opération
d’échange permet de ensuite de le positionner.

64
 Résultats
Nous pouvons déjà observer qu’il nous faut le même nombre de comparaisons que pour le tri à
bulle. Donc, cet algorithme aura la même performance en O(n2 ). Néanmoins, comme il effectue
beaucoup moins d’échanges, cet algorithme sera généralement plus rapide.

5.3.3 Le tri par insertion

Le tri par insertion utilise une approche complètement différente. Il maintient une sous-liste
triée au début de la liste d’élément. Chaque nouvel élément est ajouté de façon à ce que la sous-liste
reste toujours triée.
Pour le comprendre, examinons l’exemple donné sur la Fig. ??. La partie grisée représente la
partie de la liste déjà triée. Nous commençons par créer un liste d’un élément (12), qui est forcément
déjà triée. Les éléments sont ensuite ajoutés un par un en décalant les éléments plus grand jusqu’à
leur position finale. Par exemple, pour ajouter 7 à la sous-liste [12,23,45,75,91], [23,45,75,91] sont
séquentiellement décalés et 7 inséré à la place laissée libre.
Notons que ce tri implique de faire souvent un décalage. C’est une opération qui reste coûteuse
mais beaucoup moins qu’une inversion (une affectation à la place de trois).

F IGURE 5.10 – Exemple de tri par insertion sur une liste d’entier.

Nous obtenons le code suivant :

1 function tableau = t r i I n s e r t i o n ( tableau )
2 niteration = length ( tableau ) ;
3 for index = 1: n i t e r a t i o n
4 valeur = tableau ( index ) ;
5 p o s i t i o n = index ;
6 w h i l e p o s i t i o n > 1 && t a b l e a u ( p o s i t i o n −1) > v a l e u r
7 t a b l e a u ( p o s i t i o n ) = t a b l e a u ( p o s i t i o n −1) ;
8 position = position − 1;
9 end
10 tableau ( position ) = valeur ;

65
 11 end
12 end

Résultats
On peut observer, notamment dans l’implémentation proposée du tri par insertion, que cet
algorithme nécessite n − 1 passage pour trier une liste de n éléments. Il commence par l’index 1 et
continue jusqu’à avoir intégré tous les éléments dans la liste trié. La ligne 8 de l’implémentation
effectue l’opération de décalage.
Il faut faire au maximum n − 1 comparaison pour trouver la place d’un élément, on retrouve une
performance en O(n2 ). Néanmoins, comme l’opération de décalage nécessite moins d’opération
que l’échange, ce tri ce comporte globalement mieux que les tri précédents.

5.3.4 Le tri shell

Le tri shell est une amélioration du tri par insertion. Pour cela, il va permettre de diminuer les
décalages. Il découpe la liste initiale en plusieurs sous-listes qui seront chacune triée par un tri
par insertion. Bien évidemment, Le choix de ses sous-listes est la clé de la performance de cet
algorithme. Pour cela, à la place d’utiliser des morceaux continus de la liste originale, les liste sont
formés d’éléments espacés d’un pas i (appelé le gap).
Dans l’exemple de la Fig. 5.12, la liste a 9 éléments et le gap est fixé à 3. On commence donc
par construire 3 sous-listes, qui seront chacune triées par un tri par insertion. Cette opération permet
de rapprocher chaque élément de sa position finale. Une dernière étape de tri par insertion permet
d’obtenir la liste triée avec beaucoup moins de décalage que ceux nécessaires sur la liste initiale.

F IGURE 5.11 – Exemple de tri shell sur une liste d’entier.

Cet exemple utilise deux étapes. Mais la performance de l’algorithme peut être encore améliorée
en effectuant plusieurs passes avec un gap de plus en plus petit. Ainsi, après le gap de 3, nous
aurions pu faire une étape supplémentaire avec un gap de 2 avant le tri par insertion final. C’est
cette méthode qui est utilisée dans le code suivant avec une suite de valeur de gap prédéfinie :
1 function tableau = triShell ( tableau )
2 gap = [ 7 0 1 , 3 0 1 , 1 3 2 , 5 7 , 2 3 , 1 0 , 4 , 1 ] ;
3 f o r i = gap
4 i f i < length ( tableau ) /2
5 for j = 1: i
6 tableau = triInsertionGap ( tableau , j , i ) ;

66
 7 end
8 end
9 end
10 end
11

12 f u n c t i o n t a b l e a u = t r i I n s e r t i o n G a p ( t a b l e a u , d e b u t , gap )
13 niteration = length ( tableau ) ;
14 f o r i n d e x = d e b u t + gap : gap : n i t e r a t i o n
15 valeur = tableau ( index ) ;
16 p o s i t i o n = index ;
17 w h i l e ( p o s i t i o n −gap ) > 0 && l i s t e ( p o s i t i o n −gap ) > v a l e u r
18 t a b l e a u ( p o s i t i o n ) = t a b l e a u ( p o s i t i o n −gap ) ;
19 p o s i t i o n = p o s i t i o n − gap ;
20 end
21 tableau ( position ) = valeur ;
22 end
23 end

Résultats
Bien que cet algorithme utilise le tri par insertion à chaque étape, ses performances sont
grandement améliorées. En effet, cela limite le nombre de comparaisons et de décalage à effectuer.
On peut montrer que les performances de cet algorithme dépendent fortement de la valeur des
gaps choisis. Le tableau suivant résume quelques une des lois qui ont été étudiées depuis sa création
(n est la taille du problème et k le numéro de l’itération) :
Formule du gap gaps Performance (pire cas) Année
E 2nk E n2 , E 4n , ... , 1 O(n2 )




1959
n n
O(n3/2 )



2E 2k+1 + 1 2E 4 , ... , 3 , 1 1960
k
2 −1 ...,63,31,15,7,3,1 3/2
O(n ) 1963
Nombre successifs ...,12,9,8,6,4,3,2,1 O(nlog2 n) 1971
de la forme 2 3 p q

1 puis 4k + 3 · 2k−1 + 1 ...,281,77,23,8,1 O(n4/3 ) 1986

Expérimentale 701, 301, 132, 57, 23, 10, 4, 1 Inconnue 2001

5.3.5 Le tri fusion

Le tri fusion est un exemple d’une approche dite "diviser pour conquérir". La liste d’élément
initiale va être découpée en deux puis redécoupée par des appels récursifs à la fonction de tri. Cette
opération sera effectuée jusqu’à obtenir des listes de un élément. Les listes seront ensuite fusionnées
deux à deux en une seule liste triée.
Un exemple de cet algorithme est donné sur la Fig. 5.12. La liste de 9 éléments est tout d’abord
découpée en 9 listes de 1 éléments. Nous pouvons observer que la taille n’a pas besoin d’être paire
au moment de la découpe. Les deux sous-listes résultantes peuvent avoir une taille qui diffère de 1.
Dans la seconde étape, la liste est reconstruite en triant à chaque fois les données.
On obtient ainsi l’implémentation possible suivante :
1 function l i s t e = triMerge ( l i s t e )
2 if length ( l i s t e ) > 1
3 milieu = floor ( length ( l i s t e ) /2) ;
4 l i s t e H a u t e = l i s t e ( m i l i e u + 1 : end ) ;
5 listeBasse = l i s t e (1: milieu ) ;

67
 F IGURE 5.12 – Exemple de tri fusion sur une liste d’entier.

6 listeHaute = triMerge ( listeHaute ) ;

7 listeBasse = triMerge ( listeBasse ) ;
8 i =1; j =1; k =1;
9 w h i l e i <= l e n g t h ( l i s t e H a u t e ) && j <= l e n g t h ( l i s t e B a s s e )
10 if listeHaute ( i ) < listeBasse ( j )
11 liste (k) = listeHaute ( i ) ;
12 i = i + 1;
13 else
14 liste (k) = listeBasse ( j ) ;
15 j = j + 1;
16 end
17 k = k + 1;
18 end
19 w h i l e i <= l e n g t h ( l i s t e H a u t e )
20 liste (k) = listeHaute ( i ) ;
21 i = i + 1;
22 k = k + 1;
23 end
24 w h i l e j <= l e n g t h ( l i s t e B a s s e )
25 liste (k) = listeBasse ( j ) ;
26 j = j + 1;
27 k = k + 1;
28 end

68
 29 end
30 end

La première partie de cette implémentation s’occupe du découpage (ligne 2 à 7) en appelant

récursivement la fonction triMerge. Les appels s’arrêtent lorsque la liste obtenue est de taille 1
(ligne 2).
A partir de ce point (ligne 8 à 29) les lignes sont re-fusionnées dans la liste qui a été passée lors
de l’appel de la fonction.
Nous pouvons observer que, à chaque appel de la fonction triMerge, le tableau est recopié dans
les deux variables listeHaute et listeBasse. Ainsi le surcoût mémoire de cet algorithme est loin
d’être négligeable et peut s’avérer problématique pour de grandes tailles de problèmes.

Résultats
L’analyse des performance du tri fusion se déroule en deux parties : le découpage et la fusion.
Nous avons déjà vu, dans la recherche dichotomique, que la performance de la division successive
en deux est en O(log(n)). La deuxième étape est la fusion. Chaque élément des deux sous-listes
va être placé dans la liste fusionnée. Il n’y a qu’un passage dans tous les éléments, il faut donc n
opération pour une liste de taille n. En combinant les deux, chaque division appelle une fusion. La
performance du tri fusion est donc en O(n log(n)).

5.3.6 Le tri rapide (Quick-sort)

Pour finir, nous allons étudier l’algorithme de tri implémenté dans Matlab. L’algorithme de tri
rapide utilise la même approche que le tri fusion (diviser pour régner) sans le surcoût mémoire de
ce dernier. Pour cela nous allons nous baser sur une valeur pivot pour découper le problème initial.
Les valeurs en dessous de la valeur pivot seront placé avant dans la liste, et les valeurs au dessus de
la valeur pivot seront placées après dans la liste. Nous aurons, après cette étape, deux partitions
(avant et après la valeur pivot) sur lesquelles nous pouvons ré-appliquer l’algorithme.
La performance de l’algorithme dépend fortement du choix des valeurs pivots. Dans ce cours,
nous allons prendre à chaque fois la première valeur de la liste. Un exemple est donné sur la Fig.
5.13. Dans la première étape, 45 est choisit comme valeur pivot. Les nombres en dessous de 45
sont placés dans la partie basse et les nombres au dessus sont placés au dessus. Les deux sous-listes
résultantes n’étant pas triée, l’algorithme de tri rapide est appliqué à quelqu’une d’entre elle.

F IGURE 5.13 – Exemple de tri rapide sur une liste d’entier.

69
 Nous pouvons proposer l’implémentation suivante :

1 function l i s t e = triRapide ( l i s t e )
2 l i s t e = sousTriRapide ( l i s t e , 1 , length ( l i s t e ) ) ;
3 end
4

5 f u n c t i o n l i s t e = sousTriRapide ( l i s t e , debut , f i n )
6 i f debut < f i n
7 [ l i s t e , p i v o t P o s i t i o n ] = p a r t i t i o n ( l i s t e , debut , f i n ) ;
8 l i s t e = s o u s T r i R a p i d e ( l i s t e , d e b u t , p i v o t P o s i t i o n −1) ;
9 l i s t e = s o u s T r i R a p i d e ( l i s t e , p i v o t P o s i t i o n +1 , f i n ) ;
10 end
11 end
12

13 f u n c t i o n [ l i s t e , p i v o t P o s i t i o n ] = p a r t i t i o n ( l i s t e , debut , f i n )
14 liste = liste ;
15 pivot = l i s t e ( debut ) ;
16 iGauche = debut + 1 ;
17 iDroite = fin ;
18

19 fini = false ;
20 w h i l e f i n i == f a l s e
21 w h i l e i G a u c h e <= i D r o i t e && l i s t e ( i G a u c h e ) <= p i v o t
22 iGauche = iGauche + 1 ;
23 end
24 w h i l e i G a u c h e <= i D r o i t e && l i s t e ( i D r o i t e ) >= p i v o t
25 iDroite = iDroite − 1;
26 end
27 i f iGauche > i D r o i t e
28 fini = true ;
29 else
30 temp = l i s t e ( i G a u c h e ) ;
31 l i s t e ( iGauche ) = l i s t e ( i D r o i t e ) ;
32 l i s t e ( i D r o i t e ) = temp ;
33 end
34 end
35 temp = l i s t e ( d e b u t ) ;
36 l i s t e ( debut ) = l i s t e ( i D r o i t e ) ;
37 l i s t e ( i D r o i t e ) = temp ;
38 pivotPosition = iDroite ;
39 end

Cette implémentation commence par placer le pivot à sa position pour ensuite placer les
éléments dans leur partition en utilisant la fonction partition. Le tri rapide est alors récursivement
appelé sur les deux sous-partitions.
Le tri des éléments se fait avec deux indices iDroite et iGauche qui parcourent la liste en sens
opposés. Les éléments pointés par iGauche et iDroite sont inversé quand les éléments pointés sont
respectivement supérieur et inférieur à la valeur pivot. Finalement, le pivot est placé à sa position
finale.

70
 Résultats
La performance du tri rapide dépend fortement des pivots utilisés. Si par un choix adéquat le
pivot se positionne toujours au milieu de chaque partition, alors nous avons la performance en
O(n log(n)) habituelle.
Dans le pire cas, c’est à dire si le pivot se retrouve toujours en première ou en dernière position,
alors la performance est en O(n2 ). Par exemple pour une liste déjà trié, si le pivot est choisit comme
la première valeur, comme dans notre exemple, alors nous serons dans le pire cas.

5.4 Conclusion
A partir de deux problèmes simples (comment trouver un élément dans une liste et comment
trier une liste), nous avons vu que de nombreuses approches différentes sont possibles. Elles ont
des performances différentes et une complexité différente.
Nous avons également remarqué que même les algorithmes les plus performants montrés ici ne
seront pas toujours les plus adaptés aux problèmes étudiés car les étapes de préparation peuvent se
révéler plus longues que le calcul lui même.

71
 Introduction au calcul numérique (Ref. [7])

6.1 Introduction
Après avoir étudié la recherche et le tri, deux problèmes clés de l’algorithmie, nous allons
introduire le calcul numérique. Pour cela, nous allons partir d’un problème physique (le circuit
RLC) pour lequel la solution exacte est connu et facilement calculable. A partir des équations
différentielles obtenues par une application de la loi des mailles, nous dériverons un programme
informatique qui nous permettra d’obtenir directement la solution.
Avant d’aborder ce problème, nous allons débuter ce chapitre en rappelant brièvement les
notions de précision et de stabilité. Nous retrouverons les problèmes de représentation des nombres
que nous avons abordé dans la première partie de ce cours.
Dans la seconde partie, nous traiterons l’exemple du circuit RLC. Pour cela, nous aborderons et
utiliserons une des méthodes les plus simple : les différences finies.

6.2 Erreur, stabilité et precision

Nous avons vu, dans le chapitre 2, que les variables utilisés dans les programmes informatiques
ont un type. Elles sont enregistrées dans la mémoire de l’ordinateur sur un certain nombres
d’octets (1 octet = 8 bits). Pour les nombres, nous pouvons distinguer deux grandes familles de
représentation :
— La virgule-fixe : ce sont uniquement les entiers en Matlab.
— La virgule-flottante : ce sont les types double ou single en Matlab.
Une variable entière a toujours une représentation exacte. Il n’y a pas d’erreur d’approximation.
Tous les calcules effectués sur des nombres entiers sont également exact si leur résultat reste dans
les bornes (par exemple entre 0 et 255 pour un entier sur 1 octet) et si le résultat d’une division est
le nombre entier résultant (le reste est ignoré).
Le cas des variables en virgule flottante est plus complexe. Nous avions vu qu’un nombre v de
type double ou single pouvait se mettre sous la forme :
v = s × M × 2e−E (6.1)
s est le signe, M est la mantisse (nombre entier), E le bias de l’exposant et e est l’exposant. Dans le
cas du type E vaut 127, la mantisse est sur 23 bits et l’exposant sur 8 bits.
Voici deux exemples de nombres de type single :

Nombre Signe Exposant Bit caché Mantisse

4 0 10000001 1 00000000000000000000000
10−8 0 01100100 1 01010111100110001110111
Si l’on effectue l’opération 4 + 10−8 avec des variables de type single en Matlab, l’opérateur le
plus faible (10−8 ) doit tout d’abord être avec le même exposant e que 4. 4 a pour exposant 129 (22 )
et 10−8 a pour exposant 100 (2−27 ). Cette opération se fait en décalant les bits de la mantisse vers
la droite. Dans ce cas, il faut décaler 29 fois, et la valeur 10−8 sera remplacée par 0.
Le nombre le plus petit qui peut être ajouté à 1 et en modifier la valeur est appelé la précision
de la machine εm . Pour des variables de type single, εm ≈ 6 · 10−8 . εm est uniquement relatif à
la mantisse et ne dépend que du nombre de bit qui la compose. Ce n’est pas le nombre minimal
représentable qui lui dépend du nombre de bit de l’exposant (1.2 · 10−38 pour le type single).

73
 Ainsi, toutes les opérations arithmétiques sur deux nombres en virgule flottante induirons une
erreur appelé erreur d’arrondi. Ces erreur s’accumulent
√ avec le nombre d’opérations. Malheureuse-
ment, l’erreur résultante ne sera souvent pas en Nεm (hypothèses d’erreurs totalement aléatoires)
pour principalement deux raisons :
— Il arrive fréquemment que toutes les erreurs soit dans la même direction, dans ce cas l’erreur
résultante sera en Nεm ;
— Certaines opérations augmentent drastiquement l’erreur d’arrondi.
Pour ce dernier cas, prenons l’exemple des équations du second degré de la forme ax2 + bx + c = 0.
Elles admettent deux solutions :
√
−b + b2 − 4ac
x= (6.2)
2a
L’erreur d’arrondi prend toute son importance lorsque b2 >> 4ac où toutes les solutions seront
assimilées à 0.
L’erreur d’arrondi est relative à la représentation des nombres dans un ordinateur. Elle ne peut
pas être évitée. Considérons un deuxième exemple, le nombre d’or φ :
√
5−1
φ= (6.3)
2
On peut montrer que nous pouvons obtenir la puissance n de φ (φ n ) par la relation de récurrence :

φ n = φ n−2 − φ n−1 (6.4)

La tableau suivant donne les valeurs obtenues pour un calcul en 32 bits :

n Valeur attendue (φ n ) Valeur obtenue (φ n = φ n−2 − φ n−1 )
2 0.38197 0.38197
3 0.23607 0.23607
4 0.1459 0.1459
5 0.09017 0.09017
6 0.055728 0.055728
7 0.034442 0.034442
8 0.021286 0.021286
9 0.013156 0.013156
10 0.0081306 0.0081297
11 0.005025 0.0050265
12 0.0031056 0.0031033
13 0.0019194 0.0019232
14 0.0011862 0.0011801
15 0.00073314 0.00074315
16 0.0004531 0.0004369
17 0.00028003 0.00030625
18 0.00017307 0.00013065
19 0.00010696 0.0001756
20 6.6107e-05 -4.4942e-05
La valeur trouvé s’éloigne de plus en plus de la valeur attendue, et nous ne pouvons rapidement
plus utiliser la fonction de récurrence pour calculer la valeur de φ n . Nous dirons que la fonction
diverge, dans ce cas elle est instable. Ce problème est, ici aussi, inhérent à la précision de la
machine.
Le dernier type d’erreur que nous allons aborder dans cette partie est l’l’erreur de troncation.
Il arrive couramment que les algorithmes numériques utilisent des approximations discrètes de

74
 fonctions continues. Par exemple, pour faire le calcul d’une intégrale, nous ferons la somme de la
valeur de la fonction en un nombre fini de point. Plus nous prendrons de points, plus le résultat sera
précis. C’est un paramètre ajustable de cet algorithme. La différence entre le résultat attendu et le
résultat obtenu est appelé l’l’erreur de troncation. Cette erreur sera toujours présente, même sur un
ordinateur qui représente les nombres avec une précision infinie.
Pour conclure cette partie, l’erreur d’arrondi est intrinsèque à l’utilisation d’un ordinateur. Nous
ne pouvons que nous assurer qu’elle n’induira pas d’instabilités ou de divergences. L’erreur de
troncation est sous la responsabilité du programmeur qui doit apprendre à la contrôler.

6.3 Exemple : circuit RLC

Pour finir ce cours, nous allons étudier un cas pratique que nous allons résoudre numériquement :
le circuit RLC (Fig. 6.1. Il s’agit de l’association en série d’une résistance, d’une capacité et d’une
bobine.

F IGURE 6.1 – Un circuit RLC.

L’équation du circuit s’écrit :

di 1
Z
Ve = Ri + L + idt (6.5)
dt C
i est le courant qui traverse le circuit (en Ampères), t est le temps (n secondes), Ve la tension fournie
par le générateur (en Volts), R la résistance (en Ohm), C la valeur de la capacité (en Farad) et L la
valeur de l’inductance (en Henry).
Nous souhaitons résoudre numériquement l’équation 6.5. Tout d’abord, nous pouvons remarquer
qu’elle est également facilement soluble analytiquement pour un échelon de tension (Ve = 0 pour
t < 0, Ve = E constant pour t > 0. Les solutions du problème dépendent du paramètre ∆ défini
comme :
L
∆ = R2 − 4 (6.6)
C
√ q
On définit également les deux variables ω0 = 1/ LC et η = 2 CL . Nous pouvons dès lors mettre
R

le discriminant sous la forme :

∆ = ω02 (η 2 − 1) (6.7)

Nous avons alors trois régimes :

— η < 1 : Régime oscillatoire amorti

E p
Uc (t) = E − p · e−ηω0 ·t · sin(ω0 1 − η 2t − φ ) (6.8)
1 − η2

75
 Avec :
p !
1 − η2
φ = −arctg (6.9)
η

— η = 1 : Régime critique

Uc (t) = E − E · e−ω0 ·t (1 + ω0 · t) (6.10)

— η > 1 : Régime apériodique

E
∗ p2 e p1t − p1 e p2t


Uc (t) = E + (6.11)
p1 − p2

Avec : p 2
p1 = −ηω0 + ω0 (η − 1) (6.12)
p 2
p2 = −ηω0 − ω0 (η − 1) (6.13)
(6.14)
Ce qui pourrait nous donner l’implémentation suivante en Matlab :
1 f u n c t i o n uc = r l c f _ a n a ( L , C , R , t , ve )
2 % i n i t i a l i s a t i o n des v a r i a b l e s
3 omega0 = 1 / s q r t ( L∗C ) ;
4 eta = R/ 2 ∗ s q r t ( C / L )
5 d e l t a = omega0 ^ 2 ∗ ( e t a ^2 −1) ;
6

7 % C a l c u l de l a s o l u t i o n
8 if eta < 1
9 p h i = −a t a n ( s q r t (1− e t a ^ 2 ) / e t a ) ;
10 uc = ve − ve / ( s q r t (1− e t a ^ 2 ) ) ∗ exp (− e t a ∗ omega0 ∗ t ) . ∗
s i n ( omega0 ∗ s q r t (1− e t a ^ 2 ) ∗ t −p h i ) ;
11 end
12

13 i f a b s ( e t a −1) < 1 e−8

14 uc = ve − ve ∗ exp (−omega0 ∗ t ) . ∗ ( omega0 ∗ t + 1 ) ;
15 end
16

17 if eta > 1
18 p1 = −e t a ∗ omega0 + omega0 ∗ s q r t ( e t a ^2 −1) ;
19 p2 = −e t a ∗ omega0 − omega0 ∗ s q r t ( e t a ^2 −1) ;
20 uc = ve + ve / ( p1−p2 ) ∗ ( p2 ∗ exp ( p1 ∗ t ) + p1 ∗ exp ( p2 ∗ t ) ) ;
21 end
22 end

Nous pouvons alors tester le circuit dans les trois cas précédemment définit :
1 %P a r a m e t r e s de l a s i m u l a t i o n
2 L = 11 e −3;
3 C = 2 . 2 e −6;
4 d t = 1 e −6;

76
 5 t0 = 0;
6 t n = 1 e −2;
7

8 t = t 0 + ( 1 : l e n g t h ( uc1 ) ) ∗ d t ;
9 %S i m u l a t i o n
10 R = s q r t (4∗L /C) / 5 ;
11 uc1 = r l c f ( L , C , R , d t , t 0 , t n ) ;
12 uc1_a= r l c f _ a n a (L , C, R, t , 5 ) ;
13

14 R = s q r t (4∗L /C) ;
15 uc2 = r l c f ( L , C , R , d t , t 0 , t n ) ;
16 uc2_a= r l c f _ a n a (L , C, R, t , 5 ) ;
17

18 R = s q r t (4∗L /C) ∗ 5 ;
19 uc3 = r l c f ( L , C , R , d t , t 0 , t n ) ;
20 uc3_a= r l c f _ a n a (L , C, R, t , 5 ) ;
21

22 %A f f i c h a g e
23 figure ;
24 h o l d on ;
25 p l o t ( t , uc1_a )
26 p l o t ( t , uc2_a )
27 p l o t ( t , uc3_a )
28 l e g e n d ({ ’ Amorti ’ , ’ C r i t i q u e ’ , ’ A p e r i o d i q u e ’ })
29 x l a b e l ( ’ temps ’ )
30 y l a b e l ( ’ Uc ’ )

Le résultat est donné Fig. 6.2.

F IGURE 6.2 – Résultat de la simulation.

77
6.3.1 Calcul de l’intégrale
L’intégration numérique, également appelé quadrature, regroupe une vaste famille d’algorithme.
Nous allons en étudier quelques un, les plus simples, qui sont plus des références historiques que
des méthodes modernes.
Le problème initial est de calculer la valeur I de l’intégrale d’une fonction f entre deux bornes
a et b pouvant s’écrire sous la forme :
Z b
I= f (t)dt (6.15)
a

Les méthodes que nous allons étudier supposent que l’on connaît la valeur de la fonction f en un
certain nombre de points t0 ,t1 ,t2 , ...tn régulièrement espacés d’un pas h. On peut alors écrire, avec
i = 0, 1, ..., n :

ti = t0 + i × h (6.16)

Nous allons chercher à intégrer f par intervalle. Pour cela, nous pouvons simplement supposer que
la fonction est constante sur l’intervalle [t1 ,t2 [. C’est l’approximation des rectangles :
Z t2
f (t)dt ≈ h f (t1 ) + O(h2 f 0 ) (6.17)
t1

La notation en O() signifie que l’erreur sur l’estimation de la valeur de l’intégrale est de l’ordre de
grandeur du coefficient h2 multiplié par la dérivé première de f . Cette formule peut être modifiée
en supposant que la fonction est linéaire sur l’intervalle [t1 ,t2 [. Cette approximation est appelée
l’approximation par trapèze :
Z t2  
f (t1 ) + f (t2 )
f (t)dt ≈ h + O(h3 f 00 ) (6.18)
t1 2

Cette nouvelle équation fait intervenir deux point au lieu d’un seul. Nous pouvons continuer à
augmenter l’ordre de 1, c’est la formule de Simpson :
Z t3  
1 4 1
f (t)dt ≈ h f (t1 ) + f (t2 ) + f (t3 ) + O(h5 f (3) ) (6.19)
t1 3 3 3

Il est bien évident qu’il existe des formules pour des ordres encore supérieur. Néanmoins, nous
nous arrêterons ici dans le cadre de ce cours.
A partir des formules que nous avons définie pour un intervalle [t1 ,t2 [, nous devons étendre
l’intégration aux n intervalles entre t0 et tn . Nous aurons ainsi, pour la forme étendue des trapèzes :

(tn − t1 )3 00
Z tn    
1 1
f (t)dt ≈ h f (t0 ) + f (t1 ) + f (t2 ) + ... + f (tn−1 ) + f (tn ) + O f (6.20)
t0 2 2 n2

Comme les bornes d’intégration sont fixes, dans la suite nous ne marquerons dans le O() que
la variation en n. C’est à dire uniquement la variation dépendante du choix que nous faisons en
découpant l’intervalle. Nous pouvons bien entendu étendre également la formule de Simpson :
Z tn 
1 4 2 4 2
f (t)dt ≈ h f (t0 ) + f (t1 ) + f (t2 ) + f (t3 ) + f (t4 )... (6.21)
t0 3 3 3 3 3
  
2 4 1 1
+ f (tn−2 ) + f (tn−1 ) + f (tn ) + O 4 (6.22)
3 3 3+ n

78
6.3.2 Calcul de la dérivée
Pour calculer la dérivée d’une fonction, nous allons commencer par partir de la définition :
f (t + h) − f (t)
f 0 (t) = lim (6.23)
h→0 h
Il suffit donc de reprendre la formulation précédente :
f (ti ) − f (t − i − 1)
f 0 (ti ) ≈ (6.24)
h
Cette formule a deux sources d’erreurs principales : l’erreur de troncation et une erreur d’approxi-
mation.
Afin de déterminer l’erreur de troncation, nous pouvons utiliser un développement de la fonction
f en série de Taylor. Nous avons alors :
1 1
f (t + h) = f (t) + h f 0 (t) + h2 f 00 (t) + h3 f 000 (t) + ... (6.25)
2 6
On a donc :
f (t + h) − f (t) 1
= f 0 (t) + h f 00 (t)+ (6.26)
h 2
Ce qui nous donne l’erreur de troncation et :

et ≈ |h f 00 (t)| (6.27)

Pour ce qui est de l’erreur d’approximation ea , nous avons vu que l’addition de deux nombre
était sensible à la précision de la machine εm . Pour des fonctions assez simples (ou l’erreur sur la
fonction est de l’ordre εm ), on peut montrer que :
f 00 (t)
ea ≈ εm | | (6.28)
h
On trouve alors un optimum pour la valeur de h qui minimise ces deux erreurs :
s
εm f √
h≈ ≈ εm xc (6.29)
f 00

où xc = sqrt f (t)/ f 00 (t) est l’échelle caractéristique de la fonction. En l’absence de ces informations,
on suppose que xc = x.

6.3.3 Retour au circuit RLC

Nous allons appliquer ces méthodes à partir de l’équation définie précédemment :
di(t) 1
Z
Ve (t) = Ri(t) + L + i(t)dt (6.30)
dt C
On définit un intervalle de temps ∆t et un instant t0 de début de la simulation. Afin de simplifier
l’écriture, nous écrirons :

i(n) = i(t0 + n × ∆t) (6.31)

Nous pouvons alors écrire la dérivation comme :

di(n) i(n) − i(n − 1)
= (6.32)
dt ∆t

79
 Nous écrirons l’intégrale comme :
Z t0 +n∆t
S(n) = i(t)dt (6.33)
t0

Si nous appliquons la formule des triangles étendues, nous aurons :

 
1 1
S(n) = ∆t i(0) + i(1) + ... + i(n − 1) + i(n) (6.34)
2 2
Comme nous avons également :
 
1 1
S(n − 1) = ∆t i(0) + i(1) + ... + i(n − 2) + i(n − 1) (6.35)
2 2
Nous en déduisons :
 
1 1
S(n) = S(n − 1) + ∆t i(n − 1) + i(n) (6.36)
2 2
En posant S(0) = 0, nous retrouvons bien la formule habituelle des triangles :
 
1 1
S(1) = S(0) + ∆t i(0) + i(1) (6.37)
2 2
 
1
= ∆t (i(0) + i(1)) (6.38)
2
Nous avons donc une nouvelle forme pour l’équation 6.30 :
i(n) − i(n − 1) 1
Ve (n) = Ri(n) + L + S(n) (6.39)
∆t C  
i(n) − i(n − 1) 1 1 1
= Ri(n) + L + S(n − 1) + ∆t i(n − 1) + i(n) (6.40)
∆t C 2 2
Que nous pouvons mettre sous la forme :
   
1 L ∆t 1
i(n) = Ve (n) + − i(n − 1) − S(n − 1) (6.41)
R + ∆tL + 2C
∆t ∆t 2C C
Nous obtenons donc deux relations de récurrences (les équations 6.36 et 6.41). Elles permettent, a
partir de la connaissance de i(n − 1) et S(n − 1) de calculer les nouvelles valeurs de i(n) et S(n).
Ce qui pourrait nous donner l’implémentation suivante en Matlab :
1 f u n c t i o n uc = r l c f ( L , C , R , d t , t 0 , t n )
2 % i n i t i a l i s a t i o n des v a l i a b l e s
3 n i t e r = f l o o r ( ( t n −t 0 ) / d t ) ;
4 s = zeros (1 , n i t e r ) ;
5 i = zeros (1 , n i t e r ) ;
6 a l p h a = 1 / ( R+L / d t + d t / 2 ∗C ) ;
7 b e t a = ( L / d t −d t / 2 ∗C ) ;
8 % b o u c l e de r e c u r r e n c e
9 for n = 2: n i t e r
10 i ( n ) = a l p h a ∗ ( ve ( n∗ d t ) + b e t a ∗ i ( n −1) − 1 / C∗ s ( n −1) ) ;
11 s ( n ) = s ( n −1) + d t ∗ ( 1 / 2 ∗ i ( n −1) + 1 / 2 ∗ i ( n ) ) ;
12 end
13 uc = s ∗ 1 /C ;
14 end

80
 R
Cette fonction calcule les valeurs de i et de i itérativement en faisant croître le temps. A partir de
ces valeur, elle retourne la tension au borne de la capacité Uc . Nous avons besoin de définir une
tension Ve , par exemple un échelon :
1 function r e s = ve ( t )
2 if t > 2 e−3
3 res = 5;
4 else
5 res = 0;
6 end
7 end

Nous pouvons alors tester le circuit dans les trois cas précédemment définit :
1 %P a r a m e t r e s de l a s i m u l a t i o n
2 L = 11 e −3;
3 C = 2 . 2 e −6;
4 d t = 1 e −6;
5 t0 = 0;
6 t n = 1 e −2;
7

8 %S i m u l a t i o n
9 R = s q r t (4∗L /C) / 5 ;
10 uc1 = r l c f ( L , C , R , d t , t 0 , t n ) ;
11

12 R = s q r t (4∗L /C) ;
13 uc2 = r l c f ( L , C , R , d t , t 0 , t n ) ;
14

15 R = s q r t (4∗L /C) ∗ 5 ;
16 uc3 = r l c f ( L , C , R , d t , t 0 , t n ) ;
17

18 %A f f i c h a g e
19 t = t 0 + ( 1 : l e n g t h ( uc1 ) ) ∗ d t ;
20 figure ;
21 h o l d on ;
22 p l o t ( t , uc1 )
23 p l o t ( t , uc2 )
24 p l o t ( t , uc3 )
25 l e g e n d ({ ’ Amorti ’ , ’ C r i t i q u e ’ , ’ A p e r i o d i q u e ’ })
26 x l a b e l ( ’ temps ’ )
27 y l a b e l ( ’ Uc ’ )

Le résultat est donné Fig. 6.4. Nous retrouvons bien les même courbes qu’en utilisant les formules
analytiques (Fig. 6.2).
Nous n’avons pas discuté de la valeur de ∆t. Elle est directement liée à l’erreur de troncation.
Par exemple, la Fig. montre exactement la même figure avec ∆t = 10−4 et ∆t = 10−3 . Nous pouvons
observer que le résultat devient de plus en plus faux. Il faut que ∆t soit petit devant les phénomènes
étudiés (par exemple les oscillations dans ce cas) et qu’il permette d’effectuer le calcul dans un
temps raisonnable. Dans cet exemple simple, nous pouvons le prendre aussi petit que nous le
voulons, le temps de calcul restera faible. C’est moins vrai pour calculer la météo en france.
Nous pouvons nous demander l’intérêt de cet approche puisque nous connaissons déjà la

81
 F IGURE 6.3 – Résultat de la simulation.

F IGURE 6.4 – Résultat de la simulation avec ∆t = 10−4 (gauche) et ∆t = 10−3 (droite).

solution analytique du problème. Prenons un exemple différent en remplaçant l’échelon sur Ve

par Ve (t) = cos(ωt) ∗ log(ωt)/ωt. La solution de l’équation différentielle devient beaucoup plus
difficile à trouver alors que la solution numérique est triviale (Fig. 6.5).
Nous pouvons alors nous demander l’intérêt des solutions analytiques puisque la solution
numérique permet de résoudre plus simplement une partie des problèmes posés. Elle permet de :
— Comprendre les phénomènes observés ;
— Avoir une solution qui n’est pas sensibles aux erreurs numériques ;
— Servir de point de comparaison pour les approches numériques.
Ce sont donc deux approches complémentaires.

6.3.4 Calcul du sinus

Lors de l’exemple précédent, nous avons directement utilisé les fonction cosinus et logarithme
sans nous occuper de la manière dont ces fonctions sont calculées. Un micro-processeur n’est pas

82
 F IGURE 6.5 – Résultat d’une autre simulation.

capable de les produire directement ce résultat. Ce sont des algorithmes.

Une des méthode pour évaluer la valeur d’une fonction est d’utiliser son développement en
série entière. En effet, les fonctions analytiques peuvent s’écrire sous la forme au voisinage de x0 :
∞
f (x) = ∑ ak (x − x0 )k (6.42)
k=0

Ces fonctions sont assez facile à évaluer.

Par exemple, le développement en série entière du sinus en 0 est donné par la fonction suivante :

∞
(−1)k
sin(x) = ∑ (2k + 1)! x2k+1 (6.43)
k=0

Son rayon de convergence est infini, c’est à dire que nous pouvons déterminer la valeur du sinus
avec cette fonction quelque soit la valeur de x. Malheureusement, nous ne pouvons pas effectuer
une somme jusqu’à l’infini et nous devrons nous contenter de la somme partielle :
n
(−1)k 2k+1
sin(x) = ∑ x (6.44)
k=0 (2k + 1)!

Il faut alors savoir pour quelle valeur de n nous pouvons arrêter de calculer la série. Ce choix est
lié à l’erreur de troncation. Cette erreur sera toujours présente, quelque soit n. Un critère usuel est
lorsque le nouvel élément ajouté est inférieur à la précision souhaitée par le programmeur.
Nous pouvons utiliser, par exemple, l’implémentation suivante pour déterminer la valeur du
sinus :
1 function res=taylorSin (x , n)
2 k = 0;
3 xk = x;

83
 4 x2 = x.^2;
5 fact = 1;
6 signe= 1;
7 r e s = s i g n e / f a c t ∗ xk ;
8 for k = 1: n
9 s i g n e = −1.∗ s i g n e ;
10 xk = xk . ∗ x2 ;
11 f a c t = f a c t ∗ (2∗ k ) ∗(2∗ k +1) ;
12 r e s = r e s + s i g n e / f a c t . ∗ xk ;
13 end
14 end

Nous pouvons remarquer que nous n’effectuons bien évidemment pas la puissance de k à chaque
étape et que nous utilisons ici une procédure itérative en multipliant par (x − x0 ) le résultat précédent
(ligne 9).
Notons que, pour que ces formules soient utile, il faut que la somme partielle converge rapide-
ment vers la valeur de sin(x). Ce n’est as le cas pour des grandes valeurs de x. La série de Taylor ne
converge pas, dans le cas du sinus, avant que k  |x| (Fig. 6.6.

F IGURE 6.6 – Approximation du sinus.

Nous pouvons, bien évidemment, utiliser une des propriétés du sinus et replacer x dans l’inter-
valle [0, 2π[ pour améliorer ces résultats.

6.3.5 Conclusion
Dans cette partie du document, nous avons présenté un algorithme de résolution d’une équation
différentielle représentant un oscillateur RLC. Cet algorithme relativement simple (différence
finie) nous a permit de calculer la réponse du système dans des cas complexes ou une résolution
analytique est compliquée voir impossible.
Ainsi nous avons abordé :
— Les sources d’erreurs de calculs ;
— Le calcul des dérivés et des intégrales ;

84
— L’approximation de fonctions.

85
 Bibliographie

[1] Fabien Baillon et Jean-Louis Dirion, Initiation à Matlab, pour la résolution de problèmes
numériques, Mines Albi
[2] Documentation Matlab
[3] Stuart Madnick, Little Man Computer
[4] Sylvain Montagny, Architecture des ordinateurs, Université de Savoie
[5] George T. Heineman, Gary Pollice, Stanley Selkow, Algorithms in a Nutshell, 2nd edition
[6] Brad Miller, David Ranum, Problem Solving with Algorithms and Data Structures
[7] W.H. Press, S.A. Teukolsky, W.T. Vetterling, B.P. Flannery, Numerical Recipes in C

87