L'Aiguille de Buffon en Python

Ne converge pas vers Pi...

L’auteur de ce sujet a trouvé une solution à son problème.
Auteur du sujet

Bonjour,

J’essaye de coder une simulation de l’Aiguille de Buffon Où l’on lance des aiguilles sur un plancher. Selon la loi des grands nombre la fréquence des aiguilles tombant entre deux lattes de plancher divisé par le nombre de lancer total est censé tendre vers 2/Pi

Seulement dans ma simulation , et ce même jusqu’à 10000 aiguilles, la valeur de Pi semble converger vers 2,7 .. o_O

J’ai essayé de modifier le nombre de planches, la taille des aiguilles mais rien n’y fait.

Voilà mon code:


from tkinter import *
from random import *
from math import *

root = Tk()
root.title('Aiguille de Buffon')

can = Canvas(root,width=1000,height=1000,bd=3) 

def verifierIntersection(x1aig,y1aig,x2aig,y2aig,x1seg,y1seg,x2seg,y2seg,dessiner):
    #On cherche les équations paramétriques des droites engandrées par les vecteurs
    #aiguille
    xaig = x1aig
    yaig = y1aig
    xvectaig = x2aig-x1aig
    yvectaig= y2aig - y1aig
    # On note le paramètre t
    #autresegment
    xseg = x1seg
    yseg = y1seg
    xvectseg = x2seg - x1seg
    yvectseg = y2seg - y1seg
    # On note le paramètre b
    # Si le déterminant est non nul:
    if (xvectaig*yvectseg-xvectseg*yvectaig != 0):
        # On cherche les paramètres t et b codant le point d'intersection des droites
        b = -(-xvectaig*yaig+xvectaig*yseg+yvectaig*xaig-yvectaig*xseg)/(xvectaig*yvectseg-yvectaig*xvectseg)
        t = -(xaig*yvectseg-xseg*yvectseg-xvectseg*yaig+xvectseg*yseg)/(xvectaig*yvectseg-yvectaig*xvectseg)
        if (0 <= b and 1>=b and 0 <= t and 1>=t):
            #L'intersection est sur les deux segments
            
            if dessiner == "oui":
                
                xint = xaig + t*xvectaig
                yint = yaig + t*yvectaig
                can.create_oval(xint,yint,xint+4,yint+4,fill="red")
            
            
            return True 
        else:
            return False
    else:
        return False



def afficherPlancher():
    
    can.create_line(200,300,800,300)
    can.create_line(200,400,800,400)
    can.create_line(200,500,800,500)
    can.create_line(200,600,800,600)
    can.create_line(200,700,800,700)

def afficherBords():
    can.create_line(200,200,800,200)
    can.create_line(200,200,200,800)
    can.create_line(800,800,800,200)
    can.create_line(800,800,200,800)


def chuteAiguille():
    #Vérifie si l'aiguille sort du cadre
    
    incorrect = True
    while incorrect == True:
    # On part d'un vecteur de norme 200 
    # le vecteur vérifie alors sqrt(x^2+y^2) = 200
    # Les x et y vérifiant cette conditions forment l'équation d'un cercle de rxaigon 200.
    # On isole un des paramètres pour trouver l'autre et voilà qu'on a généré les vecteurs aléatoirements
        x3 = randrange(-100,101,1)
        #change aléatoirement la coordonées y en négative ou positive
        pos = randrange(0,2,1)
        if pos == 0:
            y3 = sqrt(10000-(x3**2))
        else: 
            y3 = -sqrt(10000-(x3**2))
    # Il ne reste plus qu'à générer des points aléatoirement qui seront l'origine de ces vecteurs.
    # Une mauvaise restriction limite la probabilité de tomber dans les coins, ce qui empêche l'expérience
        x1 = randrange(200,800,1)
        y1 = randrange(200,800,1)
        #Coordonnées du point d'arrivée
        x2 = x1 + x3
        y2 = y1 + y3
    
        #Il faut vérifier si l'aiguille ne sorte pas du cadre
        if verifierIntersection(x1, y1, x2, y2, 200, 200, 800, 200,"non") == False and verifierIntersection(x1, y1, x2, y2, 200, 200, 200, 800,"non") == False and verifierIntersection(x1, y1, x2, y2, 800, 800, 800, 200,"non") == False and verifierIntersection(x1, y1, x2, y2, 800, 800, 200, 800,"non") == False:
            incorrect = False
            can.create_line(x1,y1,x2,y2,fill="blue")
            if verifierIntersection(x1,y1,x2,y2,200,300,800,300,"oui") == True or verifierIntersection(x1,y1,x2,y2,200,400,800,400,"oui") == True or verifierIntersection(x1,y1,x2,y2,200,500,800,500,"oui") == True or verifierIntersection(x1,y1,x2,y2,200,600,800,600,"oui") == True or verifierIntersection(x1,y1,x2,y2,200,700,800,700,"oui") == True:
                return True
            else:
                return False

