最优化大作业

最优化方法大作业---------用优化算法求解函数最值问题摘要最优化(optimization)是应用数学的重要研究领域.它是研究在给定约束之下如何寻求某些因素(的量)，以使某一(或某些)指标达到最优的一些学科的总称。最优化问题一般包括最小化问题和最大化问题，而最大化问题可以通过简单的转化使之成最最小化问题。最小化问题分为两类，即约束最小化和无约束最小化问题。在此报告中，前两个问题属于无约束最小化问题的求解，报告中分别使用了“牛顿法”和“共轭梯度法”。后两个问题属于有约束最小化问题的求解，报告中分别用“外点法”和“内点法”求解。虽然命名不一样，其实质都是构造“惩罚函数”或者“障碍函数”，通过拉格朗日乘子法将有约束问题转化为无约束问题进行求解。再此报告中，“外点法”和“内点法”分别用了直接求导和调用“牛顿法”来求解无约束优化问题。在此实验中，用“共轭梯度法”对“牛顿法”所解函数进行求解时出现错误，报告中另取一函数用“共轭梯度法”求解得到正确的结果。此实验中所有的函数其理论值都是显见的，分析计算结果可知程序正确，所求结果误差处于可接受范围内。报告中对所用到的四种方法在其使用以前都有理论说明，对“外点法”中惩罚函数和“内点法”中障碍函数的选择也有相应的说明，另外，对此次试验中的收获也在报告的三部分给出。本报告中所用程序代码一律用MATLAB编写。【关键字】函数最优化牛顿法共轭梯度法内点法外点法MATLAB一，问题描述1，分别用共轭梯度发法和牛顿法来求解一下优化问题441432243221102510minxxxxxxxxxf2，分别用外点法和内点发求解一下优化问题01.min212231xxtsxx二、问题求解1.1用牛顿法求解441432243221102510minxxxxxxxxxf1.1.1问题分析：取步长为1而沿着牛顿方向迭代的方法称为牛顿法，在牛顿法中，初始点的取值随意，在以后的每次迭代中，kkkkxfxfxx121，直到终止条件成立时停止。1.1.2问题求解注：本程序编程语言为MATLAB，终止条件为162110xf，初始取值为[10101010]M文件（求解函数）如下：functions=newton1(f,c,eps)%c是初值，eps为允许误差值ifnargin==2eps=1.0e-16;elseifnargin1error('')%returnendsymsx1x2x3x4x=[x1x2x3x4].';grad=jacobian(f).';hesse=jacobian(grad);a=grad;b=hesse;i=1;gradk=subs(a,[x1x2x3x4],[c(1)c(2)c(3)c(4)]);hessek=subs(b,[x1x2x3x4],[c(1)c(2)c(3)c(4)]);pk=-1*(hessek\gradk);x=tihuan(c);whilenorm(gradk)=epsx=x+pk;gradk=subs(a,[x1x2x3x4],[x(1)x(2)x(3)x(4)]);hessek=subs(b,[x1x2x3x4],[x(1)x(2)x(3)x(4)]);pk=-hessek\gradk;i=i+1;enddisp('thetimesofiterationis:')disp(i)disp('Thegradis:')disp(gradk)disp('andtheresultis:')x=x.';disp(x)return“tihuan”子函数:functionx=tihuan(x)x(1)=x(1);x(2)=x(2);x(3)=x(3);x(4)=x(4);end调用方式如下：symsx1x2x3x4f=(x1+10*x2)^2+5*(x3-x4)^2+(x2-2*x3)^4+10*(x1-x4)^4;c=[10101010]';%初始值newton1(f,c,eps);1.1.3计算结果如下：由上述结果可知，当迭代次数达到47次时满足终止条件，此时x为1.0e-005*[-0.11110.01110.00950.0095],显然，此题的理论解为[0000],分析上述结果，与理论解的误差处于可接受范围之内。求解完成。1.2用共轭梯度法求解函数441432243221102510minxxxxxxxxxf用共轭梯度法求解上述函数的程序代码如下：1.2.1问题分析：取00xfp，当搜索到1kx时，共轭方向2,...,1,0,11nkpxfpkkkk，此时，1kp与kpA共轭，用kAp右乘上式得kTkkkkkkAppApxfApp11，由01kTkApp得2,...,1,01nkAppApxfkTpkTkk，若不满足条件，进行下一次迭代。1.1.2问题求解注：程序所用语言为MATLAB，精度为1610epssymsx1x2x3x4t0t1f=(x1+10*x2)^2+5*(x3-x4)^2+(x2-2*x3)^4+10*(x1-x4)^4;c=[10;10;10;10];grad1=diff(f,x1);grad2=diff(f,x2);grad3=diff(f,x3);grad4=diff(f,x4);grad=[grad1;grad2;grad3;grad4];a=grad;i=1;n=40;gradk=subs(a,[x1x2x3x4],[c(1)c(2)c(3)c(4)]);x=tihuan(c);p0=0;whilenorm(gradk)=epsp0=-gradk;y=x;x=x+t0*p0;k=0;gradk=subs(a,[x1x2x3x4],[x(1)x(2)x(3)x(4)]);w=solve(gradk(1)+gradk(2)+gradk(3)+gradk(4));t0=real(w);t0=eval(t0);t0=t0(1);x=y+t0*p0;gradk=subs(a,[x1x2x3x4],[x(1)x(2)x(3)x(4)]);whilenorm(gradk)=epsifk+1~=ngradk2=subs(a,[x1x2x3x4],[x(1)x(2)x(3)x(4)]);gradk1=subs(a,[x1x2x3x4],[y(1)y(2)y(3)y(4)]);lamda=norm(gradk2).^2/norm(gradk1).^2;p0=-gradk2+lamda*p0;k=k+1;elsek=0;p0=-subs(a,[x1x2x3x4],[x(1)x(2)x(3)x(4)]);endcleary;y=x;x=x+t1*p0;gradk=subs(f,[x1x2x3x4],[x(1)x(2)x(3)x(4)]);m=solve(gradk);t1=real(m);t1=eval(t1(1));x=x+t1*p0;x=eval(x);clearm;cleart1;symst1gradk=subs(a,[x1x2x3x4],[x(1)x(2)x(3)x(4)]);enddisp(x.')return;enddisp(x.')此程序为一初步程序，假设初值为[10;10;10;10]，则第一次运算得t0=0.0011,lamda=0.9291,迭代后的x=NaN。现用共轭梯度法求解另一函数222125minxxxf对上述程序稍加改动来求解本题的代码如下：注：程序所用语言为MATLAB，精度为1610epsfunctions=gongegrad2(f,c,eps)%c是初值，eps为允许误差值ifnargin==2%eps=1.0e-16;elseifnargin1error('')returnendticsymsx1x2t0t1grad1=diff(f,x1);grad2=diff(f,x2);grad=[grad1;grad2];a=grad;i=1;n=40;gradk=subs(a,[x1x2],[c(1)c(2)]);x=tihuan(c);p0=0;whilenorm(gradk)=epsp0=-gradk;y=x;x=x+t0*p0;k=0;gradk=subs(f,[x1x2],[x(1)x(2)]);w=solve(gradk);t0=real(w);t0=eval(t0);t0=t0(1);x=y+t0*p0;gradk=subs(a,[x1x2],[x(1)x(2)]);whilenorm(gradk)=epsifk+1~=ngradk2=subs(a,[x1x2],[x(1)x(2)]);gradk1=subs(a,[x1x2],[y(1)y(2)]);lamda=norm(gradk2)^2/norm(gradk1)^2;p0=-gradk2+lamda*p0;k=k+1;elsek=0;p0=-subs(a,[x1x2],[x(1)x(2)]);endcleary;y=x;x=x+t1*p0;gradk=subs(f,[x1x2],[x(1)x(2)]);m=solve(gradk);t1=real(m);t1=eval(t1(1));x=y+t1*p0;clearm;cleart1;symst1gradk=subs(a,[x1x2],[x(1)x(2)]);enddisp(sprintf('thelastpointwewantis[%f%f]',x(1),x(2)));disp(sprintf('thetimesusedtorecursionis%f',k));disp(sprintf('thefunctionvalueis%f',x(1)^2+25*x(2)^2));tocreturn;enddisp(sprintf('thelastpointwewantis[%f%f]',x(1),x(2)));disp(sprintf('thetimesusedtorecursionis%f',k));disp(sprintf('thefunctionvalueis%f',x(1)^2+25*x(2)^2));toc“tihuan”子函数为：functionx=tihuan(x)[v,g]=size(x);fori=1:vx(i)=x(i);end程序调用方式为：clearallclcsymsx1x2t0t1f=x1^2+25*x2^2;c=[2;2];%初值gongegrad2(f,c,eps)程序结果如下：由上述结果知，用共轭梯度法对上述函数求解需要迭代三次得到最优解0，此时x为[00];符合理论分析的结果，求解完成。2.1用外点法法求解函数01.min212231xxtsxx2.1.1问题分析外点法为序列无约束最优化方法，其基本思想为将条件函数吸收到目标函数中进行求解。其在数学上的直观理解是拉格朗日乘子法：miixgMxfMXT12,0minmin;min，MXT;min为总代价，xf为价格，miixgM12,0min为罚款。即在经济学上总代价为价格和罚款的和。此时mixgxgxgxgiiii~1000,0min22，当，当称miixgMxfMXT1;为增广目标函数，通常取0002xgxgxgxgiiii当当2.1.2问题求解两种方法求解程序如下：2.1.2.1注：程序所用语言为MATLAB，终止条件为1610xgi，程序中无约束优化部分通过求导实现。M文件如下：ticclc%c是初值，eps为允许误差值ifnargin==1eps=1.0e-16;elseifnargin1error