7.8. Correction des exercices¶

%matplotlib inline
import numpy as np
from IPython.display import HTML,display,IFrame,Video,Markdown

7.8.1. Test de colinéarité de 2 vecteurs¶

\[\mathbf{U} \in \mathcal{R}^n , \mathbf{V} \in \mathcal{R}^n \mbox{ si } \left|\mathbf{U}.\mathbf{V}\right| = \left\Vert\mathbf{U}\right\Vert \left\Vert\mathbf{V}\right\Vert \mbox{ alors } \mathbf{U} \mbox{ et } \mathbf{V} \mbox{ sont colinéaires} \]

Ecrire une fonction qui teste la colinéarité

    colineaire(U,V)

Attention à la précision

7.8.1.1. Algorithme: Test de colinéarité de 2 vecteurs¶

attention: calcul avec une précision $\epsilon$

\[ \left|\, \left|\mathbf{U}.\mathbf{V}\right| - \left\Vert\mathbf{U}\right\Vert \left\Vert\mathbf{V}\right\Vert \,\right| < \epsilon\]

Algorithme colineaire(U,V)
    n = dimension U et V
    normU = 0.
    normV = 0.
    scalUV = 0.
    pour i de 0 a n-1
        normU  = normU  + U[i]*U[i]
        normV  = normV  + V[i]*V[i]
        scalUV = scalUV + U[i]*V[i]
    fin pour
    normU = sqrt(normU)
    normV = sqrt(normV)
    test  = abs(scalUV) - normU*normV
    si abs(test)<epsilon
        retour VRAI
    sinon
        retour FAUX

7.8.1.2. Programme Python¶

import numpy as np
eps = 1.e-12
def colineaire(U, V):
    """ 
    test si U et V sont colinéaires 
    """
    n = U.size
    normU = normV = scalUV = 0.0
    for i in range(n):
        normU = normU + U[i]**2
        normV = normV + V[i]**2
        scalUV = scalUV + U[i] * V[i]
    normU = np.sqrt(normU)
    normV = np.sqrt(normV)
    test  = np.abs(scalUV) - normU * normV
    return np.abs(test) < eps
#
U = np.random.rand(3)
V = np.random.rand(3)
print("U =", U, "V =", V)
print(colineaire(U, V))
V = 2. * U
print("U =", U, "V =", V)
print(colineaire(U, V))

U = [0.92591445 0.73594867 0.2124627 ] V = [0.45236574 0.56627863 0.66608434]
False
U = [0.92591445 0.73594867 0.2124627 ] V = [1.8518289  1.47189734 0.42492541]
True

7.8.2. Normalisation d’un vecteur X¶

\[\mathbf{U} \in \mathcal{R}^n \mbox{ calcul } \mathbf{U} \leftarrow \frac{\mathbf{U}}{\left\Vert\mathbf{U} \right\Vert}\]

Ecrire une fonction python qui normalise un vecteur X

Xn = normalisation(X)

attention au cas X=0

7.8.2.1. Algorithme: normalisation d’un vecteur X¶

Algorithme normalisation(U)
    n dimension de U
    normU = 0.
    pour i de 0 a n-1
        normU = normU + U[i]**2
    if normU > eps
        U = U / normU
    retour

7.8.2.2. Programme Python: normalisation d’un vecteur X¶

from numpy import *
eps = 1.e-12
def normalisation(U):
    """ normalisation du vecteur U: resultat dans U """
    n = size(U)
    normU = 0.
    for i in range(n):
        normU = normU + U[i]**2
    normU = sqrt(normU)
    if normU > eps:
        for i in range(n):
            U[i] = U[i] / normU
    return
#
U = random.rand(3)
print("U=",U)
normalisation(U)
print("U normalise=",U)

U= [0.75684132 0.17318733 0.29545636]
U normalise= [0.91106578 0.20847837 0.35566264]

7.8.3. Droite des moindres carrées¶

Soient $n$ couples de points expérimentaux: $(X_i,Y_i)_{i=0,n-1}$, on cherche à lisser ces points par la droite des moindres carrés: $y = a x +b $

7.8.3.1. Principe¶

on détermine la droite $y = a x + b$ qui minimise l’erreur $Err$ entre les points de mesure et les points de la droite:

\[ Err(a,b) = \sum_{i=0}^{n-1} \Delta_i ^2 \]

où $\Delta_i ^2 = (a X_i +b - Y_i)^2$ représente le carré de l’erreur entre un point de mesure $(X_i,Y_i)$ et le point correspondant de la droite $(X_i, a X_i + b)$

Les valeurs de $a$ et $b$ sont déterminées en minimisant l’erreur $Err(a,b)$ par rapport à $a$ et $b$, i.e. tel que $\frac{\partial Err}{\partial a} = 0 $ et $\frac{\partial Err}{\partial b} = 0$

On obtient ainsi le système linéaire de 2 équations:

\[ a \overline{X} + b = \overline{Y} \mbox{ avec } \overline{X}=\frac{1}{n}\sum_{i=0}^{n-1}{X_i} \mbox{ et } \overline{Y}=\frac{1}{n}\sum_{i=0}^{n-1}{Y_i}\]

\[ a \overline{X^2} + b \overline{X} = \overline{XY} \mbox{ avec } \overline{X^2}=\frac{1}{n}\sum_{i=0}^{n-1}{X_i^2} \mbox{ et } \overline{XY}=\frac{1}{n}\sum_{i=0}^{n-1}{X_i Y_i}\]