can.pack()

afficherPlancher()
afficherBords()

lancers = 0
collisions = 0
for i in range(10):
    lancers = i
    if chuteAiguille() == True:
        collisions += 1

print(lancers+1, collisions)

 
can.create_text(400,100,text="PI = " + str((lancers/collisions)*2),width=700,font=("Arial",20))


mainloop()

Et voilà une petite image: Simualtion avec 100 aiguilles

Si quelqu’un pourrait m’orienter vers une solution ? Merci d’avance :D

Édité par jerkoco

+0 -0

Salut,

Tu parles d’une norme de 200 dans tes commentaires, on est d’accord qu’il est plutôt question de 100 ?

Sinon je suis justement dubitatif sur le calcul du vecteur aléatoire. Je comprends que tu ne veuilles pas faire intervenir pi pour calculer un angle aléatoire, mais je ne suis pas sûr que le tirage depuis une projection sur l’axe des abscisses donne la même répartition et je me dis que ça peut alors être source d’erreur (à vérifier).

Par contre j’ai pas le courage de relire la fonction verifierIntersection.

je ne suis pas sûr que le tirage depuis une projection sur l’axe des abscisses donne la même répartition

C’est bien le problème, faire un tirage uniforme en xx revient à faire un tirage avec une densité sinusoïdale de l’angle.

Édité par adri1

I don’t mind that you think slowly, but I do mind that you are publishing faster. — W. Pauli

+0 -0

Oui, c’est ce dont je voulais m’assurer, ne réussissant pas vraiment à bien visualiser la chose (mais je sentais que ça n’était pas clair). Et je confirme donc que ça ne correspond pas du tout (c’est très concentré autour des « petits » angles).

Auteur du sujet

Salut,

Tu parles d’une norme de 200 dans tes commentaires, on est d’accord qu’il est plutôt question de 100 ?

Sinon je suis justement dubitatif sur le calcul du vecteur aléatoire. Je comprends que tu ne veuilles pas faire intervenir pi pour calculer un angle aléatoire, mais je ne suis pas sûr que le tirage depuis une projection sur l’axe des abscisses donne la même répartition et je me dis que ça peut alors être source d’erreur (à vérifier).

Par contre j’ai pas le courage de relire la fonction verifierIntersection.

entwanne

Dans le code partagé , il est bien question de norme de 100, Les commentaires datent de la version antérieure du code.

Vous avez pointé le soucis je crois, le problème de répartition viens bien du tirage uniforme en x. Ducoup je ne vois plus trop comment générer mes vecteurs.

+0 -0

Ducoup je ne vois plus trop comment générer mes vecteurs.

Choisis un angle au hasard plutôt qu’un xx.

Aussi, ton code n’est pas lisible du tout donc je suis pas sûr de que tu as fait, mais soit sûr d’avoir un système qui se comporte de façon périodique spatialement (pour mimer un plancher infiniment grand) pour ne pas avoir d’effets de bord qui viendrait perturber la chose. En particulier, il me parait étrange que tu n’as aucune aiguille sur les bords de ta boîte.

