TP 1- Conditionnement et Décomposition LU

Analyse Numérique Matricielle

Polytech SU, 2023-2024

import numpy as np
from scipy.linalg import hilbert
import matplotlib.pyplot as plt

Exercice 1

1. Construire la matrice $ H^{[5]} $ en utilisant la bibliothèque scipy et la fonction scipy.linalg.hilbert ou avec des boucles for.

On considère la matrice de Hilbert $ H^{[n]} \in M_n(\mathbb{R}) $ définie pour tout $n \geq 1$ par $$ H^{[n]}_{ij}=\frac{1}{i+j-1},\, 1\leq i,j\leq n. $$

1. Construire la matrice $ H^{[5]} $ en utilisant la bibliothèque scipy et la fonction scipy.linalg.hilbert ou avec des boucles for.

# 1. Construire la matrice H^{[5]}

1. Calculer $ Cond(H^{[n]}) $ pour la norme $ \Vert\cdot\Vert_2 $ et pour$ 1\leq n\leq 10 $ à l'aide de la fonction numpy.linalg.cond. Visualisez le résultat en échelle semilogarithmique en $ y $ (à l'aide de la fonction matplotlib.pyplot.semilogy). La matrice $ H^{[n]} $ est-elle bien ou mal conditionnée selon vous ?

# 2. Calculer Cond(H^{[n]}) pour la norme \Vert\cdot\Vert_2 et visualiser

1. Répéter le point précédent avec $ Id+H^{[n]} $ au lieu de $ H^{[n]} $. Que remarquez vous ?

# 3. Répéter avec Id+H^{[n]}

Exercice 2

On considère la matrice de Vandermonde $$A=\begin{pmatrix} y_1^{n-1}&\dots&y_1^2&y_1&1\\ \vdots&\vdots&y_2^2&y_2&1\\ \vdots&\vdots&\vdots&\vdots&\vdots\\ y_n^{n-1}&\dots&y_n^2&y_n&1 \end{pmatrix},$$ où les points $y_1,\dots,y_n$ sont équirépartis sur $[0,1]$. Cette matrice peut être construite en Python à l'aide des fonctions suivantes :

y = np.linspace(0, 1, n)

A = np.vander(y, increasing=True)

On souhaite résoudre le système linéaire $Ax=b$ où $$b=\begin{pmatrix} 1+y_1^2\\\vdots\\1+y_n^2 \end{pmatrix}.$$ La solution exacte est connue et vaut $x=(0,\ldots,0,1,0,1)^T$.

1. Résoudre le système linéaire en utilisant la fonction numpy.linalg.solve. On note $x_c$ la solution obtenue. Evaluer l'erreur relative $\epsilon_{n} = \frac{\| x_{c} - x\|}{\|x\|}$ pour $n=4$.

### 1. Résoudre le système linéaire pour \( n = 4 \) et évaluer l'erreur relative :

import numpy as np

1. La variable d'arrondi en Python, notée $\mathbf{eps}$, peut être obtenue via numpy.finfo(float).eps. Une borne supérieure de l'erreur relative $\epsilon_{n}$ est donnée par $$ \eta_{n} = \text{Cond}(A) \mathbf{eps} . $$ Comparer l'erreur relative $\epsilon_{n}$ calculée précédemment avec la borne obtenue $\eta_{n}$ pour $n=4$.

### 2. Comparer l'erreur relative avec la borne pour \( n = 4 \):

1. Répéter les points $1.$ et $2.$ pour $n=4,6,8,\ldots,20$. Visualiser l'erreur $\epsilon_n$, la borne supérieure $\eta_n$ ainsi que le résidu normalisé $r_n = \frac{\| b - Ax_{c} \|}{\|b\|}$ en fonction de $n$. Pour cela, tracer deux graphes, un avec une échelle logarithmique à l'aide de la fonction matplotlib.pyplot.loglog et l'autre en échelle semilogarithmique à l'aide de la fonction matplotlib.pyplot.semilogy. Quel type de croissance de l'erreur $\epsilon_n$ observe-t-on ? Est-ce que le résidu $r_n$ est un bon indicateur de l'erreur $\epsilon_n$ ? Et $\eta_n$ ?

### 3. Répéter pour \( n=4,6,8,\ldots,20 \) et visualiser les erreurs:

1. (a) Construire la matrice $$A=\begin{pmatrix} 2&-1&&&\\-1&2&\ddots&&\\&\ddots&\ddots&\ddots&\\&&\ddots&\ddots&-1\\ &&&-1&2 \end{pmatrix}.$$ pour $n\geq 1$ en utilisant la bibliothèque numpy.

(b) Répéter le point $3.$ pour cette matrice avec $b=(2,2,\ldots,2)^T$ pour $n=5,10,\ldots,100$. Commenter les résultats obtenus.

### 4. Matrice et visualisation pour une nouvelle matrice

#(a) Construction de la matrice :

# (b) Répéter le point 3 :