Transcript TP2
L2-mathsCours Informatique appliquée aux mathématiques
TP2 A rendre le 15 octobre 2014 avec le TP3
Trouver les zéros d’une fonction
1
Comment déclarer une fonction en matlab
On veut déclarer la fonction x2 − 1.
• Definir f= inline(’x^2-1’). Calculer la valeur de f en x = 3 en utilisant
la fonction feval de matlab.
• Deuxième possibilité, que produit la commande >> h= @(x)(x^2-1) ? Calculer également la valeur de h en x = 3.
• Créer dans un fichier nommé g.m la ’fonction’
function y=g (x)
y=x^2-1.
Que se passe t-il si l’on utilise la fonction feval de matlab telle quelle ? Comment transformer le script g en une fonction au moyen du signe @ ? Que
produit feval(@g,3) ?
• Dans chacun des cas précédents, on veut prendre la valeur de la fonction en
un certains nombres de points donnés dans le vecteur x = −2 : 0.1 : 2. Ecrire
le résultat y = f (x) au moyen d’une boucle.
• Dans le premier cas, que donne la commande >>y=f(x) ? Comment modifier
la fonction pour supprimer le message d’erreur ? Même question dans les 2
autres cas.
• Tracer dans la figure (1) la fonction x2 − 1 sur l’intervalle (−2, 2).
2
Utiliser la fonction fzero de matlab
• On calcule la racine carrée positive de a = 1. Lire le help de matlab pour
utiliser la fonction fzero. Que signifie fzero(f,0) et fzero(f,[0 2 ]) ? Que
produit chacune de ces commandes ? Même question pour g et h.
• On calcule maintenant la racine carrée de a, pour un certain nombre de valeurs
de a rangées dans le vecteur Aa=1:0.1:10. Modifier le script donnant g en gp
qui prend en entrée x et a. Ecrire une boucle qui donne un zero z(j) de gp
pour chaque valeur Aa(j) de a. Tracer dans la figure (2) les valeurs de z en
fonction de a.
• Chargez la figure fig.png avec la commande imread dans un tableau A, puis
affichez la avec la commande image dans la figure (3). Avec la fonction textttginput mesurez sur l ?image les abcisses des zéros des 5 courbes présentes.
• Placez ces points dans la figure (2).
3
Un problème concret : trouver les fréquences de
vibration d’un diapason
(d’après http://supernovae.in2p3.fr/~llg/ En première approximation, on
peut modéliser un diapason comme l’assemblage de deux poutres, chacune de ces
poutres possédant une extrémité fixe (x = 0) et une extrémité libre (x = L). Les
vibrations y(x, t) transverses le long de chaque barre vérifient l’équation suivante en
fonction du temps
∂ 4y
ρ ∂ 2y
=
−
.
∂x4
Eκ2 ∂t2
Les conditions aux limites sont
(
∂2y
y = 0 et ∂x
à l’extrémité fixe x = 0,
2 = 0
∂y
∂3y
= 0 et ∂x3 = 0 à l’extrémité libre x = L.
∂x
On peut montrer que la fréquence propre ν est déterminée par l’équation non linéaire
suivante
s
π
Eκ2 2
cos(πb) cosh(πb) + 1 = 0, ν = 2
b,
2`
ρ
Où E est le module d ?Young du matériau, ρ sa masse volumique, ` la longueur
du diapason et κ une constante caractéristique de la section du diapason. Il y a
une infinité de valeurs de b, ordonnées de façon croissante, et donc une infinité de
fréquences propres νi , ordonnées de façon croissante.
• Tracer la courbe f (b) = cos(πb) cosh(πb) + 1 sur (0, 2) et trouver un encadrement de b1 . Utiliser alors fzero pour le calculer.
• Trouver tous les bi entre 0 et 10.
4
Ecriture d’un code de dichotomie
Nous cherchons à calculer un zéro x¯ de f .
• En supposant que vous connaissez un intervalle [a, b] pour lequel f (a)f (b) < 0,
écrire une fonction dicho dont les variables d’entrée sont la fonction f , les deux
bornes de l’intervalle a et b, et la tolérance tol. Les variables de sortie seront
[x, y] où x est l’approximation cherchée et y la valeur de la fonction en ce
point. Le script comportera des lignes de commentaire sur le mode
% BISECTION computes a root of a scalar equation
% [x,y]=Bisection(f,a,b,tol) finds a root x of the scalar function
% f in the interval [a,b] up to a tolerance tol. y is the
% function value at the solution
Pour être complet, le script comportera un test de vérification de ce que
f (a)f (b) < 0.
• Modifier la fonction pour compter le nombre d’itérations effectuées. Inclure en
variable d’entrée un nombre maximal d’itérations.
• Appliquer l’algorithme à la fonction x + exp(x). Attention, tenir compte de la
première partie du TP. Tracer la courbe de convergence.
• Appliquer l’algorithme à la fonction de la partie 3, calculant les modes propres
du diapason. Pour l’application numérique On prendra un diapason en fer
(ρ = 7874kg/m3 , E = 196GP a) de section circulaire (r = 5mm et κ = r/2 =
2, 5mm) et de longueur ` = 12, 6cm.
5
Ecriture d’un code de point fixe
On considère la méthode du point fixe :
xk+1 = G(xk ), où G(x) = x ⇐⇒ F (x) = 0
• Utiliser l’en-tête suivant pour écrire un script :
function x = fixedPoint(g,x0,tol,maxiter)
% FIXEDPOINT solves g(x) = x by fixed-point iteration
%
x = fixedPoint(g,x0,tol,maxiter) solves the equation g(x) = x by fixed-point iteration
%
with initial guess x0. The function terminates when a solution is found within toleran
%
tol, or when the maximum number of iterations maxiter is reached.
On utilisera le critère d’arrêt
αk
kxk − xk−1 k ≤ tol,
1 − αk
où αk = max{1/2, kxk − xk−1 k/kxk−1 − xk−2 k} est une approximation du
facteur de contractivité.
• Tester votre programme à l’aide de l’exemple F (x) = x, utiliser la fonction
Gi (x) = x − γF (x). Comment faut-il choisir γ pour que l’itération converge
vers la solution de F (x) = 0 ? Comparer la vitesse de convergence avec celle
de la méthode de dichotomie, en faisant une esquisse de l’erreur kxk − x∗ k par
rapport à k.
• Appliquer l’algorithme à la fonction de la partie 3, calculant les modes propres
du diapason.
6
Ecriture d’un code de Newton
Mêmes questions avec Newton. Comparer toutes les méthodes sur l’exemple du
diapason.