Aussi je serais toi je ferais la simulation de façon indépendante de la représentation finale. Et plutôt que calculer la position exacte de la collision, il suffit de vérifier si les deux extrémitée d’une aiguille ont des ordonnées sur la même planche ou pas. Une seule condition ultra simple à la place d’un calcul indigeste…

Édité par adri1

I don’t mind that you think slowly, but I do mind that you are publishing faster. — W. Pauli

+1 -0

Dans l’expérience telle qu’on pouvait la voir au Palais de la Découverte,l’aiguille est attirée par un électroaimant tournant à l’horizontale, tandis qu’un tapis roulant de lames métalliques défile en dessous. Quand l’électroaimant est désactivé, l’aiguille tombe sur le tapis avec une orientation angulaire aléatoire.
Un compteur affiche la proportion de cas où l’aiguille chevauche deux lames.
Le résultat tend vers l’inverse de Pi : 0,318309886183790…

Il se faut s’entraider, c’est la loi de la nature. (Jean de La Fontaine, l’âne et le chien)

+0 -1
Auteur du sujet

Ducoup je ne vois plus trop comment générer mes vecteurs.

Choisis un angle au hasard plutôt qu’un xx.

Aussi, ton code n’est pas lisible du tout donc je suis pas sûr de que tu as fait, mais soit sûr d’avoir un système qui se comporte de façon périodique spatialement (pour mimer un plancher infiniment grand) pour ne pas avoir d’effets de bord qui viendrait perturber la chose. En particulier, il me parait étrange que tu n’as aucune aiguille sur les bords de ta boîte.

Aussi je serais toi je ferais la simulation de façon indépendante de la représentation finale. Et plutôt que calculer la position exacte de la collision, il suffit de vérifier si les deux extrémitée d’une aiguille ont des ordonnées sur la même planche ou pas. Une seule condition ultra simple à la place d’un calcul indigeste…

adri1

Pour le calcul des intersections, je l’ai fait pour le challenge. Sinon merci de votre aide, j’ai bien pris en compte les effets de bords (bien que je n’ai pas encore trop saisi pourquoi il fallait répartir de cette façon le plancher) , et j’ai créé les vecteurs à partir d’un angle aléatoire. Le résultat semble bien tendre vers Pi.

Voilà le code (un peu sale)

from tkinter import *
from random import *
from math import *

root = Tk()
root.title('Aiguille de Buffon')

can = Canvas(root,width=1000,height=1000,bd=3) 

def verifierIntersection(x1aig,y1aig,x2aig,y2aig,x1seg,y1seg,x2seg,y2seg,dessiner):
    #On cherche les équations paramétriques des droites engandrées par les vecteurs
    #aiguille
    xaig = x1aig
    yaig = y1aig
    xvectaig = x2aig-x1aig
    yvectaig= y2aig - y1aig
    # On note le paramètre t
    #autresegment
    xseg = x1seg
    yseg = y1seg
    xvectseg = x2seg - x1seg
    yvectseg = y2seg - y1seg
    # On note le paramètre b
    # Si le déterminant est non nul:
    if (xvectaig*yvectseg-xvectseg*yvectaig != 0):
        # On cherche les paramètres t et b codant le point d'intersection des droites
        b = -(-xvectaig*yaig+xvectaig*yseg+yvectaig*xaig-yvectaig*xseg)/(xvectaig*yvectseg-yvectaig*xvectseg)
        t = -(xaig*yvectseg-xseg*yvectseg-xvectseg*yaig+xvectseg*yseg)/(xvectaig*yvectseg-yvectaig*xvectseg)
        if (0 <= b and 1>=b and 0 <= t and 1>=t):
            #L'intersection est sur les deux segments
            
            if dessiner == "oui":
                
                xint = xaig + t*xvectaig
                yint = yaig + t*yvectaig
                can.create_oval(xint,yint,xint+4,yint+4,fill="red")
            
            
            return True 
        else:
            return False
    else:
        return False



def afficherPlancher():
    
    can.create_line(200,400,800,400)
    can.create_line(200,500,800,500)
    can.create_line(200,600,800,600)
    can.create_line(200,700,800,700)

