1.1      Retour en dimension 1

• On considère la fonction . On considère les points correspondant à une descente de gradient à pas fixe α pour la fonction . Programmer numériquement la suite   pour des valeurs de α et de  où la convergence est linéaire. Illustrer les différents types de convergence à l’aide de graphiques.
• On considère la fonction. On considère les points  correspondant à la méthode de descente de gradient à pas fixe α. Programmer la suite  lorsque et α = 1/3 .

Discuter de l’intérêt ou non de ce genre de méthode par rapport par exemple à la méthode de la section dorée.

• Appliquer la méthode du gradient à pas optimal sur la fonction   et afficher les points obtenus sur un graphique (on pourra reprendre le code de la partie 2 ci-dessus). Observer que les directions de descentes  sont tangentes aux lignes de niveaux en  et orthogonales aux lignes de niveau en   .

1.2      Méthode du gradient projeté

On considère la fonctionnelle J définie sur  par :

On appelle  le cadrant défini par :

 

• Quelle est la solution du problème de minimisation de J sans contrainte ?
• On appelle, le minimum de J sous la contrainte.
• Mettre en œuvre la méthode du gradient projeté pour résoudre le problème de minimisation

 

• Représenter les itérés par cette méthode

SOLUTION

TP 3_ROGE_POUPA

CODE

import matplotlib.pyplot as plt
import numpy as np
import math as mp

###Exercice 1###

###1###

def f(x): # définition de f
return (x**2 – (1/4)*(x**4))

def Grad_f(x): # définition du gradient de f
return (2*x – x**3)

def Grad_fixe_f(x,alpha):
#Approximation de f par la méthode du gradient à pas fixe
k = 0
kmax=1000
epsilon=10**(-4)
X = [x]
K = [k]
while np.abs(- alpha*(Grad_f(x)))>epsilon and k < kmax:
x= x – alpha*(Grad_f(x))
k = k+1
K.append(k)
X.append(x)
return (X,K)

# Affichage de la convergence de la fonction f avec la méthode du gradient à pas fixe

Liste_Alpha = [0.3,0.2,0.1,0.05,0.02,0.01]
plt.xlabel(‘Iteration’)
plt.ylabel(‘suite xk’)
plt.title(« Illustration des différents types de convergence pour x0=0.5 « )

for i in Liste_Alpha :
Liste_conv = Grad_fixe_f(0.5,i)
plt.plot(Liste_conv[1],Liste_conv[0], label = i)
plt.legend()
plt.show()

### Question 2 ###

def f_2(x): # définition de la fonction f de la question 2 nommée f_2
return (x**2/(1+x**2))

def Grad_f_2(x): # définition du gradient de la fonction
return (2*x)/((1+x**2)**2)

def Grad_fixe_f_2(x,alpha):
#Approximation de f_2 par la méthode du gradient à pas fixe
k = 0
kmax=1000
epsilon=10**(-3)
X = [x]
K = [k]
while np.abs(- alpha*(Grad_f_2(x)))>epsilon and k < kmax:
x= x – alpha*(Grad_f_2(x))
k = k+1
K.append(k)
X.append(x)
return (X,K)

def SectionDoree(a0,b0,n_max,f):
t = (1+np.sqrt(5))/2
a,b = a0,b0
for i in range(n_max):
A = a + (b-a)/(t**2)
B = a + (b-a)/t
if f(A) > f(B) :
b = B
elif f(A) < f(B):
a = A
else :
a,b = A,B
return (a+b)/2

def GradPasOptimal(X0):
epsilon=10**(-3)
kmax=100
X = [X0]
x = X0
dx = Grad_f_2(x)
k=0
K = [k]
while k < kmax and np.abs(x**2+dx**2) >=epsilon:
b = SectionDoree(x,dx,100,f_2)
x=x-b*dx
dx = Grad_f_2(x)
k+=1
K.append(k)
X.append(x)
return (X,K)

# Affichage de la convergence de la fonction f_2 par la méthode du gradient à pas fixe et par la méthode du gradient à pas optimal

GF_f2 = Grad_fixe_f_2(3,1/3)

GO_f2 = GradPasOptimal(3)

Liste_Alpha_2 = [0.5,0.4,1/3,0.2]
for i in Liste_Alpha_2 :
GF_f2 = Grad_fixe_f_2(3,i)
plt.plot(GF_f2[1],GF_f2[0], label = i)

plt.plot(GO_f2[1],GO_f2[0],label = « pas optimal »)
plt.xlabel(‘Iteration’)
plt.ylabel(‘suite xk’)
plt.title(« Comparaison de la méthode du gradient à pas optimal et à pas fixe »)
plt.legend()
plt.show()

print(« Nombre d’itérations pour la méthode du gradient à pas optimal pour g :  » +str(GO_f2[1][-1]))

###3###

a=3

