P90 第二章 11（2）
用大M法求解
min w=2x1+x2-x3-x4
s.t x1-x2+2x3-x4=2
2x1+x2-3x3+x4=6
x1+x2+x3+x4=7
xi≥0, i = 1, 2, 3, 4
用matlab求解如下：
f=[2,1,-1,-1,200,200,200];
a=[1 -1 2 -1 1 0 0;2 1 -3 1 0 1 0;1 1 1 1 0 0 1];
b=[2 6 7];
lb=zeros(7,1);
[x,fval,exitflag,output,lambda]=linprog(f,[],[],a,b,lb)
Optimization terminated.
运行结果如下：
x =
3.0000
0.0000
1.0000
3.0000
0.0000
0.0000
0.0000
fval =
2.0000
exitflag =
1
output =
iterations: 7
algorithm: 'large-scale: interior point'
cgiterations: 0
message: 'Optimization terminated.'
lambda =
ineqlin: [0x1 double]
eqlin: [3x1 double]
upper: [7x1 double]
lower: [7x1 double]
从上述运行结果可以得出：最优解为x= ，最小值约为f*=2。
P151 第三章 26
用共轭梯度算法求f(x) = (x1-1)^2+5*(x2-x1^2)^2的极小点，取初始点x0= 。
用matlab求解如下：
function mg=MG()
%
%共轭梯度法求解习题三第26题
%
clc;
clear;
n=2;
x=[2 0]';
max_k=100;
count_k=1;
trace(1,1)=x(1);
trace(2,1)=x(2);
trace(3,1)=f_fun(x);
k=0;
g1=f_dfun(x);
s=-g1;
while count_k<=max_k
if k==n
g0=f_dfun(x);
s=-g0;
k=0;
else
r_min=fminbnd(@(t) f_fun(x+t*s),-100,100);
x=x+r_min*s;
g0=g1;
g1=f_dfun(x);
if norm(g1)<10^(-6)
break;
end
m=(norm(g1)^2)/(norm(g0)^2);
s=-g1+m*s;
count_k=count_k+1;
trace(1,count_k)=x(1);
trace(2,count_k)=x(2);
trace(3,count_k)=f_fun(x);
k=k+1;
end
end
count_k
x
f=f_fun(x)
function g=f_dfun(x)
g(1,1)=20*x(1)^3-20*x(1)*x(2)+2*x(1)-2;
g(2,1)=10*x(2)-10*x(1)^2;
function f=f_fun(x)
f=5*x(1)^4-10*x(1)^2*x(2)+x(1)^2-2*x(1)+1+5*x(2)^2;
运行结果如下：
k =
59
x0 =
1.0000
1.0000
f0 =
0
从上述运行结果可以得出：最优解为x= ，最小值约为f*=0。
P151 第三章 26
用BFGS算法求f(x) = (x1-1)^2+5*(x2-x1^2)^2的极小点，取初始点x0= 。
用matlab求解如下：
syms x1 x2;
f=(x1-1)^2+5*(x2-x1^2)^2;
v=[x1,x2];
df=jacobian(f,v);
df=df.';
x0=[2,0]';g1=subs(df,{x1,x2},{x0(1,1),x0(2,1)});k=0;
H0=[1,0;0,1];
while(norm(g1)>0.001 & k<300)
if k==0
p=-H0*g1;
else
vk=sk/(sk'*yk)-(H0*yk)/(yk'*H0*yk);
w1=(yk'*H0*yk)*vk*vk';
H1=H0-(H0*yk*yk'*H0)/(yk'*H0*yk)+(sk*sk')/(sk'*yk)+w1;
p=-H1*g1;
H0=H1;
end
x00=x0;
result=Usearch1(f,x1,x2,df,x0,p);
arf=result(1);
x0=x0+arf*p;
g0=g1;
g1=subs(df,{x1,x2},{x0(1,1),x0(2,1)});
p0=p;
yk=g1-g0; sk=x0-x00;
k=k+1;
end;
k
x0
f0=subs(f,{x1,x2},{x0(1,1),x0(2,1)})
function result=Usearch1(f,x1,x2,df,x0,p)
mu=0.001;sgma=0.99;a=0;b=inf;arf=1;
pk=p;
x3=x0;
x4=x3+ar
f*pk;
f1=subs(f,{x1,x2},{x3(1,1),x3(2,1)});
f2=subs(f,{x1,x2},{x4(1,1),x4(2,1)});
gk1=subs(df,{x1,x2},{x3(1,1),x3(2,1)});
gk2=subs(df,{x1,x2},{x4(1,1),x4(2,1)});
while (f1-f2<=-mu*arf*gk1'*pk)
b=arf;arf=(a+arf)/2;x4=x3+arf*pk;f2=subs(f,{x1,x2},{x4(1,1),x4(2,1)});gk2=subs(df,{x1,x2},{x4(1,1),x4(2,1)});
end;
while (1>0)
if(gk2'*pk
while (f1-f2<=-mu*arf*gk1'*pk)
b=arf;arf=(a+arf)/2;x4=x3+arf*pk;f2=subs(f,{x1,x2},{x4(1,1),x4(2,1)});gk2=subs(df,{x1,x2},{x4(1,1),x4(2,1)});
end;
else
break;
end
end;
result=[arf];
运行结果如下：
k =
19
x0 =
1.0000
1.0000
f0 =
4.4505e-013
从上述运行结果可以得出：最优解为x= ，最小值约为f*=0。
P229 第四章 10
用SUMT内点法求解
(1) min f(x) = x1+x2
s.t g1(x) = -x1^2+x2≥0
g2(x) = x1≥0
用matlab求解如下：
function sumt=SUMT(x0,e0,max_k0)
%SUMT内点法
global x s M
x=x0;
e=e0;
max_k=max_k0;
trace(1,1)=x(1);
trace(2,1)=x(2);
M=100;
c=3;
e_FR=10^-10;
max_FR=200;
for k=0:max_k
x=FR(x,e_FR,max_FR);
trace(1,k+2)=x(1);
trace(2,k+2)=x(2);
if f_pfun(x)
end
M=c*M;
end
x
f=f_fun(x)
k
trace
function f=FR(x0,e,max_k)
global x s;
count_k=1;
k=0;
x=x0;
g1=f_dfun(x);
s=-g1;
while count_k<=max_k
if k==n
g0=f_dfun(x);
s=-g0;
k=0;
else
r_min=fminbnd(@(t) f_fun(x+t*s),-100,100);
x=x+r_min*s;%
g0=g1;
g1=f_dfun(x)
if norm(g1)
end
m=(norm(g1)^2)/(norm(g0)^2);
s=-g1+m*s;
count_k=count_k+1;
k=k+1;
end
end
count_k=count_k-1
f=x
注意：以下三个函数体在运行计算习题四第8(1)题时使用
%*****************************************************
function f=f_fun(x)
global M;
f=x(1)+x(2)+M*f_pfun(x);%题8(1)
function g=f_dfun(x)
global M;
a=-x(1)^2+x(2);
b=x(1);
if a>=0&b>=0
g(1,1)=1;
g(2,1)=1;
elseif a<0&b>=0
g(1,1)=1+2*M*(-x(1)^2+x(2))*(-2*x(1));
g(2,1)=1+2*M*(-x(1)^2+x(2));
elseif a>=0&b<0
g(1,1)=1+2*M*x(1);
g(2,1)=1;
else
g(1,1)=1+2*M*(-x(1)^2+x(2))*(-2*x(1))+2*M*x(1);
g(2,1)=1+2*M*(-x(1)^2+x(2));
end%题8(1)
function p=f_pfun(x)
%罚函数
a=-x(1)^2+x(2);
b=x(1);
if a>=0&b>=0
p=0;
elseif a>=0&b<0
p=b^2;
elseif a<0&b>=0
p=a^2;
else
p=a^2+b^2;
end%题8(1)
运行结果如下：
x =
1.0e-004 *
-0.2058
-0.2058
f =
-2.0576e-005
k =
5
trace =
100.0000 -0.0050 -0.0017 -0.0006 -0.0002 -0.0001 -0.0000
3.0000 -0.0050 -0.0017 -0.0006 -0.0002 -0.0001 -0.0000
从上述运行结果可以得出：最优解为x= ，最小值约为f*=0。
(2) min f(x) = x1^2+x2^2
S.t 2x1+x2
-2≤0
-x2+1≤0
用matlab求解如下：
主程序与（1）相同；
注意：以下三个函数体在运行计算习题四第8(2)题时使用
%*****************************************************
function g=f_dfun(x)
global M;
a=-2*x(1)-x(2)+2;
b=x(2)-1;
if a>=0&b>=0
g(1,1)=2*x(1);
g(2,1)=2*x(2);
elseif a<0&b>=0
g(1,1)=2*x(1)+2*M*(-2*x(1)-x(2)+2)*(-2);
g(2,1)=2*x(2)+2*M*(-2*x(1)-x(2)+2)*(-1);
elseif a>=0&b<0
g(1,1)=2*x(1);
g(2,1)=2*x(2)+2*M*(x(2)-1);
else
g(1,1)=2*x(1)+2*M*(-2*x(1)-x(2)+2)*(-2);
g(2,1)=2*x(2)+2*M*(-2*x(1)-x(2)+2)*(-1)+2*M*(x(2)-1);
end%题8(2)
%*****************************************************
function f=f_fun(x)
global M;
f=x(1)^2+x(2)^2+M*f_pfun(x);%题8(2)
function p=f_pfun(x)
a=-2*x(1)-x(2)+2;
b=x(2)-1;
if a>=0&b>=0
p=0;
elseif a>=0&b<0
p=b^2;
elseif a<0&b>=0
p=a^2;
else
p=a^2+b^2;
end%题8(2)
运行结果如下：
x =
-0.0000
1.0000
f =
1.0000
k =
6
trace =
100.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000
3.0000 0.9901 0.9967 0.9989 0.9996 0.9999 1.0000 1.0000
从上述运行结果可以得出：最优解为x= ，最小值约为f*=1。
P232 第四章 25
用梯度投影法求解下列线性约束优化问题
min f(x) = x1^2+x1*x2+2*x2^2+2*x3^2+2*x2*x3+4*x1+6*x2+12*x3
s.t x1 + x2 + x3 ≤6
﹣x1﹣x2 +2x3 ≥2
xi≥0 , i = 1, 2, 3
取x1= ,
用matlab求解如下：
f=‘x1^2+x1*x2+2*x2^2+2*x3^2+2*x2*x3+4*x1+6*x2+12*x3’;
a=[1 1 1;1 1 -2];
b=[6;-2];
l=zeros(3,1);
x0=[1 1 3];
[x,fval,exitflag,output,lambda,grad,hessian]=fmincon(f,x0,a,b,[],[],l,[])
运行结果如下：
x =
0 0 1
fval =
14
exitflag =
1
output =
iterations: 2
funcCount: 14
stepsize: 1
algorithm: 'medium-scale: SQP, Quasi-Newton, line-search'
firstorderopt: 0
cgiterations: []
message: [1x144 char]
lambda =
lower: [3x1 double]
upper: [3x1 double]
eqlin: [0x1 double]
eqnonlin: [0x1 double]
ineqlin: [2x1 double]
ineqnonlin: [0x1 double]
grad =
4.0000
8.0000
16.0000
hessian =
1.1146 0.6771 0.6042
0.6771 3.3646 2.4792
0.6042 2.4792 3.4583
从上述运行结果可以得出：最优解为x= ，最小值约为f*=14。
P234 第四章 34（1）
用乘子法求解
max f(x) = 10x1+4.4x2^2 +2x3
s.t x1+4x2+5x3≤32
x1+3x2+2x3≤29
x3^2/2+x2^2≥3
x1≥2, x2≥0, x3≥0
用matlab求解如下：
function [x,minf] = ymh434(l)
format long;
syms x1 x2 x3
f=10*(x1)+4.4*(x2)^2+2*(x3);
h=[-x1-4*x2-5*x3+32,-x1-3*x2-2*x3+29,((x3)^2)/2+(x2)^2-3,x1-2,x2,x3];
x0=[2 2 0];
v=[1 0 0 0 0 0];
M=2;
alpha=2;
gama=0.25;
var=[x1 x2 x3];
eps=1.0e-4;
if nargin == 8
eps = 1.0e-4;
end
m1 = transpose(x0);
m2 = inf;
while l
FE = 0;
u=subs(h,{
x1,x2,x3},[m1(1),m1(2),m1(3)]);
for i=1:length(h)
if (v(i)+M*u(i))<=0
FE = FE +((v(i)+M*u(i))^2-(v(i))^2);
else
FE=FE+(v(i))^2;
end
end
SumF = f + (1/(2*M))*FE;
[m2,minf] = minNT(SumF,transpose(m1),var,eps);
Hm2 =subs(h,{x1,x2,x3},[m2(1),m2(2),m2(3)]);
Hm1 =subs(h,{x1,x2,x3},[m1(1),m1(2),m1(3)]);
Hx2 = Funval(h,var,x2);
Hx1 = Funval(h,var,x1);
if norm(Hx2) < eps
x = x2;
break;
else
if Hx2/Hx1 >= gama
M = alpha*M;
x1 = x2;
else
v = v - M*transpose(Hx2);
x1 = x2;
end
end
end
minf = Funval(f,var,x);
format short;
运行结果如下：
x =
19.75 0 2.45
minf =
-202.42
maxf1 =
202.42
从上述运行结果可以得出：最优解为x= ，最大值约为f*=202.42, 最小值约为f*=-202.42。
P235 第四章 35（2）
用序列二次规划法求解
min f(x)=-5x1-5x2-4x3-x1x3-6x4-5x5/(1+x5)-8x6/(1+x6)-10(1-2e^-x7+e^-2x7)
s.t g1(x)=2x4+x5+0.8x6+x7-5=0
g2(x)=x2^2+x3^2+x5^2+x6^2-5=0
g3(x)=x1+x2+x3+x4+x5+x6+x7≤10
g4(x)=x1+x2+x3+x4≤5
g5(x)=x1+x3+x5+x6^2-x7^2-5≤0
xi≥0, i=1,2,…,7
用matlab求解如下：
function f=objfun35(x)
f=-5*x(1)-5*x(2)-4*x(3)-x(1)*x(3)-6*x(4)-5*x(5)/(1+x(5))-8*x(6)/(1+x(6))-10*(1-2*exp(-x(7))+exp(-2*x(7)));
function [c,ceq] = confun35(x)
ceq=[x(2)^2+x(3)^2+x(5)^2+x(6)^2-5];
c=[x(1)+x(3)+x(5)+x(6)^2-x(7)^2-5];
clc
clear
x0=[1,1,1,1,1,1,1];
aeq=[0,0,0,2,1,0.8,1];
beq=5;
a=[1,1,1,1,1,1,1;1,1,1,1,0,0,0];
b=[10,5]';
lb=[0,0,0,0,0,0,0];
ub=[];
[x,fval,exitflag,output,lambda,grad,hessian]=fmincon('objfun35',x0,a,b,aeq,beq,lb,ub,'confun35',options)
[c,ceq]=confun35(x)
运行结果如下：
Iter F-count f(x) constraint Step-size derivative optimality Procedure
0 8 -31.4958 1 Infeasible start point
1 17 -41.7175 0.8642 1 -6.3 1.62
2 26 -43.0309 0.5358 1 -1.09 0.924
3 35 -44.5199 0.8348 1 -1.4 0.438
4 44 -44.4543 0.04621 1 0.126 0.247
5 53 -44.4636 0.004354 1 -0.00456 0.127
6 62 -44.4697 0.005324 1 -0.0007 0.00717
7 71 -44.4687 1.631e-005 1 0.000984 0.000523
8 80 -44.4687 6.735e-009 1 3.03e-006 7.25e-006
9 89 -44.4687 4.308e-013 1 1.25e-009 4.11e-007
x = 3.2418 0.0000 1.6342 0.1240 0.8896 1.2402 2.8702
fval = -44.4687
……
hessian =
0.5738 0.3772 -0.2369 0.1036 -0.1250 -0.0757 -0.1739
0.3772 0.6898 0.4769 0.0641 -0.0816 -0.0879 -0.0803
-0.2369 0.4769 1.4
590 0.0628 0.0657 0.1229 -0.0055
0.1036 0.0641 0.0628 0.8811 0.0721 0.0120 0.0699
-0.1250 -0.0816 0.0657 0.0721 1.7848 -0.2129 -0.0244
-0.0757 -0.0879 0.1229 0.0120 -0.2129 1.4472 -0.1761
-0.1739 -0.0803 -0.0055 0.0699 -0.0244 -0.1761 1.0268
c = -5.9342
ceq = 4.3077e-013
从上述运行结果可以得出：
最优解为x= ，最小值约为 f*=-44.4687。