function [x,val,k]=revisenm(fun,gfun, Hess, x0)?
%功能:用修正牛頓法求解無約束問題:min f(x)?
%輸入:x0是初始點,fun, gfun,Hess分別是求目標函數值,梯度,Hesse矩陣的函數
%輸出:x,val分別是近似最優點和最優值,k是迭代次數
n=length(x0);maxk=150;?
rho=0.55;sigma=0.4;tau=0.0;
k=0; epsilon=1e-5;
while (k<maxk)?
? ? ?gk=feval(gfun,x0);%計算梯度
? ? ?muk=norm (gk)^(1+tau);?
? ? ?Gk=feval(Hess,x0);%計算Hesse矩陣
? ? ?Ak=Gk+muk*eye (n);?
? ? ?dk=-Ak\gk;%解方程組Gk*dk=-gk,計算搜索方向
? ? ?if(norm(gk)<epsilon),break; end %檢驗終止準則
? ? ?m=0;mk=0;?
? ? ?while (m<20) %用Armij0搜索求步長
? ? ? ? ?if(feval(fun, x0+rho^m*dk)<feval (fun, x0)+sigma*rho^m*gk'*dk)?
? ? ? ? ? ? ?mk=m; break;?
? ? ? ? ?end?
? ? ? ? ?m=m+1;?
? ? ?end?
? ? ?x0=x0+rho^mk*dk;?
? ? ?k=k+1;?
end?
x=x0;?
val=feval (fun,x);
function f=fun(x)
f=100*(x(1)^2-x(2))^2+(x(1)-1)^2;
function g=gfun(x)
g=[400*x(1)*(x(1)^2-x(2))+2*(x(1)-1),-200*(x(1)^2-x(2))]';
function He=Hess(x)?
n=length(x);
He=zeros(n,n);
He=[1200*x(1)^2-400*x(2)+2, -400*x(1);
? ? ?-400*x(1), ? ? ? ? ? ? ? ?200 ?];