Les programmes de 2nde, 1ère et Tle du bac 2021 sont constellés d’exemples d’algorithmes. Ceux-ci permettent une réflexion sur le vocabulaire mathématique et les repères historiques eux aussi au programme
Cet article peut être librement diffusé à l’identique dans la limite d’une utilisation non commerciale suivant la licence CC-by-nc-nd
Les nouveaux programmes le précisent bien, « le langage choisi est Python ». Il est donc logique de décrire les algorithmes sous forme de scripts Python. Par souci de cohérence avec SNT et NSI, ce sera le plus souvent Python3 qui sera choisi ici. De plus, le programme ne fait plus aucune allusion à des notions comme l’entrée textuelle de données ou l’affichage ; on leur substitue celle de fonction :
- Au lieu d’entrer des données au clavier, on choisit les antécédents à fournir à la fonction ;
- Au lieu d’afficher des résultats, on calcule l’image des antécédents par la fonction.
Ainsi, on remplacera presque systématiquement « algorithme » par « fonction » dans cet article. Au fait, l’article est placé sous licence CC-by-SA et les scripts Python sous licence Gnu/Gpl ce qui veut dire qu’on a le droit de les copier-coller (d’ailleurs un lien de téléchargement est placé sous chaque script Python).
Les notations mathématiques sont listées à la fin du programme mais présentées comme transversales. Il en est de même pour les repères historiques, qui apparaissent conformément au rapport Torossian-Villani. Pour chaque algorithme, on essayera autant que possible de situer chronologiquement les problèmes étudiés, et de relier autant que possible les objets Python utilisés, aux objets mathématiques qu’ils permettent de définir.
Seconde
Voir aussi cet article.
Réels
En Python, chaque nombre réel est représenté par une approximation par une fraction dyadique, c’est-à-dire dont le dénominateur est une puissance de 2. Ces nombres s’appellent des flottants et
- En entrant
type(3.0)
on apprend que 3 est un flottant (quotient de 3 par 1 qui est bien une puissance de 2) ; - pour convertir 3 (qui est en réalité un entier) en flottant, on lui applique la fonction de conversion en flottant, avec
float(3)
.
Comme certains nombres décimaux ne sont pas dyadiques, l’égalité entre expressions faisant intervenir des nombres décimaux (comme 0,3 et 0,1×3 par exemple) n’est pas assurée. Il est donc prudent, au prix d’une complexité accrue des scripts, de remplacer les tests d’égalité entre réels a et b, par des tests du type (a-b)²<10-12. Cela concerne notamment les tests de colinéarité et d’orthogonalité de vecteurs et le test d’indépendance entre événements.
Déterminer par balayage un encadrement de √2 d’amplitude inférieure ou égale à 10−n
La question est qu’est-ce que c’est qu’un encadrement ?. On propose de définir cette notion comme un intervalle, ou comme la donnée d’une approximation décimale par défaut et une autre approximation, par excès.
Puisque le résultat dépend de l’entier n, on gagne à programmer cela comme une fonction. Celle-ci renverra un couple de réels (de flottants, dans le vocabulaire de Python) :
- def balayage(n):
- x = 1
- while x**2<2:
- x += 10**(-n)
- return [x-10**(-n),x]
Quelle que soit la valeur de n ∈ N, le réel (en fait, décimal) 10−n est strictement positif, donc la suite arithmétique parcourue par cet algorithme est strictement croissante et tend vers l’infini. Elle dépassera tout seuil strictement positif fixé à l’avance au bout d’un temps fini. Cette remarque est dûe à Archimède qui fut le premier à ne pas considérer cette propriété comme évidente.
Entiers
Déterminer si un entier naturel a est multiple d’un entier naturel b
Les élèves ont vu au cycle 4 des critères de divisibilité par 2, par 3, par 5, par 9 et par 10. Mais pour savoir si un entier est divisible par 7, on ne connaît pas de meilleure méthode que de faire la division euclidienne par 7 et regarder s’il y a un reste. Or le reste dans la division euclidienne de a par b se note a%b en Python. La fonction suivante répond donc au problème :
- def est_multiple(a,b):
- return a%b==0
La division euclidienne tire son nom d’Euclide qui, dans son livre VII, a décrit un algorithme de division par soustractions successives, défini les nombres premiers (voir plus bas) et appliqué son algorithme de division à des tests de divisibilité, et de primalité. Si Euclide avait programmé en Python, il aurait probablement plutôt écrit quelque chose comme :
- def est_multiple(a,b):
- while a>=0:
- a -= b
- a += b
- return a==0
La ligne 4 vient de ce que la sortie de la boucle se fait trop tard, lorsque a est devenu négatif, et il est nécessaire de compenser ce pas en arrière de trop en avançant à nouveau, une dernière fois.
Pour des entiers a et b donnés, déterminer le plus grand multiple de b inférieur ou égal à a
Remarque : Les entiers sont supposés naturels, et l’entier a est supposé non nul.
Traduction en robotique : Un robot, appelé Tortue, avance par pas de b, mais ne soit pas dépasser la distance a. Où se trouve-t-il lorsqu’il est allé le plus loin possible ?
L’algorithme ci-dessus a pour effet que le robot est allé trop loin lorsqu’on sort de la boucle. Il lui est donc nécessaire de reculer après avoir quitté la boucle, afin que la condition « multiple de b inférieur ou égal à a » soit respectée.
Transposé dans le champ numérique et en Python, on obtient ceci :
- def pgm(b,a):
- multiple_de_b = 0
- while multiple_de_b<=a:
- multiple_de_b += b
- multiple_de_b -= b
- return multiple_de_b
La suite des multiples de b est arithmétique de raison b. Ce fait est utilisé dans l’algorithme. L’explication de la ligne 5 est similaire à l’algorithme précédent, et l’ajout d’une variable compteur de boucle (initialisée à 0, incrémentée de 1 à chaque passage dans la boucle) permet de transformer cet algorithme, en un algorithme de division euclidienne (calcul du quotient).
Pour cet exemple, il est difficile de faire plus court que
- def pgm(b,a):
- return int(a/b)*b
Remarque : Comme évoqué dans l’onglet précédent, l’existence et l’unicité de ce multiple de b, a été découverte par Archimède, ce qui amène à dire que N est archimédien. R l’est également.
Déterminer si un entier naturel est premier
- def est_premier(n):
- for d in range(2,n):
- if est_multiple(n,d):
- return False
- return n>1
La ligne 4 quitte la boucle en renvoyant la valeur False dès qu’un diviseur de n est trouvé (car s’il existe un tel diviseur, n est composé, et pas premier). La ligne 5 n’est donc atteinte que si aucun diviseur n’a été trouvé, ce qui veut dire que n est premier. Renvoyer True aurait alors été inapproprié dans le cas particulier n=1, mais dans tous les autres cas n>1 vaut True comme souhaité.
Cet algorithme provient lui aussi du livre VII d’Euclide évoqué plus haut.
Calcul
Déterminer la première puissance d’un nombre positif donné supérieure ou inférieure à une valeur donnée
Si le nombre est supérieur à 1, on applique l’algorithme suivant :
- def pppss(nombre_positif,seuil):
- assert nombre_positif>1
- puissance = 1
- while puissance<=seuil:
- puissance *= nombre_positif
- return puissance
La ligne 2 est un garde-fou, destiné à ne pas boucler s’il n’y a pas lieu de boucler : Si le nombre est supérieur à 1, la suite géométrique de ses puissances tend vers l’infini et finira par dépasser le seuil. Sinon ce n’est pas garanti.
Une autre manière de dire cela est de considérer que la fonction pppss (« plus petite puissance supérieure au seuil ») n’est pas définie pour les nombres inférieurs à 1.
Si le nombre est inférieur à 1, cela a plus de sens de chercher la plus petite puissance inférieure au seuil. On s’assure au passage que le nombre est bel et bien positif :
- def pppis(nombre_positif,seuil):
- assert 0<nombre_positif<1
- puissance = 1
- while puissance>=seuil:
- puissance *= nombre_positif
- return puissance
Vecteurs
Le type tuple de Python n’est pas au programme. Cependant, il est question de « calculer les coordonnées d’un vecteur » et comme il est d’usage de grouper les coordonnées des points (et des vecteurs parfois) dans des couples de nombres, on se permettra ici de représenter les points et les vecteurs du plan par des couples de flottants. Aucun exemple d’algorithme n’est proposé sur les vecteurs, mais on va traiter ici des capacités attendues, en les traitant en Python. Cela a un double avantage :
- Meilleure cohérence entre les scripts Python ci-dessous et les notations mathématiques du cours ;
- Les scripts ci-dessous facilitent grandement la rédaction de ceux portant sur l’alignement et les équations de droite.
Pour se simplifier la vie, on suppose que les fonctions du module math ont été importées par from math import *
. En fait c’est surtout la fonction hypot qui sera utilisée.
Vecteur
Pour calculer les coordonnées du vecteur allant de A à B, on fait ainsi :
- def vecteur(A,B):
- (xA,yA) = A
- (xB,yB) = B
- return (xB-xA,yB-yA)
Au début du 19e siècle, Jean-Robert Argand notait par un simple trait (sans pointe de flèche) surmontant « AB », le vecteur allant de A à B. Il ne notait pas (3,2) le vecteur de coordonnées 3 et 2, mais il le notait 3+2√(-1). Giusto Bellavitis (avec ses « équipollents ») puis Hermann Grassmann (créateur de la notion de déterminant, entre autres) ont introduit vers le milieu du 19e siècle la notion de vecteur et le calcul vectoriel.
Colinéarité
Si un vecteur u est un couple de nombres, l’abscisse de u se note u[0] en Python, et l’ordonnée de u se note u[1]. Ce qui permet d’écrire assez facilement un test de colinéarité :
- def sont_colinéaires(u,v):
- return u[0]*v[1]==u[1]*v[0]
Cette fonction renvoie True si u et v sont colinéaires, False sinon.
Remarque : Le phénomène évoqué dans cet article a pour conséquence que ce test échoue parfois, en raison des erreurs d’approximation. on propose plus bas une variante donnant de meilleurs résultats, et basée sur la notion de déterminant.
Plutôt que toutes ces écritures avec crochets, on préférera ici nommer explicitement les coordonnées des points et des vecteurs comme on l’a déjà fait pour le vecteur entre deux points ci-dessus.
Opérations
Pour calculer la somme de deux vecteurs, par exemple, on fait ainsi :
- def somme(u,v):
- (xu,yu) = u
- (xv,yv) = v
- return (xu+xv,yu+yv)
La somme de deux vecteurs est un vecteur, représenté ici par un couple de flottants. C’est un peu ce à quoi on s’attend lorsqu’on additionne deux vecteurs : à obtenir un vecteur.
Pour multiplier un vecteur par un réel, on peut faire ainsi :
- def produit_par_un_flottant(v,k):
- (xv,yv) = v
- return (xv*k,yv*k)
Norme et distance
Pour calculer la norme d’un vecteur, on applique à ses coordonnées la fonction hypot du module math. Comme son nom l’indique, elle utilise le théorème de Pythagore pour calculer l’hypoténuse d’un triangle à partir des côtés de l’angle droit. C’est exactement ce qu’il faut pour calculer une norme :
- def norme(u):
- (xu,yu) = u
- return hypot(xu,yu)
On objectera qu’il ne s’agit pas en cours de maths, de former des experts en Python. Je contre-objecterai que les fonctionnalités de Python qui facilitent le cours de maths valent la peine qu’on s’y attarde un peu pour aller plus vite après.
Avec ça c’est très court de calculer la distance entre deux points A et B dans un repère orthonormé :
- def distance(A,B):
- return norme(vecteur(A,B))
Par exemple distance( (2,1) , (3,-2) )
permet de calculer la distance entre les points de coordonnées (2,1) et (3,-2).
Milieu
Pour calculer les coordonnées du milieu, on peut créer une fonction moyenne qui rend un peu plus naturel le calcul des coordonnées du milieu : Ce sont la moyenne des abscisses et la moyenne des ordonnées :
- def moyenne(x,y):
- return (x+y)/2
- def milieu(A,B):
- (xA,yA)=A
- (xB,yB)=B
- return (moyenne(xA,xB),moyenne(yA,yB))
Déterminant
Pour calculer le déterminant de deux vecteurs, on peut modifier le test de colinéarité ci-dessus (remplacer l’égalité par une soustraction) ou faire
- def déterminant(u,v):
- (xu,yu)=u
- (xv,yv)=v
- return xu*yv-xv*yu
Inversement, on peut utiliser le déterminant pour améliorer le test de colinéarité : Deux vecteurs sont colinéaires uniquement lorsque leur déterminant vaut zéro. Mais comme l’égalité entre deux flottants n’est pas toujours exacte, on se propose de remplacer le test de nullité du déterminant par un test portant sur son carré : On peut estimer que si le carré du déterminant est inférieur à 10-12 c’est que le déterminant est pratiquement égal à zéro. Voici le test de colinéarité amélioré :
- def sont_colinéaires(u,v):
- return déterminant(u,v)**2 < 1e-12
Droites
Étudier l’alignement de trois points dans le plan
Ce qui a été fait dans l’onglet précédent sur les vecteurs permet de vite résoudre ce problème en Python :
- def sont_alignés(A,B,C):
- return sont_colinéaires(vecteur(A,B),vecteur(A,C))
Il ne s’agit guère que d’une traduction en Python de la définition du cours.
Déterminer une équation de droite passant par deux points donnés
La première question à poser est « qu’est-ce que c’est qu’une équation de droite ? ». Voici ce que dit wikipedia, sur la notion d’équation cartésienne :
En géométrie analytique, les solutions d’une équation E d’inconnues x et y peuvent être interprétées comme un ensemble de points M(x, y) du plan affine, rapporté à un repère cartésien. Quand ces points forment une courbe, on dit que E est une équation cartésienne de cette courbe.
Toujours selon wikipedia, une équation est « une relation contenant une ou plusieurs variables ». On demande donc une fonction Python qui détermine une relation entre deux variables x et y. Mais au fait qu’est-ce que c’est qu’une relation ? Selon Gottlob Frege (seconde moitié du 19e siècle), une relation est assimilable à un prédicat, c’est-à-dire une fonction qui, au couple (x,y), associe la valeur True ou False selon que la relation est vérifiée, ou pas.
La fonction Python suivante répond donc au problème :
- def équation_cartésienne(A,B):
- def test(x,y):
- return sont_alignés(A,B,(x,y))
- return test
Voici un exemple d’utilisation :
- A = (2,3)
- B = (5,8)
- eqAB = équation_cartésienne(A,B)
- eqAB(13,21)
pour savoir si le point de coordonnées 13 et 21 appartient ou non à la droite (AB).
L’adjectif cartésienne vient de René Descartes, qui, en ramenant des problèmes de lieux géométriques à de l’algèbre, a été le premier à mettre l’accent sur cette notion d’équation de courbe.
Le problème avec cette version est qu’elle ne permet pas de connaître l’équation. Voici alors une variante qui, au lieu de l’équation, renvoie sa représentation littérale, à partir des coordonnées d’un vecteur directeur :
- def équation_droite(A,B):
- (a,b)=vecteur(A,B)
- return "{}×x+({})×y = {}".format(-b,a,-b*A[0]+a*A[1])
Pour préserver une certaine concision, on utilise l’écriture formatée de Python. Chaque paire d’accolades est à considérer comme un trou à combler en y plaçant une expression. La liste des expressions est entre parenthèses après le mot format.
Pour afficher une équation de la droite passant par les points de coordonnées (2,3) et (5,8), il suffit d’écrire
équation_droite( (2,3) , (5,8) )
En géométrie analytique, une équation cartésienne de droite est de la forme a×x+b×y=c, l’égalité étant vraie si et seulement si le point de coordonnées (x,y) est sur la droite. Mais en statistique, on appelle parfois équation de régression non pas une équation au sens mathématique, mais la fonction affine représentée par la droite. C’est cette fonction affine qui est utile pour résoudre des problèmes d’interpolation ou d’extrapolation. En présupposant que les points A et B n’ont pas la même abscisse, voici comment on peut fabriquer cette fonction affine :
- def équation_réduite(A,B):
- (xA,yA)=A
- (xB,yB)=B
- (Dx,Dy) = vecteur(A,B)
- pente = Dy/Dx
- ordonnée_origine = yA-pente*xA
- def RegEq(x):
- return pente*x+ordonnée_origine
- return RegEq
Pour trouver la valeur que prend en 8 la fonction affine qui prend la valeur 3 en 5 et la valeur 8 en 13, on peut faire
- A = (5,3)
- B = (13,8)
- équation_réduite(A,B)(8)
Remarques :
- Le vocabulaire ordonnée à l’origine et pente est celui du cours, le mot pente n’a de sens que dans un repère orthonormé (c’est la tangente de l’angle que fait la droite avec l’axe des abscisses), dans le cas affine, le vocabulaire coefficient_directeur est plus approprié.
- L’application de régression affine de la calculatrice est plus pratique que Python pour résoudre ce genre de problème : On entre la liste des deux abscisses et la liste des deux ordonnées, puis on fait une régression affine (« linéaire » sur Ti) qui donne l’équation réduite de la droite et permet même de dessiner celle-ci.
Fonctions
Se constituer un répertoire de fonctions de référence
Fonctions carré et cube
- def carré(x):
- return x*x
- def cube(x):
- return x**3
Fonction inverse
La fonction suivante ne convient pas vraiment :
- def inverse(x):
- return 1/x
En effet si on essaye de calculer inverse(0)
on a le message suivant :
ZeroDivisionError: division by zero
En effet la fonction inverse n’est pas définie sur R mais seulement sur ]-∞ ;0[∪]0 ;+∞[. Pour restreindre son domaine de définition, on va donc supposer (« assert ») que x n’est pas égal à 0 :
- def inverse(x):
- assert x!=0
- return 1/x
Avec cette fonction cohérente avec le cours de maths, une tentative de calcul de l’inverse de 0 donne toujours un message d’erreur mais cette fois-ci c’est :
assert x!=0
AssertionError
C’est logique : Si, à un moment, on a besoin de calculer l’inverse d’un réel, ce réel est censé être non nul, et tenter de calculer l’image par une fonction, d’un nombre qui n’est pas dans le domaine de définition de cette fonction, est assimilé à un bug.
Fonction racine
On suppose les fonctions du module math chargées par from math import *
; c’est le cas en particulier de sqrt (square root est la traduction anglaise de « racine carrée »). Or si on essaye de calculer la racine carrée de -1 à l’aide de cette fonction, on a le message d’erreur [1] suivant :
ValueError: math domain error
Autrement dit, une erreur de domaine de définition. Alors on propose [2]
- def racine(x):
- assert x>=0
- return sqrt(x)
qui, avec le calcul de racine(-1)
donne le message
assert x>=0
AssertionError
qui en dit plus sur la nature de l’erreur (seuls les réels positifs ont une racine).
Pour une fonction dont le tableau de variations est donné, algorithmes d’approximation numérique d’un extremum (balayage, dichotomie)
La méthode de Nelder-Mead utilisée pour chercher un extrémum, s’apparente, en dimension 1, plus à de la trichotomie, qu’à de la dichotomie. Il est donc un peu étonnant qu’on propose un algorithme de dichotomie pour chercher un extrémum, alors que cette méthode est plus adaptée à la recherche d’une solution. Certes, pour une fonction f dérivable, la recherche d’un extrémum de f revient à résoudre l’équation f’(x)=0, mais la dérivabilité n’est pas au programme de 2nde.
En tout cas, si on connaît le tableau de variation d’une fonction, il peut arriver que l’on sache qu’elle décroît puis croît et il reste à déterminer (une valeur approchée de) son minimum. Par dichotomie on peut faire ainsi :
- def minimum(f,a,b):
- m = f((a+b)/2)
- while b-a>1e-6:
- if f(a)<f(b):
- b = m
- else:
- a = m
- m = (a+b)/2
- return f(m)
Explication : f est une fonction qui décroît puis croît sur l’intervalle [a,b] (que l’on va resserrer au cours de l’exécution de l’algorithme). On s’arrête de boucler lorsque l’intervalle est devenu suffisamment petit pour que le minimum soit proche à la fois de f(a) et de f(b) (ligne 3). Dans la boucle, de deux choses l’une :
- ou bien f(a) est plus petit que f(b) et alors f(m) est un meilleur candidat que f(b) pour le minimum de f ;
- ou bien c’est le contraire, c’est f(b) qui est le plus petit et f(m) est un meilleur candidat que f(a) pour le minimum de f.
Dans le premier cas on remplace donc b par m, dans le second on remplace a par m. Dans les deux cas, l’intervalle voit sa longueur divisée par deux et l’algorithme converge rapidement. Par exemple
minimum(carré,-1,2)
donne
3.552713678800501e-15
alors que
minimum(cos,0,4)
donne
-0.9999999999999987
C’est à la fois rapide et précis.
La méthode par balayage est plus simple à programmer : On approche la courbe représentant f par un polygone et on cherche la plus petite ordonnée des sommets de ce polygone. Pour cela, comme on sait que la fonction décroît puis croît, on s’arrête de boucler dès que l’extrémité gauche du segment est devenue plus basse que l’extrémité droite. Si les coordonnées de ces sommets sont (x1,f(x1)) et (x2,f(x2)) et que x2=x1+dx (avec dx petit mais constant), on obtient la fonction suivante :
- def minimum(f,a,b):
- dx = 0.001
- x1 = a
- x2 = a+dx
- while f(x1)>f(x2):
- x1 += dx
- x2 += dx
- return f(x1)
Avec cette fonction, minimum(carré,-1,2)
donne 7.765831018255432e-31
alors que minimum(cos,0,4)
donne -0.9999999170344522
Algorithme de calcul approché de longueur d’une portion de courbe représentative de fonction
On peut se contenter de modifier l’algorithme par balayage ci-dessus : On additionne les longueurs des segments (x1,x2). Comme x2=x1+dx, la longueur du segment s’obtient par hypot(dx,f(x+dx)-f(x)) et on va, au pris d’une perte d’efficacité de l’algorithme, se passer des x1 et x2 et n’avoir que x comme variable :
- def longueur_courbe(f,a,b):
- dx = 0.000001
- x = a
- longueur_parcourue = 0
- while x<b:
- longueur_parcourue += hypot(dx,f(x+dx)-f(x))
- x += dx
- return longueur_parcourue
Remarque : La ligne 6 peut s’écrire, en appelant dy la différence f(x+dx)-f(x), sous la forme ds²=dx²+dy². C’est en affinant cette forme qu’Albert Einstein a créé la théorie de la relativité générale, en se basant sur la théorie des géodésiques de Minkovski. En bref, alors qu’on croit parler de fonctions, on fait de la géométrie :-D
Probas
On formalise la notion de loi de probabilité en s’appuyant sur le langage des ensembles
Le mot ensemble se traduit en anglais par set. Et Python possède une fonction également appelée set qui convertit en ensemble tout ce qui peut être converti en ensemble. Par exemple
- l’ensemble vide se contruit en ne donnant aucun argument à la fonction :
set()
- set convertit un couple en paire ou en singleton avec
set((2,3))
ouset((2,2))
- set peut construire un ensemble à partir d’une liste (cas le plus fréquent) avec
set([2,3,5,7])
- set peut même convertir un range en un ensemble. Par exemple pour l’univers de probabilité obtenu en lançant un dé, on peut écrire
U = set(range(1,7))
Dans la suite on suppose qu’on a calculé l’univers U correspondant à l’expérience probabiliste de lancer d’un dé icosaédrique, avec
- >>> U = set(range(1,21))
- >>> U
- {1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20}
On considère également l’événement « le résultat est pair » noté A :
- >>> A = set(range(2,21,2))
- >>> A
- {2, 4, 6, 8, 10, 12, 14, 16, 18, 20}
- >>> B = set([2,3,5,7,11,13,17,19])
- >>> B
- {2, 3, 5, 7, 11, 13, 17, 19}
et l’événement B est évidemment « le résultat est premier ».
Réunion
Pour calculer l’événement « le résultat est pair ou premier » on utilise la fonction union :
- >>> A.union(B)
- {2, 3, 4, 5, 6, 7, 8, 10, 11, 12, 13, 14, 16, 17, 18, 19, 20}
Intersection
Pour calculer l’événement » le résultat est à la fois pair et premier, on utilise la fonction intersection :
- >>> A.intersection(B)
- {2}
Complémentaire
Pour calculer l’événement « le résultat est impair » ou l’événement « le résultat est composé » on soustrait l’événement A ou l’événement B à l’univers U :
- >>> U.difference(A)
- {1, 3, 5, 7, 9, 11, 13, 15, 17, 19}
- >>> U.difference(B)
- {1, 4, 6, 8, 9, 10, 12, 14, 15, 16, 18, 20}
Ces opérations ont été définies par George Boole au milieu du 19e siècle. La théorie des ensembles remonte au début du 20e siècle avec Georg Cantor et Ernst Zermelo. C’est Émile Borel qui semble avoir eu l’idée de décrire les événements par des ensembles, ce qui a permis à Andreï Kolmogorov d’axiomatiser les probabilités en 1933.
Probabilité d’un événement
Dans le cas équiprobable, on définit la probabilité comme une fonction qui à chaque événement associe un nombre [3] compris entre 0 et 1.
- def P(E):
- global U
- return len(E)/len(U)
Avec cette fonction, on peut calculer des probabilités dans la console :
- >>> P(A)
- 0.5
- >>> P(B)
- 0.4
- >>> P(U)
- 1.0
et vérifier la formule du cours :
- >>> P(A.union(B))+P(A.intersection(B))
- 0.9
- >>> P(A)+P(B)
- 0.9
- >>> P(A.union(B))+P(A.intersection(B)) == P(A)+P(B)
- True
Loi de probabilité
Une loi de probabilité est la donnée, pour chaque événement élémentaire (n’ayant qu’un seul élément dans le langage des ensembles) de sa probabilité. On peut donc la modéliser comme une fonction qui, à chaque face du dé (par exemple ; on peut aussi penser aux deux faces d’une pièce) associe la probabilité d’obtenir cette face. Autrement dit, à tout élément x de l’univers U (x∈U) , on associe la probabilité de l’événement set((x))
(un singleton).
Or en Python, une fonction sur un ensemble fini se modélise par un dictionnaire, c’est-à-dire une liste de couples du type (x,y) où y est la l’image de x par la loi de probabilité. Par exemple, avec un dé à 6 faces non équilibré (pipé), en supposant que les probabilités des 6 faces soient
1 | 2 | 3 | 4 | 5 | 6 |
0,16 | 0,16 | 0,16 | 0,16 | 0,16 | 0,2 |
(le 6 a plus de chances d’arriver que les autres faces)
on peut écrire
loi = { 1: 0.16, 2: 0.16, 3: 0.16, 4: 0.16, 5: 0.16, 6: 0.2 }
et la fonction permettant de calculer la probabilité d’un événement avec U = set(range(1,7))
est
- def P(E):
- global U
- return sum(loi[x] for x in E)
Autrement dit, la probabilité de l’événement E est la somme des probabilités des événements élémentaires constituant E.
En définissant l’événement « le résultat est pair » par résultat_pair = set(range(2,7,2))
sa probabilité P(résultat_pair)
est 0,52 : On a plus de chances d’avoir un résultat pair, qu’un résultat impair, parce que 6 est un nombre pair et il est plus probable que les autres résultats.
Stats
Pour des données réelles ou issues d’une simulation, lire et comprendre une fonction écrite en Python renvoyant la moyenne m, l’écart type s, et la proportion d’éléments appartenant à [m-2s,m+2s]
Les données issues d’une simulation peuvent être fabriquées l’une après l’autre dans une boucle (c’est ce qu’on fera en sections technologiques dans un autre article), mais pour homogénéiser les approches « données réelles » et « simulation », il convient de se poser la question « sur quoi ces fonctions doivent-elles porter ? ». Le type de Python le plus adapté à des statistiques est la liste (tableau dans les autres langages de programmation) notée entre crochets. On propose donc de définir des fonctions moyenne, écart_type et proportion_2_sigmas s’appliquant à des listes de nombres et renvoyant des nombres (flottants) [4].
Comme les listes prises en entrée sont des listes de nombres, il est possible de calculer
- la somme des éléments de la liste (avec la fonction sum de Python)
- le nombre d’éléments de la liste (avec la fonction len, abrégée de length, de Python)
- la liste des carrés des éléments de la liste (avec la fonction map de Python, en « mappant » la fonction carré définie précédemment, sur la liste).
On arrive alors à
- def moyenne(liste):
- return sum(liste)/len(liste)
- def variance(liste):
- liste_des_carrés = list(map(carré,liste))
- return moyenne(liste_des_carrés)-carré(moyenne(liste))
- def écart_type(liste):
- return racine(variance(liste))
Explication de la fonction variance : map(carré,liste)
ne construit pas la liste des carrés, mais un objet de type map. Il est donc nécessaire de convertir cet objet en liste (avec la fonction list) pour pouvoir se livrer à des calculs statistiques sur cet objet.
Remarque sur la définition de la variance : On connaît deux algorithmes classiques de calcul de la variance :
- Moyenne des carrés des écarts à la moyenne (choisi en général comme définition) ;
- La moyenne des carrés moins le carré de la moyenne (formule de Huygens, choisie en général comme propriété caractéristique).
Le mot « variance » ne figurant pas au programme de 2nde, on peut effectuer le choix que l’on veut, et définir la variance par la formule de Huygens. L’essentiel étant d’aller rapidement à la notion d’écart-type, qui est la racine carrée de la variance.
On remarque que les fonctions carré et racine définies dans l’onglet sur les fonctions, trouvent assez rapidement leur utilité.
Pour calculer la probabilité d’appartenance à un intervalle à 2σ près, on peut utiliser les fonctions précédentes ainsi :
- def proportion_2_sigmas(liste):
- m = moyenne(liste)
- s = écart_type(liste)
- nombre_de_succès = 0
- for x in liste:
- if m-2*s<x<m+2*s:
- nombre_de_succès += 1
- return nombre_de_succès/len(liste)
Remarque : Les listes en compréhension permettent de faire plus court, mais elles n’apparaissent qu’au programme de 1ère générale. On a donc préféré faire comme en 1ère technologique, en créant une variable nombre_de_succès initialisée à 0 et incrémentée à chaque succès. La fréquence de succès est le quotient du nombre de succès par l’effectif total.
Lire et comprendre une fonction Python renvoyant le nombre ou la fréquence de succès dans un échantillon de taille n pour une expérience aléatoire à deux issues
La première question à se poser est « qu’est-ce qu’un échantillon ? ». Comme en 1ère seules les listes sont considérées, on va donc modéliser ici un échantillon par une liste. Et par souci de concision on va la décrire, osons-le, en compréhension :
- def échantillon_2_issues(n,p):
- return [random()<p for k in range(n)]
Le fait qu’on construit une liste se voit aux crochets, qui dénotent les listes en Python. La partie for k in range(n)
garantit qu’il y a n valeurs dans la liste (elle est de longueur n, dans le jargon de Python). Mais que sont ces valeurs ? random()<p
est un booléen, égal à True
si random() est inférieur à p et False
sinon. Cette technique sera utilisée dans le programme de 1ère technologique. La probabilité d’avoir True
est p justement. Et on est bien obligé de gérer ce second paramètre p puisque l’hypothèse d’équiprobabilité entre les deux issues ne figure pas dans le programme.
Il y a d’autres moyens d’obtenir un échantillon pour une expérience aléatoire à deux issues, par exemple lancer une punaise 1000 fois et écrire True si elle est tombée la pointe en haut, False sinon. Mais cela risque d’être vécu comme fastidieux par les élèves. La question essentielle ici est de savoir ce qu’on fait de l’échantillon une fois qu’il est constitué :
- def nombre_succès(n,p):
- return échantillon_2_issues(n,p).count(True)
- def fréquence_succès(n,p):
- return nombre_succès(n,p)/n
Pour avoir le nombre de succès, on demande à Python de compter les occurences de True dans l’échantillon. Pour avoir la fréquence des succès, on divise le nombre de succès par n qui se trouve être la taille de l’échantillon.
La théorie du pile ou face biaisé est de Jacques Bernoulli qui est aussi l’auteur d’autres algorithmes que l’on verra plus bas. Il cherchait avant tout à appliquer les mathématiques aux sciences sociales et le nombre de succès dans un échantillon est présent, entre autres, dans la théorie des suffrages.
Observer la loi des grands nombres à l’aide d’une simulation sur Python ou tableur
La loi des grands nombres dit, en substance, que si n grandit, la fréquence de succès se rapproche de p, en fluctuant à cause du hasard, mais de moins en moins. On peut la vérifier en affichant (avec print, horreur !) la fréquence de succès pour n de plus en plus grand :
- for n in range(1,101):
- print(fréquence_succès(n,1/6))
(ici, la fréquence de « 6 » lorsqu’on lance un dé équilibré à 6 faces, se rapproche de 0,1666666... mais pas très vite).
Simuler N échantillons de taille n d’une expérience aléatoire à deux issues. Si p est la probabilité d’une issue et f sa fréquence observée dans un échantillon, calculer la proportion des cas où l’écart entre p et f est inférieur ou égal à 1/√(N)
Simuler N échantillons, c’est répéter N fois la simulation d’un échantillon. Pour cela, on boucle N fois, et dans la boucle on fait appel à la fonction fréquence_succès(n,p) définie plus haut, en incrémentant la variable normal uniquement si l’écart entre f et p est inférieur ou égal à 1/√(N) :
- def répétitions(N,n,p):
- normal = 0
- for k in range(N):
- f = fréquence_succès(n,p)
- if p-1/racine(N)<f<p+1/racine(N):
- normal += 1
- return normal/N
Première générale
Voir aussi cet article.
logique
Les élèves apprennent en situation à utiliser les quantificateurs
En Python, ∀ se note all
et ∃ se note any
. Ce sont des fonctions qui s’appliquent à une liste (ou un ensemble) et renvoient True ou False : Précéder une équation d’un quantificateur sur son inconnue en fait donc une proposition logique, comme dans le cours de maths. L’idée est de Gottlob Frege dans la seconde moitié du 19e siècle.
Au fait, le mot-clé in
correspond à ∈. C’est une relation (d’appartenance) entre éléments et ensembles mais en Python elle peut aussi s’appliquer entre un élément et une liste voire un range.
Listes
Vers le milieu du 19e siècle, George Boole, avec ses lois de la pensée, a commencé à assimiler
- l’intersection avec la conjonction logique
- la réunion avec la disjonction (inclusive)
- le complémentaire avec la négation logique
En bref, il a assimilé les concepts avec les ensembles. Par exemple, entre le concept de nombre pair et l’ensemble des nombres pairs, Boole ne fait pas réellement la différence. Plus tard c’est Gottlob Frege qui a formalisé la logique des prédicats avec les quantificateurs et l’appartenance. Un prédicat est une fonction à valeurs booléennes. Par exemple la propriété « être pair » se modélise chez Frege par la fonction qui à tout nombre divisible par 2, associe True, et à tout nombre non divisible par 2, associe False. Il est alors apparu que des ensembles peuvent facilement se construire comme images réciproques de True par un prédicat. Par exemple l’ensemble des nombres pairs peut être défini comme « l’ensemble des entiers naturels n tels que n modulo 2 vaut 0 ».
Il a fallu attendre le début du 20e siècle pour qu’Ernst Zermelo se rende compte que l’existence d’un ensemble regroupant les objets vérifiant une propriété n’est pas si évidente que cela : C’est l’axiome de compréhension. Python permet de construire des ensembles par compréhension, mais aussi d’autres collections comme ... les listes.
Suites
Une première question à poser, avant d’étudier les algorithmes sur les suites, est celle de la définition d’une suite : Qu’est-ce que c’est qu’une suite ?
Comme le programme insiste sur la notion de fonction, on peut proposer cette définition : Une suite est une fonction définie sur N. Plus précisément, si l’ensemble d’arrivée est un ensemble de nombres (typiquement R) la suite est dite numérique.
Se pose alors le problème de la représentation d’une suite en Python : La machine sur laquelle tourne Python ne disposant que d’une quantité de mémoire finie, on ne peut y caser une infinité de nombres u(0)=u0, u(1)=u1, etc.
Le langage Haskell gère les suites qu’il considère comme des listes infinies. Pas Python. On a alors deux solutions :
- On opte pour l’infini potentiel en représentant les suites par des itérateurs
- On reste dans le domaine fini en ne calculant qu’un nombre fini de termes de la suite.
C’est le second choix qui sera fait ici, et les termes en question seront stockés dans une liste Python. Celle-ci sera typiquement construite par compréhension pour une suite explicite, et par ajouts successifs (avec append) pour une suite récurrente.
Calcul de termes d’une suite, de sommes de termes, de seuil
Termes
Un exemple de suite définie explicitement est celle des carrés dits « parfaits » :
[carré(n) for n in range(20)]
s’évalue à
[0, 1, 4, 9, 16, 25, 36, 49, 64, 81, 100, 121, 144, 169, 196, 225, 256, 289, 324, 361]
Si on affecte une variable u à cette liste, le premier terme u0 sera obtenu par u[0]
, le second terme u1 par u[1]
etc. Et la somme des termes peut s’obtenir en appliquant la fonction sum à cette liste.
La suite des carrés est un exemple de nombres figurés. Les nombre figuré sont surtout connus parce que Pythagore leur accordait des propriétés magiques.
Un exemple de suite récurrente est présent dans l’étude de la tour d’Hanoï, jeu à un joueur inventé par Édouard Lucas vers la fin du 19e siècle :
Pour résoudre la tour d’Hanoï à zéro disque, on effectue zéro mouvement : h0=0. Pour résoudre la tour d’Hanoï à n+1 disques, on effectue
- hn mouvements pour transporter tous les disques sauf le plus grand vers le piquet intermédiaire ;
- 1 mouvement pour placer le plus grand disque sur le piquet final ;
- hn mouvements pour transporter les disques restants depuis le piquet intermédiaire vers le piquet final.
Donc hn+1=hn+1+hn=2×hn+1 :
- def mouvements_Hanoi(n):
- h = [0]
- for n in range(1,n):
- h.append(2*h[-1]+1)
- return h
Explication de la ligne 4 : Le dernier élément de la liste h se note h[-1]
. On ajoute (« append ») 1 de plus que son double (relation de récurrence) à la fin de la liste h.
Somme de termes
Pour additionner les n premiers termes de la suite u représentée par une liste de n nombres, on peut faire
- S = 0
- for terme in u:
- S += terme
ou bien
sum(u)
calcul de seuil
Par exemple, ion se demande combien il faut de disques à la tour d’Hanoï pour nécessiter au moins 1000 mouvements :
- def disques_Hanoi(nmvts):
- h = 0
- ndisques = 0
- while h<1000:
- h = 2*h+1
- ndisques += 1
- return ndisques
Calcul de factorielle
Au fait c’est quoi une factorielle ? On propose de définir la factorielle de n comme produit des n premiers entiers non nuls. La question étant de savoir comment ça se calcule. On propose l’algorithme de sommation, dans lequel ici on remplace les additions par des multiplications. Algorithme trouvé difficile par plusieurs étudiants de BTS. Voici l’algorithme décrit en Sofus :
Et en Python ça donne (par exemple à l’aide de SofusPy) :
- def factorielle(n):
- assert n==int(n) and n>=0
- p = 1
- for k in range(1,n+1):
- p *= k
- return p
Explication de la ligne 2 : La fonction factorielle n’est définie que pour les entiers naturels, on restreint donc son calcul aux entiers naturels.
Remarque : Une fois importées les fonctions du module math, l’expression suivante donne le même résultat :
[factorial(n) for n in range(10)]
Voici les factorielles sur OEIS
Liste des premiers termes d’une suite : suite de Syracuse, suite de Fibonacci
Suite de Collatz
La suite de Collatz consiste à itérer une transformation appelée ici mélange car elle mélange les entiers naturels :
- def est_pair(n):
- return n%2==0
- def mélange(n):
- if est_pair(n):
- return n//2
- else:
- return 3*n+1
- def Collatz(u0,n):
- u = [u0]
- for k in range(n):
- u.append(mélange(u[-1]))
- return u
Le premier argument de cette fonction est le premier terme de la suite (un entier naturel au choix) et le second argument est le nombre de termes à calculer.
Fibonacci
Voici la suite de Fibonacci sur OEIS.
Et en Python :
- def Fibonacci(n):
- assert n==int(n) and n>=0
- if n==0:
- return [1]
- else:
- F = [1,1]
- for k in range(2,n+1):
- F.append(F[-1]+F[-2])
- return F
On peut voir d’autres exemples de suites dans ces regards croisés sur le sujet.
Équations
Second degré
Résoudre une équation, c’est donner l’ensemble de ses solutions. Il est donc logique de définir une fonction qui, aux réels a, b, et c, associe l’ensemble des solutions de l’équation a×x²+b×x+c=0. Cette fonction est remarquablement simple à calculer si on utilise les propriétés des ensembles (chaque élément n’apparaît qu’une fois et on ne tient pas compte de l’orde des éléments [7]) :
- def solutions(a,b,c):
- """résout l'équation a*x**2+b*x+c==0"""
- S = set()
- Delta = b**2-4*a*c
- if Delta>=0:
- r = racine(Delta)
- S.add((-b-r)/(2*a))
- S.add((-b+r)/(2*a))
- return S
Par exemple pour résoudre l’équation x²-5x+6=0 on fait solutions(1,-5,6)
qui affiche
{2.0, 3.0}
: il s’agit bien d’un ensemble de solutions.
Remarques : La fonction racine a été définie dans l’onglet sur les fonctions de référence dans la partie sur la classe de 2nde. Si le discriminant est négatif, on a set()
qui désigne l’ensemble vide.
Méthode de Newton, en se limitant à des cas favorables
On va ici considérer qu’un cas est favorable si
- on connaît la dérivée de f ;
- f’(a) est petit, si a est la solution cherchée ;
- on connaît une estimation de la solution cherchée.
Par exemple, si on veut, non plus résoudre l’équation x²-5x+6=0 (c’est-à-dire chercher toutes ses solutions) mais seulement estimer une bonne valeur approchée d’une de ses solutions, on pose f(x)=x²-5x+6, puis on calcule f’(x)=2x-5 et on définit
- def f(x):
- return x**2-5*x+6
- def fprime(x):
- return 2*x-5
- def zpn(f,fprime,a):
- for k in range(10):
- a = a-f(a)/fprime(a)
- return a
(zpn est une abréviation de « zéro par Newton »)
avec zpn(f,fprime,4)
on trouve la solution 3.0 égale à 3, à la précision de la machine près, mais de type float puisque tous les calculs ont été effectués sur des réels.
Remarque historique : La méthode de Newton a été inventée par Joseph Raphson mais avant qu’il ait le temps de la publier, l’immense Isaac Newton l’a réinventée et publiée. Pauvre Raphson, qui n’a pas laissé son nom à sa méthode ! Il est d’ailleurs d’usage de parler plutôt de méthode de Newton-Raphson pour rendre hommage à cet inconnu.
Dérivée
Écrire la liste des coefficients directeurs des sécantes pour un pas donné
Plutôt qu’écrire, on préférera ici produire (renvoyer comme le fait toute fonction) la liste. Le pas est donné directement dans le corps de la fonction, et nommé dx. Cela permet de modifier légèrement le script de calcul de longueur de courbe vu en Seconde (au lieu de l’hypoténuse on calcule la tangente de l’angle) :
- def lcds(f,a,b):
- dx = 0.000001
- x = a
- liste_coeff = []
- while x<b:
- liste_coeff.append((f(x+dx)-f(x))/dx)
- x += dx
- return liste_coeff
Que signifient les lettres lcds ? liste des coefficients directeurs des sécantes.
Détermination d’une valeur approchée de e à l’aide de la suite (1+1/n)n
Au début du 18e siècle, Jacques Bernoulli a commencé à s’intéresser à l’application des maths aux finances. Cela l’a amené à se poser cette question : De combien un capital augmente-t-il lorsqu’on lui applique des augmentations de pourcentages infinitésimaux équivalents en tout à 100 pourcents ?
Autrement dit, si on augmente n fois une quantité d’un pourcentage de 100/n, l’augmentation totale revient à une multiplication par e. Or augmenter de 100/n pourcents revient à appliquer un coefficient multiplicateur égal à 1+1/n et dans le langage des suites, le résultat de Jacques Bernoulli revient à dire que la suite explicite de terme général égal à (1+1/n)n tend vers e [8].
La suite étant définie explicitement, une liste en compréhension s’impose :
- def approx_e(n):
- return [(1+1/k)**k for k in range(1,n+1)]
Construction de l’exponentielle par la méthode d’Euler
Jacques Bernoulli était le père de Daniel Bernoulli. Lorsque celui-ci habitait encore à Bâle, il jouait parfois avec le fils du pasteur, dont Jacques a rapidement remarqué le talent en mathématiques. Cet enfant admiré par toute la famille Bernoulli s’appelait Leonhard Euler. Calculateur acharné, il a trouvé un moyen de résoudre graphiquement des équations différentielles, en remplaçant le coefficient directeur de la sécante par celui de la tangente. Avec l’équation différentielle y’=y cela donne le remplacement de ex par (ex+dx-ex)/dx avec dx petit. On a alors ex+dx-ex≈ ex×dx soit ex+dx≈ex+ex×dx=ex(1+dx) : Les valeurs prises par l’exponentielle sur une subdivision d’un intervalle du type [1,b] à pas constant sont approchées par les termes d’une suite géométrique de raison 1+dx.
Incidemment ce résultat généralise celui sur le calcul approché de e, puisque cela revient à dire que si le nombre n de pas tend vers l’infini, la quantité (1+x/n)n tend vers ex.
La suite géométrique peut se calculer par une compréhension :
- def exp_euler(n,x):
- r = 1+x/n # raison de la suite
- y = 1 # exp(0)==1
- return [y*r**k for k in range(n+1)]
Vecteurs
On peut faire suite au cours de 2nde (vecteurs représentés par des couples de flottants) avec le produit scalaire :
- def produit_scalaire(u,v):
- (xu,yu) = u
- (xv,yv) = v
- return xu*xv+yu*yv
Cette fonction (qui associe un nombre à deux vecteurs) peut être utilisée pour calculer le cosinus de l’angle entre deux vecteurs non nuls :
- def cosinus(u,v):
- assert u!=(0,0) and v!=(0,0)
- return produit_scalaire(u,v)/norme(u)/norme(v)
- def angle(u,v):
- return degrees(acos(cosinus(u,v)))
Au passage on voit comment, à partir d’un cosinus, Python peut donner l’angle en degrés correspondant : La fonction acos permet de retrouver l’angle en radians, et en la composant par la fonction degrees, on convertit cet angle en degrés.
On peut ensuite créer une fonction Python pour tester si deux vecteurs (nuls ou non) sont orthogonaux :
- def sont_orthogonaux(u,v):
- return produit_scalaire(u,v)==0
Cercles
Le produit scalaire permet de refaire avec les équations de cercle ce qui a été vu en 2nde sur les équations de droite :
équation comme prédicat
L’équation d’un cercle de diamètre donné peut se calculer ainsi :
- def équation_cercle_diamètre(A,B):
- def test(M):
- return sont_orthogonaux(vecteur(A,M),vecteur(B,M))
- return test
représentation textuelle de l’équation
Par centre et rayon :
- def équation_cercle(centre,rayon):
- (Cx,Cy) = centre
- R = rayon
- return "x²+y²+({})×x+({})×y = {}".format(-2*Cx,-2*Cy,R**2-Cx**2-Cy**2)
Approximation de π par la méthode d’Archimède
On peut montrer que lorsque x est petit (en radians), sin(x) et tan(x) sont proches de x et de plus encadrent x. Alors si on découpe un cercle de rayon 1 en un grand nombre n d’arcs de longueur 2π/n, le demi-périmètre π est compris entre celui n×sin(2π/n)/2 du polygone régulier inscrit, et celui n×tan(2π/n)/2 du polygone régulier circonscrit au cercle. De plus, lorsque n est grand, ces deux nombres sont proches l’un de l’autre et a fortiori de π.
Vu comme ça cet algorithme semble tourner en rond puisque pour calculer sin(2π/n) et tan(2π/n) on s’attend à avoir besoin d’une valeur approchée de π. Mais la trigonométrie permet de calculer sin(x/2) et tan(x/2) à partir de sin(x) et tan(x), en n’utilisant que les 4 opérations et la racine carrée (calcul de moyennes arithmétique et géométrique). En partant comme l’a fait Archimède d’hexagones inscrit et circonscrit au cercle unité, on construit cette suite (a est le périmètre de l’hexagone inscrit, b celui de l’hexagone circonscrit) :
- def archimede(n):
- a = 3
- b = 6/racine(3)
- p = [moyenne(a,b)]
- for k in range(n):
- b *= a/p[-1]
- a = racine(a*b)
- p.append(moyenne(a,b))
- return p
Cette fonction ne renvoie pas la suite des encadrements de π obtenus par Archimède mais celle des approximations de π. Pour retrouver le résultat d’Archimède il faut faire archimede(4)
.
Les fonctions racine et moyenne ont été définies en Seconde.
Probas
Python permet de résumer en quelques lignes la théorie de Bayes :
- def proba_sachant_que(B,A):
- assert len(B)>0
- return len(A.intersection(B))/len(B)
(la ligne 2 rappelle que la notion de probabilité conditionnelle [9] n’a aucun sens si l’événement conditionnant est impossible)
Avec les événements A et B vus en 2nde, on découvre avec proba_sachant_que(B,A)
et proba_sachant_que(A,B)
que PB(A) et PA(B) ne sont pas égales, la première valant une chance sur 10 et la seconde, une chance sur 8. D’ailleurs sont_indépendants(A,B)
répond False une fois définie cette fonction :
- def sont_indépendants(A,B):
- return P(A.intersection(B))==P(A)*P(B)
Variables aléatoires
Le programme propose de modéliser une variable aléatoire par une fonction de U dans R où U est un univers fini. Par exemple, je viens d’inventer un jeu [10] :
Chaque joueur à son tour lance un dé icosaédrique, et remporte un nombre de points égal au nombre de diviseurs du résultat du lancer de dé. Le premier joueur ayant dépassé 45 points a gagné le jeu.
La variable aléatoire « nombre de diviseurs de la face du dé icosaédrique » est donc une variable aléatoire entière. On commence par définir (uniquement pour cet exemple particulier) une fonction définie sur N et donnant le nombre nd de diviseurs d’un entier non nul :
- def nd(n):
- return len([d for d in range(1,n+1) if n%d==0])
La liste définie en compréhension est celle des diviseurs d de n comme on l’a vu en 2nde. Sa longueur est donc le nombre de diviseurs.
Comme U est fini, on modélisera mieux la variable aléatoire par un dictionnaire de Python, qui est la structure la mieux adaptée pour modéliser les fonctions sur des ensembles finis (voir l’exemple de loi de probabilité vu en 2nde). Ce dictionnaire sera décrit lui aussi en compréhension, et noté X comme il est d’usage avec les variables aléatoires :
X = { k:nd(k) for k in U}
Pour calculer la loi de probabilité (encore un dictionnaire) de X, on va commencer par construire l’événement X=k (un ensemble cette fois-ci). Et c’est là que s’exprime toute la puissance des compréhensions : L’événement « la face a 6 diviseurs » est formé de tous les entiers entre 1 et 20 ayant 6 diviseurs, et cet ensemble se calcule en résolvant l’équation X(x)=6 sur U, c’est-à-dire en regroupant en un ensemble tous les x de U tels que X(x)=6. En Python cela s’écrit
- def égal(X,k):
- return {x for x in X.keys() if X[x]==k}
U a été retraduit en X.keys()
pour que cette fonction reste calculable au cas où on aurait changé le nom de l’univers (Ω au lieu de U par exemple). Ceci permet de rapidement construire la loi de la variable aléatoire X, en associant à chacune de ses valeurs k (qui sont listées dans X.values()
) la probabilité que X égale k :
- def loi(X):
- return {k:P(égal(X,k)) for k in X.values()}
la loi du nombre de diviseurs du dé icosaédrique est
{1: 0.05, 2: 0.4, 3: 0.1, 4: 0.25, 5: 0.05, 6: 0.15}
Algorithme renvoyant l’espérance, la variance ou l‘écart type d’une variable aléatoire
Espérance
L’espérance se définit comme la somme des produits du type k×P(X=k) :
- def E(X):
- return sum({k*P(égal(X,k)) for k in X.values()})
En calculant E(X) on découvre qu’en moyenne un dé icosaédrique a 3,3 diviseurs. Qu’on se le dise !
Variance
Quant à la variance, elle est égale à 2,31, obtenu en évaluant cette fonction (calculée avec la formule de Huygens) :
- def V(X):
- X2 = { k:X[k]**2 for k in X.keys() }
- return E(X2)-E(X)**2
Écart-type
Ce qui donne un écart-type d’environ 1,52 d’après l’évaluation de cette fonction :
- def sigma(X):
- return racine(V(X))
Stats
Méthode de Monte-Carlo : estimation de l’aire sous la parabole, estimation du nombre π
Stanislaw Ulam a déjà été cité plus haut, à propos de la suite de Collatz. Travaillant sur le projet Manhattan, qui, comme son nom l’indique, se déroulait à Los Alamos, il a vite découvert l’efficacité de la simulation par le hasard pour résoudre des équations différentielles : C’est la méthode de Monte-Carlo, ainsi appelée parce que pendant la guerre aux USA Monte-Carlo passait pour la capitale des casinos.
Comme nous sommes en temps de paix, on va se contenter ici, de jouer au paintball. Mais c’est un paintball très spécial, puisque la couleur de la peinture ne se décide qu’au moment où elle touche la cible : Le point de peinture devient rouge s’il est proche de (0,0) et bleu sinon. L’estimation de π/4 se fera en comptant les points rouges du nuage, ou plus précisément leur fréquence.
Et l’intérêt essentiel de cette modélisation étant la possibilité de voir les points de peinture, on va utiliser le module kandinsky de la NumWorks. On définit les couleurs rouge et bleu et on définit deux fonctions dessinant chacune un point dans sa couleur, à partir de x et y qui sont réels compris entre 0 et 1 :
- from kandinsky import *
- rouge = color(255,0,0)
- bleu = color(0,0,255)
- def point_rouge(x,y):
- set_pixel(int(200*x),int(200-200*y),rouge)
- def point_bleu(x,y):
- set_pixel(int(200*x),int(200-200*y),bleu)
Ensuite pour jouer au paintball sur l’écran de la Numworks, on boucle N fois (typiquement N vaut quelques milliers) et dans la boucle, on peint en rouge et on incrémente S si le point est rouge et on peint en bleu sans toucher à S si le point est bleu :
- S = 0
- for k in range(N):
- (x,y) = (random(),random())
- if hypot(x,y)<1:
- S += 1
- point_rouge(x,y)
- else:
- point_bleu(x,y)
Si S a été initialisée à 0, à la fin, S vaut le nombre de points rouges (« succès ») et son quotient par N donne une estimation de π/4. Voici un exemple avec un nuage de 10000 points donnant une valeur approchée de 0,7816 pour π/4 soit 3,1264 pour π :
Pour estimer l’aire sous la parabole, on change juste le test, le point étant en rouge quand on est sous la parabole soit lorsque y < x**2
:
Le calcul intégral permet de savoir que le tiers du carré est colorié en rouge et 0,3293 n’est pas éloigné de 1/3.
Pour programmer avec le module kandinsky, il n’est point besoin de posséder une Numworks, il suffit d’utiliser son émulateur en ligne (ou celui du smartphone) ou directement tester les scripts sur cette page
Fréquence d’apparition des lettres d’un texte donné, en français, en anglais
Le texte à analyser est ici :
Pour l’analyser avec Python, il faut donc déjà ouvrir ce fichier. Pour cela on propose la syntaxe with ... as ... qui ouvre le fichier (placé dans le même dossier que le script Python) sous le nom de fichier. En envoyant à ce fichier une commande de lecture (read), on obtient une variable texte qui est une liste de lettres. On fabrique alors un dictionnaire de Python (le tableau d’effectifs) en associant à chaque lettre du texte, le nombre d’occurences de cette lettre.
- with open("dream.txt") as fichier:
- texte = list(fichier.read().lower())
- lettres = {lettre for lettre in texte if lettre.isalpha()}
- tableau_effectifs = { lettre:texte.count(lettre) for lettre in lettres }
- print(tableau_effectifs)
La fonction lower() met toutes les lettres en minuscules. L’objet lettres n’est pas une liste des lettres du texte, mais un ensemble. Cela évite d’écrire plusieurs fois la même lettre. Et on n’a placé dans l’ensemble des lettres, que les vraies lettres (isalpha() filtre ce qui est lettre de l’alphabet). Enfin le tableau des effectifs est un dictionnaire défini en compréhension, où, à chaque lettre, on associe le nombre d’occurences de cette lettre dans le texte.
Fabriquer un tel dictionnaire est plus facile avec la fonction Counter du module collections (cliquer sur « les sacs de Python ») :
- from collections import *
- with open("dream.txt") as fichier:
- texte = list(fichier.read().lower())
- print(Counter(texte))
Et si, au lieu des effectifs, on préfère les fréquences, il suffit de diviser chaque effectif par l’effectif total pour connaître la fréquence correspondante. Pour calculer l’effectif total, on effectue une sommation sur tableau_effectifs.values() qui est la liste des effectifs :
- effectif_total = sum(tableau_effectifs.values())
- for lettre,effectif in tableau_effectifs.items():
- print(lettre,effectif/effectif_total)
Remarque : C’est Al-Kindi qui, durant le 9e siècle, a, pour des besoins de cryptanalyse, eu l’idée de comparer la fréquence d’apparition des lettres dans un texte crypté selon le chiffre de César avec celle des lettres du message crypté. Pour replacer cet épisode dans un contexte historique, voir ici. Un TP sur le sujet a déjà été fait en 2nde il y a plus de 7 ans...
Spé en Terminale
Voir aussi cet article
Ensembles
Encore une fois, on cherche à éviter toute confusion en se limitant aux listes et pas d’autres collections. Au fait qu’est-ce que c’est qu’une collection en Python ? Disons que c’est un objet capable de collecter ou récolter. Autrement dit, il contient d’autres objets (nombres par exemple) et on peut ajouter (append pour une liste, add pour un ensemble) un ou plusieurs objets à la collection (le collectionneur peut augmenter sa collection, en chinant...) ou aussi en enlever.
Comme la partie sur les probabilités évoque les puissances cartésiennes d’ensembles et que la partie sur la combinatoire évoque la génération d’ensembles, on va ici faire le choix (comme en 1ère) de manipuler autant que possible les ensembles de Python, et autant que possible les compréhensions.
Par exemple pour la paire contenant 0 et 1, on écrira set([0,1])
ou set(range(2))
;
et pour le produit cartésien de A et B (ensemble de couples (x,y) tels que x ∈ A et y ∈ B) on écrira
{(x,y) for x in A for y in B}
Combinatoire
Pour un entier n donné, génération de la liste des coefficients (nk) à l’aide de la relation de Pascal
On va comprendre l’énoncé de la manière suivante : n est l’argument passé à la fonction, qui renverra une liste représentant la ligne numéro n du triangle de Pascal. Et on va utiliser la récursivité pour ne pas allonger à l’excès le script :
- def coeff_binomiaux(n):
- if n==0:
- return [1]
- else:
- l = coeff_binomiaux(n-1)[:]+[0]
- return [1]+[l[k]+l[k+1] for k in range(n)]
Mais si on précède le tout d’un from itertools import *
, on peut engendrer les combinaisons (ici, d’une liste à n éléments entiers) et demander à Python de les compter avec la fonction len :
- def combinaisons(n,k):
- return list(combinations(list(range(n)),k))
- def coeff_binomiaux(n):
- return [len(combinaisons(n,k)) for k in range(n+1)]
Attention toutefois à l’inefficacité de cet algorithme, inutilisable si n est grand. L’algorithme présenté dans l’article sur la section technologique est plus pratique.
L’algorithme du bac Pondichéry 2012 peut être adapté aux deux exemples suivants :
Génération des permutations d’un ensemble fini, ou tirage aléatoire d’une permutation
Pour ce qui est de la génération (via un générateur) des permutations d’une liste [11], le sujet a été traité par Richard Gomez dans cet article.
Avec itertools c’est assez rapide ; par exemple list(permutations(list(range(4))))
calcule la liste de toutes les permutations de 4 éléments.
Pour engendrer une permutation aléatoire, outre l’algorithme du sujet de bac Pondichéry et le choix au hasard d’un élément de la liste que l’on vient de fabriquer ci-dessus, il y a la fonction shuffle du module random. Après from random import *
on peut faire
- L = list(range(4))
- shuffle(L)
- L
pour voir apparaître une permutation aléatoire.
Génération des parties à 2, 3 éléments d’un ensemble fini
Une partie à 2 ou 3 éléments de E est une combinaison de E à 2 ou 3 éléments. itertools permet donc d’en fabriquer avec
- list(combinations(E,2))
- list(combinations(E,3))
Comme itertools donne des couples ou des triplets au lieu d’ensembles, on peut obtenir le même résultat sous forme de compréhensions :
- {(a,b) for a in E for b in E if a!=b}
- {(a,b,c) for a in E for b in E for c in E if a!=b!=c!=a}
Suites
Recherche de seuils
La recherche de seuils, pour une suite géométrique, a de fait été déjà vue en 2nde sous la forme de deux algorithmes :
- L’un cherchant quand on passe au-dessus du seuil ; il s’applique à une suite de raison supérieure à 1 ;
- L’autre cherchant quand on passe en-dessous du seuil ; il s’applique à une suite de raison inférieure à 1.
On peut regrouper les deux en un seul, en testant la raison dans la fonction :
- def mu(terme,raison,seuil):
- assert terme>0 and raison>0 and raison!=1
- n = 0
- if raison>1:
- while terme<seuil:
- n += 1
- terme *= raison
- if raison<1:
- while terme>seuil:
- n += 1
- terme *= raison
- return n
La fonction μ de Kleene, fondamentale en théorie de la calculabilité, est définie par la recherche d’un seuil (par exemple le seuil 1 dans la suite de Collatz).
Recherche de valeurs approchées de π, e, √2,(1+√5)/2, ln(2),etc
π
Un premier algorithme de calcul de π a été vu en 1ère, il s’agit de la méthode d’Archimède. En voici un autre, plus récent de 18 siècles puisque dû à François Viète :
- def viete(n):
- c,p = 0,1
- for k in range(n):
- c += 1
- c /= 2
- c = c**0.5
- p *= c
- return 2/p
En fait l’algorithme de Viète calcule une approximation de π² mais l’extraction d’une racine (par exemple avec la méthode de Heron) permet ensuite de calculer une approximation de π.
e
On a vu en 1ère que lorsqu’il était enfant, Euler avait impressionné Jacques Bernoulli par son talent en mathématiques. Il ne pouvait donc ignorer que la suite (1+1/n)n tend vers e, quoiqu’assez lentement. Parmi ses innombrables recherches, il a développé l’expression (1+1/n)n et en a déduit une autre suite tendant vers e : Celle de la somme des inverses des factorielles :
- def approx_e(n):
- return sum(1/factorial(k) for k in range(n+1))
L’algorithme a été traité en TP. Il converge si vite qu’il permet de calculer e avec beaucoup de décimales, avec le module decimal de Python.
√2
On a vu en 2nde un algorithme permettant de calculer la racine de 2 : La méthode de Heron (qu’on a généralisée en 1ère avec la méthode de Newton). Une autre méthode est basée sur le développement de √2 en fractions continues, déjà traité au bac S. L’idée remonte à Brahmagupta qui cherchait à résoudre des équations diophantiennes du type y²=2x²+k : Si (x,y) est une solution proche de l’asymptote d’équation y=x√2, la fraction y/x est proche de √2 ; or Brahmagupta a trouvé une méthode permettant, à partir d’une solution, d’en obtenir une autre, et ainsi de construire une suite tendant vers √2.
Par exemple, si r est une approximation rationnelle de √2, il en est de même pour (r+2)/(r+1) [12]. En appelant x et y le numérateur et le dénominateur de r, on a l’algorithme suivant :
- def approx_racine_2(n):
- (x,y) = (1,1)
- for k in range(n):
- (x,y) = (x+2*y,x+y)
- return x/y
L’image d’une fraction par l’homographie (r+2)/(r+1) étant une fraction, on a une suite de fractions tendant vers √2.
φ
On calcule le nombre d’Or comme fraction de deux nombres de Fibonacci successifs :
- def approx_phi(n):
- x = 1
- for k in range(n):
- x = 1/x
- x += 1
- return x
ln(2)
Pour calculer une valeur approchée de ln(2) il y a plusieurs méthodes vues ci-après (calcul intégral, calcul de primitive par la méthode d’Euler, Brouncker, résolution de l’équation ex=2). En voici une autre, basée sur les séries : ln(2) est la somme alternée des inverses des entiers (Newton) :
- def approx_ln_2(n):
- S = 0
- for k in range(1,n+1):
- S -= (-1)**k/k
- return S
Fonctions
Pour calculer des valeurs approchées de ln(2) on va résoudre l’équation f(x)=0 où f est la fonction qui, à x, associe ex-2 qui s’annule en ln(2) :
- def f(x):
- return exp(x)-2
Méthode de dichotomie
Dans les sujets de bac, le test porte sur le signe de f(a)×f(m) (m étant la moyenne des bornes a et b de l’intervalle). Cela permet de savoir si f(a) et f(m) sont, ou non, du même signe, sans présumer du sens de variation de f. Mais c’est un peu compliqué alors ici on va écrire un algorithme de dichotomie pour une fonction croissante :
- def zpd(f,a,b):
- assert a<b and f(a)<0 and f(b)>0
- while b-a>1e-6:
- m=(a+b)/2
- if f(m)<0:
- a = m
- else:
- b = m
- return (a+b)/2
Le nom zpd signifie zéro par dichotomie puisque la valeur ln(2) qui annule f est un « zéro » de f. Un zpd(f,0,1)
donne 0.6931471824645996
comme valeur approchée de ln(2).
Méthode de Newton, méthode de la sécante
- def zpn(f,fprime,a):
- while abs(f(a))>1e-12:
- a = a-f(a)/fprime(a)
- return a
zpn(f,exp,0)
donne 0.6931471805600254
. La fonction zpn a été nommée ainsi (zéro par Newton) en hommage à Isaac Newton qui a publié la méthode en 1671. L’idée est la suivante : Si on ne sait pas calculer l’abscisse du point d’intersection de la courbe avec l’axe des abscisses, on remplace la courbe par sa tangente ce qui, si celle-ci est suffisamment pentue, donne une meilleure approximation de la solution ; puis on recommence ...
L’inconvénient de la méthode de Newton est qu’elle nécessite de connaître (savoir calculer) la dérivée fprime de f. Avec l’exemple choisi on sait que la dérivée de ex-2 est ex et de fait dans l’appel à zpn on a remplacé fprime par exp. Mais si on essaye de calculer numériquement la dérivée, on risque de perdre en précision.
D’où l’idée de la méthode de la sécante : On remplace la courbe par sa sécante (il faut donc 2 abscisses au départ) :
- def zps(f,a,b):
- assert a<b and f(a)<0 and f(b)>0
- while abs(b-a)>1e-12:
- (a,b) = (b,b-(b-a)/(f(b)-f(a))*f(b))
- return (a+b)/2
Avec zps(f,0,1)
on obtient 0.6931471805599523
qui est très proche de la valeur obtenue par la dichotomie.
Algorithme de Briggs pour le calcul du logarithme
L’algorithme, déjà programmé en pseudo-code avec CaRMetal (en Tle STI2D), s’écrit ainsi en Python :
- def ln(x):
- assert x>0
- n = 0
- while not(0.999<x<1.001):
- x **= 0.5
- n += 1
- return (x-1)*2**n
Il est dû à John Briggs qui est allé le présenter à John Napier vers le début du 17e siècle. L’idée est plutôt simple : Quelle que soit la valeur de u0=x, la suite un+1=√(un) tend vers 1 donc il arrive un moment, tôt ou tard, où ln(un)≈un-1. Et comme chaque extraction de racine dans le domaine des exponentielles équivaut à une division par 2 dans le domaine des logarithmes, il suffit de compter les extractions de racines avant de soustraire 1, puis de multiplier par 2 autant de fois qu’on a extrait une racine.
Avec ln(2)
on obtient 0.6933818297000016
.
Mais on a 0.6931474097073078
avec cette variante :
- def ln(x):
- assert x>0
- n = 0
- while not(0.999999<x<1.000001):
- x **= 0.5
- n += 1
- return (x-1)*2**n
Intégrales
Résolution par la méthode d’Euler de y’=f, de y’=ay+b
Méthodes des rectangles, des milieux, des trapèzes
Méthode de Monte-Carlo
Le script vu en 1ère permet de calculer ln(2) par la méthode de Monte-Carlo. Mais comme il est prévu pour calculer des intégrales entre 0 et 1, il faut changer de variable : Au lieu d’intégrer la fonction inverse entre 1 et 2, on va intégrer entre 0 et 1 la fonction qui, à x, associe 1/(x+1). Depuis Alonzo Church, on note cette fonction λ x. 1/(x+1) et, en Python, c’est lambda x: 1/(x+1)
. C’est donc l’expression suivante :
monte_carlo(lambda x: 1/(x+1), 10000)
qui produit l’estimation suivante de ln(2) :
Dans ce cas, l’estimation 0,6887 est trop petite mais en répétant l’expérience, on trouvera typiquement aussi souvent plus que 0,69, que moins comme ici. Pour voir la théorie de cette approximation avec la loi binomiale, voir cette excellente vidéo d’Olivier Roizès.
Algorithme de Brouncker pour le calcul de ln(2)
Dans les années 1660, Isaac Newton a repris la série 1-1/2+1/3-1/4 etc convergeant vers ln(2) mais il a regroupé les termes 2 à 2 :
- 1-1/2 = 1/1/2
- 1/3-1/4 = 1/3/4
- 1/5-1/6 = 1/5/6
- 1/7 - 1/8 = 1/7/8
- etc
Or la suite de terme général 1/(2k+1)/(2k+2) tend vers 0 plus vite que la suite 1/(n+1) et on peut calculer ln(2) avec une meilleure précision ou plus de rapidité. En 1668, William Brouncker, un autre contemporain de Newton, a publié la série ci-dessus et cette méthode pour calculer l’intégrale de 1/x entre 1 et 2 :
- def brouncker(n):
- S = 0
- for k in range(n):
- S += 1/(2*k+1)/(2*k+2)
- return S
Avec brouncker(100)
on obtient 0.6906534304818241
comme approximation de ln(2).
Le problème pour Newton c’est qu’il n’avait rien publié là-dessus avant Brouncker ...
Probas
Dans le cadre d’une résolution de problème modélisé par une variable binomiale X, calculer numériquement une probabilité du type P(X=k),P(X≤k), P(k≤X≤k′), en s’aidant au besoin d’un algorithme ; chercher un intervalle I pour lequel la probabilité P(X∈I) est inférieure à une valeur donnée α, ou supérieure à 1−α
Pour calculer plus vite mais de façon moins précise [13] les coefficients binomiaux, on peut utiliser les factorielles obtenues par from math import *
. Le nombre de manières de choisir p objets (succès) parmi n (nombre total de tentatives) est
- def C(n,p):
- return factorial(n)/factorial(p)/factorial(n-p)
(« C » comme « combinaisons »)
Ceci permet de calculer la loi binomiale, comme liste des probabilités que X=k, indexées par les valeurs successives de k :
- def loi_binomiale(n,p):
- return [C(n,k)*p**k*(1-p)**(n-k) for k in range(n+1)]
Par exemple, la probabilité d’avoir 3 dés donnant 6 si on lance 10 dés (cubiques) est loi_binomiale(10,1/6)[3]
:
La probabilité d’avoir moins de 3 dés (soit, au maximum 2 dés) donnant 6 s’obtient avec sum(loi_binomiale(10,1/6)[:3])
:
Et la probabilité d’avoir entre 1 et 3 dés donnant 6 s’obtient avec sum(loi_binomiale(10,1/6)[1:4])
:
Galton
Simulation d’une marche aléatoire
en dimension 1
Une marche aléatoire en dimension 1, c’est la trajectoire d’une tortue qui, à chaque étape, choisit au hasard si elle va avancer ou reculer. Quelque chose comme
- for n in range(20):
- if random()<0.5:
- backward(10)
- else:
- forward(10)
(on suppose les fonctions des modules turtle et random chargées)
Mais comme reculer de 10 pixels revient à avancer de -10 pixels, on simule la même expérience avec
- for n in range(20):
- forward(choice([-10,10]))
(on choisit au hasard une valeur dans la liste comprenant -10 et 10)
Et finalement puisque l’abscisse de la tortue est la somme des déplacements successifs, on peut connaître celle-ci en faisant
sum([choice([-10,10]) for n in range(20)])
On peut encore mettre la longueur 10 en facteur et obtenir
10*sum([choice([-1,1]) for n in range(20)])
ou
10*sum([randrange(-1,2,2) for n in range(20)])
ou encore
sum([randrange(-10,11,20) for n in range(20)])
Ces simplifications allant vers l’abstrait ont pour but essentiel de donner des scripts plus rapides et ainsi de pouvoir répéter plus souvent l’expérience aléatoire, faisant ainsi des statistiques sur des échantillons de taille suffisamment importante pour atténuer l’effet de la fluctuation d’échantillonnage.
Ceci dit, à fonction affine près, sum([randrange(-1,2,2) for n in range(20)])
est la somme de variables de Bernoulli indépendantes et la marche aléatoire en dimension 1 suit à fonction affine près une loi binomiale, ce qui explique ce qui suivra.
en dimension 2
Maintenant on fait avancer la tortue d’une distance positive (par exemple 4) mais à chaque fois on l’oriente au hasard vers un des points cardinaux :
- for n in range(200):
- setheading(randrange(0,271,90))
- forward(4)
Simulation de la planche de Galton
Voici un script Python simulant une planche de Galton à N billes, où 1<N<159.
Si on fait abstraction du mouvement vertical de la bille, sa position finale dépend uniquement du nombre de fois qu’elle est allée à gauche et du nombre de fois qu’elle est allée à droite : La composante horizontale du mouvement de la bille est une marche aléatoire à une dimension.
Cela peut se faire avec la tortue de Python :
Par exemple, en reprenant les dimensions dessinées par Sir Francis Galton lui-même en 1889, on compte 13 clous, dont 6 de chaque côté du clou de départ :
On va donc 6 fois de suite faire avancer ou reculer la tortue au hasard, en partant de l’abscisse 6 :
- planche = [0]*13
- for bille in range(100):
- planche[6+sum([randrange(-1,2,2) for k in range(6)])] += 1
- print(planche)
Variables aléatoires
Simulation d’un échantillon d’une variable aléatoire
Plus généralement, si on connaît uniquement la loi d’une variable aléatoire (et pas son mode de fabrication), on peut simuler la variable aléatoire à l’aide de sa fonction de répartition, qui dans le cas du nombre de diviseurs d’un dé icosaédrique, est, rappelons-le :
{1: 0.05, 2: 0.45, 3: 0.55, 4: 0.8, 5: 0.85, 6: 1.0}
Comment faire pour utiliser cette fonction de répartition ? On applique son inverse à une variable aléatoire uniforme random()
:
- def nda():
- r = random()
- if r<0.05:
- return 1
- elif r<0.45:
- return 2
- elif r<0.55:
- return 3
- elif r<0.8:
- return 4
- elif r<0.85:
- return 5
- else:
- return 6
La fonction nda() (« nombre de diviseurs aléatoire ») a pour but de produire une réalisation de la variable aléatoire. Pour simuler un échantillon de taille 20, il suffit alors d’évaluer la compréhension
[ nda() for k in range(20) ]
Simuler N échantillons de taille n d’une variable aléatoire d’espérance µ et d’écart type σ. Calculer l’écart type s de la série des moyennes des échantillons observés, à comparer à σ/√n. Calculer la proportion des échantillons pour lesquels l’écart entre la moyenne et µ est inférieur ou égal à ks, ou à kσ/√n, pour k = 1,2,3
Pour simuler un échantillon de taille n=100 par exemple, il suffit de remplacer ci-dessus 20 par 100. Pour simuler N échantillons de taille 100 on peut utiliser cette fonction :
- def simul_echant(N):
- return [ [ nda() for k in range(100) ] for n in range(N) ]
Et pour obtenir la liste des moyennes :
- def simul_moyennes(N):
- return [ moyenne([ nda() for k in range(100) ]) for n in range(N) ]
Comme √100=10 et que l’écart-type du nombre de diviseurs d’un dé icosaédrique est environ 1,52, le nombre avec lequel comparer l’écart-type des moyennes est 0,152. Or on obtient quelque chose comme 0.15568434699730013
avec écart_type(simul_moyennes(1000))
.
Comme l’espérance du nombre de diviseurs du dé icosaédrique est 3,3 on va demander à Python de compter les échantillons dont la moyenne est entre 3,3-k×0,152 et 3,3+k×0,152. Pour cela on construit par compréhension une liste des tests réalisés sur les 1000 échantillons. C’est une liste de 1000 booléens, dans laquelle Python va compter les occurences de True
ce qui permet d’avoir les fréquences :
- def propfluct(k):
- liste = [ 3.3-k*0.152<moyenne([nda() for k in range(100)])<3.3+k*0.152 for n in range(1000) ]
- return liste.count(True)/1000
Avec [ propfluct(k) for k in [1,2,3] ]
on obtient [0.686, 0.956, 0.995]
qui est proche des résultats sur les variables aléatoires normales.
Calculer la probabilité de (| Sn − pn | > √n), où Sn est une variable aléatoire qui suit une loi binomiale ℬ(n,p). Comparer avec l’inégalité de Bienaymé-Tchebychev
L’événement contraire de (| Sn − pn | > √n) est (| Sn − pn | ≤ √n) soit (pn-√n≤Sn≤pn+√n). Comme il y a peu de chances que les bornes de cet intervalle soient entières, on les tronque avec la fonction int ce qui permet d’utiliser ce qui a été vu auparavant, sur la loi binomiale :
- def Bienaymé(n,p):
- bI = int(p*n-racine(n))+1
- bS = int(p*n+racine(n))
- return 1-sum( loi_binomiale(n,p)[bI:bS+1] )
On constate que le résultat est souvent nettement plus petit que la majoration donnée par Jules Bienaymé et Pafnouti Tchebychev à la fin du 19e siècle.
Maths expertes
Complexes
Pour calculer avec les nombres complexes, il suffit d’importer le module cmath au lieu de math. Alors le complexe i se note 1j, la fonction abs calcule le module, la fonction phase calcule l’argument d’un complexe.
PGCD
Le livre VII des éléments d’Euclide a été évoqué en 2nde. Euclide se posait un problème d’origine géométrique : Comment paver une salle de 39m sur 24m en utilisant des carreaux les plus grands possibles ? La constatation importante d’Euclide est que la taille des carreaux ne change pas si on enlève un carré (en l’occurence de 24m de côté) le plus grand possible. Autrement dit, le pgcd de 39 et 24 est le même que le pgcd de 39-24=15 et 24. Le test d’arrêt de l’algorithme est le fait qu’il reste un carré (donc juste avant le rectangle était formé de deux carrés). En Python cela donne
- def pgcd(a,b):
- while a!=b:
- if a>=b:
- a -= b
- else:
- b -= a
- return a
Mais Euclide est aussi connu pour avoir utilisé abondamment la division avec reste (dite euclidienne) et avec celle-ci on peut calculer le pgcd plus rapidement :
- def pgcd(a,b):
- while b>0:
- (a,b) = (b,a%b)
- return a
C’est ce nouvel algorithme qu’on appelle maintenant d’Euclide alors que ce n’est pas ainsi qu’Euclide opérait. Aujourd’hui il est plus rapide que l’algorithme originel d’Euclide (parce qu’il boucle moins souvent) mais autrefois les ordinateurs mettaient beaucoup plus de temps à effectuer des divisions que des soustractions et l’algorithme originel d’Euclide était efficace.
En étudiant des problèmes de transvasement de liquides, Claude-Gaspard Bachet de Méziriac a découvert qu’il est possible, grâce à l’algorithme d’Euclide, d’obtenir le pgcd de a et b par combinaison linéaire de a et b à coefficients entiers (les coefficients de Bézout). Mais pourquoi ce travail est-il attribué à Étienne Bézout ? Peut-être parce que ce dernier était un grand mathématicien du siècle des lumières, alors que Bachet n’était qu’un amateur de mathématiques récréatives ?
Toujours est-il que l’algorithme d’Euclide étendu est en ligne sur le site de la calculatrice Numworks. Et en voici une version récursive extraite d’un cours de l’université de Lyon :
- def bezout(a,b):
- if b==0:
- return (a,1,0)
- (d,x,y) = bezout(b,a%b)
- (x,y) = (y,x-a//b*y)
- assert a*x+b*y==d
- return (d,x,y)
Commentaire de Bernard Parisse :
L’article ne donne que la version récursive. Je pense que la version itérative mériterait d’y être, avec un return de la liste [u,v,r]
[14] :
- def bezout(a,b):
- l1=[1,0,a]
- l2=[0,1,b]
- while l2[2]!=0:
- q=l1[2]//l2[2]
- l1,l2=l2,[l1[0]-q*l2[0],l1[1]-q*l2[1],l1[2]-q*l2[2]]
- return l1
En syntaxe Python dans Xcas, la version itérative de Bezout est plus simple, parce que les listes sont considérées comme des vecteurs pour les opérations d’addition et multiplication, cela donne :
- def bezout(a,b):
- l1=[1,0,a]
- l2=[0,1,b]
- while l2[2]!=0:
- q=l1[2]//l2[2]
- l1,l2=l2,l1-q*l2
- return l1
J’en profite pour signaler que ces scripts fonctionnent en général aussi dans Xcas et Xcas pour Firefox, qui permet d’échanger facilement des sessions et de les exécuter dans un navigateur sur divers appareils. Par exemple, pour Bezout écrit en syntaxe Python ou en syntaxe Xcas (langage naturel) voir ici
Crible
Eratosthène est connu pour avoir été le premier à estimer le rayon de la Terre par des méthodes mathématiques, mais aussi pour avoir trouvé un moyen rapide de calculer, non pas un nombre premier, mais tous (du moins inférieurs à un seuil N) :
- def crible(N):
- """nombres premiers jusque N"""
- liste = list(range(2,N+1))
- for d in range(2,int(sqrt(N))+1):
- liste = [x for x in liste if x<=d or x%d!=0]
- return liste
Comme on le voit, la fonction associe, à tout entier N, une liste d’entiers (premiers inférieurs ou égaux à N). Il semble donc naturel d’utiliser les listes en compréhension : À chaque étape, on crée une nouvelle liste en n’y gardant que les entiers inférieurs à d (on sait qu’ils sont premiers et les entiers non divisibles par d. Autrement dit, on élimine (d’où le nom de crible) les entiers divisibles par d.
Une variante ne nécessitant pas de calculer sur des listes a été proposée par Tony Hoare puis répandue par Donald Knuth. En voici une version manipulable en ligne.
Commentaire de Bernard Parisse :
Cet algorithme est assez différent de l’algorithme classique tel que décrit par exemple sur wikipedia (Crible d’Ératosthène), qui est plus simple car il ne necessite de faire que des additions (il se fait donc facilement avec papier-crayon), alors que celui présenté ici fait des divisions. Or l’algorithme plus simple est aussi plus efficace. En voici une implémentation :
- def premiers(n):
- prem=list(range(0,n))
- k=2
- while k*k<n :
- for j in range(2*k,n,k):
- prem[j]=0
- k += 1
- while prem[k]==0: k += 1
- return [k for k in prem if k!=0]
Pour calculer premiers(10**6)
, le temps d’exécution est de 0,4 secondes sur ma machine contre 15 secondes pour celui avec des divisions. On peut améliorer un peu celui avec les divisions en ne testant les d que s’ils sont dans liste,
- from math import sqrt
- def crible(N):
- """nombres premiers jusque N"""
- liste = list(range(2,N+1))
- for d in range(2,int(sqrt(N))+1):
- if d in liste:
- liste = [x for x in liste if x<=d or x%d!=0]
- return liste
(le rajout est la ligne 6)
avec l=crible(10**6)
on passe ainsi de 15 secondes à 6 secondes, mais on est encore un ordre de grandeur au-dessus du crible classique.
La raison de cette inefficacité vient du nombre d’opérations à effectuer : pour la fonction premiers, le nombre d’itérations de la boucle intérieure est de n*(1/2+1/3+...+1/int(sqrt(n))) donc d’ordre n*ln(sqrt(n)), essentiellement linéaire en n [15]. Pour le deuxième (non optimisé), la liste en compréhension est en fait une boucle exécutée √(N) fois sur les éléments de liste, dont le cardinal est supérieur ou égal au nombre de nombres premiers inférieurs ou égaux à N soit asymptotiquement au moins √(N)*N/ln(N)). On a donc au moins √(N)/ln(N)² fois plus d’itérations qu’avec le crible classique, qui plus est avec des divisions au lieu d’additions. En fait, l’algorithme crible revient essentiellement à tester la primalité de tous les entiers de 2 à N individuellement (il est même un peu moins efficace, car on reteste des nombres dont on sait déjà qu’ils sont premiers).
L’analyse de complexité peut probablement être abordée en maths expertes, si on admet que le nombre de premiers inférieurs à N est asymptotiquement équivalent à N/ln(N).
Ceci illustre aussi à mon sens qu’il faut mettre en garde les élèves (on est en maths expertes ici) sur l’utilisation de fonctionnalités avancées d’un langage de programmation qui rendent le code plus compact mais masquent l’existence en interne de boucles. Ceci tant que les fondamentaux ne sont pas acquis.
facteurs
Comme certains facteurs premiers peuvent apparaître plusieurs fois, la liste est une structure de données adaptée à ce problème. On va donc créer une liste puis, pour chaque diviseur d possible, regarder si n est divisible par d. Si oui, on effectue la division et on ajoute (append) d à la liste des diviseurs. Sinon, on passe au diviseur suivant :
- def facteurs_premiers(n):
- liste = []
- d = 2
- while d<=n:
- if n%d==0:
- liste.append(d)
- n //= d
- else:
- d += 1
- return liste
Exercice : Prouver que les diviseurs stockés dans la liste sont tous premiers.
arithmétique
Les activités de cryptographie présentées dans cet article peuvent aisément se transposer à Python, avec ces fonctions :
- ord() transforme une lettre en un entier (par exemple
ord('e')
donne 101) - chr() au contraire transforme un entier (inférieur à 128) en un mot d’une lettre (par exemple
chr(101)
donne ’e’)
Les équations de Pell-Fermat données au bac S ont déjà été résolues en Python.
La recherche de triplets pythagoriciens par double balayage peut se faire ainsi :
- for x in range(1,200):
- for y in range(1,x):
- z = racine(x**2+y**2)
- if z==int(z):
- print(x,y,int(z))
Graphes
La connexité des graphes et la notion de graphe eulérien ont été étudiées avec CaRMetal.
Les chaînes de Markov ont été étudiées dans cet article et simulées en Python ici. Note : Andrei Markov a inventé ce concept de parcours aléatoire d’un graphe orienté pour essayer de caractériser le style d’un auteur plus finement que la simple fréquence d’apparition des lettres (vu en 2nde : Celle-ci caractérise plutôt la langue que l’auteur, selon Al-Kindi). C’est donc un concept linguistique qui est à l’origine du succès de Larry Page...
L’algorithme de Page justement, est programmé en Python dans ce manuel de SNT.
Le modèle des urnes d’Ehrenfest a été simulé ici avant d’être programmé en Python sur le site de la Numworks (simulable en ligne).
Nicolas Patrois a programmé l’algorithme de Dijkstra sur Numworks.
C’est Arthur Cayley qui a le premier pensé à représenter un graphe (notion dont il est l’auteur) par une matrice (notion dont il est également l’auteur), dite d’adjacence. Les puissances de cette matrice permettent d’évaluer le nombre de chemins de longueur donnée. Cette application en ligne permet de vérifier sur des exemples les applications de ce concept à la connexité du graphe et à la notion de distance entre sommets, utile en SNT.
Maths complémentaires
Suites
Recherche de seuils
Voir l’onglet suites de la spé 1ère et de la Tle spé.
Pour une suite récurrente un+1=f(un), calcul des termes successifs
Voir l’onglet suites de 1ère spé
Recherche de valeurs approchées de constantes mathématiques, par exemple π, ln 2, √2
Voir l’onglet suites de Tle spé.
Fonctions
Méthodes de recherche de valeurs approchées d’une solution d’équation du type f(x)=k:balayage, dichotomie, méthode de Newton
Le balayage a été vu en 2nde ; la dichotomie et Newton dans l’onglet fonctions de Tle spé
Algorithme de Briggs pour le calcul de logarithmes
Voir l’onglet fonctions de Tle spé.
Sur des exemples, résolution approchée d’une équation différentielle par la méthode d’Euler
Voir l’onglet intégrales de Tle spé.
Intégration
Méthode des rectangles, des trapèzes
Voir l’onglet intégrales de Tle spé.
Méthode de Monte-Carlo pour un calcul d’aire
Voir l’onglet intégrales de Tle spé.
Probas
Simulation d’une variable de Bernoulli ou d’un lancer de dé (ou d’une variable uniforme sur un ensemble fini) à partir d’une variable aléatoire uniforme sur [0 ;1]
Après un from random import *
, la fonction random()
simule une variable aléatoire uniforme sur [0,1]. La probabilité qu’elle soit inférieure à un seuil p est ∫0pdx=p. Autrement dit l’événement X = (random()<p)
a la probabilité p de se réaliser (être égal à True
) et 1-p de ne pas se réaliser (être égal à False
).
Or la fonction int()
de Python (transformation en entier), transforme True
en 1 et False
en 0. On peut donc créer une variable aléatoire de Bernoulli en appliquant la fonction int()
à l’événement X.
Le résultat est présenté dans l’onglet thèmes, dans la partie sur l’échantillonnage.
Pour simuler un lancer de dé, on peut faire
- def dé():
- return int(1+6*random())
Simulation du comportement de la somme de n variables aléatoires indépendantes et de même loi
Voir la partie sur l’échantillonnage de l’onglet thèmes pour un exemple.