该方程组，则称之为牛顿拉格朗日方法.因此，根据牛顿法的性质，该 方法具有局部二次收敛性质 下面写出牛顿-拉格朗日方法的详细算法步骤 算法12.1（牛顿-拉格朗日方法） 步0选取x0∈Rm,0∈R,p,Y∈(0,1)，0≤e≤1.令k=0. 步1计算7L(x,)川的值.若VL(xk,)川≤E,停算.否则， 转步2. 步2解方程组(12.6)得pk=（(dk,V) 步3若 I川VL(ck+dk,k+)川2≤(1-Y)川7L(ck,)l川2， (12.7) 则置=1，转步5；否则，转步4. Back Close
6/84 JJ II J I Back Close Têß|, K°Éè⁄Ó-.ÇKFê{. œd, ä‚⁄Ó{5ü, T ê{‰k¤‹g¬Ò5ü. e°—⁄Ó-.ÇKFê{ç[é{⁄½. é{ 12.1 (⁄Ó-.ÇKFê{) ⁄ 0 ¿ x0 ∈ R n , µ0 ∈ R l , ρ, γ ∈ (0, 1), 0 ≤ ε 1. - k := 0. ⁄ 1 Oé k∇L(xk, µk)k ä. e k∇L(xk, µk)k ≤ ε, é. ƒK, =⁄ 2. ⁄ 2 )êß| (12.6) pk = (dk, νk). ⁄ 3 e k∇L(xk + dk, µk + νk)k 2 ≤ (1 − γ)k∇L(xk, µk)k 2 , (12.7) Kò αk := 1, =⁄ 5; ƒK, =⁄ 4
步4令mk是使下面的不等式成立的最小非负整数m: VL(ck +pmdk,uk+pmvk2<(1-Ypm)L(xk,uk)2,(12.8) 置ak=pmk. 步5令xk+1=xk十Qkdk,k+1=k十QVk.置k=k十1，转 步1. S12.1.2牛顿-拉格朗日法的Matlab程序 本小节通过一个具体的例子来介绍算法12.1（牛顿-拉格朗日方 法)的Matlab实现. Back Close
7/84 JJ II J I Back Close ⁄ 4 - mk ¥¶e°ÿ ™§·ÅöKÍ m : k∇L(xk + ρ mdk, µk + ρ mνk)k 2 ≤ (1 − γρm)k∇L(xk, µk)k 2 , (12.8) ò αk = ρ mk . ⁄ 5 - xk+1 = xk + αkdk, µk+1 = µk + αkνk. ò k := k + 1, = ⁄ 1. §12.1.2 ⁄Ó-.ÇKF{ Matlab ßS !œLòá‰N~f50 é{ 12.1 (⁄Ó-.ÇKFê {) Matlab ¢y.
例12.1用算法12.1编程计算下列最优化问题的极小点 min f(c)=ei2243425- ++1只 1 s.t. 晚+吃+喝+x+喝-10=0， T2C3-5x4x5=0, x3+x号+1=0. 该问题有最优解x*=(-1.71,1.59,1.82,-0.763,-0.763)T,最优值f(x*)= 0.0539. 解编制Matlab程序如下 function [x,mu,val,mh,k]=newtlagr(x0,mu0) %功能：用牛顿-拉格朗日法求解约束优化问题： %minf(x),s.t.h_i(x)=0,=1,.,1. %输入：x0是初始点，mu0是乘子向量的初始值 %输出：x,m山分别是近似最优点及相应的乘子， %val是最优值，mh是约束函数的模，k是迭代次数. Back maxk=500;%最大迭代次数 Close
8/84 JJ II J I Back Close ~ 12.1 ^é{ 12.1 ?ßOéeÅ`zØK4: min f(x) = ex1x2x3x4x5 − 1 2 (x 3 1 + x 3 2 + 1)2 , s.t. x2 1 + x 2 2 + x 2 3 + x 2 4 + x 2 5 − 10 = 0, x2x3 − 5x4x5 = 0, x 3 1 + x 3 2 + 1 = 0. TØKkÅ`) x ∗ = (−1.71, 1.59, 1.82, −0.763, −0.763)T , Å`ä f(x ∗ ) = 0.0539. ) ?õ Matlab ßSXe function [x,mu,val,mh,k]=newtlagr(x0,mu0) % ıU:^⁄Ó-.ÇKF{¶)Â`zØK: % min f(x), s.t. h_i(x)=0, i=1,., l. %—\:x0¥–©:, mu0¥¶fï˛–©ä %——: x,mu©O¥CqÅ`:9ÉA¶f, %val¥Å`ä,mh¥ÂºÍ ,k¥SìgÍ. maxk=500; %ÅåSìgÍ
n=length(x0);1=length(mu0); rho=0.5;gamma=0.4; x=x0;mu=mu0; k=0;epsilon=1e-8; while(k<maxk) dl=d1a(x,mu);%计算乘子函数的梯度 if(norm(dl)<epsilon),break;end%检验终止准则 N=N1(x,mu);%计算拉格朗日矩阵 dz=-N\dl; %解方程组得搜索方向 dx=dz(1:n);du=dz(n+1:n+l) m=0;mk=0; while(m<20) %Armijo搜索 if (norm(dla(x+rho m*dx,mu+rho m*du))2<=(1-gamma*rho"m)*norm(dl)2) mk=m;break; end m=m+1; end x=x+rho'mk*dx;mu=mu+rho mk*du; k=k+1; end Back Close
9/84 JJ II J I Back Close n=length(x0); l=length(mu0); rho=0.5; gamma=0.4; x=x0; mu=mu0; k=0; epsilon=1e-8; while(k<maxk) dl=dla(x,mu); %Oé¶fºÍF› if(norm(dl)<epsilon), break; end %u™éOK N=N1(x,mu); % Oé.ÇKF› dz=-N\dl; %)êß||¢êï dx=dz(1:n); du=dz(n+1:n+l); m=0; mk=0; while(m<20) % Armijo|¢ if(norm(dla(x+rho^m*dx,mu+rho^m*du))^2<=(1-gamma*rho^m)*norm(dl)^2) mk=m; break; end m=m+1; end x=x+rho^mk*dx; mu=mu+rho^mk*du; k=k+1; end
val=f1(x); mh=norm(h1(x),inf); %%%%%拉格朗日函数L(x,mu)%%%%%%%% function l=la(x,mu) f=f1(x);%调用目标函数文件 h=h1(x);%调用约束函数文件 1=f-mu'*h;%计算乘子函数 %%%%%拉格朗日函数的梯度%%%%%% function dl=dla(x,mu) df=df1(x);%调用目标函数梯度文件 h=h1(x);%调用约束函数文件 dh=dh1(x);%调用约束函数Jacobi:矩阵文件 dl=[df-dh'*mu;-h];%计算乘子函数梯度文件 %%%%%拉格朗日函数的Hesse阵%%%%%% function d21=d2la(x,mu) d2f=d2f1(x);%调用目标函数Hesse阵文件 [d2h1,d2h2,d2h3]=d2h(x);%调用约束函数二阶导数文件 d21=d2f-mu(1)*d2h1-mu(2)*d2h2-mu(3)*d2h3;%计算乘子函数的Hesse阵 %%%%%系数矩阵N(x,mu)%%%%%% function N=N1(x,mu)%计算拉格朗日矩阵 Back Close
10/84 JJ II J I Back Close val=f1(x); mh=norm(h1(x),inf); %%%%%%%%% .ÇKFºÍ L(x,mu) %%%%%%%%%%%%% function l=la(x,mu) f=f1(x); %N^8IºÍ©á h=h1(x); %N^ºͩá l=f-mu’*h; % Oé¶fºÍ %%%%%%%%% .ÇKFºÍF› %%%%%%%%%%%%% function dl=dla(x,mu) df=df1(x); %N^8IºÍF›©á h=h1(x); %N^ºͩá dh=dh1(x); %N^ºÍJacobi› ©á dl=[df-dh’*mu; -h]; %Oé¶fºÍF›©á %%%%%%%%% .ÇKFºÍHesse %%%%%%%%%%% function d2l=d2la(x,mu) d2f=d2f1(x); %N^8IºÍHesse ©á [d2h1,d2h2,d2h3]=d2h(x); %N^ÂºÍ Í©á d2l=d2f-mu(1)*d2h1-mu(2)*d2h2-mu(3)*d2h3; %Oé¶fºÍHesse %%%%%%%%% XÍ› N(x,mu) %%%%%%%%%%% function N=N1(x,mu) %Oé.ÇKF›