(完整word版)优化设计Matlab编程作业
(完整word版)优化设计Matlab编程作业
优化设计无约束优化min f(x)= 21x +22x -21x 2x -41x初选x0=[1,1]程序：Step 1: Write an M-file objfun1.m.function f1=objfun1(x)f1=x(1)^2+2*x(2)^2-2*x(1)*x(2)-4*x(1);Step 2: Invoke one of the unconstrained optimization routinesx0=[1,1];>> options = optimset('LargeScale','off');>> [x,fval,exitflag,output] = fminunc(@objfun1,x0,options)运行结果：x =4.0000 2.0000fval =-8.0000exitflag =1output =iterations: 3funcCount: 12stepsize: 1firstorderopt: 2.3842e-007algorithm: 'medium-scale: Quasi-Newton line search'message: [1x85 char]非线性有约束优化1. Min f(x)=321x +22x +21x -32x +5Subject to: 1g (x)=1x +2x +18≤02g (x)=51x -32x -25≤03g (x)=131x -4122x 0≤4g (x)=14≤1x 130≤5g (x)=2≤2x 57≤初选x0=[10,10]Step 1: Write an M-file objfun2.mfunction f2=objfun2(x)f2=3*x(1)^2+x(2)^2+2*x(1)-3*x(2)+5;Step 2: Write an M-file confun1.m for the constraints.function [c,ceq]=confun1(x)% Nonlinear inequality constraintsc=[x(1)+x(2)+18;5*x(1)-3*x(2)-25;13*x(1)-41*x(2)^2;14-x(1);x(1)-130;2-x(2);x(2)-57];% Nonlinear inequality constraintsceq=[];Step 3: Invoke constrained optimization routinex0=[10,10]; % Make a starting guess at the solution>> options = optimset('LargeScale','off');>> [x, fval] = ...fmincon(@objfun2,x0,[],[],[],[],[],[],@confun1,options)运行结果：x =3.6755 -7.0744fval =124.14952. min f （x ）=222154x x +s.t. 0632)(211≤-+=x x x g01)(212≥+=x x x g 初选x0=[1,1]Step 1: Write an M-file objfun3.mfunction f=objfun3(x)f=4*x(1)^2+5*x(2)^2Step 2: Write an M-file confun3.m for the constraints.function [c,ceq]=confun3(x)%Nonlinear inequality constraintsc=[2*x(1)+3*x(2)-6;-x(1)*x(2)-1];% Nonlinear equality constraintsceq[];Step 3: Invoke constrained optimization routinex0=[1,1];% Make a starting guess at the solution>> options = optimset('LargeScale','off');>> [x, fval] = ...fmincon(@objfun,x0,[],[],[],[],[],[],@confun,options)运行结果：Optimization terminated: no feasible solution found. Magnitude of searchdirection less than 2*options.TolX but constraints are not satisfied.x =1 1fval =-13实例：螺栓连接的优化设计图示为一压气机气缸与缸盖连接的示意图。
大连理工优化方法大作业MATLAB编程
function [x,dk,k]=fjqx(x,s) flag=0;a=0;b=0;k=0;d=1;while(flag==0)[p,q]=getpq(x,d,s);if (p<0)b=d;d=(d+a)/2;endif(p>=0)&&(q>=0)dk=d;x=x+d*s;flag=1;endk=k+1;if(p>=0)&&(q<0)a=d;d=min{2*d,(d+b)/2};endend%定义求函数值的函数fun，当输入为x0=（x1，x2）时，输出为f function f=fun(x)f=(x(2)-x(1)^2)^2+(1-x(1))^2;function gf=gfun(x)gf=[-4*x(1)*(x(2)-x(1)^2)+2*(x(1)-1),2*(x(2)-x(1)^2)]; function [p,q]=getpq(x,d,s)p=fun(x)-fun(x+d*s)+0.20*d*gfun(x)*s';q=gfun(x+d*s)*s'-0.60*gfun(x)*s';结果：x=[0,1];s=[-1,1];[x,dk,k]=fjqx(x,s)x =-0.0000 1.0000dk =1.1102e-016k =54function f= fun( X )%所求问题目标函数f=X(1)^2-2*X(1)*X(2)+2*X(2)^2+X(3)^2+ X(4)^2-X(2)*X(3)+2*X(1)+3*X(2)-X(3);endfunction g= gfun( X )%所求问题目标函数梯度g=[2*X(1)-2*X(2)+2,-2*X(1)+4*X(2)-X(3)+3,2*X(3)-X(2)-1,2*X(4)];endfunction [ x,val,k ] = frcg( fun,gfun,x0 )%功能：用FR共轭梯度法求无约束问题最小值%输入：x0是初始点，fun和gfun分别是目标函数和梯度%输出：x、val分别是最优点和最优值，k是迭代次数maxk=5000;%最大迭代次数rho=0.5;sigma=0.4;k=0;eps=10e-6;n=length(x0);while(k<maxk)g=feval(gfun,x0);%计算梯度itern=k-(n+1)*floor(k/(n+1));itern=itern+1;%计算搜索方向if(itern==1)d=-g;elsebeta=(g*g')/(g0*g0');d=-g+beta*d0;gd=g'*d;if(gd>=0.0)d=-g;endendif(norm(g)<eps)break;endm=0;mk=0;while(m<20)if(feval(fun,x0+rho^m*d)<feval(fun,x0)+sigma*rho^m*g'*d) mk=m;break;endm=m+1;endx0=x0+rho^mk*d;val=feval(fun,x0);g0=g;d0=d;k=k+1;endx=x0;val=feval(fun,x0);end结果：>> x0=[0,0,0,0];>> [ x,val,k ] = frcg( 'fun','gfun',x0 )x =-4.0000 -3.0000 -1.0000 0val =-8.0000k =21或者function [x,f,k]=second(x)k=0;dk=dfun(x);g0=gfun(x);s=-g0;x=x+dk*s;g1=gfun(x);while(norm(g1)>=0.02)if(k==3)k=0;g0=gfun(x);s=-g0;x=x+dk*s;g1=gfun(x);else if(k<3)u=((norm(g1))^2)/(norm(g0)^2); s=-g1+u*s;k=k+1;g0=g1;dk=dfun(x);x=x+dk*s;g1=gfun(x);endendf=fun(x);endfunction f=fun(x)f=x(1)^2-2*x(1)*x(2)+2*x(2)^2+x(3)^2+x(4)^2-x(2)*x(3)+2*x(1)+3*x(2)-x(3); function gf=gfun(x)gf=[2*x(1)-2*x(2)+2,-2*x(1)+4*x(2)-x(3)+3,2*x(3)-x(2)-1,2*x(4)];function [p,q]=con(x,d)ss=-gfun(x);p=fun(x)-fun(x+d*ss)+0.2*d*gfun(x)*(ss)';q=gfun(x+d*ss)*(ss)'-0.6*gfun(x)*(ss)';function dk=dfun(x)flag=0;a=0;d=1;while(flag==0)[p,q]=con(x,d);if (p<0)b=d;d=(d+a)/2;endif(p>=0)&&(q>=0)dk=d;flag=1;endif(p>=0)&&(q<0)a=d;d=min{2*d,(d+b)/2};endEnd结果：x=[0,0,0,0];>> [x,f,k]=second(x)x =-4.0147 -3.0132 -1.0090 0 f = -7.9999k = 1function [f,x,k]=third_1(x)k=0;g=gfun(x);while(norm(g)>=0.001)s=-g;dk=dfun(x,s);x=x+dk*s;k=k+1;g=gfun(x);f=fun(x);endfunction f=fun(x)f=x(1)+2*x(2)^2+exp(x(1)^2+x(2)^2);function gf=gfun(x)gf=[1+2*x(1)*exp(x(1)^2+x(2)^2),4*x(2)+2*x(2)*(x(1)^2+x(2)^2)];function [j_1,j_2]=con(x,d,s)j_1=fun(x)-fun(x+d*s)+0.1*d*gfun(x)*(s)'; j_2=gfun(x+d*s)*(s)'-0.5*gfun(x)*(s)'; function dk=dfun(x,s)%获取步长flag=0;a=0;d=1;while(flag==0)[p,q]=con(x,d,s);if (p<0)b=d;d=(d+a)/2;endif(p>=0)&&(q>=0)dk=d;flag=1;endif(p>=0)&&(q<0)a=d;d=min{2*d,(d+b)/2};endend结果：x=[0,1];[f,x,k]=third_1(x)f =0.7729x = -0.4196 0.0001k =8（1）程序：function [f,x,k]=third_2(x)k=0;H=inv(ggfun(x));g=gfun(x);while(norm(g)>=0.001)s=(-H*g')';dk=dfun(x,s);x=x+dk*s;k=k+1;g=gfun(x);f=fun(x);endfunction f=fun(x)f=x(1)+2*x(2)^2+exp(x(1)^2+x(2)^2); function gf=gfun(x)gf=[1+2*x(1)*exp(x(1)^2+x(2)^2),4*x(2)+2*x(2)*(x(1)^2+x(2)^2)]; function ggf=ggfun(x)ggf=[(4*x(1)^2+2)*exp(x(1)^2+x(2)^2),4*x(1)*x(2)*exp(x(1)^2+x(2)^2);4*x(1)*x(2)*exp(x(1)^2+x(2)^2),4+(4*x(2)^2+2)*exp(x(1)^2+x(2)^2)]; function [j_1,j_2]=con(x,d,s)j_1=fun(x)-fun(x+d*s)+0.1*d*gfun(x)*(s)';j_2=gfun(x+d*s)*(s)'-0.5*gfun(x)*(s)';function dk=dfun(x,s)% 步长获取flag=0;a=0;d=1;b=10000;while(flag==0)[p,q]=con(x,d,s);if (p<0)b=d;d=(d+a)/2;endif(p>=0)&&(q>=0)dk=d;flag=1;endif(p>=0)&&(q<0)a=d;if 2*d>=(d+b)/2d=(d+b)/2;else d=2*d;endendEnd结果：x=[0,1];[f,x,k]=third_2(x)f =0.7729x = -0.4193 0.0001k =8（2）程序：function [f,x,k]=third_3(x) k=0;X=cell(2);g=cell(2);X{1}=x;H=eye(2);g{1}=gfun(X{1});s=(-H*g{1}')';dk=dfun(X{1},s);X{2}=X{1}+dk*s;g{2}=gfun(X{2});while(norm(g{2})>=0.001)dx=X{2}-X{1};dg=g{2}-g{1};v=dx/(dx*dg')-(H*dg')'/(dg*H*dg'); h1=H*dg'*dg*H/(dg*H*dg');h2=dx'*dx/(dx*dx');h3=dg*H*dg'*v'*v;H=H-h1+h2+h3;k=k+1;X{1}=X{2};g{1}=gfun(X{1});s=(-H*g{1}')';dk=dfun(X{1},s);X{2}=X{1}+dk*s;g{2}=gfun(X{2});norm(g{2});f=fun(x);x=X{2};endfunction f=fun(x)f=x(1)+2*x(2)^2+exp(x(1)^2+x(2)^2);function gf=gfun(x)gf=[1+2*x(1)*exp(x(1)^2+x(2)^2),4*x(2)+2*x(2)*(x(1)^2+x(2)^2)];function ggf=ggfun(x)ggf=[(4*x(1)^2+2)*exp(x(1)^2+x(2)^2),4*x(1)*x(2)*exp(x(1)^2+x(2)^2);4*x(1)*x(2)* exp(x(1)^2+x(2)^2),4+(4*x(2)^2+2)*exp(x(1)^2+x(2)^2);function [p,q]=con(x,d,s)p=fun(x)-fun(x+d*s)+0.1*d*gfun(x)*(s)';q=gfun(x+d*s)*(s)'-0.5*gfun(x)*(s)';function dk=dfun(x,s)flag=0;a=0;d=1;b=10000;while(flag==0)[p,q]=con(x,d,s);if (p<0)b=d;d=(d+a)/2;endif(p>=0)&&(q>=0) dk=d;flag=1;endif(p>=0)&&(q<0)a=d;if 2*d>=(d+b)/2d=(d+b)/2;else d=2*d;endendend结果：x=[0,1];[f,x,k]=third_3(x)f =0.7729x = -0.4195 0.0000 k=6function callqpactH=[2 0; 0 2];c=[-2 -5]';Ae=[ ]; be=[ ];Ai=[1 -2; -1 -2; -1 2;1 0;0 1];bi=[-2 -6 -2 0 0]';x0=[0 0]';[x,lambda,exitflag,output]=qpact(H,c,Ae,be,Ai,bi,x0) function [x,lamk,exitflag,output]=qpact(H,c,Ae,be,Ai,bi,x0) epsilon=1.0e-9; err=1.0e-6;k=0; x=x0; n=length(x); kmax=1.0e3;ne=length(be); ni=length(bi); lamk=zeros(ne+ni,1); index=ones(ni,1);for (i=1:ni)if(Ai(i,:)*x>bi(i)+epsilon), index(i)=0; endendwhile(k<=kmax)Aee=[];if(ne>0), Aee=Ae; endfor(j=1:ni)if(index(j)>0), Aee=[Aee; Ai(j,:)]; end endgk=H*x+c;[m1,n1] = size(Aee);[dk,lamk]=qsubp(H,gk,Aee,zeros(m1,1)); if(norm(dk)<=err)y=0.0;if(length(lamk)>ne)[y,jk]=min(lamk(ne+1:length(lamk))); endif(y>=0)exitflag=0;elseexitflag=1;for(i=1:ni)if(index(i)&(ne+sum(index(1:i)))==jk) index(i)=0; break;endendendk=k+1;elseexitflag=1;alpha=1.0; tm=1.0;for(i=1:ni)if((index(i)==0)&(Ai(i,:)*dk<0)) tm1=(bi(i)-Ai(i,:)*x)/(Ai(i,:)*dk); if(tm1<tm)tm=tm1; ti=i;endendendalpha=min(alpha,tm);x=x+alpha*dk;if(tm<1), index(ti)=1; end endif(exitflag==0), break; endk=k+1;endoutput.fval=0.5*x'*H*x+c'*x; output.iter=k;function [x,lambda]=qsubp(H,c,Ae,be) ginvH=pinv(H);[m,n]=size(Ae);if(m>0)rb=Ae*ginvH*c + be;lambda=pinv(Ae*ginvH*Ae')*rb;x=ginvH*(Ae'*lambda-c);elsex=-ginvH*c;lambda=0;end结果>>callqpactx =1.40001.7000lambda =0.8000exitflag =output =fval: -6.4500iter: 7function [x,mu,lambda,output]=multphr(fun,hf,gf,dfun,dhf,dgf,x0)%功能: 用乘子法解一般约束问题: min f(x), s.t. h(x)=0, g(x).=0%输入: x0是初始点, fun, dfun分别是目标函数及其梯度；% hf, dhf分别是等式约束（向量）函数及其Jacobi矩阵的转置；% gf, dgf分别是不等式约束（向量）函数及其Jacobi矩阵的转置；%输出: x是近似最优点，mu, lambda分别是相应于等式约束和不等式约束的乘子向量; % output是结构变量, 输出近似极小值f, 迭代次数, 内迭代次数等maxk=500;c=2.0;eta=2.0;theta=0.8;k=0;ink=0;epsilon=0.00001;x=x0;he=feval(hf,x);gi=feval(gf,x);n=length(x);l=length(he);m=length(gi);mu=zeros(l,1);lambda=zeros(m,1);btak=10;btaold=10;while(btak>epsilon&&k<maxk)%调用BFGS算法程序求解无约束子问题[x,ival,ik]=bfgs('mpsi','dmpsi',x0,fun,hf,gf,dfun,dhf,dgf,mu,lambda,c);ink=ink+ik;he=feval(hf,x);gi=feval(gf,x);btak=0;for i=1:lbtak=btak+he(i)^2;end%更新乘子向量for i=1:mtemp=min(gi(i),lambda(i)/c);btak=btak+temp^2;endbtak=sqrt(btak);if btak>epsilonif k>=2&&btak>theta*btaoldc=eta*c;endfor i=1:lmu(i)=mu(i)-c*he(i);endfor i=1:mlambda(i)=max(0,lambda(i)-c*gi(i));endk=k+1;btaold=btak;x0=x;endendf=feval(fun,x);output.fval=f;output.iter=k;%增广拉格朗日函数function psi=mpsi(x,fun,hf,gf,dfun,dhf,dgf,mu,lambda,c) f=feval(fun,x);he=feval(hf,x);gi=feval(gf,x);l=length(he);m=length(gi);psi=f;s1=0;for i=1:lpsi=psi-he(i)*mu(i);s1=s1+he(i)^2;endpsi=psi+0.5*c*s1;s2=0;for i=1:ms3=max(0,lambda(i)-c*gi(i));s2=s2+s3^2-lambda(i)^2;endpsi=psi+s2/(2*c);%不等式约束函数文件g1.mfunction gi=g1(x)gi=10*x(1)-x(1)^2+10*x(2)-x(2)^2-34;%目标函数的梯度文件df1.mfunction g=df1(x)g=[4, -2*x(2)]';%等式约束（向量）函数的Jacobi矩阵（转置）文件dh1.m function dhe=dh1(x)dhe=[-2*x(1), -2*x(2)]'%不等式约束（向量）函数的Jacobi矩阵（转置）文件dg1.m function dgi=dg1(x)dgi=[10-2*x(1), 10-2*x(2)]';function [x,val,k]=bfgs(fun,gfun,x0,varargin)maxk=500;rho=0.55;sigma=0.4;epsilon=0.00001;k=0;n=length(x0);Bk=eye(n);while(k<maxk)gk=feval(gfun,x0,varargin{:});if(norm(gk)<epsilon)break;enddk=-Bk\gk;m=0;mk=0;while(m<20)newf=feval(fun,x0+rho^m*dk,varargin{:});oldf=feval(fun,x0,varargin{:});if(newf<oldf+sigma*rho^m*gk'*dk)mk=m;break;endm=m+1;endx=x0+rho^mk*dk;sk=x-x0;yk=feval(gfun,x,varargin{:})-gk;if(yk'*sk>0)Bk=Bk-(Bk*sk*sk'*Bk)/(sk'*Bk*sk)+(yk*yk')/(yk'*sk);endk=k+1;x0=x;endval=feval(fun,x0,varargin{:});结果x=[2 2]';[x,mu,lambda,output]=multphr('fun','hf','gf1','df','dh','dg',x0) x =1.00134.8987mu =0.7701lambda =0.9434output =fval: -31.9923iter: 4f=[3,1,1];A=[2,1,1;1,-1,-1];b=[2;-1];lb=[0,0,0];x=linprog(f,A,b,zeros(3),[0,0,0]',lb)结果：Optimization terminated.x =0.00000.50000.5000。
(完整word版)MATLAB程序设计教程(第二版)课后答案
MATLAB第二版课后答案unit3—8 unit3实验指导1、 n=input（'请输入一个三位数：'）;a=fix（n/100)；b=fix（（n-a*100)/10）；c=n-a＊100-b＊10；d=c＊100+b*10+a2（1）n=input('请输入成绩’）；switch ncase num2cell（90：100)p='A'；case num2cell（80:89)p='B'；case num2cell(70：79）p=’C';case num2cell（60:69)p='D'；otherwisep='E';endprice=p(2）n=input(’请输入成绩');if n〉=90&n〈=100p='A’；elseif n>=80&n<=89p='B';elseif n〉=70&n<=79p=’C’;elseif n>=60＆n<=69p='D'；elsep='E'；endprice=p(3）tryn;catchprice='erroe'end3n=［1,5,56，4,3,476，45，6，3，76，45,6,4,3，6，4,23,76,908，6]; a=n(1);b=n(1)；for m=2:20if n（m）>aa=n（m）；elseif n（m)<bb=n（m）;endendmax=amin=b法2n=[1,5，56,4，3，476,45，6,3,76,45，6，4，3,6，4，23，76,908,6］;min=min(n)max=max(n）4b=［—3.0:0.1:3.0];for n=1：61a=b(n）；y(n)=（exp（0.3＊a)-exp(—0。
(完整word版)matlab经典习题及解答
第1章 MATLAB概论1.1与其他计算机语言相比较，MATLAB语言突出的特点是什么？MATLAB具有功能强大、使用方便、输入简捷、库函数丰富、开放性强等特点.1.2 MATLAB系统由那些部分组成?MATLAB系统主要由开发环境、MATLAB数学函数库、MATLAB语言、图形功能和应用程序接口五个部分组成。
1.4 MATLAB操作桌面有几个窗口?如何使某个窗口脱离桌面成为独立窗口?又如何将脱离出去的窗口重新放置到桌面上?在MATLAB操作桌面上有五个窗口,在每个窗口的右上角有两个小按钮，一个是关闭窗口的Close按钮，一个是可以使窗口成为独立窗口的Undock按钮，点击Undock按钮就可以使该窗口脱离桌面成为独立窗口,在独立窗口的view菜单中选择Dock ……菜单项就可以将独立的窗口重新防止的桌面上。
1.5 如何启动M文件编辑/调试器？在操作桌面上选择“建立新文件”或“打开文件"操作时，M文件编辑/调试器将被启动.在命令窗口中键入edit命令时也可以启动M文件编辑/调试器.1.6 存储在工作空间中的数组能编辑吗？如何操作？存储在工作空间的数组可以通过数组编辑器进行编辑：在工作空间浏览器中双击要编辑的数组名打开数组编辑器，再选中要修改的数据单元，输入修改内容即可。
1。
7 命令历史窗口除了可以观察前面键入的命令外,还有什么用途?命令历史窗口除了用于查询以前键入的命令外，还可以直接执行命令历史窗口中选定的内容、将选定的内容拷贝到剪贴板中、将选定内容直接拷贝到M文件中。
1。
8 如何设置当前目录和搜索路径,在当前目录上的文件和在搜索路径上的文件有什么区别？当前目录可以在当前目录浏览器窗口左上方的输入栏中设置，搜索路径可以通过选择操作桌面的file菜单中的Set Path 菜单项来完成。
在没有特别说明的情况下，只有当前目录和搜索路径上的函数和文件能够被MATLAB运行和调用，如果在当前目录上有与搜索路径上相同文件名的文件时则优先执行当前目录上的文件，如果没有特别说明，数据文件将存储在当前目录上。
Matlab优化设计作业答案
长江大学机械工程学院机械优化设计大作业班级2012年5月31-361、⎩⎨⎧=+-=++⋅+-++=22422min 321321321232221x x x x x x t s x x x x x x f解： function f=fun1(x)f=x(1)^2+2*x(2)^2+x(3)^2-2*x(1)*x(2)+x(3)clearclcx0=[1;1;1];A=[];b=[];Aeq=[1 1 1;2 -1 1];beq=[4;2];lb=[];ub=[];[x,fval,exitflag]=fmincon('fun1',x0,A,b,Aeq,beq,lb,ub)x =1.90911.95450.1364fval =3.9773exitflag =12、221212min ()2130f x x x x s t x =+-+⋅-≤解： function f=fun2(x)f=x(1)^2+x(2)^2-2*x(1)+1clearclcx0=[1;1];A=[0 -1];b=-3;Aeq=[];beq=[];lb=[];ub=[];[x,fval,exitflag]=fmincon('fun2',x0,A,b,Aeq,beq,lb,ub)x =13fval =9exitflag =13、 2212min ()(2)f X x x =-+112222312..()0()0()10s t g X x g X x g X x x =-≤=-≤=-+-≤解：function f=fun3(x)f=(x(1)-2)^2+x(2)^2function [g,ceq]=mycon3(x)g=[-x(1);-x(2);-x(1)^2+x(2)^2-1];ceq=[];clearclcx0=[0;0];A=[];b=[];Aeq=[];beq=[];lb=[];ub=[];[x,fval]=fmincon('fun3',x0,A,b,Aeq,beq,lb,ub,'mycon3')x =2.0000fval =2.2204e-0164、⎩⎨⎧≤--≤+--⋅++++=01005.1)12424(min 21212122122211x x x x x x t s x x x x x e f x解： function f=fun4(x)f=exp(x(1))*(4*x(1)^2+2*x(2)^2+4*x(1)*x(2)+2*x(2)+1);function [g,ceq]=mycon4(x)g(1)=1.5+x(1)*x(2)-x(1)-x(2);g(2)=-x(1)*x(2)-10;ceq=[];clearclcx0=[0;0];A=[];b=[];Aeq=[];beq=[];lb=[];ub=[];[x,fval]=fmincon('fun4',x0,A,b,Aeq,beq,lb,ub,'mycon4')x =1.1825-1.7398fval =3.06085 喜糖问题：需要购买甲乙两种喜糖，喜糖甲10元/斤，喜糖乙20元/斤。
(完整word版)Matlab求解优化问题
预备知识：M 文件简介在MATLAB 中，用户可以利用Edtior （编辑器）建立M 文件，然后在命令窗口中的“>>”提示符下键入M 文件的主文件名，回车执行.MATLAB 的M 文件有两类：命令文件和函数文件。
将原本要在MATLAB 环境下直接输入的语句，放在一个以 .m 为后缀的文件中，这一文件就称为命令文件；函数文件由五部分组成：函数定义行、H1行、函数帮助文本、函数体、注释，MATLAB 的内部函数都是由函数文件定义的。
1.11 优化(最值、数学规划)在数学上，优化问题包括最值问题和数学规划问题等，后者又包括线性规划、整数规划（含0-1规划）、二次规划等.在MATLAB 中，求解最值问题的命令主要有：fminbnd (f,x1,x2) 求一元函数f 在区间[x1,x2]上的最小值点[x,fval]=fminbnd(f,x1,x2) 求一元函数f 在区间[x1,x2]上的最小值点和最小值 fminsearch (’f’,x0) 求多元函数f 在点x0附近的最小值点[x,fval]=fminsearch(’f’,x0) 求多元函数f 在点x0附近的最小值点和最小值例1.11.1 求函数23)(2++=x x x f 在区间]5,5[-上的最小值点和最小值. >> [x,fval]=fminbnd('x^2+3*x+2',-5,5) x =-1.5000 fval =-0.2500例1.11.2 求函数21212122),(x x x x x x f ++=在点)1,1(附近的最小值点和最小值. >> [x,fval]= fminsearch('x(1)*x(2)+2/x(1)+2/x(2)',[1 1]) x =1.2599 1.2599 fval =4.7622在MATLAB 中，求解数学规划问题的命令主要有：（1）线性规划⎪⎪⎩⎪⎪⎨⎧≤≤=⋅≤=ub x lb beq x Aeq bAx t s x c z T ..min命令：[x,fval]=linprog (c,A,b,Aeq,beq,lb,ub)在上述命令中，当某些参数空缺时，可用[]代替或省略，下同。
(完整word版)matlab(7,4)汉明码和(7,4)循环码的编程设计
二、创新实验设计创新实验一：（7,4）汉明码的编码与译码实现1、实验目的实现（7,4）汉明码的编码与译码，通过这次实验不但加深了对汉明码编码和译码原理了解，而且对线性分组码有所了解。
2、实验原理线性分组码的构造方法比较简单、理论较为成熟，应用比较广泛。
汉明码是一种能够纠正一个错码的效率比较高的线性分组码，下面以（7,4）码为例就汉明码的编码与译码分别进行介绍：（1）编码原理一般来说，若汉明码长为n ，信息位数为k ，则监督位数r=n-k 。
若希望用r 个监督位构造出r 个监督关系式来指示一位错码的n 种可能位置，则要求21r n -≥或211rk r -≥++ (1)设汉明码（n,k ）中k=4，为了纠正一位错码，由式（1）可知，要求监督位数r ≥3。
若取r=3,则n=k+r=7。
这样就构成了（7,4）码。
用6543210a a a a a a a 来表示这7个码元，用123s s s 的值表示3个监督关系式中的校正子，则123s s s 的值与错误码元位置的对应关系可以规定如表1所列。
表2.1 校正子和错码位置的关系则由表1可得监督关系式：16542s a a a a =⊕⊕⊕()226531s a a a a =⊕⊕⊕()3 36430s a a a a =⊕⊕⊕()4 在发送端编码时，信息位6543a a a a 的值决定于输入信号，因此它们是随机的。
监督位2a 、1a 、a 应根据信息位的取值按监督关系来确定，为使所编的码中无错码，则123,,S S S 等于0，即65426531643000(5)0a a a a a a a a a a a a ⊕⊕⊕=⎧⎪⊕⊕⊕=⎨⎪⊕⊕⊕=⎩方程组（5）可等效成如下矩阵形式6543210111010001101010010110010a a a a a a a ⎡⎤⎢⎥⎢⎥⎢⎥⎡⎤⎡⎤⎢⎥⎢⎥⎢⎥=⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦⎣⎦⎢⎥⎢⎥⎢⎥⎣⎦（6）式(6)可简化为0T T HA =，H 为监督矩阵，则由式(6)可得到监督矩阵11101001101010=[P I ] (7)1011001r H ⎡⎤⎢⎥=⎢⎥⎢⎥⎣⎦因为生成矩阵'=[I Q]=[I ]k k G P ，所以由(7)得生成矩阵G 如下：[]k 10001110100110[']00101010001011k G I Q I P ⎡⎤⎢⎥⎢⎥===⎢⎥⎢⎥⎣⎦然后利用信息位和生成矩阵G 相乘产生整个码组，即有[][]65432106543=(8)A a a a a a a a a a a a G=其中A 为整个码组矩阵，6543a a a a 是信息位。
最优化方法matlab作业
实用最优化方法——matlab编程作业初值为[-1;1]其中g0、g1分别为不同x值下得导数，f0、f1为函数值MATLAB程序：x0=[-1;1];s0=[1;1];c1=0.1;c2=0.5;a=0;b=inf;d=1;n=0;x1=x0+d*s0;g0=[-400*(x0(2)-x0(1)^2)*x0(1)-2*(1-x0(1));200*(x0(2)-x0(1) ^2)];g1=[-400*(x1(2)-x1(1)^2)*x1(1)-2*(1-x1(1));200*(x1(2)-x1(1) ^2)];f1=100*(x1(2)-x1(1)^2)^2+(1-x1(1))^2;f0=100*(x0(2)-x0(1)^2)^2+(1-x0(1))^2;while((f0-f1<-c1*d*g0'*s0)||(g1'*s0<c2*g0'*s0))if ((f0-f1)<(-c1*d*g0'*s0))b=d;d=(d+a)/2;x1=x0+d*s0;g0=[-400*(x0(2)-x0(1)^2)*x0(1)-2*(1-x0(1));200*(x0 (2)-x0(1)^2)];g1=[-400*(x1(2)-x1(1)^2)*x1(1)-2*(1-x1(1));200*(x1 (2)-x1(1)^2)];f1=100*(x1(2)-x1(1)^2)^2+(1-x1(1))^2;f0=100*(x0(2)-x0(1)^2)^2+(1-x0(1))^2;elseif (((g1')*s0)<(c2*(g0')*s0))a=d;if(2*d<=(d+b)/2)d=2*d;elsed=(d+b)/2;endx1=x0+d*s0;g0=[-400*(x0(2)-x0(1)^2)*x0(1)-2*(1-x0(1));200*(x0(2) -x0(1)^2)];g1=[-400*(x1(2)-x1(1)^2)*x1(1)-2*(1-x1(1));200*(x1(2 )-x1(1)^2)];f1=100*(x1(2)-x1(1)^2)^2+(1-x1(1))^2;f0=100*(x0(2)-x0(1)^2)^2+(1-x0(1))^2;endx1df1=100*(x1(2)-x1(1)^2)^2+(1-x1(1))^2 计算结果：最优点：x1 = -0.99611.0039步长： d = 0.0039最优解：f1 = 3.9981目标函数：function f = fun2( x )f=x(1)^2-2*x(1)*x(2)+2*x(2)^2+x(3)^2-x(2)*x(3)+2*x(1)+3*x(2 )-x(3);end目标函数梯度：function g = gfun2( x )g=[2 -2 0;-2 4 -1;0 -1 2]*x+[2;3;-1];end源代码：（以初值为为（0；0；0））x0=[0;0;0]; %初始值eps=1.0e-5; %精度g0=gfun2(x0);s0=-g0;n=0;syms d1;while norm(g0)>epsif n<3g=gfun2(x0+d1*s0);d= double(solve(s0'*g));x1=x0+d*s0;g1=gfun2(x1);if norm(g1)<epsn=n+1;x0=x1;breakelses0=-g1+(norm(g1)^2/norm(g0)^2)*s0;x0=x1;g0=g1;endelseif n==3x0=x1;g0=gfun2(x0);s0=-g0;n=0;endn=n+1;x0nfun2(x0)计算结果：最优点：x0 = -4 -3 -1 迭代次数： n = 3 最优值： ans = -8（1）最速下降法：目标函数：function f= fun3_1(x )f=x(1)+2*x(2)^2+exp(x(1)^2+x(2)^2);end目标函数的梯度：function g= gfun3_1(x)g=[1+2*x(1)*exp(x(1)^2+x(2)^2);4*x(2)+2*x(2)*exp(x(1)^2+x(2 )^2)];end源代码（初值为（1；0））：x0=[1;0];%初始值eps=1.0e-5;%精度n=0;g0=gfun3(x0);syms d1;while norm(g0)>=epss0=-g0;g=gfun3(x0+d1*s0);d= double(solve(s0'*g));x1=x0+d*s0;g1=gfun3(x1);if( norm(g1)<eps)n=n+1;x0=x1;break;elsex0=x1;g0=gfun3(x0);endn=n+1;endf0=fun3(x0)x0n计算结果：最优值： f0 = 0.7729 最优点：x0 = -0.4194 0迭代次数：n = 1（2）牛顿法目标函数：function f= fun3_1(x )f=x(1)+2*x(2)^2+exp(x(1)^2+x(2)^2);end目标函数梯度：function g= gfun3_1(x)g=[1+2*x(1)*exp(x(1)^2+x(2)^2);4*x(2)+2*x(2)*exp(x(1)^2+x(2 )^2)];end目标函数的Hesse阵：function g2 = g2fun3(x)g2=[2*exp(x(1)^2+x(2)^2)+4*x(1)^2*exp(x(1)^2+x(2)^2),4*x(1) *x(2)*exp(x(1)^2+x(2)^2)4*x(1)*x(2)*exp(x(1)^2+x(2)^2),4+2*exp(x(1)^2+x(2)^2)+4*x(2 )^2*exp(x(1)^2+x(2)^2)];end源代码（初值为（1；0））：x0=[1;0];%初始值eps=1.0e-5;%精度n=0;g0=gfun3(x0);g20=g2fun3(x0);while norm(g0)>=epsd=-g20\g0;x1=x0+d;g1=gfun3(x1);if( norm(g1)<eps)n=n+1;x0=x1;break;elsex0=x1;g0=gfun3(x0);endn=n+1;endf0=fun3(x0)x0n计算结果：最优值：f0 = 0.7729最优点：x0 = -0.4194迭代次数：n = 63（3）利用BGFS法：目标函数：function f= fun3_1(x )f=x(1)+2*x(2)^2+exp(x(1)^2+x(2)^2);end目标函数梯度function g= gfun3_1(x)g=[1+2*x(1)*exp(x(1)^2+x(2)^2);4*x(2)+2*x(2)*exp(x(1)^2+x(2 )^2)];end源代码（初值为（1；0））：x0=[1;0];%初始值eps=1.0e-5;%精度n=0;g0=gfun3(x0);syms d1;h0=eye(2);while norm(g0)>=epss0=-h0*g0;g=gfun3(x0+d1*s0);d= double(solve(s0'*g));x1=x0+d*s0;g1=gfun3(x1);if( norm(g1)<eps)n=n+1;x0=x1;break;elseh0=h0-(h0*(g1-g0)*(g1-g0)'*h0)/((g1-g0)'*h0*(g1-g0))+...((x1-x0)*(x1-x0)')/((x1-x0)'*(g1-g0))+...((g1-g0)'*h0*(g1-g0))*((x1-x0)*(x1-x0)');x0=x1;g0=gfun3(x0);endn=n+1;endf0=fun3(x0)x0n计算结果：最优值：f0 = 0.7729最优点：x0 = -0.4194迭代次数：n = 1题四：求解子程序：function [x,lambda]=qsubp(H,c,Ae,be) ginvH=pinv(H);[m,n]=size(Ae);if(m>0)rb=Ae*ginvH*c + be;lambda=pinv(Ae*ginvH*Ae')*rb;x=ginvH*(Ae'*lambda-c);elsex=-ginvH*c;lambda=0;endend源代码：H=[2 0;0 2];%目标函数的hesse阵c=[-2,-4];Ae=[0 0];be=[0;0];Ai=[1/2 0;-1 3];bi=[1;2];x0=[-1;0];%初始值内部点eps=1.0e-9; %µ±ax-b=epsʱµ±×öax-b=0´¦Àíerr=1.0e-6;k=0;x0(1)=x0(1)+x0(2);x=x0;n=length(x);max=1.0e3;ne=length(be);ni=length(bi);lamk=zeros(ne+ni,1);index=ones(ni,1);for i=(1:ni)if(bi(i)-Ai(i,:)*x>eps)index(i)=0;。
优化方法MATLAB编程——大连理工大学
优化方法上机大作业学院：姓名：学号：指导老师：肖现涛第一题源程序如下：function zy_x = di1ti(x)%di1ti是用来求解优化作业第一题的函数。
x0=x; yimuxulong=0.000001;g0=g(x0);s0=-g0;A=2*ones(100,100);k=0;while k<100lanmed=-(g0)'*s0/(s0'*A*s0);x=x0+lanmed*s0;g=g(x);k=k+1;if norm(g)<yimuxulongzy_x=x;fprintf('After %d iterations,obtain the optimal solution.\n \n The optimal solution is \n %f.\n\nThe optimal "x" is "ans".',k,f(x) )break;endmiu=norm(g)^2/norm(g0)^2;s=-g+miu*s0;g0=g; s0=s;x0=x;endfunction f=f(x)f=(x'*ones(100,1))^2-x'*ones(100,1);function g=g(x)g=(2*x'*ones(100,1))*ones(100,1)-ones(100,1);代入x0,运行结果如下：>> x=zeros(100,1);>> di1ti(x)After 1 iterations,obtain the optimal solution.The optimal solution is-0.250000.The optimal "x" is "ans".ans =0.005*ones(100,1).第二题1.最速下降法。
机械优化设计MATLAB程序
机械优化设计MATLA程序(总11页)-CAL-FENGHAI.-(YICAI)-Company One 1■CAL•本页仅作为文档封面.使用请直接删除机械优化设计作业2•用二次插值法求函数刃)=(/ + 1卜极小值，精度e=o在MATLAB的M文件编辑器中编写的M文件，如下：f=inline('(t+l)*(t-2)A2,/,t,)a=0;b=3;epsil on二；tl=a;fl=f(tl);t3=b;f3=f(t3); t2=*(tl+t3);f2=f(t2); cl=(f3-fl)/(t3-tl); C2=((f2-fl)/(t2-tl)-cl)/(t2-t3); t4=*(tl+t3-cl/c2);f4=f(t4); k=0;while(abs(t4-t2)>=epsilon)if t2<t4if f2>f4fl=f2;tl=t2;t2=t4;f2=f4;elsef3=f4;t3=t4;endelseif f2>f4f3=f2;t3=t2;t2=t4;f2=f4;elsefl=f4;t2=t4;endendcl=(f3-fl)/(t3-tl);C2=((f2-fl)/(t2-tl)-cl)/(t2-t3); t4=*(tl+t3-cl/c2);f4=f(t4);k=k+l;end%输出最优解if f2>f4t=t4;f=f(t4);elset=t2;f=f(t2);endfprintf(l；迭代计算k=%\n',k)fprintffl/极小点坐标t=%\n',t) fprintf(l；函数值f=%\n\f)运行结果如下：迭代计算k= 7极小点坐标t= 2函数值f=2.用黄金分割法求函数池)=疔-亦T的极小值，精度在MATLAB的M文件编辑器中编写的M 文件，如下：f=inline(,tA(2/3)-(t A2+l)A(l/3)l/t,);a=0;b=3;epsil on二；tl=*(b-a);fl=f(tl);t2=a+*(b-a);f2=f(t2);k=l;while abs(b-a)>=epsilonif fl<f2b=t2;t2=tl;f2=fl;tl=*(b-a);fl=f(tl);elsea=tl;tl=t2;fl=f2; t2=a+*(b-a);f2=f(t2);endt=*(b+a);k=k+l;fo=f(t);endfprintf(l；迭代次数k=% \n',k)fprintf(l,'迭代区间一左端a=%\n'/a)fprintf(l；试点 1 坐标值tl=%\n'/tl)fprintf(l；函数值fl=%\n',f(tl))fprintf(l/迭代区间一右端b=%\n',b)fprintf(l；试点 2 坐标值t2=%\n',t2)fprintf(l；函数值f2=%\n',f(t2))fprintf(l；区间中点t=%\n',t)fprintffl；函数值fO=%\n\f(t))运行结果如下：迭代次数k= 13迭代区间一左端护试点1坐标值tl=函数值fl=迭代区间一右端b=试点2坐标值12=函数值f2=区间中点匸函数值fo=山黄金分割法在初始区间［0, 3］求得的极小值点为tr极小值为。
大连理工优化方法大作业MATLAB编程
fun ctio n [x,dk,k]=fjqx(x,s) flag=0;a=0;b=0;k=0;d=1;while (flag==0)[p,q]=getpq(x,d,s);if (P<0)b=d;d=(d+a)/2;endif(p>=0) &&( q>=0)dk=d;x=x+d*s;flag=1;endk=k+1;if (p>=0)&&(q<0)a=d;d=min{2*d,(d+b)/2};endend%定义求函数值的函数 fun ，当输入为 x0= (x1 ， x2 )时，输出为 f function f=fun(x)f=(x(2)-x(1)A2)A2+(1-x(1)F2;function gf=gfun(x)gf=[-4*x(1)*(x (2) -x(1)A2)+2*(x(1)-1),2*(x(2)-x(1)A2)];function [p,q]=getpq(x,d,s)p=fun(x)-fun(x+d*s)+0.20*d*gfun(x)*s';q=gfun(x+d*s)*s'-0.60*gfun(x)*s';结果：x=[0,1];s=[-1,1];[x,dk,k]=fjqx(x,s)x =-0.0000 1.0000dk =1.1102e-016k =54取初始= (0.0. 0,0)r^'l用兵柜梯皮法求解下面无约東优化问题：min f (x) = x孑—2x^X2 十2x孑 + x孑H-爲—X2天3 十 2xj + 3|X2 —*3,其中步长g的选取可利用习題1戎精确一维披索.注:通过比习题验证共範梯度法求辉门无二次西数极小点至多需要“次迭代.fun ctio n f= fun( X )%所求问题目标函数f=X(1)A2-2*X(1)*X (2)+2*X(2)A2+X(3)A2+ X(4) A2-X( 2)*X(3)+2*X(1)+3*X(2)-X(3);end function g= gfun( X )%所求问题目标函数梯度g=[2*X(1)-2*X(2)+2,-2*X(1)+4*X(2)-X(3)+3,2*X (3) -X (2)-1,2*X(4)];end function [ x,val,k ] = frcg( fun,gfun,xO )%功能：用FR共轭梯度法求无约束问题最小值%输入：x0是初始点，fun和gfun分别是目标函数和梯度%输出：x、val分别是最优点和最优值，k是迭代次数maxk=5000; %最大迭代次数rho=0.5;sigma=0.4;k=0;eps=10e-6;n=length(x0);while (k<maxk)g=feval(gfun,x0); % 计算梯度 itern=k-(n+1)*floor(k/(n+1));itern=itern+1;%计算搜索方向if (itern==1)d=-g;elsebeta=(g*g')/(g0*g0');d=-g+beta*d0;gd=g'*d;if (gd>=0.0)d=-g;endendif (norm(g)<eps)break ;endm=0;mk=0;while (m<20)if(feval(fu n,xO+rhoAm*d)<feval(fu n,xO)+sigma*rhoAm*g'*d) mk=m; break ;endm=m+1;endx0=x0+rho A mk*d;val=feval(fun,x0);g0=g;d0=d;k=k+1;endx=x0;val=feval(fun,x0);end结果：>> x0=[0,0,0,0];>> [ x,val,k ] = frcg( 'fun','gfun',x0 ) x =-4.0000 -3.0000 -1.0000 0val =-8.0000k =或者function [x,f,k]=second(x)k=0;dk=dfun(x);g0=gfun(x);s=-g0;x=x+dk*s;g1=gfun(x);while (norm(g1)>=0.02)if (k==3)k=0;g0=gfun(x);s=-g0;x=x+dk*s;g1=gfun(x);else if (k<3)u=(( norm(g1))A2)/( norm(gO)A2); s=-g1+u*s;k=k+1;g0=g1;dk=dfun(x);x=x+dk*s;g1=gfun(x);endendf=fun(x);endfunction f=fun(x)f=x(1F2-2*x(1)*x (2)+2*x (2)A2+x(3)A2+x(4)A2-x (2) *x (3)+2*x(1)+3*x(2)-x(3); function gf=gfun(x)gf=[2*x(1)-2*x(2)+2,-2*x(1)+4*x(2)-x(3)+3,2*x(3)-x(2)-1,2*x(4)];function [p,q]=con(x,d)ss=-gfun(x);p=fun(x)-fun(x+d*ss)+0.2*d*gfun(x)*(ss)';q=gfun(x+d*ss)*(ss)'-0.6*gfun(x)*(ss)';function dk=dfun(x)flag=0;a=0;d=1;while (flag==0)[p,q]=con(x,d);if (p<0)b=d;d=(d+a)/2;endif (p>=0)&&(q>=0)dk=d;flag=1;endif (p>=0)&&(q<0)a=d;d=min{2*d,(d+b)/2};endEnd结果： x=[0,0,0,0];>> [x,f,k]=second(x)x =-4.0147 -3.0132-1.0090 0 f = -7.9999k = 1取初始点3 = （0」）二考虑下面无约東优化问题：min f（x）二冷 + 2x2 + exp（xf + 天孑），其中歩长Qk的选取可別用习题1或精确一维搜索•搜索方向为一HNW ♦取垃=b•取皿=R2f防）］"9耳丈啟为BFG5公式亠通过此习题体会上述三种算法的收敛速度.fun ctio n [f,x,k]=third_1(x) k=0;g=gfu n(x);while (norm(g)>=0.001) s=-g;dk=dfu n( x,s);x=x+dk*s;k=k+1;g=gfu n(x);f=fun( x);endfun ctio n f=fun(x)f=x(1)+2*x(2)A2+exp(x(1)A2+x(2)A2);fun ctio n gf=gfu n(x)gf=[1+2*x(1)*exp(x(1)A2+x(2)A2),4*x(2)+2*x(2)*(x(1)A2+x(2)A2)]; function[j_1,j_2]=con(x,d,s)j_1=fun(x)-fun(x+d*s)+0.1*d*gfun(x)*(s)'; j_2=gfun(x+d*s)*(s)'-0.5*gfun(x)*(s)'; function dk=dfun(x,s) % 获取步长 flag=0;a=0;d=1;while (flag==0)[p,q]=con(x,d,s);if (p<0)b=d;d=(d+a)/2;endif (p>=0)&&(q>=0)dk=d;flag=1;endif (p>=0)&&(q<0)a=d;d=min{2*d,(d+b)/2}; end结果：x=[0,1];[f,x,k]=third_1(x)f =0.7729x = -0.4196 0.0001k =8(1 ) 程序：function [f,x,k]=third_2(x)k=0;H=inv(ggfun(x));g=gfun(x);while (norm(g)>=0.001)s=(-H*g')';dk=dfun(x,s);x=x+dk*s;k=k+1;g=gfun(x);f=fun(x);endfunction f=fun(x)f=x(1)+2*x(2)A2+exp(x(1F2+x(2)A2);function gf=gfun(x) gf=[1+2*x(1)*exp(x(1F2+x(2)A2),4*x(2)+2*x(2)*(x(1F2+x(2)A2)]; function ggf=ggfun(x)ggf=[(4*x(1)A2+2)*exp(x(1)A2+x (2) A2),4*x(1)*x (2) *exp(x(1)A2+x(2)A2);4*x(1)*x(2)*exp(x(1)A2+x(2)A2),4+(4*x(2)A2+2)*exp(x(1)A2+x(2)A2)];function [j_1,j_2]=con(x,d,s)j_1=fun(x)-fun(x+d*s)+0.1*d*gfun(x)*(s)';j_2=gfun(x+d*s)*(s)'-0.5*gfun(x)*(s)'; function dk=dfun(x,s) % 步长获取flag=0;a=0;d=1;b=10000;while (flag==0)[p,q]=con(x,d,s);if (p<0)b=d;d=(d+a)/2;endif(p>=0)&&(q>=0)dk=d;flag=1;endif (p>=0)&&(q<0)a=d;if 2*d>=(d+b)/2d=(d+b)/2;endendEnd结果：x=[0,1];[f,x,k]=third_2(x)f =0.7729x = -0.4193 0.0001k =8(2) 程序：function [f,x,k]=third_3(x) k=0;X=cell(2);g=cell(2);X{1}=x;H=eye(2);g{1}=gfun(X{1});s=(-H*g{1}')';dk=dfun(X{1},s);X{2}=X{1}+dk*s;g{2}=gfun(X{2});while (norm(g{2})>=0.001)dg=g{2}-g{1};v=dx/(dx*dg')-(H*dg')'/(dg*H*dg');h1=H*dg'*dg*H/(dg*H*dg');h2=dx'*dx/(dx*dx');h3=dg*H*dg'*v'*v;H=H-h1+h2+h3;k=k+1;X{1}=X{2};g{1}=gfun(X{1});s=(-H*g{1}')';dk=dfun(X{1},s);X{2}=X{1}+dk*s;g{2}=gfun(X{2});norm(g{2});f=fun(x);x=X{2};endfunction f=fun(x)f=x(1)+2*x(2)A2+exp(x(1F2+x(2)A2);function gf=gfun(x)gf=[1+2*x(1)*exp(x(1)A2+x(2)A2),4*x(2)+2*x(2)*(x(1)A2+x(2)A2)];function ggf=ggfun(x)ggf=[(4*x(1)A2+2)*exp(x(1)A2+x(2)A2),4*x(1)*x(2)*exp(x(1)A2+x(2)A2);4*x(1)*x(2)* exp(x(1)A2+x(2)A2),4+(4*x(2)A2+2)*exp(x(1)A2+x(2)A2);function [p,q]=con(x,d,s)p=fun(x)-fun(x+d*s)+0.1*d*gfun(x)*(s)';q=gfun(x+d*s)*(s)'-0.5*gfun(x)*(s)';function dk=dfun(x,s)flag=0;a=0;d=1;b=10000;while (flag==0)[p,q]=con(x,d,s);if (p<0)b=d;d=(d+a)/2;if (p>=0)&&(q>=0)dk=d;flag=1;endif (p>=0)&&(q<0)a=d;if 2*d>=(d+b)/2d=(d+b)/2;else d=2*d;endendend结果：x=[0,1];[f,x,k]=third_3(x)f =0.7729x = -0.41950.0000 k=6*U 用有效集法求解下面勺勺二次规划问题:(XI 一 I)2 + (x 2 一 2.5)2 X1 - 2X2 + 2 > 0-Xi — 2>(2 + 6 > 0-Xi + 2X2 + 2 > 0xi,x 2 > 0function callqpactH=[2 0; 0 2];c=[-2 -5]';Ae=[ ]; be=[];Ai=[1 -2; -1 -2; -1 2;1 0;0 1];bi=[-2 -6 -2 0 0]';x0=[0 0]';[x,lambda,exitflag,output]=qpact(H,c,Ae,be,Ai,bi,xO)fun ctio n [x,lamk,exitflag,output]=qpact(H,c,Ae,be,Ai,bi,x0) epsilo n=1.0e-9; err=1.0e-6;k=0; x=x0; n=len gth(x); kmax=1.0e3;n e=le ngth(be); ni=le ngth(bi); lamk=zeros( ne+n i,1); in dex=ones(n i,1);for (i=1:ni)if(Ai(i,:)*x>bi(i)+epsil on), i ndex(i)=0; end while (k<=kmax)mmSi.Aee=[];if (ne>0), Aee=Ae; endfor (j=1:ni)if (index(j)>0), Aee=[Aee; Ai(j,:)]; end endgk=H*x+c;[m1,n1] = size(Aee);[dk,lamk]=qsubp(H,gk,Aee,zeros(m1,1)); if (norm(dk)<=err)y=0.0;if (length(lamk)>ne)[y,jk]=min(lamk(ne+1:length(lamk))); endif (y>=0)exitflag=0;elseexitflag=1;for (i=1:ni)if (index(i)&(ne+sum(index(1:i)))==jk) index(i)=0; break ;endendk=k+1;elseexitflag=1;alpha=1.0; tm=1.0;for (i=1:ni)if ((index(i)==0)&(Ai(i,:)*dk<0))tm1=(bi(i)-Ai(i,:)*x)/(Ai(i,:)*dk);if (tm1<tm)tm=tm1; ti=i;endendendalpha=min(alpha,tm);x=x+alpha*dk;if (tm<1), index(ti)=1; endendif (exitflag==0), break ; endk=k+1;endoutput.fval=0.5*x'*H*x+c'*x;output.iter=k;function [x,lambda]=qsubp(H,c,Ae,be) ginvH=pinv(H); [m,n]=size(Ae);if (m>0)rb=Ae*ginvH*c + be;lambda=pinv(Ae*ginvH*Ae')*rb; x=ginvH*(Ae'*lambda-c);elsex=-ginvH*c;lambda=0;end结果>>callqpactx =1.40001.7000lambda =0.8000exitflag =output =fval: -6.4500iter: 7function [x,mu,lambda,output]=multphr(fu n, hf,gf,dfu n, dhf,dgf,xO)%功能：用乘子法解一般约束问题：min f(x), s.t. h(x)=0, g(x).=0%输入：x0是初始点，fun, dfun分别是目标函数及其梯度；% hf, dhf分别是等式约束(向量)函数及其 Jacobi矩阵的转置；% gf, dgf分别是不等式约束(向量)函数及其 Jacobi矩阵的转置；%输出：x是近似最优点，mu, lambda分别是相应于等式约束和不等式约束的乘子向量% output是结构变量，输出近似极小值f,迭代次数，内迭代次数等maxk=500;c=2.0;eta=2.0;theta=0.8;k=0;i nk=0;epsilo n=0.00001;x=xO;he=feval(hf,x);gi=feval(gf,x);n=len gth(x);l=le ngth(he);m=le ngth(gi);mu=zeros(l,1);lambda=zeros(m,1);btak=10;btaold=10;while (btak>epsilon&&k<maxk)%调用BFGS算法程序求解无约束子问题[x,ival,ik]=bfgs( 'mpsi' ,'dmpsi' ,x0,fun,hf,gf,dfun,dhf,dgf,mu,lambda,c);ink=ink+ik;he=feval(hf,x);gi=feval(gf,x);btak=0;for i=1:lbtak=btak+he(y2;end% 更新乘子向量for i=1:mtemp=min(gi(i),lambda(i)/c);btak=btak+temp A2;endbtak=sqrt(btak);if btak>epsilonif k>=2&&btak>theta*btaoldc=eta*c;endfor i=1:lmu(i)=mu(i)-c*he(i);endlambda(i)=max(0,lambda(i)-c*gi(i));endk=k+1;btaold=btak;x0=x;endendf=feval(fun,x);output.fval=f;output.iter=k;%增广拉格朗日函数function psi=mpsi(x,fun,hf,gf,dfun,dhf,dgf,mu,lambda,c) f=feval(fun,x);he=feval(hf,x);gi=feval(gf,x);l=length(he);m=length(gi);psi=f;s1=0;for i=1:lpsi=psi-he(i)*mu(i);s仁 s1+he(y2;psi=psi+0.5*c*s1;s2=0;for i=1:ms3=max(0,lambda(i)-c*gi(i));s2=s2+s3A2-lambda(i)A2;endpsi=psi+s2/(2*c);% 不等式约束函数文件 g1.mfunction gi=g1(x)gi=10*x(1)-x(1)A2+10*x(2)-x(2)A2-34;% 目标函数的梯度文件df1.mfunction g=df1(x)g=[4, -2*x(2)]';% 等式约束(向量)函数的Jacobi 矩阵(转置)文件 dh1.m function dhe=dh1(x)dhe=[-2*x(1), -2*x(2)]'% 不等式约束(向量)函数的Jacobi 矩阵(转置)文件 dg1.m function dgi=dg1(x)dgi=[10-2*x(1), 10-2*x(2)]';function [x,val,k]=bfgs(fun,gfun,x0,varargin) maxk=500; rho=0.55;sigma=0.4;epsilon=0.00001;k=0;n=length(x0);Bk=eye(n);while (k<maxk)gk=feval(gfun,x0,varargin{:});if (norm(gk)<epsilon)break ;enddk=-Bk\gk;m=0;mk=0;while (m<20)n ewf=feval(fu n, x0+rho A m*dk,vararg in {:});oldf=feval(fun,x0,varargin{:});if(newf<oldf+sigma*rhoAm*gk'*dk) mk=m;break ;endm=m+1;endx=x0+rhoAmk*dk;sk=x-x0;yk=feval(gfun,x,varargin{:})-gk;if (yk'*sk>0)Bk=Bk-(Bk*sk*sk'*Bk)/(sk'*Bk*sk)+(yk*yk')/(yk'*sk);endk=k+1;x0=x;endval=feval(fun,x0,varargin{:});结果x=[2 2]';[x,mu,lambda,output]=multphr( 'fun' ,'hf' ,'gf1' ,'df' ,'dh' ,'dg' ,x0) x =1.00134.8987mu =0.7701lambda =0.9434output =fval: -31.9923iter: 4利用序列二次规划方法求解习题5中的约束优化问题:min 4xi 一好一 12s.t. 25 - x? —x孑=Q10x一召 + 10旳-xj - 34 > 0 X1,X2 > 0tf=[3,1,1]；A=[2,1,1;1,-1,-1];b=[2;-1];lb=[0,0,0]； x=li nprog(f,A,b,zeros(3),[0,0,0]',lb)结果：Optimization terminated.0.00000.50000.5000。
优化设计Matlab程序
进退法步骤：1. 给定初始点(0)x ，初始步长0h ，令(1)(0)0,,0h h x x k ===2.令(4)(1),1x x h k k =+=+ 3.若(4)(1)()()f x f x <，则转4，否则转5 4.(2)(1)(1)(4)(2)(1)(1)(4),,()(),()()x x x x f x f x f x f x ====，令h =2h ，转2 5.若k =1，则转6，否则转,7 6.令h =-h ，(2)(4)(2)(4),()()x x f x f x ==，转2 7. 令(3)(2)(2)(1)(1)(4),,x x x x x x ===，停止计算，极小点包含于区间(1)(3)[,]x x 或(3)(1)[,]x x%目标函数：f%初始点：x0%初始步长：h0%精度：eps%目标函数取包含极值的区间左端点：minx%目标函数取包含极值的区间右端点：maxxfunction [minx,maxx] = minJT(f,x0,h0,eps)format long ;if nargin == 3eps = 1.0e-6;endx1 = x0;k = 0;h = h0;while 1x4 = x1 + h;k = k+1;f4 = subs(f, findsym(f),x4); ! subs : Symbolic substitution in symbolic expression or matrixf1 = subs(f, findsym(f),x1); ! findsym : Determine variables in symbolic expression or matrixif f4 < f1x2 = x1;x1 = x4;f2 = f1;f1 = f4;h = 2*h;elseif k==1h = -h;x2 = x4;f2 = f4;elsex3 = x2;x2 = x1;x1 = x4;break;endendendminx = min(x1,x3); maxx = x1+x3 - minx; format short;syms t;f=t^4-t^2-2*t+5;[x1,x2]=minJT(f,0,0.1)黄金分割法：% 目标函数：f% 极值区间左端点：a% 极值区间右端点：b% 精度：eps% 目标函数取最小值时的自变量值：x % 目标函数的最小值：minffunction [x,minf] = minHJ(f,a,b,eps) format long;if nargin == 3eps = 1.0e-6;endl = a + 0.382*(b-a);u = a + 0.618*(b-a);k=1;tol = b-a;while tol>eps && k<100000fl = subs(f , findsym(f), l);fu = subs(f , findsym(f), u);if fl > fua = l;l = u;u = a + 0.618*(b - a);elseb = u;u = l;l = a + 0.382*(b-a);endk = k+1;tol = abs(b - a);endif k == 100000disp('找不到最小值');x = NaN;minf = NaN;return;endx = (a+b)/2;minf = subs(f, findsym(f),x);format short;抛物线法：% 目标函数：f% 极值区间左端点：a% 极值区间右端点：b% 精度：eps% 目标函数取最小值时的自变量值：x% 目标函数的最小值：minffunction [x,minf] = minPWX(f,a,b,eps)format long;if nargin == 3eps = 1.0e-6;endt0 = (a+b)/2;k = 0;tol = 1;while tol>epsfa = subs(f,findsym(f),a);fb = subs(f,findsym(f),b);ft0 = subs(f,findsym(f),t0);tu = fa*(b^2 - t0^2)+fb*(t0^2 - a^2)+ft0*(a^2 - b^2);td = fa*(b - t0)+fb*(t0 - a)+ft0*(a - b);t1 = tu/2/td;ft1 = subs(f,findsym(f),t1);tol = abs(t1 - t0);if ft1 <= ft0if t1<= t0b = t0;t0 = t1;elsea = t0;t0 = t1;endk = k+1;elseif t1<= t0a = t1;elseb = t1;endk = k+1;endendx = t1;minf = subs(f,findsym(f),x); format short;一维牛顿法：% 目标函数：f% 初始点：x0% 精度：eps% 目标函数取最小值时的自变量值：x % 目标函数的最小值：minffunction [x,minf] = minNewton(f,x0,eps) format long;if nargin == 2eps = 1.0e-6;enddf = diff(f);d2f = diff(df);k = 0;tol = 1;while tol>epsdfx = subs(df,findsym(df),x0);d2fx=subs(d2f,findsym(d2f),x0;x1=x0-dfx/d2fx;k = k + 1;tol = abs(dfx);x0 = x1;endx = x1;minf = subs(f,findsym(f),x);format short;最速下降法：%: 目标函数：f%: 初始点：x0%: 自变量向量：var%: 精度：eps%: 目标函数取最小值时的自变量值：x%: 目标函数的最小值：minffunction [x,minf] = minFD(f,x0,var,eps)format long;if nargin == 3eps = 1.0e-6;endsyms l;tol = 1;gradf = - jacobian(f,var);while tol>epsv = Funval(gradf,var,x0);!Objective function value of the current point tol = norm(v); ! Vector and matrix normsy = x0 + l*v;yf = Funval(f,var,y);[a,b] = minJT(yf,0,0.1); %初始点0，步长为0.1xm = minHJ(yf,a,b);x1 = x0 + xm*v;x0 = x1;endx = x1;minf = Funval(f,var,x);format short;用共轭梯度法求无约束问题minf(x),n x R ∈的算法步骤如下：1. 给定初始点(0)x ，及精度ε2. 若(0)()f x ε∇≤，停止，极小点为(0)x ，否则转3.3. 取(0)(0)()d f x =-∇，且置k=04. 用一维搜索法求k α，使得()()()()()min ()k k k k k f x d f x d αα+=+令(1)()()k k k x x d α+=+，转55. 若(1)()k f x ε+∇≤，停止，极小点为(1)k x +，否则转66. 若k+1=n ，令(0)()n x x =转3，否则转77. 令2(1)(1)(1)()2()()(),()k k k k k k k f x d f x d f x λλ+++∇=-∇+=∇，置1k k =+，转4程序举例：function [x,minf] = minGETD(f,x0,var,eps) format long ;if nargin == 3eps = 1.0e-6;endx0 = transpose(x0);n = length(var);syms l ;gradf = jacobian(f,var);v0 = Funval(gradf,var,x0);d = -transpose(v0);k = 0;while 1v = Funval(gradf,var,x0);tol = norm(v);if tol<=epsx = x0;break;endy = x0 + l*d;yf = Funval(f,var,y);[a,b] = minJT(yf,0,0.1);xm = minHJ(yf,a,b);x1 = x0 + xm*d;vk = Funval(gradf,var,x1);tol = norm(vk);if tol<=epsx = x1;break;endif k+1==nx0 = x1;continue;elselamda = dot(vk,vk)/dot(v,v);d = -transpose(vk) + lamda*d;k = k+1;x0 = x1;endendminf = Funval(f,var,x);format short;用牛顿法求无约束问题，步骤：1. 给定初始点(0)x ，及精度ε2. 若(0)(f x ε∇≤，停止，极小点为(0)x ，否则转3.3. 计算12()()k f x -⎡⎤∇⎣⎦，令2()1()[()]();k k k d f x f x -=-∇∇ 4. 令(1)()(),1k k k x x d k k +=+=+，转2程序举例：function [x,minf] = minNT(f,x0,var,eps) format long ;if nargin == 3eps = 1.0e-6;endtol = 1;x0 = transpose(x0);gradf = jacobian(f,var);jacf = jacobian(gradf,var);while tol>epsv = Funval(gradf,var,x0); tol = norm(v);pv = Funval(jacf,var,x0);p = -inv(pv)*transpose(v);p = double(p);x1 = x0 + p;x0 = x1;endx = x1;minf = Funval(f,var,x);format short;syms t s;f=(t-4)^2+(s+2)^2+1;[x, mf]=minNewton(f, [0 0],[t s])。
优化设计matlab
第5章 优化问题5.1 线性规划问题线性规划问题是目标函数和约束条件均为线性函数的问题，MATLAB6.0解决的线性规划问题的标准形式为：min n R x x f ∈' sub.to ：b x A ≤⋅beq x Aeq =⋅ ub x lb ≤≤其中f 、x 、b 、beq 、lb 、ub 为向量，A 、Aeq 为矩阵。
其它形式的线性规划问题都可经过适当变换化为此标准形式。
在MA TLAB6.0版中，线性规划问题（Linear Programming ）已用函数linprog 取代了MATLAB5.x 版中的lp 函数。
当然，由于版本的向下兼容性，一般说来，低版本中的函数在6.0版中仍可使用。
函数 linprog格式 x = linprog(f,A,b) %求min f ' *x sub.to b x A ≤⋅线性规划的最优解。
x = linprog(f,A,b,Aeq,beq) %等式约束beq x Aeq =⋅，若没有不等式约束b x A ≤⋅，则A=[ ]，b=[ ]。
x = linprog(f,A,b,Aeq,beq,lb,ub) %指定x 的范围ub x lb ≤≤，若没有等式约束beq x Aeq =⋅ ，则Aeq=[ ]，beq=[ ]x = linprog(f,A,b,Aeq,beq,lb,ub,x0) %设置初值x0x = linprog(f,A,b,Aeq,beq,lb,ub,x0,options) % options 为指定的优化参数 [x,fval] = linprog(…) % 返回目标函数最优值，即fval= f ' *x 。
[x,lambda,exitflag] = linprog(…) % lambda 为解x 的Lagrange 乘子。
[x, lambda,fval,exitflag] = linprog(…) % exitflag 为终止迭代的错误条件。
最优化课程设计matlab
作业：1.用薄钢板制造一体积5m 3，长度不小于4m ，无上盖的货箱，要求钢板耗量最小。
确定货箱的长x 1、宽x 2和高x 3。
试列出问题的数学模型。
解：假设钢板为单位厚度，不考虑钢板焊接或连接耗钢量及相关劳工费用。
设y 为货箱的钢板消耗量，则此问题的数学模型如下：min y= x 1 x 2+2 x 1 x 3+2 x 2 x 3s.t. x 1 x 2 x 3=5x 1>4 x 1 ,x 2 ,x 3>02.将下面的线性规划问题表示为标准型并用单纯形法求解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 解：先化为标准型：首先，将第二个约束条件两边乘以（－1），再分别对三个约束不等式添加非负的松弛变量654,,x x x ，即可化为如下标准型：3212min x x x f ---=⎪⎪⎩⎪⎪⎨⎧=≥=+++=++-=+-+6,5,4,3,2,1,06465222..632153214321i x x x x x x x x x x x x x t s i 列成表格：2 ○1 -1 1 0 0 2 2 -1 5 0 1 0 6 4 1 1 0 0 1 6 -1 -2* -1 0 0 0 0可见此表已具备1、2、3三个特点。
首先从底行中选元素-2，再在第二列三个元素中，由2/1,6/1,6/1最小者决定选第一行第二列的元素1，标以记号，迭代一次得2 1 -1 1 0 0 2 4 0 ④ 1 1 0 8 2 0 2 -1 0 1 43 0 -3* 2 0 0 4再从底行中选元素-34再迭代一次得3 1 0 5/4 1/4 0 4 1 0 1 1/4 1/4 0 2 0 0 0 -3/2 -1/2 1 0 6 0 0 11/4 3/4 0 10此时，所有的检验数均为正，停止迭代，最优解为：()0,0,0,2,4,0*=x ，最优值为：102142*-=⨯-+⨯-=Z ；3. 试用DFP 变尺度法求解下列无约束优化问题。
优化方法MATLAB编程——大连理工大学
The optimal solution is 0.000000.
The optimal "x" is "ans".
ans =
1.0000 1.0000 1.0000 1.0000 可以看出，用Newton法经过6次迭代就能求出最优解。
3. BFGS法 源程序如下： function zy_x=di2tiBFGS(x) %该函数用来解大作业第二题。 x0=x; yimuxulong=1e-5; k=0; g0=g(x0); H0=eye(4);s0=-H0*g0;
3
s0=-g0;
大连理工大学优化方法上机大作业
end end x=x0+lanmed*s0; x0=x; g0=g(x); s0=-g0; k=k+1; end end zy_x=x; zyj=f(x); fprintf('after %d iterations,obtain the optimal solution.\n\nThe optimal solution is %f.\n\n The optimal "x" is "ans".\n',k,zyj);
4
大连理工大学优化方法上机大作业
>> x=[-1.2 1 -1.2 1]'; >> di2titidu(x) after 5945 iterations,obtain the optimal solution.
The optimal solution is 0.000000.
The optimal "x" is "ans".
大连理工大学优化方法上机大作业
最优化方法的Matlab实现(公式(完整版))【范本模板】
第九章最优化方法的Matlab实现在生活和工作中，人们对于同一个问题往往会提出多个解决方案,并通过各方面的论证从中提取最佳方案。
最优化方法就是专门研究如何从多个方案中科学合理地提取出最佳方案的科学。
由于优化问题无所不在，目前最优化方法的应用和研究已经深入到了生产和科研的各个领域，如土木工程、机械工程、化学工程、运输调度、生产控制、经济规划、经济管理等,并取得了显著的经济效益和社会效益。
用最优化方法解决最优化问题的技术称为最优化技术，它包含两个方面的内容：1）建立数学模型即用数学语言来描述最优化问题。
模型中的数学关系式反映了最优化问题所要达到的目标和各种约束条件。
2）数学求解数学模型建好以后，选择合理的最优化方法进行求解。
最优化方法的发展很快，现在已经包含有多个分支,如线性规划、整数规划、非线性规划、动态规划、多目标规划等.9。
1 概述利用Matlab的优化工具箱，可以求解线性规划、非线性规划和多目标规划问题。
具体而言,包括线性、非线性最小化，最大最小化，二次规划，半无限问题，线性、非线性方程(组）的求解，线性、非线性的最小二乘问题。
另外，该工具箱还提供了线性、非线性最小化,方程求解，曲线拟合，二次规划等问题中大型课题的求解方法，为优化方法在工程中的实际应用提供了更方便快捷的途径。
9.1.1 优化工具箱中的函数优化工具箱中的函数包括下面几类：1．最小化函数表9—1 最小化函数表2．方程求解函数表9—2 方程求解函数表3．最小二乘（曲线拟合）函数表9-3 最小二乘函数表4．实用函数表9—4 实用函数表5．大型方法的演示函数表9—5 大型方法的演示函数表6．中型方法的演示函数表9—6 中型方法的演示函数表9.1.3 参数设置利用optimset函数，可以创建和编辑参数结构;利用optimget函数,可以获得opt ions优化参数。
● optimget函数功能：获得options优化参数。
val = optimget（options，’param’)val = optimget（options，’param’,default)描述:val = optimget(options,'param’) 返回优化参数options中指定的参数的值。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
优化设计
无约束优化
min f(x)= 21x +22x -21x 2x -41x
初选x0=[1,1]
程序：
Step 1: Write an M-file objfun1.m.
function f1=objfun1(x)
f1=x(1)^2+2*x(2)^2-2*x(1)*x(2)-4*x(1);
Step 2: Invoke one of the unconstrained optimization routines
x0=[1,1];
>> options = optimset('LargeScale','off');
>> [x,fval,exitflag,output] = fminunc(@objfun1,x0,options)
运行结果：
x =
4.0000 2.0000
fval =
-8.0000
exitflag =
1
output =
iterations: 3
funcCount: 12
stepsize: 1
firstorderopt: 2.3842e-007
algorithm: 'medium-scale: Quasi-Newton line search'
message: [1x85 char]
非线性有约束优化
1. Min f(x)=321x +2
2x +21x -32x +5
Subject to: 1g (x)=1x +2x +18≤0
2g (x)=51x -32x -25≤0
3g (x)=131x -412
2x 0≤
4g (x)=14≤1x 130≤
5g (x)=2≤2x 57≤
初选x0=[10,10]
Step 1: Write an M-file objfun2.m
function f2=objfun2(x)
f2=3*x(1)^2+x(2)^2+2*x(1)-3*x(2)+5;
Step 2: Write an M-file confun1.m for the constraints.
function [c,ceq]=confun1(x)
% Nonlinear inequality constraints
c=[x(1)+x(2)+18;
5*x(1)-3*x(2)-25;
13*x(1)-41*x(2)^2;
14-x(1);
x(1)-130;
2-x(2);
x(2)-57];
% Nonlinear inequality constraints
ceq=[];
Step 3: Invoke constrained optimization routine
x0=[10,10]; % Make a starting guess at the solution
>> options = optimset('LargeScale','off');
>> [x, fval] = ...
fmincon(@objfun2,x0,[],[],[],[],[],[],@confun1,options)
运行结果：
x =
3.6755 -7.0744
fval =
124.1495
2. min f （x ）=2
22154x x +
s.t. 0632)(211≤-+=x x x g
01)(212≥+=x x x g 初选x0=[1,1]
Step 1: Write an M-file objfun3.m
function f=objfun3(x)
f=4*x(1)^2+5*x(2)^2
Step 2: Write an M-file confun3.m for the constraints.
function [c,ceq]=confun3(x)
%Nonlinear inequality constraints
c=[2*x(1)+3*x(2)-6;
-x(1)*x(2)-1];
% Nonlinear equality constraints
ceq[];
Step 3: Invoke constrained optimization routine
x0=[1,1];% Make a starting guess at the solution
>> options = optimset('LargeScale','off');
>> [x, fval] = ...
fmincon(@objfun,x0,[],[],[],[],[],[],@confun,options)
运行结果：
Optimization terminated: no feasible solution found. Magnitude of search
direction less than 2*options.TolX but constraints are not satisfied.
x =
1 1
fval =
-13
实例：螺栓连接的优化设计
图示为一压气机气缸与缸盖连接的示意图。
已知D1=400mm,D2=240mm ，缸内工作压力p=8.5Mpa ，螺栓材料为45Cr ，抗拉强度Mpa 1000b =σ，屈服强度Mpa 320s =σ，拉压疲劳极限Mpa 3301-=σ，许用疲劳安全系数[7.1]S a =，取残余预紧力F F 6.1'
'=，采用铜皮石棉密封垫片，螺栓相对刚度8.0=c K 。
从安全、可靠、经济的角度来选择螺栓的个数n 和螺栓的直径d 。
解：
1. 目标函数
取螺栓组连接经济成本n C 最小为目标。
当螺栓的长度、材料和加工条件一定时，螺栓的总成本与n ，d 值成正比，故本问题优化设计的目标函数为 min 21)(x x nd C x f n ===
由此可见，设计变量为螺栓个数n 和直径d 为
n x [= 1[]x d T = T x ]2
2．约束条件
（1）强度约束条件：螺栓在脉动载荷下工作，因此螺栓组连接须满足疲劳强度条件
()()()[]a
a a S K K S ≥++-+=-min min 122σσψσψσσσσσ 其中a σ为应力幅值;min σ为最小应力；σK 为疲劳极限综合影响系数，取σK =4.4；σψ为应力折算系数，取σψ=0.23. 气缸最大载荷p D P 224
1π=
螺栓最大工作载荷n P F F +=,,0，210max 4
1d F πσ+ 螺栓最小工作载荷n p F F K F F /8.000,
-=-=σ, 21,
min 41D F πσ= 螺栓应力幅值2min
max σσσ-=a
对于普通螺纹，小径d d 85.01=，于是疲劳强度约束条件为
()[]01≤-=a a S S X g
（2）密封约束条件：考虑密封安全，螺栓间距应小于8d ，故密封约束条件为
()08400821
1
2≤-=-=x x d n D X g ππ （3）扳手工作空间约束条件：考虑扳手工作空间，螺距间距应大于2d ，故扳手工作空间约束条件为()0400221
213≤-=-
=x x n D d X g ππ （4）非约束条件 (),011≤-=x X g ()025≤-=x X g
螺栓连接的优化数学模型
综上所述，本问题的数学模型可表达如下
设计变量：[]T
x x x 21,= 目标函数：()21m in x x x f =
约束条件：()5,4,3,2,1,0..=≤u x g t s u
M 文件
function f=stud_obj(x)
f=x(1)*x(2);
global p Ksigam psai sigam_1
p=8.5*le6;
Ksigam=4.4;psai=0.23;sigam_1=330;D2=240;
p=1/4*pi*D2^2*le-6*p;
function [c,ceq]=stud_conl(x)
global p Ksigam psai sigam_1
p,x
F0=(1.6*p+p)/x(1);
d1=0.85*x(2);
A=1/4*pi*d1^2;
F1=F0-0.8*p/x(1);
sigma_max=F0/A;
sigma_min=F1/A;
sigma_1=(sigma_max-sigma_min)/2;
Sa=(2*sigma_1+(Ksigma-psai)*sigma_min)/((Ksigma+psai)*(2*sigma_a+sigm a_min));
c=[1.7-Sa;
400*pi/x(1)-8*x(2);
-400*pi/x(1)+2*x(2);
16-x(1);
-x(2)];
ceq=[];
x0=[7,20];
[x,feval]=fmincon(@stud_obj,x0,[],[],[],[],[],[],@stud_conl)
运行结果：
x =
16.00 28.84
fval =
461.39
根据实际问题的意义取整、标准化：n=16，d=30，经验证n和d的取值满足约束条件。