//
//   méthode de descente////
//  ordre 1 ou ordre 2
//  recherche linéaire: Armijo+backtracking ou Armijo+Wolfe
///////////////////////////////////
//
function z=f(X)
    x=X(1);y=X(2);
    z=10*(y-x^2)^2+(x-1)^2
endfunction
function Z=gradf(X)
     x=X(1);y=X(2);
     Z=[-20*(y-x^2)*(2*x)+2*(x-1);20*(y-x^2)]    
endfunction
function Z=hessien(X)
     x=X(1);y=X(2);
     Z=[-40*(y-x^2)+80*x^2+2, -40*x;-40*x,20];   
endfunction
//
//   méthode de descente////
//
// parametres
//
// Armijo+backtracking
alphainit=1
beta=0.001
tau=0.3
//
//// Armijo+Wolfe
m1=0.001;m2=0.9;  
lambda=1.5;
//
//
X0=[0;1];
Niter=evstr(x_dialog('nombre d''iterations','20'));
Ndir=evstr(x_dialog('direction de descente (1: gradient, 2: hessien)','1'));
rec=evstr(x_dialog('recherche linéaire (1: Armijo+backtrack, 2: Armijo+Wolfe)','1'));
X=X0;
Xtot=X;
for i=1:Niter
if (Ndir==1) then
    d=-gradf(X)  // direction de descente du hessien//
else
    d=-inv(hessien(X))*gradf(X) // direction de descente du gradient//
end
alpha=alphainit;qprime0=d'*gradf(X);
if (rec==1) then    // Armijo + backtracking
  while f(X+alpha*d)>f(X)+beta*alpha*(qprime0)
  alpha=tau*alpha;
  end    
else        // Armijo +Wolfe
  tg=0;td=0;
  while (f(X+alpha*d)>f(X)+m1*alpha*qprime0)|(alpha*(d'*gradf(X+alpha*d))<m2*qprime0)
    if   f(X+alpha*d)>f(X)+m1*alpha*(qprime0) then
     td=alpha;  
    else 
     tg=alpha;   
    end
   if td==0 then 
       alpha=lambda*alpha
   else
       alpha=(tg+td)/2
   end           
  end     
end  
X=X+alpha*d;Xtot=[Xtot,X];      
end
disp('valeur obtenue:');disp(X)
//
///////////////// tracé//////////
Nt=100;
x=linspace(-1,2,Nt);
y=linspace(-1,2,Nt);
z=zeros(Nt,Nt);
for i=1:Nt
    for j=1:Nt
       z(i,j)=f([x(i),y(j)]) 
    end
end
//xset('window',0)
//clf()
//plot3d(x,y,z)

xset('window',1)
clf()
contour(x,y,z,[0:0.1:1,2:10])
plot2d(Xtot(1,:),Xtot(2,:))