dont la solution s’écrit:

\[ a = \frac{\overline{XY}- \overline{X}\ \overline{Y}}{\overline{X^2}- (\overline{X})^2} \mbox{ et } b = \overline{Y} - a \overline{X} \]

7.8.3.2. Algorithme: moindre carré¶

initialiser n 
X,Y tableaux de n reels
Xm=0.; XXm=0.
Ym=0.; XYm=0.
pour i de 0 a n-1
    Xm  = Xm  + X[i]
    XXm = XXm + X[i]*X[i]
    Ym  = Ym  + Y[i]
    XYm = XYm + X[i]*Y[i]
fin pour
Xm = Xm / n; XXm = XXm / n 
Ym = Ym / n; XYm = XYm / n
a = (XYm - Xm*Ym)/(XXm - Xm*Xm)
b = Ym - a * Xm

7.8.3.3. Programme Python « moindre carré »¶

X = np.array([1., 2., 3., 4.])
Y = np.array([1.01, 1.99, 3.1, 4.05])
n = X.size
Xm = XXm = 0 
Ym = XYm = 0
for i in range(n):
    Xm = Xm + X[i]
    Ym = Ym + Y[i]
    XXm = XXm + X[i]**2
    XYm = XYm + X[i] * Y[i]
Xm = Xm / n; XXm = XXm / n
Ym = Ym / n; XYm = XYm / n
a = (XYm-Xm*Ym) / (XXm-Xm**2)
b = Ym - a * Xm
print("droite des moindres carres y = ", a, " x + ", b)

droite des moindres carres y =  1.023000000000001  x +  -0.020000000000003126

Visualisation de l’execution du programme « moindre carré »

vous pouvez copier l’exemple de code python sur le site pour l’exécuter

display(Markdown("**Visualisation de l'execution sur le site pythontutor**"))
Video("VIDEO_COURS/pythonlive_exo1.mp4", embed=True,width=700, height=300)

Visualisation de l'execution sur le site pythontutor

Remarque quel résultat obtient-on si on remplace la première ligne par:

X=array([1,2,3,4])

7.8.4. Produit matrice-matrice¶

$\mathbf{C} = \mathbf{A}\ \mathbf{B} $ pour $\mathbf{A} \in \mathcal{R}^n\times\mathcal{R}^p$ , $\mathbf{B} \in \mathcal{R}^p\times\mathcal{R}^q$ et $\mathbf{C} \in \mathcal{R}^n\times\mathcal{R}^q$ $$ C_{i,j} = \sum_{k=0}^{p-1} A_{i,k} B_{k,j}$$

7.8.5. Algorithme: « produit matrice-matrice »¶

initialise n,p,q (dimension)
A matrice dimension n,p
B matrice dimension p,q
C matrice dimension n,q

pour i de 0 a n-1
  pour j de 0 a q-1
    C[i,j]=0.
    pour k de 0 a p-1
        C[i,j]=C[i,j]+A[i,k]*B[k,j]
    fin pour k
  fin pour j
fin pour i
affiche C

7.8.5.1. Programme Python: « produit matrice-matrice »¶

n, p, q = 2, 2, 3 
A = np.reshape(np.arange(n*p), (n, p))
B = np.reshape(np.arange(p*q), (p, q))
C = np.zeros((n, q))
for i in range(n):
    for j in range(q):
        for k in range(p):
            C[i, j] = C[i, j] + A[i, k] * B[k, j]
print("A =",A)
print("B =",B)
print("C = A.B",C)

A = [[0 1]
 [2 3]]
B = [[0 1 2]
 [3 4 5]]
C = A.B [[ 3.  4.  5.]
 [ 9. 14. 19.]]

7.8.6. Calcul matrice de Vandermonde¶

\[\begin{split} V = \left[ \begin{array}{ccccc} 1 & x_0 & x_0^2 & \ldots & x_0^{n} \\ 1 & x_1 & x_1^2 & \ldots & x_1^{n}\\ \ldots & \ldots & \ldots & \ldots & \ldots \\ 1 & x_n & x_n^2 & \ldots & x_n^{n}\\ \end{array} \right] \end{split}\]

Intervient dans le problème de la détermination du polynôme d’interpolation de degré $n$ passant par $n+1$ points $\{x_i,y_i\}_{i=0,n}$

7.8.6.1. Algorithme: vandermonde¶

initialiser n
X vecteur dimension n+1
V matrice dimension (n+1,n+1)
pour i de 0 a n
    V[i,0]=1.0
    pour j de 1 a n
        V[i,j] = V[i,j]*X[i]
    fin pour j
fin pour i
print V

7.8.6.2. Programme Python: vandermonde¶

n = 4
X = np.arange(1, n+2)
print(X)
V = np.zeros((n+1, n+1))
for i in range(n+1):
    V[i, 0] = 1.
    for j in range(1, n+1):
        V[i, j] = V[i, j-1] * X[i]
print(V)

[1 2 3 4 5]
[[  1.   1.   1.   1.   1.]
 [  1.   2.   4.   8.  16.]
 [  1.   3.   9.  27.  81.]
 [  1.   4.  16.  64. 256.]
 [  1.   5.  25. 125. 625.]]

7.8.7. FIN des exercices¶

INPROS INtroduction à la PROgrammation Scientifique en Python

Correction des exercices

Contenu