• Gradient
• Gradient et extremum relatif

Gradient et extremums relatifs:
       une approche numrique

Pierre Lantagne (avril 2001)

Collge de Maisonneuve

plantag@edu.cmaisonneuve.qc.ca

http://math.cmaisonneuve.qc.ca/plantagne

Cette feuille a pour but de prsenter une approche numrique dans la recherche d'un extremum relatif. Nous allons utiliser le gradient d'une fonction afin de btir, point par point, un chemin vers un maximum ou vers un minimum relatif.

Gradient

Rendons d'abord disponibles les macro-commandes des extensions linalg, plottools et plots.

>	with ( linalg ): with ( plottools ): with ( plots ):

Warning, the protected names norm and trace have been redefined and unprotected

Warning, the name changecoords has been redefined

Warning, the previous binding of the name arrow has been removed and it now has an assigned value

Considrons la fonction g dfinie par .

>	g:=(x,y)->x*exp(-(x-1)^2-(y+0.25)^2)+(y)*exp(-x^2-y^2)+3;

Traons la surface d'quation sur le pav [-2; 3,5] x [-2,5; 2].

>	Surface:=plot3d([x,y,g(x,y)],x=-2..3.5,y=-2.5..2,grid=[30,30]): display(Surface,labels=[x,y,z],        lightmodel=light2,style=hidden,shading=zhue,color=navy,        axes=framed,font=[TIMES,ROMAN,8],        orientation=[-160,60]);

Cette surface possde des extremums relatifs. Traons alors des courbes de niveau pour mieux estimer les valeurs x et  y qui permettent d'atteindre ces extremums.

>	Courbes_niveau:=contourplot(g(x,y),x=-2..3,y=-2.5..2,contours=15,numpoints=1300,color=orange): Courbes_niveau;

>

Les courbes de niveau nous montrent que le maximum est atteint avec (x,y) au voisinage du point (1,4; -0,3) et que le minimum est atteint avec (x,y) dans le voisinage du point (-0,2; -0,7).

Visualisons le fait que le gradient en un point est toujours perpendiculaire la courbe de niveau passant par ce point. Obtenons d'abord un champ de vecteurs gradient. La macro-commande gradplot de l'extension plots sera utilise.

>	gp:=gradplot(g(x,y),x=-2..3.5,y=-2.5..2,grid=[25,25],arrows=slim,color=navy): gp;

>

Superposons ensuite ce champ de vecteurs au trac des courbes de niveau.

>	display([Courbes_niveau,gp]);

>

On observe effectivement que les flches sont toujours perpendiculaire aux courbes de niveau. Slectionnez le graphique prcdent et cliquez sur la plus grosse loupe (ou utilisez le raccourcis-clavier ctrl-n) pour agrandir le graphique. Cette illustration nous indique la direction avec laquelle la fonction possde le plus grand taux de variation. Le sens de ces vecteurs confirme qu'il se prsente effectivement un maximum relatif (flches rentrantes) dans le voisinage du point (1,4; -0,3) et un minimum relatif (flches sortantes) dans le voisinage du point (-0,2; -0,7).

 

Pour rechercher les couples (x,y) avec lesquels la fonction g prsentera un extremum relatif, appliquons d'abord le test d'extremum local. En fait, nous nous contenterons de trouver les points critiques et sur la base du trac de la surface, nous identifierons la nature de chacun de ces points.

>	gx:=D[1](g): gy:=D[2](g): Systeme:={gx(x,y),gy(x,y)};

>	Sol_symboliques:=solve({gx(x,y)=0,gy(x,y)=0,x<infinity},{x,y});

>

Il fallait s'y attendre, l'valuateur n'a trouv aucune solution ce systme. Plus souvent qu'autrement, il faut plutt rsoudre numriquement de tels systmes d'quations. Rsolvons d'abord numriquement ce systme au voisinage de (1,4; -0,3) avec la macro-commande fsolve.

>	xy_max:=fsolve({gx(x,y)=0,gy(x,y)=0},{x,y},{x=0..1.4,y=-0.3..0});

>

Rsolvons ensuite numriquement le systme au voisinage de (-0,2; -0,7).

>	xy_min:=fsolve({gx(x,y)=0,gy(x,y)=0},{x,y},{x=-0.2..0,y=-0.7..-0.5});

>

Contrlons si chaque rsultat annule effectivement le gradient.

>	Gradient:=[gx(x,y),gy(x,y)];

>	evalf(subs(xy_max,Gradient)); evalf(subs(xy_min,Gradient));

Les valeurs des composantes sont presque nulles. Obtenons les valeurs de ces extremums relatifs.

>	z_max:=evalf(subs(xy_max,g(x,y))); z_min:=evalf(subs(xy_min,g(x,y)));

>

Gradient et extremum relatif

