罚函数法罚函数法是能够处理一般的约束优化问题:min()()0,1,2,()0,1,2,,iifxhxikgxjm的一类方法。其基本思想是将约束优化问题卑微无约束问题来求解。罚函数是由目标函数和约束函数的某种组合得到的函数,对于等式约束的优化问题min()()0,1,2,ifxhxik,可以定义如下的罚函数:21()()()kiiFxfxchx将约束优化问题转化为无约束优化问题;对于不等式约束的优化问题min()()0,1,2,,ifxgxjm可以定义如下的罚函数:11()()()mjjFxfxCgx对于同时存在等式约束和不等式约束的优化问题,可以去上面两个罚函数的组合。当然罚函数还有其他的取法,但是构造罚函数的思想都是一样的,即使得在可行点罚函数等于原来的目标函数值,在不可行点罚函数等于一个很大的数。外点罚函数法1.算法原理外点罚函数法是通过一系列罚因子{}ic,求罚函数的极小值来逼近原约束问题的最有点。之所以称为外点罚函数法,是因为它是从可行域外部向约束边界逐步靠拢的。2,。算法步骤用外点罚函数法求解线性约束问题min()fxAxb的算法过程如下:1,给定初始点(0)x,罚参数列{}ic及精度0,置1k;2,构造罚函数2()()FxfxcAxb;3,用某种无约束非线性规划,以(1)kx为初始点求解min()Fx;4,设最优解为()kx,若()kx满足某种终止条件,则停止迭代输出()kx,否则令1kk,转2;罚参数列{}ic的选法:通常先选定一个初始常数1c和一个比例系数2,则其余的可表示为11iicc。终止条件可采用()Sx,其中2()SxcAxb。3算法的MATLAB实现function[x,minf]=minPF(f,x0,A,b,c1,p,var,eps)%目标函数:f;%初始点:x0;%约束矩阵:A;%约束右端向量:b;%罚参数的初始常数:c1;%罚参数的比例系数:p;%自变量向量Var;%精度:eps;%目标函数取最小值时自变量值:x;%目标函数的最小值:minf;formatlong;ifnargin==7eps=1.0e-4;endk=0;FE=0;fori=1:length(b)FE=FE+(var*transpose(A(1,:))-b(i))^2;endx1=transpose(x0);x2=inf;while1M=c1*p;FF=M*FE;SumF=f+FF;[x2,minf]=minNT(SumF,transpose(x1),var);ifnorm(x2-x1)=epsx=x2;break;elsec1=M;x1=x2;endendminf=Funval(f,var,x);formatshort;4.算法举例用外点罚函数法求下面的优化问题:22min(,)0.50.25,..1ftstsstts其中取10.05,2c,初始点取0(0,0)x。解:由题已知[11];1Ab,则每一步迭代需求解的罚函数为:222(,)0.50.250.05*2(1)iftststs在MATLAB命令窗口中输入:f=0.5*t^2+s^2/4;A=[11];b=1;c1=0.05;p=2;[x,minf]=minPF(f,[00],A,b,c1,p,[ts])所得结果为:x=0.33330.6666minf=0.1666对于一般的等式约束问题min()()0fxhx也可以用外点罚函数法解决:function[x,minf]=minGeneralPF(f,x0,h,c1,p,var,eps)formatlong;ifnargin==6eps=1.0e-4;endk=0;FE=0;fori=1:length(h)FE=FE+(h(i))^2;endx1=transpose(x0);x2=inf;while1M=c1*p;FF=M*FE;SumF=f+FF;[x2,minf]=minNT(SumF,transpose(x1),var);ifnorm(x2-x1)=epsx=x2;break;elsec1=M;x1=x2;endendminf=Funval(f,var,x);formatshort;例用通用罚函数法求下面的优化问题:min(,)2,..1,,0ftstssttsts其中取10.05,2cp,初始点取为0(0,0)x。外点罚函数法除了可以求含等式约束的优化问题外,还可以求混合约束的优化问题,不过其罚函数形式要修正。min()()0,1,2,()0,1,2,,iifxhxikgxjm通过引入罚函数2211()()()min(0,())kmiiijFxfxMhxMgx可以转化为无约束优化问题。function[x,minf]=minConPF(f,x0,g,h,c1,p,var,eps)formatlong;ifnargin==7eps=1.0e-4;endk=0;FE=0;fori=1:length(h)FE=FE+(h(i))^2;endx1=transpose(x0);x2=inf;while1M=c1*p;FF=M*FE;gx=Funval(g,var,x1);gF=0;fori=1:length(h)ifgx(i)0gF=gF+M*(g(i)^2);%构造罚函数endendSumF=f+FF+gF;[x2,minf]=minNT(SumF,transpose(x1),var);%用牛顿法求解无约束规划问题ifnorm(x2-x1)=epsx=x2;break;elsec1=M;x1=x2;endendminf=Funval(f,var,x);formatshort;内点罚函数法1.算法原理内点罚函数法的所有迭代过程均在可行域内进行,它每次迭代得到的点都是可行点,内点罚函数法的基本思想是用可行域内的点列来逼近最优解。2,。算法步骤用内点罚函数法求解线性约束问题min()()0,1,2,,ifxgxjm的算法过程如下:1,给定初始点(0)x,罚因子0,缩小系数1v及精度0,置1k;2,构造增广目标函数()()()FxfxBx,其中11()()miiBxgx;3,用某种无约束非线性规划,以(1)kx为初始点求解min()Fx;设最优解为()kx,若()()kBx,则终止迭代输出()kx,否则令v,1kk转2;3算法的MATLAB实现function[x,minf]=minNF(f,x0,g,u,v,var,eps)formatlong;ifnargin==6eps=1.0e-4;endk=0;FE=0;fori=1:length(g)FE=FE+1/g(i);%构造罚函数endx1=transpose(x0);x2=inf;while1FF=u*FE;SumF=f+FF;[x2,minf]=minNT(SumF,transpose(x1),var);%用牛顿法求解无约束规划Bx=Funval(FE,var,x2);ifu*Bxeps%精度判断ifnorm(x2-x1)=epsx=x2;break;elseu=v*u;%参数修正x1=x2;endelseifnorm(x2-x1)=epsx=x2;break;elseu=v*u;%参数修正x1=x2;endendendminf=Funval(f,var,x);formatshort;4算法举例2211min(,),.220tsftstsstts其中取8,0.5v,初始点取为(,)(3,1)ts。symsts;f=t^2+s^2;g=[t+s-1;2*t-s-2];[x,minf]=minNF(f,[31],g,8,0.5,[ts])x=1.00020.0001minf=1.0004乘子法1.算法原理外点罚函数法和内点罚函数法的均要求罚因子趋于无穷才能得到目标函数的最优解,但是罚因子太大的话会引起计算的困难,为克服这一缺陷就产生了乘子法。2,算法步骤用乘子法求解约束优化问题min()()0,1,2,,ifxhxjm的算法过程如下:1,给定初始点(0)x,初始乘子向量12(,,,)mvvvv,罚因子0M及放大系数1,及精度0,参数(0,1),置1k;2,构造目标函数211()()()()2mmiiiiiMFxfxhxvhx;3,用某种无约束非线性规划,以(1)kx为初始点求解min()Fx;设最优解为()kx,4,若()()khx,则终止迭代输出()kx,否则转2;5,若()(1)()()kkhxhx,则令MM,否则保持M不变,转66,令()(),1,2,,,kiiivvMhxim置1kk,转2,3算法的MATLAB实现function[x,minf]=minFactor(f,x0,h,v,M,alpha,gama,var,eps)formatlong;ifnargin==8eps=1.0e-4;endFE=0;fori=1:length(h)FE=h(i)^2;endx1=transpose(x0);x2=inf;while1FF=M*FE/2;Fh=v*h;SumF=f+FF-Fh;[x2,minf]=minNT(SumF,transpose(x1),var);Hx2=Funval(h,var,x2);Hx1=Funval(h,var,x1);ifnorm(Hx2)epsx=x2;break;elseifHx2/Hx1=gamaM=alpha*M;x1=x2;elsev=v-M*transpose(Hx2);x1=x2;endendendminf=Funval(f,var,x);formatshort;4算法举例22min(,),.21ftstsstts其中取8,3,0.7,(1),Mv,初始点取为0(0,0)x。symsts;f=t^2+s^2;h=[t+s-1;2*t-s-2];[x,minf]=minFactor(f,[00],h,2,3,0.7,[ts])x=0.20000.4000minf=0.2000