[PDF] TP7 : le théor`eme du point fixe en action sous MATLAB

View - Download TP7 : le théor`eme du point fixe en action sous MATLAB

Previous PDF

Next PDF

TP7 : le theoreme du point xe en action sous

MATLAB

Cette seance de TP7 poursuit la familiarisation avecMATLAB. Elle illustre le chapitre 4 du cours (le theoreme du point xe et ses applications en algebre lineaire).

OuvrezMATLABpour commencer.

Exercice1 (l'algorithme iteratif conduisant au calcul approche du rayon spec- tral d'une matrice carree reelle). (1) Declarez sousMATLABla matrice symetrique reelle (de taille (4;4)) sui- vante : >> A A =

0.5172 0.5473 -1.2240 0.8012

0.5473 1.3880 1.3530 -1.1120

-1.2240 1.3530 0.0364 2.8930

0.8012 -1.1120 2.8930 0.0583

En utilisant la commandeeig(A), calculez les quatre valeurs propres reelles de cette matrice reelle symetriqueA. Cette matriceAest-elle denie positive? Est-elle diagonalisable? Quel est le signe de la valeur propre de valeur absolue egale au rayon spectralr(A)de cette matrice? Quelle est la dimension (dansR4) du sous-espace propre correspondant? (2) Telechargez depuis le site la routinerayonspectral. Modiez ensuite cette routine pour realiser une routine [r,Niter] = rayonspectral1(A,X,epsilon,k); qui, etant donne un nombre maximalkd'iterations autorisees, une ma- trice carree (reelle ou complexe) de taille (N;N)A, un vecteur colonneX de longueurN, calcule de maniere approchee, dans les cas favorables, le rayon spectralrde la matriceA, l'algorithme etant initie aX, ainsi que le nombreNiterkd'iterations necessaires avant que l'erreur entre une approximation et la suivante ne vienne a passer en valeur absolue sous le seuilepsilon1. Testez cet algorithme sur la matriceAde la question 1 en prenantk= 100,epsilon=eps, et pourXrespectivement

1. Cela ne donne pasa prioriune majoration parepsilonde l'erreur commise entre le rayon

spectral et sa valeur approchee (au terme deNiteriterations), mais seulement un contr^ole de 1

2TP7 : LE THEOREME DU POINT FIXE EN ACTION SOUSMATLAB