def afficherBords():
    can.create_line(200,300,800,300,fill="red")
    can.create_line(200,300,200,750,fill="red")
    can.create_line(800,750,800,300,fill="red")
    can.create_line(800,750,200,750,fill="red")


def chuteAiguille():
    #Vérifie si l'aiguille sort du cadre
    
    incorrect = True
    while incorrect == True:
        angle = uniform(0,2*pi)
    # On part d'un vecteur de norme 200 
    # le vecteur vérifie alors sqrt(x^2+y^2) = 200
    # Les x et y vérifiant cette conditions forment l'équation d'un cercle de rxaigon 200.
    # On isole un des paramètres pour trouver l'autre et voilà qu'on a généré les vecteurs aléatoirements
        x3 = 100*cos(angle)
        #change aléatoirement la coordonées y en négative ou positive
        y3 = 100*sin(angle)
    # Il ne reste plus qu'à générer des points aléatoirement qui seront l'origine de ces vecteurs.
    # Une mauvaise restriction limite la probabilité de tomber dans les coins, ce qui empêche l'expérience
        x1 = randrange(200,800,1)
        y1 = randrange(300,750,1)
        #Coordonnées du point d'arrivée
        x2 = x1 + x3
        y2 = y1 + y3
    
        #Il faut vérifier si l'aiguille ne sorte pas du cadre
        if verifierIntersection(x1, y1, x2, y2, 200,300,800,300,"non") == False and verifierIntersection(x1, y1, x2, y2, 200,250,200,750,"non") == False and verifierIntersection(x1, y1, x2, y2, 800,750,800,250,"non") == False and verifierIntersection(x1, y1, x2, y2, 800,750,200,750,"non") == False:
            incorrect = False
            can.create_line(x1,y1,x2,y2,fill="blue")
            if verifierIntersection(x1,y1,x2,y2,200,400,800,400,"oui") == True or verifierIntersection(x1,y1,x2,y2,200,500,800,500,"oui") == True or verifierIntersection(x1,y1,x2,y2,200,600,800,600,"oui") == True or verifierIntersection(x1,y1,x2,y2,200,700,800,700,"oui") == True:
                return True
            else:
                return False

can.pack()

afficherPlancher()
afficherBords()

lancers = 0
collisions = 0
for i in range(1000):
    lancers = i
    if chuteAiguille() == True:
        collisions += 1

print(lancers+1, collisions)


can.create_text(400,100,text="PI = " + str((lancers/collisions)*2),width=700,font=("Arial",20))


mainloop()
+0 -0

Pour le calcul des intersections, je l’ai fait pour le challenge.

Ton calcul me parait extrêmement compliqué. La seule chose qu’il y a besoin de calculer, c’est son abscisse, puisque l’ordonnée est déjà connue (c’est celle de la rainure du plancher). Il me parait beaucoup plus simple de faire une bête condition (du genre int(y1 / length) != int(y2 / length)) pour vérifier si il y a collision ou non, puis calculer l’abscisse qui doit être à vue de nez x1 + (x2 - x1) * (y - y1) / (y2 - y1).

j’ai bien pris en compte les effets de bords (bien que je n’ai pas encore trop saisi pourquoi il fallait répartir de cette façon le plancher)

Honnêtement ton code n’est pas suffisament lisible pour savoir ce que tu as fais exactement, mais si il faut faire un plancher périodique (ou qui se comporte pareil), c’est juste parce que la probabilité de toucher une limite de planche va tendre vers 2/π2/\pi uniquement dans ce cas là. Si on a un plancher borné, le fait qu’il y a des bords modifie la probabilité de toucher une limite de planche.

Par ailleurs, si ton problème est résolu, merci de marquer le sujet comme tel avec le bouton dans le menu de gauche.

I don’t mind that you think slowly, but I do mind that you are publishing faster. — W. Pauli

+0 -0
Connectez-vous pour pouvoir poster un message.
Connexion

Pas encore membre ?

Créez un compte en une minute pour profiter pleinement de toutes les fonctionnalités de Zeste de Savoir. Ici, tout est gratuit et sans publicité.
Créer un compte