Bibm@th

Forum de mathématiques - Bibm@th.net

Bienvenue dans les forums du site BibM@th, des forums où on dit Bonjour (Bonsoir), Merci, S'il vous plaît...

Vous n'êtes pas identifié(e).

#1 21-05-2010 15:34:27

yoshi
Modo Ferox
Inscription : 20-11-2005
Messages : 16 946

[Python] Approximations de Pi...

Bonjour,


Approximations via le calcul de la valeur approchée de [tex]\int_0^1 \frac{4}{1+x^2}\,dx[/tex]
Justification :
          [tex]\int_0^1 \frac{4}{1+x^2}\,dx=4\big[arctan(x)\big]_0^1=4\left(\frac{\pi}{4}-0\right)=\pi[/tex]
4 méthodes au choix :
méthode des rectangles, des trapèzes, des tangentes et formule de Simpson composée.
Seul le pas vous sera demandé, exclusivement 0.1, 0.01, 0.001, 0.0001... etc.

# -*- coding: cp1252 -*-

from __future__ import division
from decimal import Decimal
import decimal as d

def titre():
    print
    print
    print
    print "                  *******************************************"
    print "                  *      +- APPROXIMATIONS de Pi -+         *"
    print "                  *******************************************"
    print
    print
    print
    print "Calcul de la valeur approchée de l'intégrale de 0 à 1 de f(x)=4/(1+x^2)"
    print "                         par :"
   
def sortie():
    print
    print
    print "   Au revoir !"

def incorrect():
    print "           ...Désolé, valeur incorrecte. Veuillez recommencer ..."

def pas_nb():
    print "             ... Ce n'est pas un nombre. Veuillez recommencer..."    

def message(cx):
    ch={0:'       ... Pour retourner au menu,',
        1:'... Pour la suite,',2:'... '}
    print
    print "       ",ch[cx],"Presser ENTREE ..."
    raw_input("")

def pas():
    ok=0
    print "\n\n"
    while not ok:
        chp=raw_input("Choisissez un pas (0.1, 0.01, 0.001...) : ")
        l=len(chp)
        try:
            p=float(chp)
            if p>0 and p<1:
                if chp=="0."+(l-3)*"0"+"1" :
                    ok=1
                else:
                    incorrect()
                    message(2)
            else:
                incorrect()
                message(2)          
        except ValueError:
            pas_nb()
            message(1)
    return p

def affiche_resultat(p,s):
    print "\n\n"
    print "-------------------------------------------------------"
    print
    print "            *****************************"
    print "   **********    Résultat du calcul      ***********"
    print "            *****************************"
    print
    print "Valeur calculée de Pi :",s, "avec un pas de",p
    print "Valeur  connue  de Pi : 3.14159265358979323846..."
    print
    print "-------------------------------------------------------\n\n"
         
def rectangles():
    # Somme des aires des rectangles de largeur x-p/2 ; x+p/2
    p=0
    p=pas()
    nt,x,s=int(1/p),p/2,0
    for k in xrange(nt):
        s+=4/(1+(x+k*p)**2)*p
    affiche_resultat(p,s)

def trapezes():
    s=0
    p=pas()
    nt,x,s=int(1/p),p/2,0
    for k in xrange(nt):
        s+=(4/(1+(k*p)**2)+4/(1+((k+1)*p)**2))/2
    s=s*p
    affiche_resultat(p,s)

def tangentes():
    s=0
    p=pas()
    nt,x=int(1/p),p/2
    for k in xrange(nt):
        s+=4/(1+(x+k*p)**2)
    s=s*p
    affiche_resultat(p,s)

def simpson():
    d.getcontext().prec = 21
    s=0
    p=pas()
    a,b,nt=0,1,int(1/p)
    for k in xrange(1,nt):
        s+=Decimal((1+k%2))*(Decimal(4)/Decimal(str((1+(a+k*p)**2))))
    s=Decimal(str(p))/Decimal(3)*(Decimal(4)/Decimal(1+a**2)+Decimal(4)/Decimal(1+b**2)+Decimal(2*s))
    affiche_resultat(p,s)
   
fin=5
while fin>0:
    titre()
    print "      1.la méthode des rectangles"
    print "      2.la méthode des trapèzes"
    print "      3.la méthode des tangentes"
    print "      4.la formule de Simpson composée"  
    print
    print " 0. Sortie du programme\n"
    chx_menu='5'
    while chx_menu=='5':
        chx_menu=raw_input(" Votre choix : ")
        try:
            fin = int(chx_menu)
            choix2 = {'1': rectangles,
                      '2': trapezes,
                      '3': tangentes,
                      '4': simpson,
                      '0': sortie}
            if not fin in [0,1,2,3,4]:
                chx_menu='5'
                incorrect()
        except ValueError:
            pas_nb()
            chx_menu='5'
    choix2[chx_menu]()

Un exemple de sortie avec la formule de Simpson :

---------------------------------------------------------------

                   *****************************
   **********          Résultat du calcul             ***********
                   *****************************

Valeur calculée de Pi : 3.14159265358979323840 avec un pas de 0.001
Valeur connue de Pi : 3.1415926535897932384626..."...

