约束优化算法：拉格朗日乘子法

1、拉格朗日乘子法约束优化问题的标准形式为：min f ( x), xRnst.gi ( x)0,i1,2,., mhj (x)0, j1,2,.,l其中f , gi , hj : RnR约束优化算法的基本思想是： 通过引入效用函数的方法将约束优化问题转换为无约束问题，再利用优化迭代过程不断地更新效用函数，以使得算法收敛。1. 罚函数法罚函数法（内点法）的主思想是：在可行域的边界上筑起一道很高的“围墙” ，当迭代点靠近边界时， 目标函数陡然增大， 以示惩罚，阻止迭代点穿越边界， 这样就可以将最优解“挡”在可行域之内了。它只适用于不等式约束：min f (x), xRnst.gi0, i1,2,.,

2、 m它的可行域为：D xRn | gi (x)0, i1,2,., m对上述约束问题，其其可行域的内点可行集D0的情况下，引入效用函数：min B( x,r )f ( x)rB (x) 、m1m其中 B(x)或 B( x)|ln(gi (x) |i 1gi ( x)i 1算法的具体步骤如下：给定控制误差0，惩罚因子的缩小系数0 c1 。步骤 1：令 k1，选定初始点x(0)D0 ，给定 r10 （一般取10）。步骤 2：以 x(k ) 为初始点，求解无约束min B( x, r )f ( x)rk B(x)m1mgi (x) | ，得最优解 x(k )其中 B(x)或 B( x)|ln(x(

3、rk )i 1gi ( x)i 1步骤 3：若 rkB( x( k) )，则 x( k) 为其近似最优解，停；否则，令rk crk , k k 1，转步骤 2.2. 拉格朗日乘子法（ 1） PH 算法：（约数为等式的情况引入）效用函数为min M ( x, u (k) ,k )f ( x)u (k) T h(x)kh(x)T h(x)判断函数为khx( k)当 k(x(k ) )时迭代停止。步骤 1：选定初始点x(0)，初始拉格朗日乘子向量u(1) ，初始罚因子1 及其放大系数 c1，控制误差0 与常数(0,1) ，令 k1。步骤 2：以 x(k 1) 为初始点，求解无约束问题：min M (

4、 x, u (k) ,k )f ( x)u (k) T h(x)kh(x)T h(x)得到无约束问题最优解x(k )步骤 3：当 h x( k)时， x( k ) 为所求的最优解，停；否则转步骤4.步骤 4：当 hx( k )/ hx( k )时，转步骤5；否则令k 1ck ，转步骤 5.步骤 5：令 u(k 1)u( k)kh(x(k ) ), kk 1，转步骤 1。（2） PHR 算法（一般约束形式的松弛变量法和指数形式法）松弛变量法：1M (u, v,)f (x)212max 0, uigi ( x)ui22lvj hj (x)j1乘子的修正公式为：l2 j 1h2j (x)v(j k

5、1)ui(k 1)判断函数为：v(j k )hj( x( k ) ), j 1,.,lmax0,u( k) igi ( x( k) ) , i 1,.,mlm(k )21/2gi ( x( k ) ), uhj2 ( x( k) )kmaxij 1i 1当k(x(k ) )时迭代停止。3. 乘子法 MATLAB程序及其作用3.1Al _ main 函数程序（ 1）：乘子法效用函数程序函数功能：将约束优化问题，根据效用函数方法，将其转变成无约束问题。function f=AL_obj(x)%拉格朗日增广函数%N_equ 等式约束个数%N_inequ不等式约束个数global r_al pena

6、N_equ N_inequ;%全局变量h_equ=0;h_inequ=0;h,g=constrains(x);%等式约束部分for i=1:N_equh_equ=h_equ+h(i)*r_al(i)+(pena/2)*h(i).2;end%不等式约束部分for i=1:N_inequh_inequ=h_inequ+(0.5/pena)*(max(0,(r_al(i)+pena*g(i).2-r_al(i).2);end%拉格朗日增广函数值f=obj(x)+h_equ+h_inequ;程序（ 2）：判断函数函数功能：判断是否符合约束条件% the compare function is the

7、stop condition function f=compare(x)global r_al pena N_equ N_inequ; h_equ=0;h_inequ=0; h,g=constrains(x); %等式部分for i=1:N_equh_equ=h_equ+h(i).2;end%不等式部分for i=1:N_inequh_inequ=h_inequ+(max(-g(i),r_al(i+N_equ)/pena).2;endf=sqrt(h_equ+h_inequ);程序（ 3） AL 算法主程序函数功能：对无约束的效用函数利用拟牛顿算法求解其最优解，更新乘子。function X,

8、FVAL=AL_main(x_al,r_al,N_equ,N_inequ)%本程序为拉格朗日乘子算法示例算法%函数输入：% x_al: 初始迭代点% r_al ：初始拉格朗日乘子% N-equ: 等式约束个数% N_inequ: 不等式约束个数%函数输出% X ：最优函数点% FVAL: 最优函数值%=程序开始 =global r_al pena N_equ N_inequ; %参数（全局变量）pena=10;%惩罚系数c_scale=2;%乘法系数乘数cta=0.5;%下降标准系数e_al=0.005;%误差控制范围max_itera=25;out_itera=1;%迭代次数%=算法迭代开始

9、=while out_itera<max_iterax_al0=x_al;r_al0=r_al;% 判断函数compareFlag=compare(x_al0);% 无约束的拟牛顿法 BFGS X,FVAL=fminunc(AL_obj,x_al0);x_al=X;%得到新迭代点% 判断停止条件if compare(x_al)<e_aldisp('we get the opt point');breakend%c判断函数下降度if compare(x_al)<cta*compareFlagpena=pena;%可以根据需要修改惩罚系数变量elsepena=mi

10、n(1000,c_scale*pena); %乘法系数最大disp('pena=2*pena');end1000% 更新拉格朗日乘子h,g=constrains(x_al);for i=1:N_equ% 等式约束部分r_al(i)=r_al(i)+pena*h(i);endfor i=1:N_inequ% 不等式约束部分r_al(i+N_equ)=max(0,(r_al(i+N_equ)+pena*g(i);endout_itera=out_itera+1;end%+迭代结束 +disp('!the iteration over!');disp('the

11、 value of the obj function');obj(x_al)disp('the value of constrains');compare(x_al)disp('the opt point');X=x_al;FVAL=obj(X);3.1.4乘子法 AL _ main 函数使用方法（1）定义目标函数及约束条件min f ( x)x1 x2x2 x3 x3 x1st.x1x2x31 0x12x22x323 0目标函数 m 文件function f=obj(x)f=-x(1)*x(2)-x(2)*x(3)-x(3)*x(1);约束函数 m 文件function h, g constrains( x)hx(1)x(2)x(3) 1;gx(1)2x(2)2x(3)23;（ 2） AL _ main 函数调用x_al=1,