Marc BUFFAT

Professeur au département de Mécanique, Lyon 1 e-mail

Blog scientifique et pédagogique utilisant des notebooks IPython et Linux

Cours Python scientifique


In [1]:
%matplotlib inline
#%autosave 300
import numpy as np
import scipy as sp
import matplotlib.pyplot as plt
from matplotlib import rcParams
rcParams['font.family'] = 'serif'
rcParams['font.size'] = 14
from IPython.core.display import HTML
from IPython.display import display
from matplotlib import animation
#from JSAnimation import IPython_display
css_file = 'style.css'
HTML(open(css_file, "r").read())
Out[1]:

Programmation et Algorithme

Marc BUFFAT, dpt mécanique, Université Claude Bernard Lyon 1

Python comme langage

Créé en 1989 par Guido van Rossum

cartoon from a webcomic of romance, sarcasm, math, and language: xkcd

Les bonnes pratiques de la programmation

les rêgles de programmation sous Python

In [2]:
import this
The Zen of Python, by Tim Peters

Beautiful is better than ugly.
Explicit is better than implicit.
Simple is better than complex.
Complex is better than complicated.
Flat is better than nested.
Sparse is better than dense.
Readability counts.
Special cases aren't special enough to break the rules.
Although practicality beats purity.
Errors should never pass silently.
Unless explicitly silenced.
In the face of ambiguity, refuse the temptation to guess.
There should be one-- and preferably only one --obvious way to do it.
Although that way may not be obvious at first unless you're Dutch.
Now is better than never.
Although never is often better than *right* now.
If the implementation is hard to explain, it's a bad idea.
If the implementation is easy to explain, it may be a good idea.
Namespaces are one honking great idea -- let's do more of those!

Un algorithme est une suite finie et non ambigüe d’opérations ou d’instructions permettant de résoudre un problème. Les algorithmes sont connus depuis l’antiquité (Euclide).

Le mot algorithme vient du nom du mathématicien perse du 9ième siècle (AJC)

  • Abu Abdullah Muhammad ibn Musa al-Khwarizmi (photo).

L’algorithmique correspond à la phase préparatoire avant une quelconque programmation. Elle permet de décrire un problème sous une forme que l’on peut ensuite programmer sur un ordinateur et ceci dans un langage naturel, indépendant d’un langage de programmation.

algorithme numérique suite finie et non ambiguë d’opérations ou d’instructions sur des nombres permettant de résoudre un problème.

Et il n’est pas nécessaire d’avoir un ordinateur pour exécuter un algorithme (machine de Turing inventé en 1936 avant l’ordinateur)!

ordinateur = machine qui exécute (rapidement) des algorithmes

Algorithme itératif

Problème d’Euclide: pavage d’une piece

Soit une pièce rectangulaire de coté a et b (en cm) que l’on souhaites paver avec des petits carreaux carrés hxh.

Quelle est la taille maximum h des carreaux ?

Algorithme

idée : découpage recursif / itératif de la pièce en carré

Algorithme PavageR(A,B)

détermine la taille maxi pour paver une piece de cotés AxB

# version recursive
si A > B alors
   retour PavageR(A-B,A)
sinon si A < B alors
   retour PavageR(A,B-A)
sinon
   retour A

Algorithme Pavage(A,B)

détermine la taille maxi pour paver une piece de cotés AxB

# version itérative
tant que A # B faire
   si A > B alors
       A = A - B
   sinon
       B = B - A
fin tant que
# on a ici A=B donc
retour A
In [3]:
def PavageR(A,B):
    '''version recursive du PGCD'''
    print("PavageR : ",A,B)
    if A > B :
        return PavageR(A-B,B)
    elif A < B :
        return PavageR(B-A,A)
    else :
        return A
def Pavage(A,B):
    '''version itérative du PGCD'''
    while A != B :
        print("Pavage : ",A,B)
        if A > B :
            A = A - B
        else :
            B = B - A
    return A
In [4]:
PavageR(15,21)
PavageR :  15 21
PavageR :  6 15
PavageR :  9 6
PavageR :  3 6
PavageR :  3 3
Out[4]:
3
In [5]:
Pavage(15,21)
Pavage :  15 21
Pavage :  15 6
Pavage :  9 6
Pavage :  3 6
Out[5]:
3

vérification

Propriétés:

  • Pavage(a,a) -> a

cas particulier

  • Pavage(15,21) - > 3