La troncature à [tex]10^{-19}[/tex] de la valeur calculée est juste !

Le module "decimal" de Python permet, à partir de la version 2.4, d'augmenter la précision d'un calcul en virgule flottante.
A noter aussi que le même calcul en C est précis, avec une déclaration "double s" également jusqu'à [tex]10^{-17}[/tex] pour p=0.001

Commentaires, réactions ?

@+


Arx Tarpeia Capitoli proxima...

Hors ligne

#2 21-05-2010 16:49:23

thadrien
Membre
Lieu : Grenoble
Inscription : 18-06-2009
Messages : 526
Site Web

Re : [Python] Approximations de Pi...

Salut,

Juste pour le fun, ajoute aussi une intégration par une méthode type monté-carlo.

Hors ligne

#3 21-05-2010 17:21:42

yoshi
Modo Ferox
Inscription : 20-11-2005
Messages : 16 946

Re : [Python] Approximations de Pi...

Re,

Je vais chercher ça, connais pas !

Merci.

@+

[EDIT]
J'ai vu...
Je comprends le "pour le fun" de thadrien, et aussi le pourquoi de "Monte Carlo" !
Nan, je vais plutôt chercher une autre méthode -sérieuse- d'approximation d'une intégrale...

Le Bac devient intéressant, Valentin n'aurait pas intérêt à le repasser.
Epreuve au Bac :
http://lyc-eparny.ac-reunion.fr/IMG/pdf/Sujet071.pdf

Dernière modification par yoshi (21-05-2010 18:08:50)


Arx Tarpeia Capitoli proxima...

Hors ligne

#4 24-05-2010 16:57:39

yoshi
Modo Ferox
Inscription : 20-11-2005
Messages : 16 946

Re : [Python] Approximations de Pi...

Re,

Petite explication sur "ma" formule de Simpson :
La formule officielle est (avec m=nb/2, nb étant le nombre total d'intervalles, obligatoirement pair) :
[tex]s=\frac{p}{3}\left(f(a)+f(b)+2\sum\limits_{k=1}^{m-1} f(x_{2k})+4\sum\limits_{k=0}^{m-1} f(x_{2k+1})\right)[/tex]
Je fais ce calcul en 2 fois :
1. Dans un premier temps avec la boucle je calcule [tex]s=\sum\limits_{k=1}^{m-1} f(x_{2k}) \;+\;2\sum\limits_{k=0}^{m-1}f(x_{2k+1})[/tex] via [tex] s= \sum\limits_{k=1}^{m-1} f(x_{2k}) \;+\;\sum\limits_{k=0}^{m-1}2f(x_{2k+1}) [/tex]

2. Dans un deuxième temps, j'apporte le correctif : [tex]s=\frac{p}{3}\left(f(a)+f(b)+2s\right)[/tex]

Certains d'entre vous m'objecteront que j'ai écrit au dessus  que les coefficients des sommes sont 2 et 4  pour k pair et impair (et non pas 1 et 2) et que je ne fais pas le distinguo entre pair et impair...
C'est exact...
Quoique pour le pair/impair ça ne l'est pas, exact :
je calcule en effet [tex] s = \sum\limits_{k=0}^{2m-1}(1+k\%2)\times f(x_k)[/tex]
Et le 1+k%2 vaut 1 si k est pair et 2 si k est impair, : cette astuce me permet d'effectuer le calcul dans la boucle d'une seule traite, sans utiliser if... else et en allant de 0 à 2m - 1 = nb-1

Toutefois, le s à la sortie de la boucle est quand même 2 fois trop petit d'où le +2s du correctif et non pas seulement le +s.
Toutes les mentions Decimal, (que je ne suis pas sûr de maîtriser encore totalement); sont là pour améliorer la représentation décimale, en
Python, des nombres à virgule, voir  à propos du module decimal :
http://docs.python.org/library/decimal.html
ou encore :
https://docs.python.org/2/library/decimal.html

@+


Arx Tarpeia Capitoli proxima...

Hors ligne

#5 11-07-2010 18:16:09

yoshi
Modo Ferox
Inscription : 20-11-2005
Messages : 16 946

Re : [Python] Approximations de Pi...

Bonsoir,


Mise à jour : une erreur s'était glissée à la dernière minute, après les tests, que j'ai rectifiée...
C'était la ligne :
print "Valeur connue de pi :........... etc...

Avec mes excuses

@+


Arx Tarpeia Capitoli proxima...

Hors ligne

Réponse rapide

Veuillez composer votre message et l'envoyer
Nom (obligatoire)

E-mail (obligatoire)

Message (obligatoire)

Programme anti-spam : Afin de lutter contre le spam, nous vous demandons de bien vouloir répondre à la question suivante. Après inscription sur le site, vous n'aurez plus à répondre à ces questions.

Quel est le résultat de l'opération suivante (donner le résultat en chiffres)?
vingt huit moins huit
Système anti-bot

Faites glisser le curseur de gauche à droite pour activer le bouton de confirmation.

Attention : Vous devez activer Javascript dans votre navigateur pour utiliser le système anti-bot.

Pied de page des forums