BFGS算法Matlab（拟牛顿）

%bfgs方法基于wolf条件求步长
%x是近似最优点，val是近似最优值，k迭代次数
function [x,val,k]=BFGS(fun,gfun,x0)
maxk=500;
rho=0.55;eps=1e-5;
k=0; n=length(x0);
I=eye(n);
wk=eye(n);%初始对称正定矩阵一般选单位矩阵
c1=0.001; c2=0.8;
while (k<maxk)
   gk=feval(gfun,x0);
    if(norm(gk)<eps)
        break;
    end
    dk=-wk*gk;
    m=0;mk=0;
     while(m<20)%用Armijo求步长
        newf=feval(fun,x0+rho^m*dk);
        oldf=feval(fun,x0);
        newgk=feval(gfun,x0+rho^m*dk);
        oldgk=feval(gfun,x0);
        if(newf<oldf+c1*rho^m*gk'*dk)%wolf准测确定步长
            if(newgk'*dk>=c2*oldgk'*dk)
                mk=m;break;
            end 
        end
        m=m+1;
     end
     %BFGS矫正
     x=x0+rho^mk*dk;
     sk=x-x0;
     yk=feval(gfun,x)-gk;
     gamak=1/(yk'*sk);
     wk=(I-gamak*yk*sk')'*wk*(I-gamak*yk*sk')+gamak*sk*(sk)';
     k=k+1;x0=x; 
end
val=feval(fun,x0);

 下面为简单例子

function f=f1(x)
f=100*(x(1)^2-x(2))^2+(x(1)-1)^2;

function g=g1(x)%梯度
g=[400*x(1)*(x(1)^2-x(2))+2*(x(1)-1),-200*(x(1)^2-x(2))]';