Traduction en Python

In [6]:
def PavageR(a,b):
    #print("PavageR:",a,b)
    if a>b:
        return PavageR(a-b,b)
    elif a<b:
        return PavageR(a,b-a)
    else:
        return a
In [7]:
def Pavage(a,b):
    while a != b :
        #print("Pavage:",a,b)
        if a > b :
            a = a - b
        else :
            b = b - a
    return a

validation

In [8]:
print(Pavage(15,21))
print(PavageR(15,21))
3
3
In [9]:
print(Pavage(15,15))
print(PavageR(15,15))
15
15
In [10]:
n1=np.random.randint(10000)
n2=np.random.randint(10000)
%timeit Pavage(n1,n2)
%timeit PavageR(n1,n2)
2.14 µs ± 1.02 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)
3.81 µs ± 7.86 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)

variables globales et locales

attention: les arguments et les variables d'une fonction sont locales

In [11]:
def ma_fonction(a):
    """a et x locales"""
    a=a-1
    x=a+2
    print("@a dans f=",id(a),a)
    print("@x dans f=",id(x),x)
    return
# programme principal
# a global
a=5
print("@a global =",id(a),a)
ma_fonction(a)
print("@a global=",id(a),a)
@a global = 10919552 5
@a dans f= 10919520 4
@x dans f= 10919584 6
@a global= 10919552 5

Calcul d’une série

Formulation mathématique: calcul de la somme d’une série

$$ S_N = x - \frac{x^2}{2} + \frac{x^3}{3} + .. = \sum_{i=1}^n (-1)^{i+1} x^i/i $$

Solutions algorithmiques

1ere solution

Algorithme Serie1(x,n)
   somme = 0
   pour i de 1 a n
      somme=somme + (-1)^(i+1)*x^i/i
   retour somme

2nde solution (optimisée)

Algorithme Serie(x,n)
   somme = 0
   coeff = x
   pour i de 1 a n
      somme = somme + coef/i
      coef  = -coef*x
   retour somme

Test et Validation

$$ \lim_{N\rightarrow\infty} S_N = log(1+x) $$

Implementation en Python

In [12]:
def Serie(x,n):
    '''calcul des n premiers termes de la série'''
    somme = 0.
    for i in range(1,n+1):
        somme = somme + (-1)**(i+1)*(x**i)/i
    return somme
def SerieO(x,n):
    '''version optimisée'''
    somme = 0.0
    coeff = -1.0
    for i in range(1,n+1):
        coeff = - coeff*x
        somme = somme + coeff/i
    return somme  
In [13]:
x=0.1
print(Serie(x,100),SerieO(x,100),np.log(1+x))
0.09531017980432488 0.09531017980432488 0.0953101798043
In [ ]:
 
In [14]:
# algorithme en Python
from numpy  import log
def serie1(x,n):
    """ calcul de la somme de n termes de la serie x-x^2/2+x^3/3- """
    somme=0.0
    for i in range(1,n+1):
        somme = somme + (-1)**(i+1)*x**i/i
    return somme
In [15]:
# algorithme en Python (version optimisé)
from numpy  import log
def serie(x,n):
    """ calcul de la somme de n termes de la serie x-x^2/2+x^3/3- """
    coef=x
    somme=0.0
    for i in range(1,n+1):
        somme = somme + coef/i
        coef  = -coef*x
    return somme
In [16]:
# utilisation et validation
x=0.2
n=10
print("Calcul de la serie pour x=",x," et n=",n)
print("somme1   = ",serie1(x,n))
print("somme    = ",serie(x,n))
print("log(1+x) = ", log(1+x))
err = serie(x,n) - log(1+x)
print("erreur: ",err)
Calcul de la serie pour x= 0.2  et n= 10
somme1   =  0.18232155522031748
somme    =  0.18232155522031748
log(1+x) =  0.182321556794
erreur:  -1.57363710951e-09

Efficacité: temps d’execution

utilisation de timeit pour calculer le temps d’exécution

ATTENTION l'important est une implémentation correcte et validée d'un lagorithme. La phase d'optimisation n'interviens qu'ensuite si nécéssaire!!!

In [17]:
%timeit serie1(x,10000)
12.6 ms ± 29.4 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)
In [18]:
%timeit serie(x,10000)
2.17 ms ± 8.7 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)

Etude du mouvement du pendule simple

Voir le notebook sur le pendule simple

FIN