Remerciez-le!

Remerciez @Admin pour avoir partagé cet document gratuitement, de la manière la plus simple, en partageant sur les réseaux sociaux.

Ecole Nationale Des Sciences Appliquées de Tétouan- ENSA 2AP-1 – 2015/2016 Init

Ecole Nationale Des Sciences Appliquées de Tétouan- ENSA 2AP-1 – 2015/2016 Initiation à MATLAB TP4 Exercice 1 Traduire en langage Matlab l’algorithme PGCD suivant : a,b entiers positifs tant-que a b faire si a b > alors a a b - sinon b b a - fin tant-que PGCD a Corrigé de l’exercice 1 Pour traduire en langage Matlab l’algorithme PGCD considéré dans cet exercice, nous avons besoin des structures de contrôle, en l’occurrence la boucle while et le test if-else-end. Nous proposons le M-file suivant pour le programme traduisant cet algorithme: clear all; close all; clc; a=input('Entrer l''entier positif a '); b=input('Entrer l''entier positif b '); a0=a; b0=b; while a~=b if a>b a=a-b; else b=b-a; 1 Ecole Nationale Des Sciences Appliquées de Tétouan- ENSA 2AP-1 – 2015/2016 end end str=['Le plus petit diviseur commun de ' num2str(a0) ' et ' num2str(b0) ' est ' num2str(a)]; disp(str); Remarques: *Le programme M-file constituant la réponse à cet exercice n’est pas unique. Chacun a sa façon personnelle de programmer. Le programme proposé est donc un exemple seulement de solution. Il existe cependant des critères à satisfaire: le programme doit être sans redondances, facile à lire et optimisé. *On a fait appel aux variables auxiliaires a0 et b0 pour stocker les valeurs initiales de a et b qui sont modifiées dans la boucle while et ne gardent pas donc leurs valeurs d’entrées à la sortie de cette boucle. Elles sont en fait écrasées dans les lignes 10 et 12. *On peut tester le programme en faisant appel à la commande gcd de Matlab qui calcule le plus grand diviseur commun de deux nombre entiers positifs. Pour cela, il suffit de taper dans notre cas (a0=18 ; b0=5). >> gcd(a0,b0) ans = 1 * Il faut faire la distinction entre la commande disp qui affiche la valeur de la variable sans écrire le nom de la variable et la commande display qui écrit le nom de la variable avant de l’afficher. Taper pour cela les deux commandes suivantes: >> display(str) str = Le plus petit diviseur commun de 18 et 25 est 1 >> disp(str) Le plus petit diviseur commun de 18 et 25 est 1 2 Ecole Nationale Des Sciences Appliquées de Tétouan- ENSA 2AP-1 – 2015/2016 Exercice 2 Ecrire la fonction vec2col.m qui transforme tout vecteur (ligne ou colonne) passé en argument en vecteur colonne. Un traitement d’erreur testera que la variable passée à la fonction est bien un vecteur et non une matrice (utiliser la commande size). Corrigé de l’exercice 2 On propose le M-file de type fonction suivant : function [y]=vect2col(x) n=size(x); n1=n(1); n2=n(2); if n1==1 y=x.'; elseif n2==1 y=x; else erreur=['la variable entrée n''est pas un vecteur ligne ou colonne']; display(erreur) end On peut le tester avec les commandes >> vect2col([1 1 1]) ans = 1 1 1 >> vect2col([1; 1; 1]) ans = 3 Ecole Nationale Des Sciences Appliquées de Tétouan- ENSA 2AP-1 – 2015/2016 1 1 1 >> vect2col([1 1; 1 2]) erreur = la variable entrée n'est pas un vecteur ligne ou colonne Exercice 3 Ecrire un script qui permet selon le choix de transformer les cordonnées cartésiennes (x,y,z) en coordonnés cylindriques (r, ,z) q ou bien en coordonnées sphériques ( , , ) r q j . Corrigé de l’exercice 3 Un exemple de script répondant à l’énoncé est donné dans la suite. Remarquer que l’on sort les coordonnées cylindriques selon l’ordre : ( ) r, ,z q et les coordonnées sphériques selon l’ordre ( ) r, , q j . function [a,b,c]=transformation(x,y,z,choix) switch choix case 'cylindrique' a=sqrt(x^2+y^2); b=atan(y/x); c=z; case 'sphérique' a=sqrt(x^2+y^2+z^2); b=atan(y/x); c=atan(sqrt(x^2+y^2)/z); otherwise display('Vous n''avez pas préciser la nature de la tansformation') end 4 Ecole Nationale Des Sciences Appliquées de Tétouan- ENSA 2AP-1 – 2015/2016 end Exemples d’utilisation : >> [a,b,c]=transformation(1,1,1,'cylindrique') a = 1.4142 b = 0.7854 c = 1 % >> [a,b,c]=transformation(2,3,4,'sphérique') a = 5.3852 b = 0.9828 c = 0.7336 Avec Matlab vous avez quatre commandes qui vous permettent de faire la transformation des coordonnées cartésiennes en coordonnées cylindriques ou bien des coordonnées cartésiennes en coordonnées sphériques, ainsi que leurs inverses. [theta,rho,z]=cart2pol(x,y,z) % transforamtion catésiennes en cylindriques [x,y,z]=pol2cart(theta,rho,z) % transforamtion cylindriques en cartésiennes [theta,phi,r]=cart2pol(x,y,z) % transforamtion catésiennes en sphériques [x,y,z]=sph2cart(theta,phi,r) % transforamtion sphériques en catésiennes 5 Ecole Nationale Des Sciences Appliquées de Tétouan- ENSA 2AP-1 – 2015/2016 Attention: * Ici tous les angles sont calculés en radians (Matlab utilise cette unité par défaut). * Matlab sort les coordonnées cylindriques selon l’ordre ( ) ,r,z q et les coordonnées sphériques selon l’ordre ( ) , ,r q j . * Dans Matlab l’angle phi est compté à partir du plan horizontal et non pas à partir de l’axe des z comme le montre la figure suivante : Définition des coordonnées sphériques selon Matlab : phi est compté à partir du plan: phi=90- Exercice 4 On considère l’algorithme suivant pour calculer π: on génère n couples ( ) { } k k x ,y de nombres aléatoires dans l’intervalle [0, 1], puis on calcule le nombre m de ceux qui se trouvent dans le premier quart du cercle unité, c'est-à-dire vérifiant la condition: 2 2 k k x y 1 + . Naturellement, π est la limite de la suite n 4m u n = . Ecrire un programme Matlab pour calculer cette suite et observer comment évolue l’erreur quand n augmente. Corrigé de l’exercice 4 On peut proposer le script suivant : 6  Notre choix Ecole Nationale Des Sciences Appliquées de Tétouan- ENSA 2AP-1 – 2015/2016 clear all; close all; clc; format long; N=input('Entrer N? ') x=rand(N,1); y=rand(N,1); z=x.^2+y.^2; m=sum(z<=1); uN=4*m/N; display(uN); erreur=(pi-uN)/pi; display(erreur); A l’exécution, on obtient : Entrer N? 2.e7 N = 20000000 uN = 3.141851800000000 erreur = -8.248886433788238e-005 La convergence est très lente avec le processus de Monte Carlo. Sur ma machine, on ne peut pas dépasser N 2.e7 = car la mémoire se sature. Ceux qui préfèrent les boucles et les tests peuvent être intéressés par la variante suivante: 7 Ecole Nationale Des Sciences Appliquées de Tétouan- ENSA 2AP-1 – 2015/2016 clear all; close all; clc; format long; N=input('Entrer N? ') m=0; tic; for k=1:N xk=rand(1); yk=rand(1); if xk^2+yk^2<=1 m=m+1; else end end uN=4*m/N; display(uN); erreur=(pi-uN)/pi; display(erreur); toc; En exécutant ce script, on obtient : Entrer N? 2.e7 N = 20000000 uN = 3.141854000000000 erreur = 8 Ecole Nationale Des Sciences Appliquées de Tétouan- ENSA 2AP-1 – 2015/2016 -8.318914608747151e-005 Elapsed time is 106.232164 seconds. On a placé tic et toc pour compter la durée d’exécution, c'est-à-dire celle que l’on mesurerait en consultant une horloge étalonnée déclenchée avant de lancer le programme et arrêtée à la fin de l’exécution du programme (A ne pas confondre avec le temps CPU que l’on peut obtenir avec la commande cputime). On a attendu dans ce cas à peu près 106 s. Avec la première version, on obtient Entrer N? 2.e7 N = 20000000 uN = 3.141189800000000 erreur = 1.282322803158388e-004 Elapsed time is 3.452282 seconds. La durée d’exécution est meilleure avec la première version du programme. Elle est seulement de 3.3% de celle qui est exigée avec la deuxième version. Cependant la deuxième version ne souffre pas de limitation au niveau de l’occupation de la mémoire car on peut dépasser facilement un nombre de tirages supérieur à N 2.e7 = , du fait que l’on ne stocke pas en mémoire les variables: on les écrase. Un programme n’est jamais unique, mais il existe des critères qui font qu’une version de script est meilleure que l’autre: la clarté, la vitesse d’exécution et l’occupation de la mémoire sont des éléments importants. L’étude de la convergence est résumée dans le tableau suivant: N 1.e3 1.e4 1.e5 1.e6 1.e7 2.e7 Erreur 1.4513e-2 4.3267e-3 1.0973e-3 2.5612e-4 -7.0494e-6 1.5446e-4 Il faut remarquer que ces erreurs varient d’un tirage à l’autre et que le tableau représente le résultat d’un tirage particulier. Par ailleurs la convergence n’est pas monotone comme le montre les deux dernières colonnes. 9 Ecole Nationale Des Sciences Appliquées de Tétouan- ENSA 2AP-1 – 2015/2016 Exercice 5 Ecrire un script qui permet de calculer par la méthode de dichotomie la racine de l’équation nonlinéaire suivante: 5sin(x) exp(x/2) 0 - = sur l’intervalle [0,1] . Le script devra permettre de calculer par l’algorithme de dichotomie une approximation à h près de la solution de l’équation f(x) = 0 dans l’intervalle [a,b], a, b et h étant saisis par l’utilisateur et la fonction f étant définie par la commande inline. Principe de la méthode de dichotomie : Soit f une fonction continue sur [a, b] telle que f(a)f(b) < 0. Nécessairement f a au moins un zéro dans ]a, b[ (ce résultat est un cas particulier du théorème des valeurs intermédiaires). Supposons pour simplifier qu’il est unique et notons le α (dans le cas où il y a plusieurs zéros, on peut localiser graphiquement, à l’aide de la commande plot, un intervalle qui n’en contient qu’un). La méthode de dichotomie consiste à diviser en deux un intervalle donné, et à choisir le sous-intervalle où uploads/Industriel/ corrige-tp4-matlab.pdf