>> X1 = [1;1;1;1]; >> X2 = [1;1;0;0]; >> X3 = [1;0;0;0]; >> X4 = rand(4,1); Quelle valeur deNitertrouvez vous dans chacun de ces cinq cas? Veriez que la valeur derobtenue alors est bien en accord avec le resultat fourni a la question 1 par la commandeeig(A)donnant les quatre valeurs propres (dans ce cas reelles) de la matriceA. Recommencez avec cette foisepsilon= 106. (3) On suppose queAest une matrice reelle et que le rayon spectralr(A)est la valeur absolue d'une valeur propre reelle simple de la matriceA. Veriez que c'est bien le cas pour la matriceAdonnee a la question 1. Modiez la routinerayonspectral1en une routineAppVPdom(≪App roximation du V ecteur P ropre associe a la valeur propre dom inante [VP,Niter] = AppVPdom(A,X,epsilon,k); qui fournisse, avec les m^emes donnees que dansrayonspectral1eninput, { le vecteurVPNiter+1de la suite initiee aX1=X=norm(X)et regie ensuite par l'equation recurrente VP k+1=AVPk norm (AVPk); k1: { le nombre d'iterationsNitereffectue dans la boucle sachant que cette boucle s'arr^ete des que, pour la premiere fois : min (norm(VP_(Niter) -VP_(Niter-1)), norm (VP_(Niter)+ VP_(Niter-1))) <= epsilon Testez cet algorithme avec la matriceAen prenantepsilon=eps, k=200et respectivementX = X1,X2,X3,X4comme a la question 2. Quelle valeur deNiterobtenez vous dans chacun de ces quatre cas? Recommencez avec cette foisepsilon= 106. Calculez les vecteurs propres (normalises) de la matriceAen utilisant : >> [V,D] = eig(A); >> V Comparez le vecteur propreY=V(:,1)ainsi obtenu (correspondant a la valeur propre de valeur absoluer(A)) et le vecteurVPobtenu via >> [VP,Niter] = AppVPdom(A,X,epsilon,k); (4) Modiez la routineAppVPdomconstruite a la question 3 en une nouvelle routine function [VP,Niter,errY] = AppVPdom1(A,X,Y,epsilon,k); qui, en plus des sortiesVPetNiter(comme pour la routineAppVPdom construite a la question 3), fournit aussi la listeerrYdes erreurs succes- sives 2 cette erreur a un facteur multiplicatif pres, ce contr^ole etant de l'ordre deepsilon=(1), ou designe le rapport entrer(A)et le module de la premiere valeur propre de module strictement inferieur ar(A). Voir pour cela la preuve de la Proposition 4.2 du cours.

2. M^eme remarque que precedemment a propos du contr^ole d'erreur : cette toleranceepsilon

contr^ole l'erreur entreVPNiteret un vrai vecteur propre normalise (pour la valeur proprer(A)) enepsilon=(1), oua ete deni dans la note 1 precedente.

TP7 : LE TH

EOREME DU POINT FIXE EN ACTION SOUSMATLAB3

min(norm(VP_k-Y), norm(VP_k+Y)), k=1,2,...,Niter lorsqueYest un vecteur colonne deR4donne eninputen plus des donnees precedentesA,X,epsilon,k. En utilisant cette routine avecepsilon= 10

6,Y=V(:,1)(vecteur propre normalise de la matriceAcorrespondant

a la valeur propre de valeur absoluer(A)) etX=X1,X2,X3,X4, comparez la rapidite de la convergence de la suite (VPk)kversYdans les quatre cas de gure (suivant la valeur du vecteur initialXdepuis lequel est lance l'algorithme). (5) Reprendre la question 4 en remplacant dans la routine l'utilisation de la norme euclidiennenormpar la norme∥ ∥1(norm(.,inf)sousMATLAB). On notera la nouvelle routine (obtenue en modiant a la marge la routine

AppVPdom1

de la question 4)AppVPdom1bis. Exercice2 (un schema simpliste pourPagerank).Cet exercice met en uvre l'approche proposee dans la section 4.1.2 des notes de cours ( ≪maquette≫simpliste du fonctionnement deGoogle). (1) Construire une routine function G=AdjacencePonderee(M); qui, etant donnee une matriceMde taille(N,N)dont les entrees sont constituees de0et de1(consideree comme lamatrice d'adjacenced'un certain graphe oriente(E,V)), calcule lamatrice d'adjacence pondereede ce m^eme graphe oriente, c'est-a-dire la matriceGdeduite de la matriceM en transformant chaque ligne deMainsi : { une ligne constituee entierement de0reste inchangee; { une ligne dans laquelle gure au moins un1est divisee par le nombre de1presents sur cette ligne. (2) Redigez une procedure function [Lequilibre,Niter] =Pagerank(M,L0,kappa,epsilon,k); qui, etant donne un graphe oriente(E,V)aNsommets, materialise par sa matrice d'incidenceM, un nombrekappastrictement entre0et1, calcule de maniere approchee le vecteur ligneLequilibre(de longueurN), unique point xe de l'application strictement contractante :

L = ((1-kappa)/N)*ones(1,N) + kappa*L*G

lorsque : {Gdesigne la matrice d'adjacence ponderee du graphe oriente(E,V) (cf.la question 1); { l'algorithme est initie avec le vecteur ligneL0(dont les entrees sont positives et de somme 1); { le nombre maximal d'iterations autorisees estk; { le nombreNiterkest le nombre d'iterations necessaires (lorsque cela est possible) jusqu'a ce que, pour la premiere fois 3, norm(L_(Niter)-L_(Niter-1)) <=epsilon

3. D'apres l'etude faite en cours,cf.la preuve du Theoreme 4.1 (du point xe), l'erreur entre

L Niteret sa limite est alors majoree parepsilon/(1-kappa).

4TP7 : LE THEOREME DU POINT FIXE EN ACTION SOUSMATLAB

(3) Generez un graphe oriente a10sommetsviala donnee de sa matrice d'adjacenceMet calculez la mesure d'equilibreLequilibrede ce graphe oriente avec le choix dekappa=0.85(le choix classiquement privilegie dans l'algorithmePagerank). Prendreepsilon= 108, puisepsilon=eps, cal- culez aussiNiter(lorsquek=100) et examinez la dependance en le choix du vecteur ligne initialL0: on pourra prendre par exemple pour ce faire comparer les resultats { une liste initiale ≪creuse≫telleL0=[1 0 0 ... 0]; { une liste initiale Affichez les resultats par exemple avecplot(Lequilibre,'d'). Exercice3 (algorithmes iteratifs de Jacobi et de Gau-Seidel).Telechargez depuis le site les deux routinesJacobi.metGaussSeidel.m, respectivement basees sur les decompositions (1) Transformez ces deux routines en des routines : function [XX,Niter] = Jacobi1(M,B,X,epsilon,k); function [XX,Niter] = GaussSeidel1(M,B,X,epsilon,k); qui, etant donnes une matriceMde taille(N,N), sans zeros sur la diagonale, et un vecteur colonneBde longueurN: { renvoient toutes les deux le message condition non remplie si la conditionr(D1E)<1 our(T1F)<1 se trouve en defaut; { generent, si cette condition se trouve remplie, l'algorithme iteratif, soit de Jacobi, soit de Gau-Seidel, initie sur le vecteur colonneX (aussi de longueurN), et visant a calculer de maniere approchee la solutionXXdu systeme de CramerMXX=B; le nombre maximal d'iterations autorisees estk(le m^eme que celui autorise pour calcu- ler en preambule le rayon spectral des matricesD1EouT1F), et l'on decide que l'algorithme iteratif s'arr^ete des que, pour la premiere fois 4, norm(XX_(Niter)-XX_(Niter-1)) <=epsilon (2) Construisez une routine [XX,Niter] = ExempleJacobi(N,B,X,epsilon,k); qui, etant donne un entier strictement positifN, un vecteurBde taille (N,1): { genere la matrice creuseM(N)de taille(N,N)dont la diagonale est constituee de3, la sur-diagonale et la sous-diagonale de-1, toutes les autres entrees de la matrice etant nulles;

4. D'apres l'etude faite en cours,cf.la preuve du Theoreme 4.1 (du point xe), l'erreur entre

XX Niteret sa limite (a savoir la solution du systeme de Cramer que l'on tente d'approcher) est alors majoree parepsilon/(1-r(A)), ouA=D1EouA=T1Fsuivant le cas (Jacobi ou

Gau-Seidel).

TP7 : LE TH

EOREME DU POINT FIXE EN ACTION SOUSMATLAB5

{ resout par la methode de Jacobi (initiee au vecteurXde taille(N,1), aveckiterations autorisees au plus et un seuil d'erreurepsilon comme a la question 1) le systeme de CramerM(N)XX=B, ou la sortieNiterkdesigne toujours le nombre d'iterations reellement effectue lors de l'algorithme de Jacobi. Faites la m^eme chose en remplacant l'algorithme de Jacobi par celui de Gau-Seidel pour construire une fonction similaire : [XX,Niter] = ExempleGaussSeidel(N,B,X,epsilon,k); (l'algorithme de Jacobi ayant cette fois ete remplace par l'algorithme de

Gau-Seidel).

(3) Generez une matrice reelleAde taille(20,20)en utilisant la routine

A = 2*(rand(20)-ones(20)/2);

Generez ensuite la matriceA*A'. L'algorithme de Gauss-Seidel converge- t-il lorsque la matriceMest la matriceM=A*A'? Pourquoi? Veriez le en generant un vecteurB=rand(20,1)et en essayant de resoudre le systeme de Cramer(AA')XX =Bde maniere iterative en utilisant la routine [XX,Niter] = GaussSeidel1(A*A',B,zeros(20,1),10^(-8),k); Comparez le resultat que vous obtenez ainsi avec celui donne par la resolution directe >> (A*A')^(-1)*B Comment faut-il choisirkpour que les deux resultats soient en coherence? Calculezeig(T1F) pour la decomposition de Gau-SeidelAA'=T-Fet expliquez pourquoi il s'averait necessaire de choisir de telles valeurs dek pour appliquer l'algorithme de Gau-Seidel. (4) Soient les trois matrices : 0 @1:75:75 :75 1:75 :75:75 11 A0 @1 22 1 1 1

2 2 11

A0 @21 1 2 2 2 11 21 A Quels algorithmes (de Jacobi ou de Gau-Seidel) sont-ils convergeant pour la resolution iterative des systemesMXX=BlorsqueMest l'une de ces trois matrices (etudiez les trois cas separement)? Exercice4 (le conditionnement des matrices et les facteurs d'amplication d'erreurs relatives dans la resolution des ystemes de Cramer perturbes). (1) En utilisant la routinesvd: >> [U,D,V] =svd(M); donnant ladecomposition en valeurs singulieresd'une matrice reelle ou complexe (cf.la section 4.4.1 du cours), ecrivez une routine function c=ConditionnementNorme2(M);

6TP7 : LE THEOREME DU POINT FIXE EN ACTION SOUSMATLAB

qui donne la valeur du conditionnement d'une matrice carree inversible Mrelativement a la norme∥ ∥2. En utilisant cette routine, calculez le conditionnement de la matrice A:=0 B

B@10 7 8 7

quotesdbs_dbs22.pdfusesText_28

[PDF] Séance de travaux pratiques n° 1

[PDF] simulations, algorithmes en probabilités et statistique(s) au - Apmep

[PDF] Loi de Bernoulli et loi binomiale, cours, première S - MathsFG - Free

[PDF] Probabilités, simulation et algorithmique (pour TI)

[PDF] Algorithmes et programmation en Pascal TD corrigés - Limuniv-mrsfr

[PDF] Notes de cours / Algo et Python

[PDF] Algorithmique et Programmation Projet : algorithme de - DI ENS

[PDF] Score ASIA

[PDF] Un algorithme de simulation pour résoudre un problème de probabilité

[PDF] Algorithme PanaMaths

[PDF] Algorithmique en classe de première avec AlgoBox - Xm1 Math

[PDF] Algorithme U prend la valeur [expression de la suite - Maths en ligne

[PDF] Rappels sur les suites - Algorithme - Lycée d Adultes

[PDF] Les tableaux - Luc Brun

[PDF] Les tableaux 1 Exercice 1 - Lipn