Implémentation et utilisation de la méthode des rectangles à gauche¶
À partir d’une fonction connue¶
import numpy as np
import matplotlib.pyplot as plt
%matplotlib inline
def rectangleGauche_f(f, a, b, m):
"""
Renvoie la valeur approchée de l'intégrale d'une fonction f
entre les bornes a et b par la méthode des rectangle à gauche
avec un intervalle découpé en m rectangles
"""
h = (b-a) / float(m)
s = 0
for i in range (m):
s += f(a + i*h)
return h * s
rectangleGauche_f(lambda x : np.sin(x), 0, np.pi, 8)
rectangleGauche_f(lambda x : np.sin(x), 0, np.pi, 16)
À partir d’un ensemble de valeurs discrètes¶
data = np.array([[0,1.5],[0.5,2],[1,2],[1.5,1.6364],[2,1.25],[2.5,0.9565]])
def rectangleGauche_l(pts):
"""
Pour des valeurs rangées en colonne
"""
h = pts[1,0] - pts[0,0] # Le pas
s = sum(pts[:-1,1]) #la somme avec la fonction sum()
return h * s
rectangleGauche_l(data)
plt.plot(data[:,0], data[:,1], "-or")
plt.grid()
plt.show()
Implémentation et utilisation de la méthode des trapèzes¶
def trapeze_f(f, a, b, m):
'''
Renvoie la valeur approchée de l'intégrale d'une fonction f
entre les bornes a et b par la méthode des trapèze avec
un intervalle découpé en m trapèzes
'''
h = (b-a) / float(m)
s = 0.5 * (f(a) + f(b))
for i in range (1,m):
s += f(a + i*h)
return h * s
trapeze_f(lambda x : np.sin(x), 0, np.pi, 8)
trapeze_f(lambda x : np.sin(x), 0, np.pi, 16)
Erreur avec m=8 : $\quad E = \dfrac{\pi^3}{768}\times sin(\xi)\quad \xi \in ]0, \pi[$
$\xi$ n'est pas connu mais on peut minorer et majorer cette erreur
$$E_{min} = \dfrac{\pi^3}{768}\times \sin(0) \simeq 0 \quad\quad E_{max} = \dfrac{\pi^3}{768}\times \sin(\frac{\pi}{2}) \simeq 0,04037$$
Erreur avec m=16 : $\quad E = \dfrac{\pi^3}{768}\times sin(\xi)\quad \xi \in ]0, \pi[$
$$E_{min} = \dfrac{\pi^3}{3072}\times \sin(0) \simeq 0 \quad\quad E_{max} = \dfrac{\pi^3}{3072}\times \sin(\frac{\pi}{2}) \simeq 0,01009$$
Pour obtenir une précision de $10^{-6}$ on a $$E\leq \dfrac{\pi^3}{12\times I^2}\times sin(\frac{\pi}{2}) \text{ donc } I^2 \geq \dfrac{\pi^3}{12\times 10^{-6}} \simeq 2583856.39 \text{ ssi } I \geq 1607$$
trapeze_f(lambda x : np.sin(x), 0, np.pi, 1607)
Valeur moyenne d'un signal¶
Soit un signal $s(x)$ (périodique ou non). Sa valeur moyenne (notée $\bar{s}$) sur un intervalle $[a ; b ]\in \mathbb{R}$ est définie par:
$$\bar{s}(a, b) = \dfrac{1}{b - a}\int_{a} ^{b} s(x)dx \quad\quad (1)$$
En physique il n'est pas rare que le signal soit défini par une acquisition en temps discret qui renvoie un nuage de points. On adoptera donc une représentation en temps discret du signal. On appelle horodate, un ensemble (fini) des mesures ($t_i, f(t_i)$) réalisées sur une période et avec un échantillonnage donnés.
La moyenne de ce signal peut être calculé avec la relation $(1)$ en adoptant l'une des méthodes (trapèze, rectangle à gauche ou à droite) explicitées dans le TD 12.
Exemple : méthode des rectangles à gauche
Dans le cas d'un signal en temps discret m(t) possédant $n$ points on note $(t_i, m(t_i))$ le i-ième point de la liste et $H$ le pas d'acquisition constant. La valeur moyenne du signal est donnée par:
$$\bar{m}(t_1, t_n) = \dfrac{H}{t_{n}-t_1}\left(\sum\limits_{k=1}^{n-1}f(t_k)\right)\quad\quad f(t_k) \text{ étant une valeur de la grandeur mesurée}$$
moyenne_rectangle = 1/(data[-1,0]-data[0,0])*rectangleGauche_l(data)
print(moyenne_rectangle)
On remarquera dans le cas de la méthode des rectangles à gauche que la valeur moyenne d'un signal en temps discret sur n points est égale à la moyenne arithmétique des $(n-1)$ premiers points
$$\bar{m}(t_1, t_n) = \dfrac{H}{t_{n}-t_1}\left(\sum\limits_{k=1}^{n-1}f(t_k)\right) = \dfrac{1}{n-1}\left(\sum\limits_{k=1}^{n-1}f(t_k)\right)$$
def moyenneArithmetique(data):
fin = len(data)
return sum(data)/fin
print("Les valeurs pour la moyenne", data[0:-1, 1])
print("Moyenne : ",moyenneArithmetique(data[0:-1, 1]))
def trapeze_l(pts):
'''
Renvoie la valeur approchée de l'intégrale d'une fonction f
entre les bornes a et b par la méthode des trapèze à partir
d'une série de données expérimentales pts de type <ndarray>
'''
fa, fb = pts[0,1], pts[-1,1] # Pour des valeurs rangées en colonne
h = pts[1,0] - pts[0,0] # Le pas
s = 0.5 * (fa + fb)
s += sum(pts[1:-1,1]) #la somme avec la fonction sum()
return h * s
trapeze_l(data)
En dehors de la méthode vue en cours, la lecture d'un fichier CSV peut également se faire avec le module csv ou avec le module numpy
Avec le module csv¶
Le plus simple est d'écrire une fonction permettant de lire les données expérimentales dans le fichier data.csv et de les placer dans un tableau.
import csv
def readCSV(filename, delim):
'''
Cette fonction renvoie les données d'un fichier au format csv dans un tableau de type list
filename : nom du fichier au format csv
delim : séparateur de champs dans le fichier csv
'''
f = open(filename, "r")
tab = []
lines = csv.reader(f, delimiter=delim)
for l in lines:
if l[0][0] != '#': #pour ignorer les lignes de commentaires dans le fichier csv
tab.append(l)
f.close()
return tab
Utilisation de la fonction readCSV
Changement du répertoire de lecture.
import os
os.chdir("/chemin/absolu/du/fichier")
readCSV("data.csv", ";")
Attention le séparateur décimal des nombres doit être un point et non une virgule. Ce symbole dépend des conventions régionales du système de numération ; communément, il est représenté par un point dans les systèmes anglo-saxons et par une virgule dans les autres systèmes.
pts = np.array(readCSV("data.csv", ";"), dtype='float')
print pts
# print "A partir d'une série de valeurs expérimentales"
print
print "Im=", trapeze_l(pts)
Lecture du fichier directement avec numpy¶
Attention pour l'argument de séparateur de champs il est nécessaire d'indiquer le paramètre delimiter = ";"
En effet la documentation indique que par défaut le deuxième paramètre est le type des données du tableau. Il est impératif de spécifier le nom d'un champs qui n'est pas celui par défaut dans la liste des paramètres.
print help(np.loadtxt)
pts = np.loadtxt("data.csv", delimiter = ";");
print "Im=", trapeze_l(pts)
La spirale de Cornu¶
def cornu(t):
x = lambda u: np.cos(u**2)
y = lambda u: np.sin(u**2)
return [trapeze_f(x, 0, t, 1000), trapeze_f(y, 0, t, 1000)]
import matplotlib.pyplot as plt
import time
%matplotlib inline
start = time.time() #étude d'une complexité temporelle
liste = np.array([cornu(t) for t in np.arange(-10,10,0.01)])
plt.plot(liste[:,0], liste[:,1])
plt.show()
print(f"Temps d'exécution pour obtenir la spirale de cornu : {time.time()-start:.3f} secondes")
def trapeze_v(f, a, b, m):
'''
Renvoie la valeur approchée de l'intégrale d'une fonction f
entre les bornes a et b par la méthode des trapèzes vectorisée
'''
s = 0.5 * (f(a) + f(b))
i = np.linspace(a, b, m + 1)
h = i[1] - i[0]
s = s + np.sum(f(i[1:-1]))
return h * s
def cornu(t):
x = lambda u: np.cos(u**2)
y = lambda u: np.sin(u**2)
return [trapeze_v(x, 0, t, 1000), trapeze_v(y, 0, t, 1000)]
start = time.time() #étude d'une complexité temporelle
liste = np.array([cornu(t) for t in np.arange(-10,10,0.01)])
plt.plot(liste[:,0], liste[:,1])
plt.show()
print(f"Temps d'exécution pour obtenir la spirale de cornu : {time.time()-start:.3f} secondes")
La complexité en temps est bien meilleur...
Une méthode basée sur une technique probabiliste : Monte-Carlo¶
def monteCarlo(f, a, b, m):
h = (b-a) / float(m)
x = (b - a) * np.random.random_sample((m,)) + a
return h*sum(f(x))
monteCarlo(lambda x: np.sin(x), 0, math.pi, 10000)