Commenons par rechercher un maximum relatif de la fonction (en supposant qu'il existe, bien sr). Rappelons que le gradient en un point indique la direction du plus grand taux de variation de la fonction. Alors, avec un point de dpart, nous allons nous dplacer lgrement, dans la direction donne par le gradient en ce point, vers un deuxime point . Si, avec le point de dpart, on n'atteint pas un maximum de la fonction, ce deuxime point nous permettra par contre de se rapprocher du point qui permet la fonction d'atteindre le maximum recherch.

Soit donc l'quation vectorielle de la droite passant par le point et ayant pour vecteur directeur. La valeur de >0 dterminera l'importance du dplacement dans le direction de . Avec une valeur >0 donne, on dduira un point particulier de cette droite. Ce point particulier, qu'on nommera , nous permettra de se rapprocher du point (x,y) permettant d'atteindre le maximum cherch.

Voyons maintenant comment appliquer ce raisonnement. Comme les calculs seront numriques, nous considrerons que le point P(x,y) permet d'atteindre un maximum relatif de la fonction si le gradient en ce point est presque gal au vecteur [0,0]. Un faon lgante de tester ce fait, c'est de vrifier si la longueur (la norme) de ce vecteur gradient est presque nulle.

Pour les besoins de la cause, le point P(-2,2) sera dsign point de dpart.

>	P[0]:=[-2,2];

Contrlons (ici, par principe) si ce point permet d'atteindre le maximum de la fonction. Utilisons une tolrance de .

>	tolerance:=0.00000001; Grad_P:=evalf(subs({x=P[0][1],y=P[0][2]},Gradient)); Grad_longueur:=norm(Grad_P,2); evalb(Grad_longueur<tolerance);

On a utilis la macro-commande norm de l'extension linalg pour calculer la longueur du gradient au point P(-2,2). De plus, rappelons que evalb commande une valuation boolenne. Le rsultat vrai pour l'ingalit Grad_longueur<tolerance signifiera que le point P(x,y) est, au niveau de tolrance admis, suffisamment prs du point permettant d'avoir le maximum cherch.

Obtenons maintenant un deuxime point avec un pas .

>	lambda:=0.1;

Obtenons le point .

>	Grad_P:=evalf(subs({x=P[0][1],y=P[0][2]},Gradient)); P[1]:=P[0]+lambda*Grad_P;

Contrlons si le point permet d'atteindre le maximum recherch ( toujours avec la tolrance de ).

>	Grad_P:=evalf(subs({x=P[1][1],y=P[1][2]},Gradient)); Grad_longueur:=norm(Grad_P,2); evalb(Grad_longueur<tolerance);

Il sera beaucoup plus commode d'utiliser une boucle de calcul pour continuer de se dplacer de point en point. Prenons soins de limiter le nombre d'itrations 1000 et rappelons les paramtres qui seront utiliss dans cette boucle.

>	k:=1000; lambda:=0.1; tolerance:=0.00000001;

>

L'algorithme de cette boucle est le suivant.

1) valuer le gradient au point .

2) Calculer la longueur (la norme) de ce vecteur gradient.

3) Contrler si cette longueur est infrieur la longueur tolre.

4) Si tel est le cas, terminer les itrations

5) Sinon, calculer un prochain point .

>	P[0]:=[-2,2];  # Point de dpart

>

for n from 0 to k do  Grad_P:=evalf(subs({x=P[n][1],y=P[n][2]},Gradient));  Grad_longueur:=norm(Grad_P,2);    if Grad_longueur<tolerance  then break  else P[n+1]:=P[n]+lambda*Grad_P;  fi; od: Grad_P:=evalf(subs({x=P[n][1],y=P[n][2]},Gradient)); norm(Grad_P,2); Nombre_iterations:=n; 'P[n]'=P[n]; n:='n':

Comparons ce point avec celui xy_max obtenu la section prcdente.

>

xy_max;

Traons, dans le plan, sans les relier, le chemin de points ( le point de dpart plus les 812 points obtenus par les itrations prcdentes).

>	Courbe:=pointplot({seq(P[n],n=0..Nombre_iterations)},color=navy): Courbe;

>

On observe trs bien sur le trac prcdent, que plus le taux de variation de la fonction g est important, plus la distance sparant les points augmente. Cela dcoule du fait que le pas de dplacement est constant et que le taux de variation de la fonction g ne l'est pas.

Illustrons ce fait en ralisant une animation du trac de ce chemin de points. Pour animer la partie significative du chemin de point, cette animation sera construite partir du 600 ime point jusqu'au 750 ime. Nous constaterons qu'avec un pas constant, plus la distance entre les points augmente, plus l'animation s'acclre.

>	Animation:=seq(pointplot([seq(P[n],n=0..k)],color=navy),k=600..750):

>	display(Animation,insequence=true);

Superposons maintenant le trac de ces points ceux du champ de vecteurs gradient et des courbes de niveau. On aura alors une ide claire du chemin parcourir dans le domaine, partir du point (-2, 2), pour se rapprocher du point (x,y) qui permet d'atteindre le maximum relatif de la fonction. Notez que le chemin obtenu est toujours, videmment, perpendiculaire toutes les courbes de niveau passant par ces points.

>	display([Courbe,Courbes_niveau,gp]);

>

Obtenons le trac de ce chemin dans l'espace en ajoutant, chacun des points calculs, une troisime composante (cote) dont la valeur sera celle de l'image par la fonction g des deux premires composantes. Relions cette fois-ci les points en prcisant l'option connect=true.

>	Courbe3d:=pointplot3d([seq([P[n][1],P[n][2],evalf(g(P[n][1],P[n][2]))],n=1..Nombre_iterations)],connect=true,thickness=3,axes=boxed,color=navy,orientation=[-160,60]): Courbe3d;

>

Superposons finalement le trac de cette courbe celui de la surface. L'effet est spectaculaire.

>	display([Surface,Courbe3d],labels=[x,y,z],        lightmodel=light2,style=hidden,shading=zhue,color=navy,        axes=framed,font=[TIMES,ROMAN,8],        orientation=[-190,50]);

>

>