最优化方法大作业
最优化方法大作业答案
1.用薄钢板制造一体积5m 3，长度不小于4m ，无上盖的货箱，要求钢板耗量最小。
确定货箱的长x 1、宽x 2和高x 3。
试列出问题的数学模型。
解：min 32312122x x x x x x z ++= s.t 5321=x x x 41≥x 0,,321≥x x x2.将下面的线性规划问题表示为标准型并用单纯形法求解max f=x 1+2x 2+x 3s ．t ．2x 1+x 2-x 3≤2 -2x 1+x 2-5x 3≥-6 4x 1+x 2+x 3≤6 x i ≥0 i=1，2，3 解：先化标准形：Min 321x x x z -+=224321=+-+x x x x 6525321=++-x x x x646321=+++x x x x列成表格：121610011460105122001112-----可见此表已具备1°，2°，3°三个特点，可采用单纯形法。
首先从底行中选元素-1，由2/2,6/2,6/4最小者决定选第一行第一列的元素2，标以记号，迭代一次得121210231040116201002121211--------再从底行中选元素-2/3，和第二列正元素1/2，迭代一次得12123230210231040116201002121211-------再从底行中选元素-3，和第二列正元素2，迭代一次得4233410120280114042001112---再迭代一次得1023021062210231010213000421021013--选取最优解：01=x 42=x 23=x3. 试用DFP 变尺度法求解下列无约束优化问题。
min f （X ）=4（x 1-5）2+（x 2-6）2取初始点X=（8，9）T ，梯度精度ε=0.01。
解：取IH=0,初始点()TX 9,8=2221)6()5(4)(-+-=x x x f⎥⎦⎤⎢⎣⎡--=∇122408)(21x x x f⎪⎪⎭⎫⎝⎛=∇624)()0(xfTx f d )6,24()()0()0(--=-∇=)0(0)0()1(dxxα+=T)69,248(00αα--=])669()5248(4min[)(min 2020)0(0)0(--+--⨯=+αααdxf 0)6()63(2)24()2458(8)(00)0(0)0(=-⨯-+-⨯--=+ααααd d xdf13077.0130170≈=α⎪⎪⎭⎫⎝⎛=⎪⎪⎭⎫ ⎝⎛--⨯+⎪⎪⎭⎫ ⎝⎛=21538.886153.462413077.098)1(x⎪⎪⎭⎫⎝⎛-=∇43077.410784.1)()1(xf进行第二次迭代：⎥⎦⎤⎢⎣⎡--=-=78463.013848.31)0()1(xxδ⎥⎦⎤⎢⎣⎡--=∇-∇=56924.110783.25)()(1)0()1(xf xf γ101011011101γγγγγδδδH HH H H TTTT-+=03172.8011=γδT86614.6321101==γγγγH T⎥⎦⎤⎢⎣⎡=61561.046249.246249.285005.911Tδδ⎥⎦⎤⎢⎣⎡==46249.240022.3940022.3940363.630110110TTHH γγγγ所以：⎪⎪⎭⎫⎝⎛--=0038.103149.003149.012695.01H⎪⎪⎭⎫⎝⎛-⨯⎪⎪⎭⎫⎝⎛---=∇-=43076.410784.10038.103149.003149.012695.0)()1(1)1(xf H d⎪⎪⎭⎫⎝⎛-=48248.428018.0令 )1(1)1()2(dx x α+=利用)()1()1(=+ααd dxdf ，求得49423.01=α，所以⎪⎪⎭⎫⎝⎛-+⎪⎪⎭⎫⎝⎛=+=21538.213848.021538.886152.449423.0)1()1()2(dxx⎪⎪⎭⎫ ⎝⎛=65因)()2(=∇xf ，于是停，)2(x 即为最优解。
最优化方法大作业-算法源程序-0.618法、抛物线法、共轭梯度法
最优化方法程序作业专业：油气储运工程班级：姓名：学号：一、开发工具该应用程序采用的开发工具是Visual studio 2005，编程语言使用的是C#。
二、程序代码（1）0.618法和抛物线法求ψ(t)=sint，初始点t0=1①主程序：using System;using System.Data;using System.Configuration;using System.Web;using System.Web.Security;using System.Web.UI;using System.Web.UI.WebControls;using System.Web.UI.WebControls.WebParts;using System.Web.UI.HtmlControls;using System.Collections;public partial class_Default : System.Web.UI.Page{JiSuan JS = new JiSuan();protected void Button1_Click(object sender, EventArgs e){double xx0 = 0.01;//步长double tt0 = 1, pp = 2, qq = 0.5;ArrayList list1 = new ArrayList();list1 = (ArrayList)JS.SuoJian(xx0, tt0, pp, qq);//调用倍增半减函数double aa = double.Parse(list1[0].ToString());double bb = double.Parse(list1[1].ToString());txtShangxian.Text = bb.ToString();//在页面上显示极小区间txtXiaxian.Text = aa.ToString();ArrayList list2 = new ArrayList();list2 = (ArrayList)JS.JiXiao1(aa, bb);//调用0.618法函数double jixiao1 = double.Parse(list2[0].ToString());double fjixiao1 = double.Parse(list2[1].ToString());txtJixiao1.Text = jixiao1.ToString();//在页面上显示极小点txtFjixiao1.Text = fjixiao1.ToString();//在页面上显示极小点处的函数值ArrayList list3 = new ArrayList();list3 = (ArrayList)JS.JiXiao2(aa, bb);//调用抛物线法函数double jixiao2 = double.Parse(list3[0].ToString());double fjixiao2 = double.Parse(list3[1].ToString());txtJixiao2.Text = jixiao2.ToString();//在页面上显示极小点txtFjixiao2.Text = fjixiao2.ToString();//在页面上显示极小点处的函数值 }}②各个子函数的程序：using System;using System.Data;using System.Configuration;using System.Web;using System.Web.Security;using System.Web.UI;using System.Web.UI.WebControls;using System.Web.UI.WebControls.WebParts;using System.Web.UI.HtmlControls;using System.Collections;///<summary>/// JiSuan 的摘要说明///</summary>public class JiSuan{public JiSuan(){//// TODO: 在此处添加构造函数逻辑//}//目标函数public double f(double z){double zz = Math.Sin(z);return zz;}//倍增半减函数public ArrayList SuoJian(double x0, double t0, double p, double q){double t1;double f0, f1;double x1, t2;double f2;double a = 0, b = 0, c = 0;double fa, fc, fb;double t3, f3;bool biaozhi1 = true;//设置是否循环的标志bool biaozhi2 = true;t1 = t0 + x0;f0 = f(t0);//调用目标函数f1 = f(t1);if (f0 > f1){while (biaozhi1){x1 = p * x0;t2 = t1 + x1;f2 = f(t2);//调用目标函数if (f1 > f2){t0 = t1;t1 = t2;f0 = f1;f1 = f2;continue;}else{a = t0;c = t1;b = t2;fa = f0;fc = f1;fb = f2;break;}}}else{t3 = t1;f3 = f1;t1 = t0;f1 = f0;t0 = t3;f0 = f3;while (biaozhi2){x1 = q * x0;t2 = t1 - x1;f2 = f(t2);//调用目标函数if (f1 > f2){t0 = t1;t1 = t2;f0 = f1;f1 = f2;continue;}else{a = t2;c = t1;b = t0;fa = f2;fc = f1;fb = f0;break;}}}ArrayList list = new ArrayList();list.Add(a);list.Add(b);return list;}//0.618法函数public ArrayList JiXiao1(double a, double b) {double jd = 0.0001;//精度double p = 0.618;double t1, t2;double f1, f2;double jixiao=0;//极小点bool biaozhi1 = true;//设置是否循环标志bool biaozhi2 = true;while (biaozhi1){t1 = a + (1 - p)*(b - a);t2 = a + p * (b - a);f1 = f(t1);//调用目标函数f2 = f(t2);while (biaozhi2){if (f1 < f2){b = t2;t2 = t1;f2 = f1;t1 = a + (1 - p)*(b - a); f1 = f(t1);//调用目标函数if (Math.Abs(b - a) > jd) {continue;}else{biaozhi1 = false;break;}}else if (f1 == f2){a = t1;b = t2;if (Math.Abs(b - a) > jd) {break;}else{biaozhi1 = false;break;}}else if (f1 > f2){a = t1;t1 = t2;f1 = f2;t2 = a + p * (b - a);f2 = f(t2);//调用目标函数if (Math.Abs(b - a) > jd){continue;}else{biaozhi1 = false;break;}}}}jixiao = (a + b) / 2;double fjixiao = f(jixiao);//调用目标函数ArrayList list = new ArrayList();list.Add(jixiao);list.Add(fjixiao);return list;}//抛物线法函数public ArrayList JiXiao2(double a, double b){double jd = 0.0001;//精度double c = a + (b - a) / 2;//将c的值设为a、b点的中点double fa, fb, fc, c1, c2;double jixiao = 0;//极小点double ta, fta;double fac;bool biaozhi1 = true;//设置是否循环标志while (biaozhi1){fa = f(a);//调用目标函数fb = f(b);fc = f(c);c1 = (fb - fa) / (b - a);c2 = ((fc - fa) / (c - a) - c1) / (c - b);if (c2 == 0){jixiao = c;break;}else{ta = 0.5 * (a + b - (c1 / c2));fta = f(ta);//调用目标函数if (Math.Abs(b - a) <= jd){jixiao = ta;break;}else{if (fc > fta){if (c > ta){b = c;fb = fc;c = ta;fc = fta;continue;}else{a = c;fa = fc;c = ta;fc = fta;continue;}}else if (fc == fta){if (c > ta){a = ta;b = c;c = (c + ta) / 2;fa = fta;fb = fc;fc = f(c);//调用目标函数continue;}else if (c == ta){fac = f((a + c) / 2);//调用目标函数if (fac < fc){c = (a + c) / 2;fc = f(c);//调用目标函数b = ta;fb = fta;continue;}else{a = (a + c) / 2;fa = f(a);//调用目标函数continue;}}else if (c < ta){a = c;c = (c + ta) / 2;b = ta;fa = fc;fb = fta;fc = f(c);//调用目标函数continue;}}else if (fc < fta){if (c > ta){a = ta;fa = fta;continue;}else{b = ta;fb = fta;continue;}}}}}double fjixiao = f(jixiao);//调用目标函数ArrayList list = new ArrayList();list.Add(jixiao);list.Add(fjixiao);return list;}}（2）共轭梯度法求函数极小点：f(x)=1.5x12+0.5x22-x1x2-2x1；初始点为（-2，4）①主程序：using System;using System.Data;using System.Configuration;using System.Web;using System.Web.Security;using System.Web.UI;using System.Web.UI.WebControls;using System.Web.UI.WebControls.WebParts;using System.Web.UI.HtmlControls;using System.Collections;public partial class_Default : System.Web.UI.Page{JiSuan JS = new JiSuan();protected void Button1_Click(object sender, EventArgs e){//共轭梯度法求极小值double[] x0 = new double[2];//定义二维数组double[] x = new double[2];x0[0] = -2;x0[1] = 4;double x00 = -2;double x01 = 4;double jd = 0.0001;//精度bool biaozhi1 = true;//设置是否循环的标志bool biaozhi2 = true;double y;//定义关于x的函数值double[] g = new double[2];//定义函数在点x处的梯度值double[] p = new double[2];double ggm, ggo=0;double jixiao;double x1=0, x2=0, fy=0;int i;for (int j = 0; j < 2; j++){x[j] = x0[j];}while (biaozhi1){i = 0;while (biaozhi2){y = JS.f(x[0], x[1]);//调用目标函数g[0] = JS.g0(x[0], x[1]);//调用梯度函数，求在点x处的梯度值 g[1] = JS.g1(x[0], x[1]);ggm = g[0] * g[0] + g[1] * g[1];if (ggm <= jd){x1 = x[0];x2 = x[1];fy = y;biaozhi1 = false;break;}else{if (i == 0){p[0] = -g[0];p[1] = -g[1];}else{p[0] = -g[0] + (ggm / ggo) * p[0];p[1] = -g[1] + (ggm / ggo) * p[1];}//调用0.618法子函数，找出极小点jixiao = JS.JiXiao1(x[0], x[1], p[0], p[1], x00, x01); x[0] = x[0] + jixiao * p[0];x[1] = x[1] + jixiao * p[1];ggo = ggm;i++;if (i >= 2){break;}else{continue;}}}}txtx1.Text = x1.ToString();txtx2.Text = x2.ToString();txty.Text = fy.ToString();}}②子函数程序：using System;using System.Data;using System.Configuration;using System.Web;using System.Web.Security;using System.Web.UI;using System.Web.UI.WebControls;using System.Web.UI.WebControls.WebParts;using System.Web.UI.HtmlControls;using System.Collections;///<summary>/// JiSuan 的摘要说明///</summary>public class JiSuan{public JiSuan(){//// TODO: 在此处添加构造函数逻辑//}//目标函数public double f(double z1,double z2){double zz = 1.5 * z1 * z1 + 0.5 * z2 * z2 - z1 * z2 - 2 * z1;return zz;}//求梯度值函数public double g0(double z1, double z2){double zz = 3 * z1 - z2 - 2;return zz;}public double g1(double z1, double z2){double zz = z2 - z1;return zz;}//倍增半减函数public ArrayList SuoJian(double x0, double t0, double p, double q, double xx0, double xx1, double p0, double p1,double x00,double x01){double t1;double f0, f1;double x1, t2;double f2;double a = 0, b = 0, c = 0;double fa, fc, fb;double t3, f3;bool biaozhi1 = true;//设置是否循环的标志bool biaozhi2 = true;double x_0t1, x_1t1, x_0t2, x_1t2;x_0t1 = xx0 + t1 * p0;x_1t1 = xx1 + t1 * p1;f0 = f(x00, x01);//调用目标函数f1 = f(x_0t1, x_1t1);if (f0 > f1){while (biaozhi1){x1 = p * x0;t2 = t1 + x1;x_0t2 = x_0t1 + t2 * p0;x_1t2 = x_1t1 + t2 * p1;f2 = f(x_0t2, x_1t2);//调用目标函数if (f1 > f2){t0 = t1;t1 = t2;f0 = f1;f1 = f2;continue;}else{a = t0;c = t1;b = t2;fa = f0;fc = f1;fb = f2;break;}}}else{t3 = t1;f3 = f1;t1 = t0;f1 = f0;t0 = t3;while (biaozhi2){x1 = q * x0;t2 = t1 - x1;x_0t2 = x_0t1 + t2 * p0;x_1t2 = x_1t1 + t2 * p1;f2 = f(x_0t2, x_1t2);//调用目标函数if (f1 > f2){t0 = t1;t1 = t2;f0 = f1;f1 = f2;continue;}else{a = t2;c = t1;b = t0;fa = f2;fc = f1;fb = f0;break;}}}ArrayList list = new ArrayList();list.Add(a);list.Add(b);return list;}//0.618法函数public double JiXiao1(double x0, double x1, double p0, double p1,double x00,double x01){double jd = 0.0001;//精度double p = 0.618;double t1, t2;double jixiao = 0;//极小点bool biaozhi1 = true;//设置是否循环标志bool biaozhi2 = true;double xx0 = 0.01;//步长double tt0 = 0, pp = 2, qq = 0.5;double x_0t1, x_1t1, x_0t2, x_1t2;ArrayList list1 = new ArrayList();//调用倍增半减函数获取极小区间list1 = (ArrayList)SuoJian(xx0, tt0, pp, qq, x0, x1, p0, p1, x00, x01);double a = double.Parse(list1[0].ToString());double b = double.Parse(list1[1].ToString());while (biaozhi1){t1 = a + (1 - p) * (b - a);t2 = a + p * (b - a);x_0t1 = x00 + t1 * p0;x_1t1 = x01 + t1 * p1;x_0t2 = x_0t1 + t2 * p0;x_1t2 = x_1t1 + t2 * p1;f1 = f(x_0t1, x_1t1);//调用目标函数f2 = f(x_0t2, x_1t2);while (biaozhi2){if (f1 < f2){b = t2;t2 = t1;f2 = f1;t1 = a + (1 - p) * (b - a);x_0t1 = x00 + t1 * p0;x_1t1 = x01 + t1 * p1;f1 = f(x_0t1, x_1t1);//调用目标函数if (Math.Abs(b - a) > jd){continue;}else{biaozhi1 = false;break;}}else if (f1 == f2){a = t1;b = t2;if (Math.Abs(b - a) > jd){break;}else{biaozhi1 = false;break;}}else if (f1 > f2){a = t1;t1 = t2;f1 = f2;t2 = a + p * (b - a);x_0t2 = x_0t1 + t2 * p0;x_1t2 = x_1t1 + t2 * p1;f2 = f(x_0t2, x_1t2);//调用目标函数if (Math.Abs(b - a) > jd){continue;}else{biaozhi1 = false;break;}}}}jixiao = (a + b) / 2;return jixiao;}}三、程序运行界面及结果（1）0.618法和抛物线法求ψ(t)=sint，初始点t0=1（2）共轭梯度法求函数极小点：f(x)=1.5x12+0.5x22-x1x2-2x1；初始点为（-2，4）。
优化方法大作业1
优化方法大作业一、Wolfe-Powell法利用MATLAB软件编写，其中初始值x0=-5；其他参数按照已知条件来取。
当b分别等于8、9、10时，均得到如下结果：而当初值x0变化时，则结果变化比较大，如将x0取-6，则计算结果如下：通过比较可以看出，b的取值对计算结果的影响较小，而初始值x0的取值则对结果影响很大。
从中也表明Wolfe-Powell准则的收敛条件比较弱，容易出现当函数还没取极小值而迭代循环已结束的情况。
具体代码见附录1.二、无约束优化1.DFP法精度ε = 0.02，确定λ使用的是一维非精确搜索算法(直接法，Wolfe-Powell准则)。
结果如下图所示：2. BFGS方法同上一种方法，精度ε = 0. 02，确定λ使用的是一维非精确搜索算法(直接法，Wolfe-Powell准则)。
结果如下图所示：比较两种方法，查看每一步的函数值，得出如下结果。
图1：DFP与BFGS法结果对比图通过图1与表1的比较，BFGS比DFP法最初收敛得更快，但是DFP法比BFGS法的最终结果更好。
DFP与BFGS法的代码分别见附录2、3。
三、Rockafellar乘子法文件myfun.m：function y=myfun(x)y=1000-x(1)^2-2*x(2)^2-x(3)^2-x(1)*x(2)-x(1)*x(3);文件mycon.m：function [c,ceq]=mycon(x)c(1)=(-x(1));c(2)=(-x(2));c(3)=(-x(3));ceq(1)=x(1)^2+x(2)^2+x(3)^2-25;ceq(2)=8*x(1)+14*x(2)+7*x(3)-56文件Q3.m：clearclcx0=[2 2 2]';[x,fval,exitflag,output]=fmincon(@myfun,x0,[],[],[],[],[],[],@mycon);附录1：%% Wolfe-Powell 准则方法clear;tic;c1=0.1;c2=0.65;a=0;b=8;syms x;f=x*x-2*x+7; %函数G0=jacobian(f,x);%求梯度x0=-6; %初始取x=-5% 由于只有一个未知数，默认S=1lambda = 1;x=x0;for k=1:1:1000f0 = eval(f); %求出函数值grad1= eval(G0); %求出梯度值figure1(k) = f0; % 用于绘出迭代过程中函数值变化x= x0+lambda;f1 = eval(f); %求出函数值grad2 = eval(G0); %求出梯度值value1 = f0 - f1 + c1* lambda * grad1 ;value2 = grad2 -c2*grad1 ;if value1 >=-1e-6 && value2 >=-1e-6disp('The variable matrix is : ');disp(x);fun=eval(f);fprintf('The minimum value of function is %f \n',fun);break; %如果满足两个条件，则退出循环。
实用最优化方法编程大作业
实用最优化方法编程大作业班级：姓名：指导老师：学号：大连理工大学2015年11月27日版本号：MATLAB 7.11.0 (R2010b)【文件名WP.m】function x = WP(f,x0,var,s0,eps)clcsyms x1 x2 ;c1=0.1;c2=0.5;b=inf;lambda=1;ifnargin==4eps=1.0e-6;endgradf=jacobian(f,var);g0=subs(gradf,var,x0);f0=subs(f,var,x0);gs=g0*s0';a=0;j=0;while j<1000new_x=x0+lambda*s0;new_f=subs(f,var,new_x);left=f0-new_f;new_g=subs(gradf,var,new_x);new_gs=new_g*s0';right=-1*c1*lambda*gs;if left<right %不满足第一个条件j=j+1;b=lambda;lambda=0.5*(lambda+a);else %满足第一个条件left2=new_gs;right2=c2*g0;if left2<right2 %不满足第二个条件j=j+1;a=lambda;lambda=min(2*lambda,0.5*(lambda+b));elsex=lambda;break;endendend在Command Window 输入：symsx1 x2WP(100*(x2-x1^2)^2+(1-x1)^2,[-1,1],[x1 x2],[1 1])程序运行后可得出结果：ans=0.0039【文件名minGETD.m】function [x,minf] = minGETD(f,x0,var,eps)clcsyms x1 x2 x3format long;n=length(x0);ifnargin==3eps=1.0e-3;endx0=x0';syms lambda;gradf=jacobian(f,var);g=subs(gradf,var,x0);s=-g';k=0;while 1tol=norm(double(g));iftol<=epsx=x0;break;endx1=x0+lambda*s;f1=subs(f,var,x1);dy1=diff(f1);lambda0=solve(dy1);x1=x0+lambda0*s;g1=subs(gradf,var,x1)tol=norm(double(g1))iftol<=epsx=x1;break;endif k+1==nx0=x1;continue;elsemiu=dot(g1,g1)/dot(g,g)s=-g1'+miu*s;k=k+1;x0=x1;endendx在Command Window 输入：symsx1 x2 x3x0=[0 0 0];var=[x1 x2 x3];f=x1^2-2*x1*x2+2*x2^2+x3^2-x2*x3+2*x1+3*x2-x3;minGETD(f,x0,var,eps)程序运行后可得出结果：x=[-236894/59711,-178563/59711,-59465/59711]可认为最终解为[-4,-3,-1]。
最优化方法与应用大作业(一)最速下降法
最优化方法与应用大作业（一）
---最速下降法部分：
1.问题描述:
用梯度下降法求解以下优化问题
min f(x)=(x1+10*x2)^2+5(x3-x4)^2+(x2-2*x3)^4+10*(x1-x4)^4
2.编程感想：
该算法需要计算Hesse矩阵，C语言在向量运算时没有Matlab方便，所以手工完成了理论计算，再输入，破坏了程序的移植性。
同时，实验表明当初始值离理想点较远且精度要求较高时，最速下降法的收敛速率极慢，迭代几乎不可能完成，这对初值的选取提出了一定限制。
3.结果分析：
编译界面（Mac os X，Xcode环境）
输入参数（设定为(0.1,0.2,0.3,0.4)）：
结果（此处列出每次迭代结果）。
明显的看到，最速下降法的收敛较慢，最终结果接近理论值（F（0，0，0，0）=0）所以该结果可以满意。
4.算法代码见下页
西安电子科技大学电子工程学院020951
李骏昊02095005。
最优化方法大作业模板
命题人：审核人：大作业学期：至学年度第学期课程：最优化方法课程代号：签到序号：使用班级：姓名：学号：题号一二三四五六七八九十总分得分一、（目标1）请从以下6种算法中任选一种，说明算法的来源、定义、基本思想和优缺点，并给出算法步骤（包含算法流程图）和例子（包含程序与运算结果）。
①禁忌搜索算法；②模拟退火算法；③遗传算法；④神经网络算法；⑤粒子群算法；⑥蚁群算法。
二、（目标1）某工厂生产甲、乙两种产品，已知生产这两种产品需要消耗三种材料A 、B 和C ，其中生产过程中材料的单位产品消耗量和总量，以及单位产品的利润如下表所示。
该如何配置安排生产计划，使得工厂所获得的利润最大？材料甲乙资源总量材料A （Kg ）3265材料B （Kg ）2140材料C （Kg ）0375单位利润（元/件）15002500-（1）要保证工厂利润的最大化，写出相应的生产计划数学模型；（2）根据对偶理论，直接写出该线性规划的对偶问题；（3）采用单纯形表法对该该线性规划问题进行求解，写出详细的计算过程；（4）采用Matlab 软件对该线性规划问题进行求解，写出完整的源程序，并给出程序运行结果；（5）讨论当材料B 的资源总量发生变化时，该线性规划问题的最优解会如何变化？课程目标目标1……题号一、二、三、四、五……分值20、25、20、20、15……得分得分三、（目标1）求解下列无约束非线性规划问题（1）采用黄金分割法求解：min 4()24f x x x =++。
初始区间为[-1.0]，精度为ε=10-4。
（要求：采用黄金分割法进行Matlab 编程求解，写出源程序，并给出运行结果，列出迭代过程的数据表格）（2）采用阻尼牛顿法求解：222121212min (,)4f x x x x x x =+-。
分别取两个初始点：x A =(1,1)T ，x B =(3,4)T 。
（要求：采用阻尼牛顿法进行Matlab 编程求解，并给出运行结果，列出迭代过程的数据表格）四、（目标1）求解下列约束非线性规划问题：22112212121212min ()23532..00f x x x x x x x x x x x s t x x =-+--+≤⎧⎪-≤⎪⎨≥⎪⎪≥⎩（1）采用罚函数法进行求解，需写出具体计算过程；（2）采用二次规划方法进行求解，需写出具体计算过程，并进行MATLAB 编程，写出源程序和运算结果；五、（目标1）（1）某商店在未来的4个月里，准备利用它的一个仓库来专门经营某种商品，仓库的最大容量为1000单位，而且该商店每月只能出卖仓库现有的货。
最优化方法课程设计-最优化大作业-用优化算法求解函数最值问题
最优化方法大作业---------用优化算法求解函数最值问题摘要最优化(optimization) 是应用数学的重要研究领域.它是研究在给定约束之下如何寻求某些因素(的量)，以使某一(或某些)指标达到最优的一些学科的总称。
最优化问题一般包括最小化问题和最大化问题，而最大化问题可以通过简单的转化使之成最最小化问题。
最小化问题分为两类，即约束最小化和无约束最小化问题。
在此报告中，前两个问题属于无约束最小化问题的求解，报告中分别使用了“牛顿法”和“共轭梯度法”。
后两个问题属于有约束最小化问题的求解，报告中分别用“外点法”和“内点法”求解。
虽然命名不一样，其实质都是构造“惩罚函数”或者“障碍函数”，通过拉格朗日乘子法将有约束问题转化为无约束问题进行求解。
再此报告中，“外点法”和“内点法”分别用了直接求导和调用“牛顿法”来求解无约束优化问题。
在此实验中，用“共轭梯度法”对“牛顿法”所解函数进行求解时出现错误，报告中另取一函数用“共轭梯度法”求解得到正确的结果。
此实验中所有的函数其理论值都是显见的，分析计算结果可知程序正确，所求结果误差处于可接受范围内。
报告中对所用到的四种方法在其使用以前都有理论说明，对“外点法”中惩罚函数和“内点法”中障碍函数的选择也有相应的说明，另外，对此次试验中的收获也在报告的三部分给出。
本报告中所用程序代码一律用MATLAB编写。
【关键字】函数最优化牛顿法共轭梯度法内点法外点法 MATLAB一，问题描述1，分别用共轭梯度发法和牛顿法来求解一下优化问题()()()()()441432243221102510min x x x x x x x x x f -+-+-++=2， 分别用外点法和内点发求解一下优化问题⎩⎨⎧≥-++01.min 212231x x t s x x二、问题求解1.1 用牛顿法求解()()()()()441432243221102510min x x x x x x x x x f -+-+-++=1.1.1问题分析：取步长为1而沿着牛顿方向迭代的方法称为牛顿法，在牛顿法中，初始点的取值随意，在以后的每次迭代中，()[]()k k k k x f x f x x ∇∇-=-+121，直到终止条件成立时停止。
最优化方法大作业
单位代码03学号《最优化方法》课程实践完成时间：2015年5月30日星期六选择题目：题目一使用优化软件，编写重要算法的程序1.第一大题：（1）学习最优流量工程问题，nonsmooth_MCFP.pdf（2）问题重述：Figure 1一个简单的网络拓扑和流量需求如Figure 1所示，网络有7 个节点，13 条弧，每条弧的容量是5 个单位. 此外有四个需求量均为4个单位的源－目的对（M=4)，具体的源节点、目的节点信息如图所示. 这里为了简单，省去了未用到的弧，此外弧上的数字表示弧的编号。
（3）极小化MAU设定变量x，为531⨯的向量，其中(53)x即为变量z。
使用linprog 函数求解极小化问题得到x。
之前确定三个约束条件。
1、Ax b≤，其中A为1353⨯的矩阵，b为131⨯的向量。
2、eq eq x b A =，其中eq A 为2853⨯的矩阵，eq b 为281⨯的向量。
3、x lb ≥，其中lb 为153⨯的向量 编程计算后得到结果如下：（4） 极小化FT 成本函数设定变量x ，为651⨯的向量，其中(53:65)x 即为变量l z 。
使用linprog 函数求解极小化问题得到x 。
之前确定三个约束条件。
1、Ax b ≤，其中A 为7865⨯的矩阵，b 为781⨯的向量。
2、eq eq x b A =，其中eq A 为2865⨯的矩阵，eq b 为281⨯的向量。
3、x lb ≥，其中lb 为165⨯的向量 编程计算后得到结果如下：2. 第二大题： 2.1. 习题5.6 2.1.1. 问题分析问题2112212()(101810)/241513q x x x x x x x =-++-+ 通过matlab 画出其等高线为：2.1.2. 最速下降法最速下降法中，取值：k k p g =-==()()k k k kk T k k T k g p g g p Gp g Ggα- x1x 2等高线-224681012(1)()k x k x k p α+=+2.1.3. 算法流程图如下图所示：2.1.4. 初始值（0,0）编程运行结构为：收敛过程曲线为：2.1.5. 初始值（-0.4,0）编程运行结构为：收敛过程曲线为：x1x 2等高线-2246810122.1.6. 初始值（10,0）编程运行结构为：收敛过程曲线为：x1x 2-2246810122.1.7. 初始值（11,0）编程运行结构为：收敛过程曲线为：x1x 2-2246810122.2. 习题5.7 2.2.1. 问题分析问题()94ln(7)f x x x =--497g x =-- 24(7)G x =- Matlab 画出在区间(7 10)的函数、一阶导数、二阶导数的变化曲线为x1x 2-22468101277.588.599.510707274767880828486xf函数变化曲线77.588.599.510-35-30-25-20-15-10-50510xg一阶导数g 变化曲线2.2.2. 牛顿法牛顿法中，取值：k k k G s g =-1k k k s xx +=+其中，如果G 不是半正定，则采用修正牛顿法(+)k k k G I s g λ=-77.588.599.51050100150200250300350400xg二阶导数G 变化曲线2.2.3.算法流程图如下图所示：2.2.4.初始值7.40编程运行结构为：收敛过程曲线为：2.2.5. 初始值7.20编程运行结构为：收敛过程曲线为：7.397.47.417.427.437.447.457.4670.24570.2570.25570.2670.265xf2.2.6. 初始值7.01编程运行结构为：收敛过程曲线为：7.17.27.37.47.57.67.770.270.470.670.87171.271.471.6xf2.2.7. 初始值7.80编程运行结构为：收敛过程曲线为：6.856.9 6.9577.057.17.157.27.257.37.35727476788082xf2.2.8. 初始值7.88编程运行结构为：收敛过程曲线为：7.17.27.37.47.57.67.77.87.97070.57171.572xf2.2.9. 分析函数在区间（7，7.8888）内是凸函数，G 恒大于零，所以单纯牛顿法保证收敛。
西电最优化上机报告(大作业)
上机报告一.最速下降法算法简述：1.在本例中，先将最速下降方向变量赋一个值，使其二范数满足大于ε的迭代条件，进入循环。
2.将函数的一阶导数化简，存在一个矩阵，将其hesse矩阵存在另一个矩阵。
依照公式求出α，进而求出下一任迭代的矩阵初值。
循环内设置一个计数功能的变量，统计迭代次数。
3.求其方向导数的二范数，进行判别，若小于ε，则跳出循环，否则将继续迭代。
4.显示最优解，终止条件，最小函数值。
心得体会：最速下降法的精髓，无疑是求梯度，然后利用梯度和hesse矩阵综合计算，求解下一个当前最优解。
但是，要求函数是严格的凸函数，结合严格凸函数的大致图像，这就给初值的选取提供了一点参考。
例如在本例中，由于含有两个变量的二次方之和，结合大致图像，想当然的，初值的选取应当在原点附近；又因为变量的二次方之和后面，还减去了变量的一次形式和一次混合积，所以初值的选取应该再向第一象限倾斜。
综合以上考量，第一次选取（1,1）作为初值，判别精度方面，取到千分位，暂定为0.001。
运行以后，结果显示迭代了25次，最优解为（3.9995，1.9996），终止条件为5.4592e-04，目标函数为-8.0000。
这个结果已经相当接近笔算结果。
整体的运行也比较流畅，运算速度也比较快。
第二次取值，决定保留判别精度不变，将初值再适当向第一象限倾斜，取（2,2）作为初值，运行后，显示只迭代了11次！最优结果显示（3.9996，1.9997），终止条件为3.6204e-04，最优解-8.0000。
可见，最优结果更接近理想值，终止条件也变小了，最关键的是，迭代次数减少至第一次的一半以下！这说明以上初选取的方向是对的！第三次再进行初值细化，判别精度仍然不变，初值取（3,3）。
结果令人兴奋，只迭代了四次！最优解已经显示为（4.0000,2.0000），终止条件为2.4952e-04，目标函数-8.0000。
第四次，判别精度不变，取初值（4,4）。
最优化理论和算法： 大作业(一)
最优化理论和算法:大作业（一）简介：这个大作业的主要目的是在Matlab下自己编写单纯形法算法来求解标准型线性规划问题：min c T xs.t.Ax=bx≥0其中b≥0,A是m×n(m≤n)的矩阵。
假设A的秩是m.特别的，A并不一定包含单位矩阵。
按照要求编写下列小程序。
程序(一)：实现单步单纯形法程序格式：function[istatus,iB,iN,xB]=simplex_step(A,b,c,iB,iN,xB)%实现一步单纯形法%输入参数：%A-(n,m)系数矩阵%ｂ-(ｍ,１)(正)右端向量%ｃ-(ｎ,１)目标函数系数%iB-(1,m)整数向量记录当前基本变量的指标数%iN-(1,n-m)整数向量记录当前非基本变量的指标数%xB-(m,1)向量代表当前基本变量的值%输出参数：%istatus-整数标记单纯形法执行状态%istatus=0正常单纯形法步完成% istatus=32问题无界% istatus=-１找到最优基本可行解%iB-(1,m)整数向量记录运行单纯形法之后的基本变量的指标数%iN-(1,n-m)整数向量记录运行单纯形法之后的非基本变量的指标数%xB-(m,1)向量记录运行单纯形法之后的基本变量的值注：该程序不考虑退化情形。
程序(二)：利用两步法中的第一步来求解一个初始基本可行解程序格式:function[istatus,iB,iN,xB]=simplex_init(A,b)%实现两步法中的第一步来求解一个初始基本可行解，通过求解下面的问题：%min y_1+...+y_m%s.t.Ax+y=b%x>=0,y>=0%A是m x n矩阵。
%输入参数：%A-(n,m)系数矩阵%ｂ-(ｍ,１)正的右端向量%输出参数：%istatus-整数标记初始化状态% istatus=1找到原问题的一个基本可行解% istatus=4问题可行域是空集% istatus=16初始化过程失败%iB-(1,m)整数向量记录运行初始化之后的基本变量的指标数（对应原问题）%iN-(1,n-m)整数向量记录运行初始化之后的非基本变量的指标数（对应原问题）%xB-(m,1)向量记录运行初始化之后的基本变量的值（对应原问题）注：为了简单化程序，若初始化过程找到的初始基本可行解包含某些人工变量y j,设置istatus=16（初始化失败）。
最优化方法大作业
发动机空燃比控制器引言：我主要从事自动化相关研究。
这里介绍我曾经接触过的发动机空燃比控制器设计中的优化问题。
发动机空燃比控制器设计中的最优化问题AFR =afm m && (1)空燃比由方程（1）定义，在发动机运行过程中如果控制AFR 稳定在14.7可以获得最好的动力性能和排放性能。
如果假设进入气缸的空气流量am &可以由相关单元检测得到，则可以通过控制进入气缸的燃油流量f m &来实现空燃比的精确控制。
由于实际发动机的燃油喷嘴并不是直接对气缸喷燃油，而是通过进气歧管喷燃油，这么做会在进气歧管壁上液化形成油膜，因此不仅是喷嘴喷出的未液化部分燃油会进入气缸，油膜蒸发部分燃油也会进入气缸，如方程（2）。
这样如何更好的喷射燃油成为了一个问题。
1110101122211ττττ⎡⎤⎡⎤-⎢⎥⎢⎥⎡⎤⎡⎤⎢⎥⎢⎥=+⎢⎥⎢⎥⎢⎥⎢⎥-⎢⎥⎢⎥⎣⎦⎣⎦⎢⎥⎢⎥-⎢⎥⎢⎥⎣⎦⎣⎦⎡⎤⎡⎤⎢⎥⎣⎦⎢⎥⎣⎦f ff v X x x u x x X x y =x && (2)其中12、,==ff fv x m x m &&=f y m &,=fi u m &这里面，表示油膜蒸发量ff m &、fvm &表示为液化部分燃油、fim &表示喷嘴喷射的燃油，在τf 、τv 、X 都已知的情况下，由现代控制理论知识，根据系统的增广状态空间模型方程（3）00000011011011114.70ττττ⎡⎤⎡⎤-⎡⎤⎢⎥⎢⎥⎡⎤⎡⎤⎢⎥=-+-⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦⎣⎦⎢⎥⎢⎥⎢⎥⎣⎦⎢⎥⎢⎥⎣⎦⎣⎦⎡⎤⎡⎤=⎢⎥⎣⎦⎣⎦ff v v a X X u +q q m y q x x x &&& (3)其中()014.7⎰taq =y -m&。
由极点配置方法，只要设计控制器方程（4），就可以使得y 无差的跟踪阶跃输入，那么y 也能较好的跟踪AFR *am /&。
最优化大作业(一)终稿
西安电子科技大学电子工程学院智能科学与技术专业最优化大作业（一）内容：1.最速下降法2.牛顿法班级：020951姓名：张普照学号：02095076日期：2011.11.24一、心得体会1.最速下降法优点；算法简单，每次迭代计算量小，内存占用小，即使从一个不好的初始点出发，往往也能收敛到局部极小点。
缺点：一个严重缺点就是收敛速度慢，迭代次数多。
特别是对于等值线（面）具有狭长深谷形状的函数，收敛速度更慢，其原因是产生锯齿现象。
2.牛顿法优点：收敛速度非常快，迭代次数相应就很少，这可以从实验结果中看出，而且具有二次收敛的优点。
缺点：当Hesse矩阵非正定时，Newton法的搜索将会失败；对初始点要求严格，一般要求比较接近或有利于接近极值点，而这在实验计算中是比较难办的；在进行某次迭代时可能求不出搜索方向；Newton方向构造困难，计算相当复杂，除了求梯度以外还需计算Hesse矩阵及其逆矩阵，占用机器内存相当大。
二、结果目标函数为：f=(x1+10*x2)^2+5*(x3-x4)^2+(x2-2*x3)^4+10*(x1-x4)^4;初始点：X0=[100,100,100,100]’;终止限:0.000001在结果中：X表示所求得的最优解，F为取最优解释的目标函数值，k为总的迭代次数。
1.最速下降法结果>> fastest_drop2()X =0.0048-0.00050.00240.0024F =1.0835e-009k=3976222.牛顿法结果>>Newton2()X =-0.00200.00020.00020.0002F = 2.0774e-010k = 34三、源代码1.最速下降法%********************主函数**************************** function fastest_drop2()X0=[100,100,100,100]';k=0;f0=fun(X0);g0=grad(X0);H0=hesse(X0);while(1)X=ls(X0,g0,H0);F=fun(X);G=grad(X);H=hesse(X);if(norm(G)<0.000001)XFkbreak;elseX0=X;f0=F;g0=G;H0=H;k=k+1;endend%待求解原函数function f=fun(x)f=(x(1,1)+10*x(2,1))^2+5*(x(3,1)-x(4,1))^2+(x(2,1)-2*x(3,1))^4+10*(x(1,1)-x(4,1))^4;end%梯度函数function g=grad(x)g=[2*x(1,1)+20*x(2,1)+40*(x(1,1)-x(4,1))^3;20*x(1,1)+200*x(2,1)+4*(x(2,1)-2*x(3,1))^3;10*x(3,1)-10*x(4,1)-8*(x(2,1)-2*x(3,1))^3;-10*x(3,1)+10*x(4,1)-40*(x(1,1)-x(4,1))^3;]; end%直线搜索function ls=ls(X,g,H)ls=X-(g'*g)/(g'*H*g)*g;end%海森矩阵function h=hesse(x)h=[ 2+120*(x(1,1)-x(4,1))^2, 20, 0,-120*(x(1,1)-x(4,1))^2; 20, 200+12*(x(2,1)-2*x(3,1))^2, -24*(x(2,1)-2*x(3,1))^2, 0; 0, -24*(x(2,1)-2*x(3,1))^2, 10+48*(x(2,1)-2*x(3,1))^2, -10; -120*(x(1,1)-x(4,1))^2, 0, -10, 10+120*(x(1,1)-x(4,1))^2];e ndend%***************最速下降法代码结束*************2.牛顿法%**********牛顿法主函数****************function Newton2()X0=[100,100,100,100]';k=0;f0=fun(X0);g0=grad(X0);while(1)H=hesse(X0);P=-inv(H)*g0;X=X0+P;F=fun(X);G=grad(X);if(norm(G)<0.000001)XFkbreak;elseX0=X;f0=F;g0=G;k=k+1;endendend%**************目标函数************************function f=fun(x)x1=x(1,1); x2=x(2,1); x3=x(3,1); x4=x(4,1);f=(x1+10*x2)^2+5*(x3-x4)^2+(x2-2*x3)^4+10*(x1-x4)^4;end%************梯度矩阵******************************* function g=grad(x)x1=x(1,1); x2=x(2,1); x3=x(3,1); x4=x(4,1);g =[ 2*x1+20*x2+40*(x1-x4)^3;20*x1+200*x2+4*(x2-2*x3)^3;10*x3-10*x4-8*(x2-2*x3)^3;-10*x3+10*x4-40*(x1-x4)^3];end%*************海森矩阵****************************** function h=hesse(x)x1=x(1,1); x2=x(2,1); x3=x(3,1); x4=x(4,1);h = [ 2+120*(x1-x4)^2, 20, 0, -120*(x1-x4)^2;20, 200+12*(x2-2*x3)^2, -24*(x2-2*x3)^2, 0;0, -24*(x2-2*x3)^2, 10+48*(x2-2*x3)^2, -10;-120*(x1-x4)^2, 0, -10, 10+120*(x1-x4)^2];end%***********牛顿法代码结束***************。
北航最优化方法大作业参考
北航最优化方法大作业参考旅行商问题是一个经典的组合优化问题，目标是找到一条最短路径，使得旅行商能够在访问所有城市后回到起始城市。
在实际应用中，旅行商问题有着广泛的应用，例如物流配送、城市规划等领域。
为了解决旅行商问题，我们可以采用启发式算法，其中一个常用的方法是遗传算法。
遗传算法是一种模拟自然进化过程的优化算法，通过模拟生物遗传的选择、交叉和变异等操作，逐步优化问题的解。
首先，我们需要对问题进行建模。
假设有N个城市，我们可以通过一个N*N的距离矩阵来表示各个城市之间的距离。
同时，我们需要定义一个染色体表示一条路径，其中每个基因表示一个城市的编号。
接下来，我们可以采用遗传算法来求解最优解。
遗传算法一般包括以下几个步骤：1.初始化种群：随机生成初始的染色体种群，每个染色体都表示一条路径。
2.适应度评价：根据染色体的路径长度来评估每个染色体的适应度，路径越短适应度越高。
3.选择操作：选择适应度较高的染色体作为父代，采用轮盘赌选择算法确定父代。
4.交叉操作：采用部分映射交叉算子对父代进行交叉操作，生成新的子代。
5.变异操作：对子代进行变异操作，以增加种群的多样性。
6.环境选择：根据适应度选择下一代种群，同时保留精英个体，避免解的丢失。
7.终止条件：当达到预设的迭代次数或者达到最优解时，终止算法。
通过以上步骤的迭代，我们可以逐步优化路径的长度，最终得到一条最短路径。
除了遗传算法，我们还可以尝试其他的优化算法，例如模拟退火算法、蚁群算法等。
这些算法在求解旅行商问题时都有一定的优势和适用性。
总结起来，旅行商问题是一个经典的组合优化问题，在北航最优化方法大作业中可以选择使用启发式算法来解决。
我们可以尝试使用遗传算法来求解最优路径，并根据实际情况选择合适的算法参数和终止条件。
通过不断地迭代和优化，我们可以得到一条最短路径，满足旅行商的需求。
以上是关于北航最优化方法大作业的参考内容，希望对你的写作有所帮助。
如果有其他疑问，欢迎继续提问。
最优化方法-习题
( x2 x1) f ( x 2) f ( x2) f ( x1)
三 、 设 f(x)=
2 1
“充分性” 由 则
x ，x
1
2
的任意性取
x = x 时，f( x )>f( x )
1
T 1 T T Qx b x c, Q Q 0 试 证 ： 共 轭 梯 度 法 的 线 性 搜 索 中 2x
3
*
T
f ( x)
x
2
1
2 x 2 2 x1 x2 4 x1 , x R
T
1、 给定问题
x
(1)
(1,1)
解：1）DFC 法 取初始对称矩阵
1 0 H1 0 1
第一次迭代：
x x x 2 x 6 x 14 x x x x 2 s.t. x 2 x 3 x 0, x 0, x 0 取初始点 x (1 ,1,0) ，用简约梯度法求其最优解
T
2) 因为 Q=
2 2 2 2 ,所以 |Q|= 2 6 =8>0 即可知 Q 是非奇异的 2 6 2 2 2 6
=8>0 ,所以 Q 是正定的,故 f(x)是正定的
f ( x1) (2,4,5)
3) 因为|2|>0,
2 d f ( x1) =(1,0,-1) 4 = -3<0 5
2 2
2
一、设 f(x)为定义在区间[a,b]上的实值函数，
2) f(x)=ln(
x1 + x x x 2 )
x1 x2 9 x3 , x1 6 x2 x3 2 , 9 x1 x2 )
西电最优化大作业
p0 f x0 ， 当 搜 索 到
xk 1
时
，
共
轭
方
向
为
pk 1 f xk 1 k pk , k 0,1,...,n 2 ，此时， pk 1 与 pk A 共轭，用 Apk 右乘上式得
T pk 1 Apk f xk 1 Apk k pk Apk
m 2 理解为拉格朗日乘子法： minT X ; M min f x M min0, g i x i 1
其中
min0, g i x 2
g i 0，i 1 ~ m
当g i x 0，
，
由
T pk 1 Apk 0
得
f x k 1 Ap k k k 0,1,...,n 2 ，若不满足条件，进行下一次迭代。 pT p Ap k
T
1.2.2 问题求解 本程序编程语言为 MATLAB，终止条件为 f x k x0 =[1 1]。 程序代码见附件conjugate.m 1.2.3 计算结果如下：
三、此次实验的收获
经过几个晚上的艰苦奋斗，努力学习，不断调试程序，最终才得以成功运行 程序并得到满意的结果。 有过山重水复疑无路的困境，但最终迎来的还是柳暗花 明又一村的喜悦。通过此次实验，我的收获主要有以下几点： 以前自己在求解函数最优化问题时都是通过求导、画图等方法，而这几个算法都 是通过不断迭代寻找最优解， 相对来说更有利于电脑编程的实现和推广，对函数 本身性质的要求也不是太高。 在最速下降法中， 沿负梯度方向函数值很快的说法容易使我们相信这一定是 最理想的搜索方向， 沿该方向搜索时收敛速度应该很快，然而在算法实现过程中 发现，梯度法的收敛速度并不快（迭代次数为 45 次） ，比其它算法收敛速度都要 慢。而共轭梯度法仅需要利用一阶导数信息，也不再要求精确的直线搜索，进而 克服了最速下降法收敛慢的缺点（迭代次数为 2 次） 。 内点法和外点法的实质都是构造 “惩罚函数” 或者 “围墙函数 （或障碍函数） ” ， 将约束函数其转化为非约束函数， 其中外点法在将函数转化为非约束函数后调用 了“牛顿法”来求解，内点法也尝试过用“牛顿法”来求解非约束函数，但由于 障碍函数为倒数形式， 导致了程序在求矩阵逆的时候产生了无穷大的量，函数无 解，所以内点法采用“直接求导”来求解非约束函数。此外，我也尝试了用求极 限的方法来求解最优点，根据手算的步骤，在求得偏导数为 0 的点后，令障碍因 子 mk→0，求得最优解，方法简单，算法易于实现。 这门课在学习过程中多以理论学习为主，如果平时不注重实践，将自己所学 运用到实际中， 就很难会有太大的收获，通过自己的努力用课堂上的算法解决了 实际问题，无疑加深了自己对算法理论的理解。
最优化方法练习题答案
练习题一1、建立优化模型应考虑哪些要素“答：决策变量、目标函数和约束条件。
2、讨论优化模型最优解的存在性、迭代算法的收敛性及停顿准则。
答：针对一般优化模型，讨论解的可行域，假设存在一()()min ()..0,1,2, 0,1,,i j f x s t g x i m h x j p≥===L L D 点，对于均有则称为优化模型最优解，最优解存在；*X D ∈X D ∀∈*()()f X f X ≤*X 迭代算法的收敛性是指迭代所得到的序列，满足，(1)(2)(),,,K X X X L L (1)()()()K K f X f X +≤则迭代法收敛；收敛的停顿准则有，，(1)()k k x x ε+-<(1)()()k k k x x xε+-<，，等等。
()()(1)()k k f x f x ε+-<()()()(1)()()k k k f x f x f x ε+-<()()k f x ε∇<练习题二1、*公司看中了例2.1中厂家所拥有的3种资源R 1、R2、和R 3，欲出价收购〔可能用于生产附加值更高的产品〕。
如果你是该公司的决策者，对这3种资源的收购报价是多少？〔该问题称为例2.1的对偶问题〕。
解：确定决策变量对3种资源报价作为本问题的决策变量。
123,,y y y 确定目标函数问题的目标很清楚——“收购价最小〞。
确定约束条件资源的报价至少应该高于原生产产品的利润，这样原厂家才可能卖。
因此有如下线性规划问题：123min 170100150w y y y =++*2、研究线性规划的对偶理论和方法〔包括对偶规划模型形式、对偶理论和对偶单纯形法〕。
答：略。
3、用单纯形法求解以下线性规划问题：〔1〕⎪⎪⎩⎪⎪⎨⎧≥≤+-≤++≤-++-=0,,43222..min32131321321321x x x x x x x x x x x t s x x x z ；〔2〕⎪⎪⎩⎪⎪⎨⎧=≥=++=+-=+-+-=)5,,2,1(052222..4min 53243232132 i x x x x x x x x x x t s x x z i 解：〔1〕引入松弛变量*4，*5，*6c j →1-11C B基b*1*2*3*4*5*60*421[1]-21000*532110100*64-101001c j -z j1-11因检验数σ2<0，故确定*2为换入非基变量，以*2的系数列的正分量对应去除常数列，最小比值所在行对应的基变量*4作为换出的基变量。
最优化方法大作业
学号《最优化方法》课程实践完成时间：2015年5月30日星期六选择题目：题目一使用优化软件，编写重要算法的程序1.第一大题：（1）学习最优流量工程问题，nonsmooth_MCFP.pdf（2）问题重述：Figure 1一个简单的网络拓扑和流量需求如Figure 1所示，网络有7 个节点，13 条弧，每条弧的容量是5 个单位. 此外有四个需求量均为4个单位的源－目的对（M=4)，具体的源节点、目的节点信息如图所示. 这里为了简单，省去了未用到的弧，此外弧上的数字表示弧的编号。
（3）极小化MAU设定变量x，为531⨯的向量，其中(53)x即为变量z。
使用linprog 函数求解极小化问题得到x。
之前确定三个约束条件。
1、Ax b⨯的矩阵，b为131⨯的向量。
≤，其中A为13532、eq eq x b A =，其中eq A 为2853⨯的矩阵，eq b 为281⨯的向量。
3、x lb ≥，其中lb 为153⨯的向量 编程计算后得到结果如下：（4） 极小化FT 成本函数设定变量x ，为651⨯的向量，其中(53:65)x 即为变量l z 。
使用linprog 函数求解极小化问题得到x 。
之前确定三个约束条件。
1、Ax b ≤，其中A 为7865⨯的矩阵，b 为781⨯的向量。
2、eq eq x b A =，其中eq A 为2865⨯的矩阵，eq b 为281⨯的向量。
3、x lb ≥，其中lb 为165⨯的向量 编程计算后得到结果如下：2. 第二大题： 2.1. 习题5.6 2.1.1. 问题分析问题2112212()(101810)/241513q x x x x x x x =-++-+ 通过matlab 画出其等高线为：2.1.2. 最速下降法最速下降法中，取值：k k p g =-==()()k k k kk T k k T k g p g g p Gp g Ggα- x1x 2等高线-224681012(1)()k x k x k p α+=+2.1.3. 算法流程图如下图所示：2.1.4. 初始值（0,0）编程运行结构为：收敛过程曲线为：2.1.5. 初始值（-0.4,0）编程运行结构为：收敛过程曲线为：x1x 2等高线-2246810122.1.6. 初始值（10,0）编程运行结构为：收敛过程曲线为：x1x 2-2246810122.1.7. 初始值（11,0）编程运行结构为：收敛过程曲线为：x1x 2-2246810122.2. 习题5.7 2.2.1. 问题分析问题()94ln(7)f x x x =--497g x =-- 24(7)G x =- Matlab 画出在区间(7 10)的函数、一阶导数、二阶导数的变化曲线为x1x 2-22468101277.588.599.510707274767880828486xf函数变化曲线77.588.599.510-35-30-25-20-15-10-50510xg一阶导数g 变化曲线2.2.2. 牛顿法牛顿法中，取值：k k k G s g =- 1k k k s xx +=+其中，如果G 不是半正定，则采用修正牛顿法(+)k k k G I s g λ=-77.588.599.51050100150200250300350400xg二阶导数G 变化曲线2.2.3.算法流程图如下图所示：2.2.4.初始值7.40编程运行结构为：收敛过程曲线为：2.2.5. 初始值7.20编程运行结构为：收敛过程曲线为：7.397.47.417.427.437.447.457.4670.24570.2570.25570.2670.265xf2.2.6. 初始值7.01编程运行结构为：收敛过程曲线为：7.17.27.37.47.57.67.770.270.470.670.87171.271.471.6xf2.2.7. 初始值7.80编程运行结构为：收敛过程曲线为：6.856.9 6.9577.057.17.157.27.257.37.35727476788082xf2.2.8. 初始值7.88编程运行结构为：收敛过程曲线为：7.17.27.37.47.57.67.77.87.97070.57171.572xf2.2.9. 分析函数在区间（7，7.8888）内是凸函数，G 恒大于零，所以单纯牛顿法保证收敛。
最优化方法课程设计-最优化大作业-用优化算法求解函数最值问题
最优化方法大作业---------用优化算法求解函数最值问题摘要最优化(optimization) 是应用数学的重要研究领域.它是研究在给定约束之下如何寻求某些因素(的量)，以使某一(或某些)指标达到最优的一些学科的总称。
最优化问题一般包括最小化问题和最大化问题，而最大化问题可以通过简单的转化使之成最最小化问题。
最小化问题分为两类，即约束最小化和无约束最小化问题。
在此报告中，前两个问题属于无约束最小化问题的求解，报告中分别使用了“牛顿法”和“共轭梯度法”。
后两个问题属于有约束最小化问题的求解，报告中分别用“外点法”和“内点法”求解。
虽然命名不一样，其实质都是构造“惩罚函数”或者“障碍函数”，通过拉格朗日乘子法将有约束问题转化为无约束问题进行求解。
再此报告中，“外点法”和“内点法”分别用了直接求导和调用“牛顿法”来求解无约束优化问题。
在此实验中，用“共轭梯度法”对“牛顿法”所解函数进行求解时出现错误，报告中另取一函数用“共轭梯度法”求解得到正确的结果。
此实验中所有的函数其理论值都是显见的，分析计算结果可知程序正确，所求结果误差处于可接受范围内。
报告中对所用到的四种方法在其使用以前都有理论说明，对“外点法”中惩罚函数和“内点法”中障碍函数的选择也有相应的说明，另外，对此次试验中的收获也在报告的三部分给出。
本报告中所用程序代码一律用MATLAB编写。
【关键字】函数最优化牛顿法共轭梯度法内点法外点法 MATLAB一，问题描述1，分别用共轭梯度发法和牛顿法来求解一下优化问题()()()()()441432243221102510min x x x x x x x x x f -+-+-++=2， 分别用外点法和内点发求解一下优化问题⎩⎨⎧≥-++01.min 212231x x t s x x二、问题求解1.1 用牛顿法求解()()()()()441432243221102510min x x x x x x x x x f -+-+-++=1.1.1问题分析：取步长为1而沿着牛顿方向迭代的方法称为牛顿法，在牛顿法中，初始点的取值随意，在以后的每次迭代中，()[]()k k k k x f x f x x ∇∇-=-+121，直到终止条件成立时停止。
最优化方法试卷与答案5套
《最优化方法》1一、填空题：1．最优化问题的数学模型一般为：____________________________，其中___________称为目标函数，___________称为约束函数，可行域D 可以表示为_____________________________，若______________________________，称*x 为问题的局部最优解，若_____________________________________，称*x 为问题的全局最优解。
2．设f(x)= 212121522x x x x x +-+，则其梯度为___________，海色矩阵___________,令,)0,1(,)2,1(T T d x ==则f(x)在x 处沿方向d 的一阶方向导数为___________，几何意义为___________________________________,二阶方向导数为___________________，几何意义为____________________________________________________________。
3．设严格凸二次规划形式为：012..222)(min 2121212221≥≥≤+--+=x x x x t s x x x x x f则其对偶规划为___________________________________________。
4．求解无约束最优化问题：n R x x f ∈),(min ，设k x 是不满足最优性条件的第k 步迭代点，则：用最速下降法求解时，搜索方向k d =___________ 用Newton 法求解时，搜索方向k d =___________ 用共轭梯度法求解时，搜索方向k d =___________________________________________________________________________。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
发动机空燃比控制器引言：我主要从事自动化相关研究。
这里介绍我曾经接触过的发动机空燃比控制器设计中的优化问题。
发动机空燃比控制器设计中的最优化问题AFR =afm m && (1)空燃比由方程（1）定义，在发动机运行过程中如果控制AFR 稳定在14.7可以获得最好的动力性能和排放性能。
如果假设进入气缸的空气流量am &可以由相关单元检测得到，则可以通过控制进入气缸的燃油流量f m &来实现空燃比的精确控制。
由于实际发动机的燃油喷嘴并不是直接对气缸喷燃油，而是通过进气歧管喷燃油，这么做会在进气歧管壁上液化形成油膜，因此不仅是喷嘴喷出的未液化部分燃油会进入气缸，油膜蒸发部分燃油也会进入气缸，如方程（2）。
这样如何更好的喷射燃油成为了一个问题。
1110101122211ττττ⎡⎤⎡⎤-⎢⎥⎢⎥⎡⎤⎡⎤⎢⎥⎢⎥=+⎢⎥⎢⎥⎢⎥⎢⎥-⎢⎥⎢⎥⎣⎦⎣⎦⎢⎥⎢⎥-⎢⎥⎢⎥⎣⎦⎣⎦⎡⎤⎡⎤⎢⎥⎣⎦⎢⎥⎣⎦f ff v X x x u x x X x y =x && (2)其中12、,==ff fv x m x m &&=f y m &,=fi u m &这里面，表示油膜蒸发量ff m &、fvm &表示为液化部分燃油、fim &表示喷嘴喷射的燃油，在τf 、τv 、X 都已知的情况下，由现代控制理论知识，根据系统的增广状态空间模型方程（3）00000011011011114.70ττττ⎡⎤⎡⎤-⎡⎤⎢⎥⎢⎥⎡⎤⎡⎤⎢⎥=-+-⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦⎣⎦⎢⎥⎢⎥⎢⎥⎣⎦⎢⎥⎢⎥⎣⎦⎣⎦⎡⎤⎡⎤=⎢⎥⎣⎦⎣⎦ff v v a X X u +q q m y q x x x &&& (3)其中()014.7⎰taq =y -m&。
由极点配置方法，只要设计控制器方程（4），就可以使得y 无差的跟踪阶跃输入，那么y 也能较好的跟踪AFR *am /&。
12-- u =K q K x (4)这里面的12、K K 确定，可由主导极点概念降维成两个参数12C ,C ，虽然都是最终稳态无差，但是目标是使得瞬态过程中y 和阶跃输入y r 的差异尽可能的小。
所以原问题转化成了一个关于参数的12C ,C 的，以()()2=-⎰W r c y y 为目标函数的优化问题。
这里面由于无法写出参数与目标函数的显示表达式，从而也就无法使用和导数有关的优化方法。
所以我选择用坐标轮换法，由于坐标轮换法，局部寻优性能好，且不需要事先知道函数的导数。
在matlab7.10.0上对其进行仿真试验（程序见附录）。
假设进入气缸的空气质量流量已知，并且已经适当量化，假设状态值已经由状态观测器测量得到，在τf = 1s 、τv = 60ms 、X = 0.5时按上述方法设计控制器。
坐标轮换法的精度10.051ε⎡⎤=⎣⎦,a ,b其中，11001010，''⎡⎤⎡⎤==⎣⎦⎣⎦a ,b ,为12'⎡⎤⎣⎦c ,c 的初次搜索区间；基于阶跃响应的系统性能指标()()2=-⎰W r c y y ，通过坐标轮换法得到最优参数1 4.0542,0.9544,2**=-=c c 得到最优极点1 4.0542+0.9544, 4.05420.92λλ**=-=--j 544j ，将最优极点带入爱克曼公式得到反馈增益矩阵9.7820,4.4379,43.3499⎡⎤-⎣⎦=K ，进一步得到控制器：()9.7820,4.437943.349914.7⎡⎤⎡⎤⎢⎥=---*-⎣⎦⎢⎥⎣⎦fffi a f fv m m mm m &&&&& （5） 画出最后的仿真跟踪效果图，如图1，可见坐标轮换法很好的完成了任务。
图1 fm &对14.7a m &的跟踪效果从而使得发动机空燃比的变化如图2。
图2 空燃比变化time(sec)time(sec)空燃比主函数tf = 1;tv = 0.06; X=0.5;k0 = 10000000;k = 0;A = [-1/tf 0;0 -1/tv];tf0 = 1;tv0 = 0.06;X0 = 0.5;A1 = [-1/tf0 0;0 -1/tv0];B1 = [X0/tf0 (1-X0)/tv0]';C = [1 1];D = 0;AA = [A1 zeros(2,1);C 0];BB = [B1 ;0];P = zblh();PJ = [-10*sqrt(P(1)*P(1)+P(2)*P(2)) (-P(1)+P(2)*1i) (-P(1)-P(2)*1i)];K = acker(AA, BB, J)K1 = [K(1) K(2)];K2 = K(3);-10*sqrt(P(1)*P(1)+P(2)*P(2))Ac = [(A-B*K1) -B*K2;C 0];Bc = [0;0;-1];Cc = [C 0];Dc = 0;hold on;sys = ss(Ac,Bc,Cc,Dc);a1 = 0:0.01:3;b1 = 1+0*a1;a2 = 3.01:0.01:4;b2 = 1+(a2-3);a3 = 4.01:0.01:6;b3 = 2+0*a3;a4 = 6.01:0.01:10;b4 = 2+(a4-6)/2;a5 = 10.01:0.01:14;b5 = 4+0*a5;a6 = 14.01:0.01:17;b6 = 4-(a6-14);a7 = 17.01:0.01:20;b7 = 1+0*a7;t0 = [a1 a2 a3 a4 a5 a6 a7];u0 = [b1 b2 b3 b4 b5 b6 b7];x0 = [0.5 0.5 -0.1769]';[y0,t0,x] = lsim(sys,u0,t0,x0);x1 = [1 0 0 ]*x';x2 = [0 1 0]*x';x3 = [0 0 1]*x';x = [x1;x2];u = -K1*x-K2*x3;l = u0./y0';hold on;plot(t0,u0,'-');plot(t0,y0,':');xlabel('time(sec)')ylabel('Output')zblu（）坐标轮换法函数function [ R ] = zblh( )m = 0;n = 10;x2 = m + 0.618*(n-m);f2 = wch(x2,0.5);x1 = m + 0.382*(n-m);f1 = wch(x1,0.5);M = 0;N = 10;y2 = M + 0.618*(N-M);F2 = wch(1,y2);y1 = M + 0.382*(N-M);F1 = wch(1,y1);flag = 1;a = 0.5;b = 0.5;k1 = 0;k2 = 0;while(abs(abs(m-n)+abs(M-N))>0.05) if flag ~= 0k1 = k1 + 1;if f1 < f2n = x2;x2 = x1;f2 = f1;x1 = m + 0.382*(n-m);f1 = wch(x1,b);else if f1 == f2;m = x1;n = x2;x2 = m + 0.618*(n-m); f2 = wch(x2,b);x1 = m + 0.382*(n-m); f1 = wch(x1,b);elsem = x1;x1 = x2;f1 = f2;x2 = m + 0.618*(n-m); f2 = wch(x2,b);endenda = (m+n)/2;flag = 0;elsek2 = k2 +1;if F1 < F2N = y2;y2 = y1;F2 = F1;y1 = M + 0.382*(N-M);F1 = wch(a,y1);else if F1 == F2;M = y1;N = y2;y2 = M + 0.618*(N-M); F2 = wch(a,y2);y1 = M + 0.382*(N-M);F1 = wch(a,y1);elseM = y1;y1 = y2;F1 = F2;y2 = M + 0.618*(N-M);F2 = wch(a,y2);endendb = (M+N)/2;flag = 1;endendk1k2R = [a b]end目标函数function [ f ] = wch(a,b)tf = 1;tv = 0.06;X=0.5;k0 = 10000;k = 0;A = [-1/tf 0;0 -1/tv];tf0 = 1;A1 = [-1/tf0 0;0 -1/tv];B = [X/tf (1-X)/tv]';B1 = [X/tf0 (1-X)/tv]';C = [1 1];D = 0;AA = [A1 zeros(2,1);C 0];BB = [B1 ;0];J = [-10*sqrt(a*a+b*b) (-a+b*1i) (-a-b*1i)]; K = acker(AA, BB, J);K1 = [K(1) K(2)];K2 = K(3);Ac = [(A-B*K1) -B*K2;C 0];Bc = [0;0;-1];Cc = [C 0];Dc = 0;t0 = 0:0.01:4;yn = heaviside(t0);sys = ss(Ac,Bc,Cc,Dc);u0 = heaviside(t0);y0 = lsim(sys,u0,t0);for t0 = 1:400;kn = (yn(t0) - y0(t0))*(yn(t0) - y0(t0)); k = k+kn;endf = k;end。