b=t2;t2=t1;f2=f1;
t1=b-0.618*（b-a）;f1=f（t1）;
else
a=t1;t1=t2;f1=f2;
t2=a+0.618*（b-a）;f2=f（t2）;
end
t=0.5*（b+a）;
k=k+1;
f0=f（t）;
end
fprintf（1,'迭代次数k=%3.0f\n',k）
fprintf（1,'迭代区间—左端a=%3.4f\n',a）
fprintf（1,'试点1坐标值t1=%3.4f\n',t1）
fprintf（1,'函数值f1=%3.4f\n',f（t1））
fprintf（1,'迭代区间—右端b=%3.4f\n',b）
fprintf（1,'试点2坐标值t2=%3.4f\n',t2）
fprintf（1,'函数值f2=%3.4f\n',f（t2））
fprintf（1,'区间中点t=%3.4f\n',t）
fprintf（1,'函数值f0=%3.4f\n',f（t））
运行结果如下：
迭代次数k=13
迭代区间—左端a=0.0000
试点1坐标值t1=0.0036
函数值f1=-0.9767
迭代区间—右端b=0.0093
试点2坐标值t2=0.0058
函数值f2=-0.9679
区间中点t=0.0047
函数值f0=-0.9721
由黄金分割法在初始区间[0，3]求得的极小值点为t=0.0047，极小值为-0.9721。
3.用牛顿法、阻尼牛顿法及变尺度法求函数
的极小点。
（1）在用牛顿法在MATLAB的M文件编辑器中编写的M文件，如下：
function[x,fx,k]=niudunfa（x0）
symsx1x2
f=（x1-2）^4+（x1-2*x2）^2;
fx=0;
v=[x1,x2];
df=jacobian（f,v）;
df=df.';
G=jacobian（df,v）;
epson=1e-12;
g1=subs（df,{x1,x2},{x0（1,1）,x0（2,1）}）;
G1=subs（G,{x1,x2},{x0（1,1）,x0（2,1）}）;
k=0;
p=-G1\g1;
x0=x0+p;
while（norm（g1）>epson）
p=-G1\g1;
x0=x0+p;
g1=subs（df,{x1,x2},{x0（1,1）,x0（2,1）}）;
G1=subs（G,{x1,x2},{x0（1,1）,x0（2,1）}）;
k=k+1;
end
x=x0;
fx=subs（f,{x1,x2},{x（1,1）,x（2,1）}）;
运行结果如下：
>> [x,fx,k]=niudunfa（[1;1]）
x =1.897
0.9483
fx =0.7647
k =23
（2）用阻尼牛顿法在MATLAB的M文件编辑器中编写的M文件，如下：
function[x,fx,k]=zuniniudunfa（x0）%阻尼牛顿法
symsx1x2
f=（x1-2）^4+（x1-2*x2）^2;
fx=0;
v=[x1,x2];
df=jacobian（f,v）;
df=df.';
G=jacobian（df,v）;
epson=1e-12;%停机原则
g1=subs（df,{x1,x2},{x0（1,1）,x0（2,1）}）;
G1=subs（G,{x1,x2},{x0（1,1）,x0（2,1）}）;
k=0;%迭代次数
p=-G1\g1;
a0=-p'*g1/（p'*G1*p）;
x0=x0+a0*p;
while（norm（a0*p）>epson）
p=-G1\g1;
a0=-p'*g1/（p'*G1*p）;
x0=x0+a0*p;
g1=subs（df,{x1,x2},{x0（1,1）,x0（2,1）}）;
G1=subs（G,{x1,x2},{x0（1,1）,x0（2,1）}）;
k=k+1;
end
x=x0;
fx=subs（f,{x1,x2},{x0（1,1）,x0（2,1）}）;
运行结果如下：
>>[x,fx,k]=zuniniudunfa（[1;1]）
x=1.897
0.9483
fx=0.7647
k=23
（3）用变尺度法在MATLAB的M文件编辑器中编写的M文件，如下：
4.用共轭梯度法求函数
的极小点
（1）用共轭梯度法在MATLAB的M文件编辑器中编写的M文件，如下：
function[y,x,k]=CG（A,b,c,x0）
%共轭梯度法解minf（x）=0.5*X'*A*X+b'x+c
eps=1e-6;%迭代停机原则
%fx=0.5*x0'.*A.*x0+b'.*x0+c;
r0=A*x0+b;
ifnorm（r0）<=eps
x=x0;
y=0.5*x'*A*x+b'*x+c;
k=0;
end
p0=-r0;
a=-r0'*p0/（p0'*A*p0）;
x1=x0+a*p0;
r1=A*x1+b;
k=0;
whilenorm（r1）>eps
beta=（r1'*r1）/（r0'*r0）;
p1=-r1+beta*p0;
alpha=-（r1'*p1）/（p1'*A*p1）;
x1=x1+alpha*p1;
r2=A*x1+b;
p0=p1;
r0=r1;
r1=r2;
k=k+1;
end
x=x1;
y=0.5*x'*A*x+b'*x+c;
运行结果如下：
[y,x,k]=CG（[3 -1;-1 1],[-2;0],0,[2;1]）
y = -1
x = 1.0000
1.0000
k = 1
（2）用变尺度法在MATLAB的M文件编辑器中编写的M文件，如下：
function[x,fx,k]=bianchidufa（A,b,c,x0）
%用变尺度法求fx=0.5*x'*A*x+b'*x+c;
epson=1e-12;
g0=A*x0+b;
G0=A;
H0=eye
（2）;
k=0;
d0=-H0*g0;
a0=-d0'*g0/（d0'*G0*d0）;
s0=a0*d0;%x（k+1）-x（k）;
y0=A*a0*d0;%g（k+1）-g（k）;
x1=x0+a0*d0;
while（norm（s0）>=epson）
switchk
case{10}
x0=x1;
g0=A*x0+b;
H0=eye
（2）;
k=0;
d0=-H0*g0;
a0=-d0'*g0/（d0'*G0*d0）;
s0=a0*d0;
x1=x0+a0*d0;
break
otherwise
g1=A*x1+b;
y0=A*a0*d0;
s0=a0*d0;
%H1=H0+s0*s0'/（s0'*y0）-H0*y0*y0'*H0/（y0'*H0*y0）;
H1=H0+（（1+y0'*H0*y0/（s0'*y0））*s0*s0'-H0*y0*s0'-s0*y0'*H0）/（s0'*y0）;
k=k+1;
d1=-H1*g1;
a1=-d1'*g1/（d1'*G0*d1）;
a0=a1;
d0=d1;
H0=H1;
s0=a0*d0;
x1=x1+a0*d0;
break
end
end
x=x1;
fx=0.5*x1'*A*x1+b'*x1+c;
运行结果如下：
》 [x,fx,k]=bianchidufa（[3 -1;-1 1],[-2;0],0,[2;1]）
H1 =
0.4031 0.2578
0.2578 0.8945
fx = -1
x =
1.0000
1.0000
fx = -1
k = 1
故函数极小点是点（1,1）
5.用鲍威尔法求函数
的极小点。
用鲍威尔法在MATLAB的M文件编辑器中编写的M文件，如下：
function[x,fx,k]=bowell（A,b,c,x0）%鲍威尔法
d01=[1;0];
d02=[0;1];
x02=[0;0];
esp=1e-12;%停机原则
k=0;%迭代次数
whilenorm（x0-x02）>=esp
k=k+1;
g01=A*x0+b;
a01=-d01'*g01/（d01'*A*d01）;
x01=x0+a01*d01;
g02=A*x01+b;
a02=-d02'*g02/（d02'*A*d02）;
x02=x01+a02*d02;
d10=x02-x0;
g10=A*x02+b;
a10=-d10'*g10/（d10'*A*d10）;
x10=x0+a01*d01;
d01=d02;
d02=d10;
x0=x10;
end
x=x0;
fx=0.5*x'*A*x+b'*x+c;
运行结果如下：
[x,fx,k]=bowell（[2 -2;-2 4],[-4;0],0,[2;1]）
fx =
-8
x =
4
2
fx =
-8
k =
3
6.用单纯形法求线性规划问题
用单纯形法在MATLAB的M文件编辑器中编写的M文件，如下：
%单纯形法matlab程序-danchunxingfa
%求解标准型线性规划:
maxc*x;s.t.A*x=b;x>=0
%本函数中的A是单纯初始表，包括:
最后一行是初始的检验数，最后一列是资源向量b
%N是初始的基变量的下标
%输出变量sol是最优解,其中松弛变量（或剩余变量）可能不为0
%输出变量val是最优目标值，kk是迭代次数
function[sol,val,kk]=danchunxingfa（A,N）
[mA,nA]=size（A）;
kk=0;%迭代次数
flag=1;
whileflag
kk=kk+1;
ifA（mA,:
）<=0%已找到最优解
flag=0;
sol=zeros（1,nA-1）;
fori=1:
mA-1
sol（N（i））=A（i,nA）;
end
val=-A（mA,nA）;
else
fori=1:
nA-1
ifA（mA,i）>0&A（1:
mA-1,i）<=0%问题有无界解
disp（'haveinfinitesolution!
'）;
flag=0;
break;
end
end
ifflag%还不是最优表，进行转轴运算
temp=0;
fori=1:
nA-1
ifA（mA,i）>temp
temp=A（mA,i）;
inb=i;%进基变量的下标
end
end
sita=zeros（1,mA-1）;
fori=1:
mA-1
ifA（i,inb）>0
sita（i）=A（i,nA）/A（i,inb）;
end
end
temp=inf;
fori=1:
mA-1
ifsita（i）>0&sita（i）temp=sita（i）;
outb=i;%出基变量下标
end
end
%以下更新N
fori=1:
mA-1
ifi==outb
N（i）=inb;
end
end
%以下进行转轴运算
A（outb,:
）=A（outb,:
）/A（outb,inb）;
fori=1:
mA
ifi~=outb
A（i,:
）=A（i,:
）-A（outb,:
）*A（i,inb）;
end
end
end
end
end;
运行结果如下：
>>A=[11104;
122.535;
1.12.2-3.34.40];N=[3;4];[sol,val,kk]=danchunxingfa（A,N）
sol=
004.00001.6667
val=
7.3333
kk=
2
所以，
7.求解线性规划问题
用单纯形法在MATLAB的M文件编辑器中编写的M文件，如下：
%单纯形法matlab程序-danchunxingfa
%求解标准型线性规划:
maxc*x;s.t.A*x=b;x>=0
%本函数中的A是单纯初始表，包括:
最后一行是初始的检验数，最后一列是资源向量b
%N是初始的基变量的下标
%输出变量sol是最优解,其中松弛变量（或剩余变量）可能不为0
%输出变量val是最优目标值，kk是迭代次数
function[sol,val,kk]=danchunxingfa（A,N）
[mA,nA]=size（A）;
kk=0;%迭代次数
flag=1;
whileflag
kk=kk+1;
ifA（mA,:
）<=0%已找到最优解
flag=0;
sol=zeros（1,nA-1）;
fori=1:
mA-1
sol（N（i））=A（i,nA）;
end
val=-A（mA,nA）;
else
fori=1:
nA-1
ifA（mA,i）>0&A（1:
mA-1,i）<=0%问题有无界解
disp（'haveinfinitesolution!
'）;
flag=0;
break;
end
end
ifflag%还不是最优表，进行转轴运算
temp=0;
fori=1:
nA-1
ifA（mA,i）>temp
temp=A（mA,i）;
inb=i;%进基变量的下标
end
end
sita=zeros（1,mA-1）;
fori=1:
mA-1
ifA（i,inb）>0
sita（i）=A（i,nA）/A（i,inb）;
end
end
temp=inf;
fori=1:
mA-1
ifsita（i）>0&sita（i）temp=sita（i）;
outb=i;%出基变量下标
end
end
%以下更新N
fori=1:
mA-1
ifi==outb
N（i）=inb;
end
end
%以下进行转轴运算
A（outb,:
）=A（outb,:
）/A（outb,inb）;
fori=1:
mA
ifi~=outb
A（i,:
）=A（i,:
）-A（outb,:
）*A（i,inb）;
end
end
end
end
end;
运行结果如下：
>>A=[94100360;
45010200;
310001300;
7120000];N=[3;4;5];[sol,val,kk]=danchunxingfa（A,N）
sol=
20248400
val=
420
kk=
3
所以，经3次转轴运算，得到的最优解为
升级会员 