def fa(X): # X estun vecteur
if a > 0:
return 1-1/(1+a*(X[0]**2)+X[1]**2)

def grad_fa(X):
return( 2*a*X[0]/(1+a*X[0]**2+X[1]**2)**2 , 2*X[1]/(1+a*X[0]**2+X[1]**2)**2)

def SectionDoree_2(a0,b0,n_max,f):
t = (1+np.sqrt(5))/2
a,b = a0,b0
A=[0,0]
B=[0,0]
for i in range(n_max):
A[0] = a[0] + (b[0]-a[0])/(t**2)
A[1] = a[0] + (b[0]-a[0])/(t**2)

B[0]= a[0] + (b[0]-a[0])/t
B[1]= a[1] + (b[1]-a[1])/t

if f(A) > f(B) :
b = B
elif f(A) < f(B):
a = A
else :
a,b = A,B

return (mp.sqrt(a[0]**2+a[1]**2)+mp.sqrt(b[0]**2+b[1]**2))/2

def GradPasOptimal_2(X0):
epsilon=10**(-3)
kmax=100
X = [X0]
x1 = [X0[0]]
y1 = [X0[1]]
[x,y] = X0
F = [fa([x,y])]

[dx,dy] = grad_fa([x,y])
Dx=[dx,dy]
k=1

while k < kmax and mp.sqrt(dx**2+dy**2) >=epsilon:

b = SectionDoree_2(X[-1],Dx,100,fa)

x,y=x-b*dx,y-b*dy
[dx,dy] = grad_fa([x,y])
k+=1
F.append(fa([x,y]))
X.append([x,y])
x1.append(x)
y1.append(y)

return (F,X,k,x1,y1,Dx)

plt.figure()
plt.xlabel(‘x’)
plt.ylabel(‘y’)

plt.title(« Tracé de la fonction fa avec mise en évidence de la convergence »)
x = np.linspace(-1.1,1.1)
y = np.linspace(-1.1,1.1)
z = np.array([[fa(np.array([x1,y1])) for x1 in x] for y1 in y])
convergence = GradPasOptimal_2([1,0])

x2 = np.linspace(-3,3,convergence[2])
y2 = np.linspace(-3,3,convergence[2])
plt.plot(convergence[3],convergence[4], »+ »)

plt.contour(x,y,z)
plt.legend([‘convergence’])
plt.show()

###Exercice 2###

###4###

# Définition de la fonction J à optimiser et de ses dérivées

def J (x,y):
return 2*(x**2)+3*x*y+2*(y**2)

def dJx(x, y):
return 4*x + 3*y

def dJy(x, y):
return 3*x + 4*y

def gradientJ(x, y):

return np.array([dJx(x, y), dJy(x, y)])

# Fonction de projection du vecteur X. Ici, la zone de projection (espace des contraintes)
# est définie par les conditions x <= x_max et y <= y_max

def proj(X):

return np.array([min(x_max, X[0]), min(y_max, X[1])])
# Pour chaque vecteur, on lui associe le point le plus proche appartenant
# à l’espace des contraintes.

# Méthode du gradient projeté. On utilise un pas fixe pour les itérations

def Gradiant_Pas_fixe_Projete(X0, pas):

n_max=1000
epsilon=1e-8
x, y = X0
i = 0

d = -gradientJ(x, y)
X1 = proj(X0 + pas*d)

while (i < n_max) and (np.abs((X1[0]-X0[0])**2+(X1[1]-X0[1])**2) >= epsilon):

d = -gradientJ(X1[0], X1[1])
X0 = X1
X1 = proj(X1 + pas*d)

return X1

X = np.array([1,0])
alpha = 0.1
x_max = -0.5
y_max = -0.5

# Tracé de la fonction et visualisation du point final

X1 = Gradiant_Pas_fixe_Projete(np.array([-1,-1]), 1e-4)

« Valeurs pour les deux variables dans le but de tracer la fonction dans le bon intervalle »
figure = plt.figure()
x1 , x2 = np.meshgrid(np.linspace(-1,1,200),np.linspace(-1,1,200))
plt.title(‘Itérés par la methode du gradiant projeté’)
plotHX = plt.contourf(x1,x2,J(x1,x2),np.linspace(0,1,20))
plt.xlabel(‘x’)
plt.ylabel(‘y’)
cbar = plt.colorbar(plotHX)
cbar.ax.set_ylabel(‘h(x,y)’)
plt.plot([X1[0]], [X1[1]], ‘ro’, label = « solution de l’optimisation »)

line_x = np.linspace(-1, x_max, 10)
line_ymax = np.array(10*[y_max])
line_y = np.linspace(-1, y_max, 10)
line_xmax = np.array(10*[x_max])

plt.plot(line_x, line_ymax, color = ‘b’, label = « Espace des contraintes : Q »)
plt.plot(line_xmax, line_y, color = ‘b’)

plt.legend()
plt.show()