syms x1 x2;
f=4*(x1)^2-4*(x1)*(x2)+3*(x2)^2+x1;%函数表达式 %f=8*(x1)^2+4*(x1)*(x2)+5*(x2)^2; X=-10:0.1:10; Y=-10:0.1:10; [X,Y]=meshgrid(X,Y);
Z=4.*X.^2-4.*X.*Y+3.*Y.^2+X; mesh(X,Y,Z) contour(X,Y,Z) x0=[10 10]';%初始点 v=[x1,x2];
mu=0.01;%最小误差
gradf=gradient(f);%函数的梯度 H=jacobian(gradf,v);
h0=subs(H,[x1;x2],x0);%在点x0处的Hessian g0=subs(gradf,[x1;x2],x0);%在点x0处的梯度值 if g0'*g0 d0=-g0;%搜索方向 alpha=-(g0'*d0)/(d0'*h0*d0);%步长 xk=x0+alpha*d0;%下一点 gk=subs(gradf,[x1;x2],xk);%梯度值 beta=gk'*gk/(g0'*g0);%求搜索方向时的系数 dk=-gk+beta*d0;%下一个方向 x0=xk;%更新点 g0=gk;%更新所在点的梯度 d0=dk;%更新方向 while g0'*g0>mu alpha=-(g0'*d0)/(d0'*h0*d0);%步长 xk=x0+alpha*d0;%下一点 gk=subs(gradf,[x1;x2],xk);%梯度值 beta=gk'*gk/(g0'*g0);%求搜索方向时的系数 dk=-gk+beta*d0;%下一搜索方向 x0=xk;%更新点 g0=gk;%更新所在点的梯度 d0=dk;%更新方向 hk=subs(H,[x1;x2],x0);%在点xk处的梯度值 h0=hk;%更新矩阵 end minf=subs(f,[x1;x2],xk)%函数的最小值 xk