MO102.
Introduction à Matlab 1
MO102. Manipulations matricielles avancées
Patrick Ciarlet
Pour s’échauffer...
Exercice 1 Soient les vecteurs colonnes et la matrice suivants
1 −5 −1 2 3 4
~u1 = 2 , ~u2 = 2 , ~u3 = −3 , A0 = 7 6 5 .
3 1 7 2 8 7
1. Entrer les données en ligne.
2. Calculer le produit A0~u1 .
3. Déterminer les commandes Matlab permettant de :
(a) calculer k~u1k2 , k~u2 k1 , k~u3 k∞ ;
(b) déterminer les dimensions de la matrice A0 , en extraire le nombre de colonnes ;
(c) calculer le déterminant et l’inverse de A0 .
4. Proposer deux méthodes permettant de résoudre le problème A0~x = ~u1 , et déterminer
les commandes Matlab associées. Quelle méthode doit-on retenir ?
Pour étudier un exemple plus pratique, on se propose de discrétiser le modèle du fil
pesant : trouver u tel que
−u′′ = f sur ]0, 1[, u(0) = u(1) = 0.
Si f est de classe C 2 ([0, 1]) (ou, ce qui est équivalent, si u est de classe C 4 ([0, 1])), on peut
vérifier que pour x ∈]0, 1[ et h tel que [x − h, x + h] ∈ [0, 1], on peut écrire
!
−u(x − h) + 2u(x) − u(x + h)
−u′′ (x) = +r(x, h), avec |r(x, h)| ≤ sup |f ′′ (x)|/12 h2 .
h2 x∈[0,1]
Ce résultat simple fournit une méthode de discrétisation et d’approximation de l’équation
du fil pesant. On parle de schéma aux différences finies. Comment procède-t-on ?
1
Pour commencer, on choisit N ∈ N, et on fixe h = . Remarquons tout de suite que
N +1
′′
pour avoir une ”bonne” approximation de u (x), il convient que h soit petit. Ceci signifie
que N est un paramètre de discrétisation qui aura vocation à devenir ”grand”, lors de la
réalisation des expériences numériques.
Comment approcher la valeur de u aux points xi = ih, pour i ∈ {0, 1, · · · , N + 1}, par des
nombres, notés (ui )0≤i≤N +1 ?
Comme on sait que u(0) = u(1) = 0, on choisira comme approximation u0 = uN +1 = 0 !
Ensuite, on définit fi = f (xi ), pour i ∈ {1, · · · , N}, et on considère les équations
−ui−1 + 2ui − ui+1
= fi , 1 ≤ i ≤ N, avec u0 = uN +1 = 0.
h2
MO102. Introduction à Matlab 2
Si on appelle ~u (resp. f~) le vecteur de RN de composantes (ui )1≤i≤N (resp. (fi )1≤i≤N ), on
peut réécrire l’ensemble des équations sous la forme vectorielle équivalente
2 −1 0 . . . 0
. .
−1 2 −1 . . ..
1 .. .. ..
A1 ~u = f~, avec A1 = 2
N ×N
h .
0 . . . 0 ∈R .
.. . . . . . . . . . −1
0 . . . 0 −1 2
Exercice 2 1. Ecrire une fonction qui construise la matrice A1 d’ordre N obtenue
par la discrétisation par différences finies du modèle du fil pesant.
2. Comparer les temps de construction ainsi que la faisabilité, selon que la structure
associée à A1 est pleine ou creuse, en fonction de N.
3. Résoudre le système linéaire A1~u = f~, lorsque f ≡ 1.
MO102. Introduction à Matlab 3
Et plus si affinités...
Lorsque l’on doit résoudre plusieurs fois des systèmes linéaires avec la même matrice,
la question se pose de savoir si l’utilisation de la commande \ est suffisamment rapide...
Si tel n’est pas le cas, il est intéressant de réaliser la factorisation de la matrice en ques-
tion a priori, puis d’utiliser cette factorisation pour résoudre très rapidement les systèmes
linéaires. Qu’en est-il pour la matrice A1 ?
Commençons par noter que la matrice A1 est symétrique définie-positive :
— Une matrice A de RN × RN est symétrique lorsque Ai,j = Aj,i , pour 1 ≤ i, j ≤ N ;
— Une matrice A de RN × RN est définie-positive lorsque (Ax, x)RN > 0, pour tout
x non-nul.
Par voie de conséquence, on peut réaliser la factorisation de Cholesky de la matrice A,
c’est-à-dire construire une matrice triangulaire inférieure L (Li,j = 0 si j > i) telle que
A = LLT . En remplaçant A par sa factorisation, on peut résoudre le système linéaire
A~u = f~
en deux temps :
— une descente : résoudre L~y = f~ ;
— une remontée : résoudre LT~u = ~y .
Ceci est très avantageux, car chacun des deux systèmes fait intervenir une matrice trian-
gulaire, pour laquelle la résolution du système linéaire associé est élémentaire...
NB. Dans le cas général, on utilisera la factorisation de Gauss (fonction lu) ; on
pourra également considérer des méthodes itératives de résolution (rubrique iterative
methods).
Exercice 3 1. Ecrire une fonction Matlab pour calculer la solution de A~u = f~ à
l’aide la méthode de Cholesky, en utilisant la fonction chol.
2. Comment procéder si on doit résoudre p systèmes linéaires, de seconds membres
f1 , · · · , fp ?
Outre la résolution de systèmes linéaires, l’autre grande catégorie de problèmes liés
aux matrices est le calcul de ses valeurs propres (et de ses vecteurs propres.) Précisons.
Soit A une matrice de CN ×N : λ ∈ C est une valeur propre de A si, et seulement si, l’une
des deux assertions ci-dessous est vérifı́ee :
— det(A − λ IN ) = 0 ;
— il existe un vecteur x non nul de RN satisfaisant à Ax = λx.
NB. Il est souvent intéressant de calculer les modes propres d’un système physique. Une
fois discrétisé (et après linéarisation éventuelle), on aboutit typiquement au problème
mentionné ci-dessus de la détermination des valeurs propres.
Exercice 4 1. En utilisant la fonction eig, calculer les valeurs propres de A0 .
2. En utilisant la fonction eig, peut-on calculer les valeurs propres de A1 , pour N =
10.000 ?
3. En utilisant la fonction eigs, calculer les 5 plus petites valeurs propres de A1 , pour
N = 100, 300, 1.000, 3.000, etc. Conclusion ?