中科大《优化设计》课程大作业之无约束优化实验报告.docx
《中科大《优化设计》课程大作业之无约束优化实验报告.docx》由会员分享，可在线阅读，更多相关《中科大《优化设计》课程大作业之无约束优化实验报告.docx（27页珍藏版）》请在冰点文库上搜索。
中科大《优化设计》课程大作业之无约束优化实验报告
无约束优化设计实验报告
力学系
型号：
联想y470
CPU：
i5-2450M
内存：
2GB
系统：
win7-64位
如下是三个目标函数（包括自定义函数）以及初值和精度选取：
1.minf（x）=x1^2+2*x2^2-2*x1*x2-4*x1
初值x0=[1,1]’；精度为：
0.1
2.minf（x,y）=x^4-2*x*x*y-2*x*y+x^2+2*y*y+4.5*x-4*y+4
初值为（-2.5,4.25）；精度为0.01
3.minf（x）=x1^2+x2^2+x3^2
初值为（3,2,1）；精度为0.01
如下是运算结果：
目标函数
无约束方法
一维搜索
所需时间
迭代次数
极值点
极值
1
最速下降法
黄金分割法
5.81
8
（3.9,1.9）
-8.0
牛顿法
2.61
8
（3.9,1.9）
-8.0
不精确法
2.96
8
（3.9,1.9）
-8.0
阻尼牛顿法
黄金分割法
1.66
1
（4.0,2.0）
-8.0
牛顿法
1.02
1
（4.0,2.0）
-8.0
不精确法
0.76
1
（4.0,2.0）
-8.0
共轭梯度法
黄金分割法
14.65
24
（3.9,2.0）
-8.0
牛顿法
7.91
24
（4.0,2.0）
-8.0
不精确法
11.28
27
（3.9,2.0）
-8.0
鲍维尔法
黄金分割法
9.53
2
（4.0,2.0）
-8.0
牛顿法
2.86
2
（4.0,2.0）
-8.0
不精确法
--
--
--
--
变尺度法
黄金分割法
8.04
2
（4.0,2.0）
-8.0
牛顿法
1.24
2
（4.0,2.0）
-8.0
不精确法
1.35
3
（4.0,2.0）
-8.0
单形替换法
无
0.02
9
（4.1,2.2）
-8.0
目标函数
无约束方法
一维搜索
所需时间
迭代次数
极值点
极值
2
最速下降法
黄金分割法
20.76
21
（1.90,3.73）
0.99
牛顿法
3.23
6
（1.93,3.82）
0.99
不精确法
10.49
27
（1.94,3.84）
0.99
阻尼牛顿法
黄金分割法
3.50
3
（-1.05,1.03）
-0.51
牛顿法
1.61
3
（-1.05,1.03）
-0.51
不精确法
1.58
4
（-1.05,1.03）
-0.51
共轭梯度法
黄金分割法
25.06
30
（-1.05,1.03）
-0.51
牛顿法
24.40
35
（-1.05,1.03）
-0.51
不精确法
19.80
33
（1.94,3.86）
0.99
鲍维尔法
黄金分割法
11.49
3
（-1.05,1.03）
-0.51
牛顿法
6.17
3
（-1.05,1.03）
-0.51
不精确法
4.75
1
（-1.50,4.25）
17.31
变尺度法
黄金分割法
3.94
3
（1.94,3.85）
0.99
牛顿法
2.39
3
（1.94,3.85）
0.99
不精确法
2.46
6
（1.94,3.84）
0.99
单形替换法
无
0.01
16
（1.92,3.81）
0.99
目标函数
无约束方法
一维搜索
所需时间
迭代次数
极值点
极值
3
最速下降法
黄金分割法
2.70
2
（0,0,0）
0
牛顿法
0.76
1
（0,0,0）
0
不精确法
0.78
1
（0,0,0）
0
阻尼牛顿法
黄金分割法
1.61
1
（0,0,0）
0
牛顿法
0.66
1
（0,0,0）
0
不精确法
0.65
1
（0,0,0）
0
共轭梯度法
黄金分割法
3.15
3
（0,0,0）
0
牛顿法
0.40
0
（0,0,0）
0
不精确法
0.44
0
（0,0,0）
0
鲍维尔法
黄金分割法
7.17
1
（0,0,0）
0
牛顿法
1.99
1
（0,0,0）
0
不精确法
4.08
0
（3.00,2.00,1.00）
14.00
变尺度法
黄金分割法
2.69
2
（0,0,0）
0
牛顿法
0.67
1
（0,0,0）
0
不精确法
0.77
1
（0,0,0）
0
单形替换法
无
0.01
24
（-0.01,0.03,-0.08）
0.01
总结及比较：
根据上面三个函数的表格可以看出：
首先，从迭代时间来看，三种一维搜索方法中黄金分割法所用时间最久，牛顿法和不精确法所用时间较少，这两种方法相比较而言牛顿法所用时间更少一些。
而六种无约束方法中，由于单形替代法不需要使用一维搜索方法，故迭代时间最少，紧接着在使用一维搜索的五种方法中以阻尼牛顿法迭代时间相对较少，共轭梯度法迭代时间最久；然后，从迭代次数来看，共轭梯度法往往需要较多的迭代次数，从而所需时间也最久；接着，从计算结果的精度来看，阻尼牛顿法的结果精度最高，而单形替换法的精度最低；最后，从编程来看，在编好一维搜索方法的情况下，最速下降法和阻尼牛顿法编程简单容易，而共轭梯度法、变尺度法和单形替代法需要两重循环实现，鲍威尔法和单形替换法则需要编程者对矩阵的操作能力有较高的要求，故编程较难。
同时，从上面的结果也可以发现，鲍威尔法在使用不精确的一维搜索方法时，对函数1无法收敛，对函数2、3收敛到错误的结果，所以鲍威尔法是依赖于精确的一维搜索过程的，而其他几个则不依赖于精确一维搜索过程。
精确的一维搜索方法通常需要花费很大的工作量，特别是当迭代点远离问题的解时，精确的求解一个一维子问题通常不是十分有效率的。
因此，只要保证目标函数值在每一步都有满意的下降，使用不是非常精确的一维搜索，就可以大大节省工作量。
在分析函数2的计算结果时，可以发现存在两个收敛结果，当然这两个结果都是极值，因为函数2是二元四次函数，存在多个极值。
故为了验证正确性，自己曾将初始点（-2.5,4.25）调成（2.5,4.25），分别代入程序中计算，计算结果都收敛于极值为0.99的这个点上。
所以，在存在多个极值点的情况下结果是和初始点的选取有关。
对于单形替换法，这种方法不需要一维搜索最佳步长，故没有一维搜索方法反复地计算最佳步长的计算时间，程序运行效率快。
但它的收敛条件不好选择，通过查找文献资料总结出以下三个收敛条件：
1.利用最坏点函数值与最好点函数值之差判别；2.利用相邻两次函数值差值的绝对值判别；3.利用各点的函数值与最好点函数值之差的均方根判别。
为了保证程序的执行可靠性，这三种常用的方法中自己选择了判别3，即：
综上所述，阻尼牛顿法是无约束方法中最有效的方法。
不仅编程简单，而且迭代次数较少，运行时间较短，结果的精度也较高。
在程序的运行方面，分别设置了可变的函数选择、无约束方法选择、一维搜索方法选择、起始点、精度这五个输入，故可以在命令窗口运行主程序main，再根据提示要求分别输入这五个参数的所需值，就可以得到运行结果。
程序如下：
1、主函数
clear;
globalk;
k=0;
disp（'1.f（x）=x1^2+2*x2^2-2*x1*x2-4*x1'）;
disp（'2.f（x,y）=x^4-2*x*x*y-2*x*y+x^2+2*y*y+4.5*x-4*y+4'）;
disp（'3.f（x）=x1^2+x2^2+x3^2'）;
while1
n0=input（'请输入上面所想选择函数的编号（1、2、3）：
'）;
ifn0==1||n0==2||n0==3
break;
end
disp（'此次输入无效.'）;
end
disp（''）;
disp（'1.最速下降法'）;
disp（'2.阻尼牛顿法'）;
disp（'3.共轭梯度法'）;
disp（'4.鲍威尔法'）;
disp（'5.变尺度法'）;
disp（'6.单纯形法'）;
while1
m0=input（'请输入上面所想选择无约束方法的编号（1、2、3、4、5、6）：
'）;
ifm0==1||m0==2||m0==3||m0==4||m0==5||m0==6
break;
end
disp（'此次输入无效.'）;
end
disp（''）;
disp（'1.黄金分割法'）;
disp（'2.牛顿法'）;
disp（'3.不精确一维搜索法'）;
while1
m1=input（'请输入上面所想选择一维搜索方法的编号（1、2、3）：
'）;
ifm1==1||m1==2||m1==3
break;
end
disp（'此次输入无效.'）;
end
disp（''）;
s=input（'请输入用空格隔开的初始值坐标向量（如：
1.12.0）：
','s'）;
x=str2num（s）;
x=x';
disp（''）;
while1
e=input（'请输入精度（建议0.1或0.01）：
'）;
ife>0
break;
end
disp（'此次输入无效.'）;
end
disp（''）;
disp（'~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~'）;
[xx,yy]=fmins（m0,m1,n0,x,e）;
fprintf（'迭代次数为：
%8.0f\n',k）;
disp（'所求极值点的坐标向量为：
'）;
fprintf（'%16.5f\n',xx）;
fprintf（'所求函数的极值为：
%16.5f\n',yy）;
2、外部多维的调用函数
function[xx,yy]=fmins（m0,m1,n0,x,e）
%UNTITLED此处显示有关此函数的摘要
%此处显示详细说明
ifm0==1
tic;[xx,yy]=zuisu（m1,n0,x,e）;toc;
elseifm0==2
tic;[xx,yy]=zuni（m1,n0,x,e）;toc;
elseifm0==3
tic;[xx,yy]=gonge（m1,n0,x,e）;toc;
elseifm0==4
tic;[xx,yy]=powell（m1,n0,x,e）;toc;
elseifm0==5
tic;[xx,yy]=bianchi（m1,n0,x,e）;toc;
elseifm0==6
tic;[xx,yy]=danxing（n0,x,e）;toc;
end
end
3、最速法
function[xx,yy]=zuisu（m1,n0,x,e）
%UNTITLED2此处显示有关此函数的摘要
%此处显示详细说明
globalk;
[g,ss]=gra（n0）;
while1
d=-double（subs（g,ss,x））;
a=fmin（m1,n0,x,d,e）;
x1=x+a*d;
ifnorm（x1-x）break;
end
x=x1;
k=k+1;
end
xx=x1;
yy=f0（n0,xx）;
end
4、阻尼法
function[xx,yy]=zuni（m1,n0,x,e）
%UNTITLED此处显示有关此函数的摘要
%此处显示详细说明
globalk;
[g,ss]=gra（n0）;
h=jacobian（g',ss）;
while1
d=-double（subs（h,ss,x）^（-1）*subs（g,ss,x））;
a=fmin（m1,n0,x,d,e）;
x1=x+a*d;
ifnorm（x1-x）break;
end
x=x1;
k=k+1;
end
xx=x1;
yy=f0（n0,xx）;
end
5、共轭梯度法
function[xx,yy]=gonge（m1,n0,x,e）
%UNTITLED2此处显示有关此函数的摘要
%此处显示详细说明
globalk;
ifn0==1
n=2;
elseifn0==2
n=4;
elseifn0==3
n=2;
end
[g,ss]=gra（n0）;
while1
kk=0;
d=-double（subs（g,ss,x））;
while1
a=fmin（m1,n0,x,d,e）;
x1=x+a*d;
gx=double（subs（g,ss,x））;
gx1=double（subs（g,ss,x1））;
ifnorm（gx1）break;
elseifkk==n
break;
end
beta=norm（gx1）^2/norm（gx）^2;
d=-gx+beta*d;
x=x1;
k=k+1;
kk=kk+1;
end
ifnorm（gx1）break;
end
x=x1;
k=k+1;
end
xx=x1;
yy=f0（n0,xx）;
end
6、鲍威尔法
function[xx,yy]=powell（m1,n0,x,e）
%UNTITLED3此处显示有关此函数的摘要
%此处显示详细说明
globalk;
ifn0==1
n=2;
elseifn0==2
n=2;
elseifn0==3
n=3;
end
d=zeros（n,n+1）;
xk=zeros（n,n+1）;
deta=zeros（n,1）;
forj=1:
n
d（j,j）=1;
end
while1
xt=x;
fori=1:
n
a=fmin（m1,n0,xt,d（:
i）,e）;
xk（:
i）=xt+a*d（:
i）;
deta（i）=f0（n0,xt）-f0（n0,xk（:
i））;
xt=xk（:
i）;
end
xt=x;
xk（:
n+1）=2*xk（:
n）-x;
ff0=f0（n0,x）;
ff2=f0（n0,xk（:
n））;
ff3=f0（n0,xk（:
n+1））;
md=max（deta）;
m=find（deta==md）;
if（ff3d（:
n+1）=xk（:
n）-xt;
a=fmin（m1,n0,xk（:
n）,d（:
n+1）,e）;
x=xk（:
n）+a*d（:
n+1）;
fori=m:
n
d（:
i）=d（:
i+1）;
end
else
ifff2x=xk（:
n）;
else
x=xk（:
n+1）;
end
end
ifnorm（x-xt）break;
end
k=k+1;
end
xx=x;
yy=f0（n0,xx）;
end
7、变尺度法
function[xx,yy]=bianchi（m1,n0,x,e）
%UNTITLED此处显示有关此函数的摘要
%此处显示详细说明
globalk;
ifn0==1
n=2;
elseifn0==2
n=2;
elseifn0==3
n=3;
end
[g,ss]=gra（n0）;
while1
gx=double（subs（g,ss,x））;
h=eye（n）;
kk=0;
while1
d=-h*gx;
a=fmin（m1,n0,x,d,e）;
xk=x+a*d;
ifnorm（xk-x）break;
end
ifkk==n
break;
end
gxk=double（subs（g,ss,xk））;
yk=gxk-gx;
sk=xk-x;
h=h+（sk*sk'）/（sk'*yk）-（（h*yk）*yk'*h）/（yk'*h*yk）;
x=xk;
gx=gxk;
kk=kk+1;
k=k+1;
end
ifnorm（xk-x）break;
end
x=xk;
k=k+1;
end
xx=xk;
yy=f0（n0,xx）;
end
8、单形替换法
function[xx,yy]=danxing（n0,x,e）
%UNTITLED此处显示有关此函数的摘要
%此处显示详细说明
globalk;
ifn0==1
n=2;
elseifn0==2
n=2;
elseifn0==3
n=3;
end
f=zeros（n+5,1）;
xk=zeros（n,n+5）;
h=2*eye（n）;
xk（:
1）=x;
fori=1:
n
xk（:
i+1）=x+h（:
i）;
end
while1
fori=1:
n+1
f（i）=f0（n0,xk（:
i））;
end
while1
f（n+2）=nan;
f（n+3）=nan;
f（n+4）=nan;
f（n+5）=nan;
fl=min（f）;
xll=find（f==fl）;
xl=xll
（1）;
fh=max（f）;
xhh=find（f==fh）;
xh=xhh
（1）;
fff=f;
fff（xh）=[];
fg=max（fff）;
fz=0;
fori=1:
n+1
fz=fz+（f（i）-f（xl））^2;
end
fz=sqrt（fz/n）;
iffzbreak;
end
xkk=xk（:
1）;
fori=1:
n
xkk=xkk+xk（:
i+1）;
end
xk（:
n+2）=（xkk-xk（:
xh））/n;
xk（:
n+3）=2*xk（:
n+2）-xk（:
xh）;
f（n+3）=f0（n0,xk（:
n+3））;
iff（n+3）xk（:
n+4）=xk（:
n+2）+2*（xk（:
n+3）-xk（:
n+2））;
f（n+4）=f0（n0,xk（:
n+4））;
iff（n+4）xk（:
xh）=xk（:
n+4）;
f（xh）=f（n+4）;
else
xk（:
xh）=xk（:
n+3）;
f（xh）=f（n+3）;
end
else
iff（n+3）xk（:
xh）=xk（:
n+3）;
f（xh）=f（n+3）;
else
iff（n+3）>=fh
xk（:
n+3）=xk（:
xh）;
end
xk（:
n+5）=xk（:
n+2）+0.5*（xk（:
n+3）-xk（:
n+2））;
f（n+5）=f0（n0,xk（:
n+5））;
iff（n+5）xk（:
xh）=xk（:
n+5）;
f（xh）=f（n+5）;
else
fori=1:
n+1
xk（:
i）=（xk（:
i）+xk（:
xl））/2;
end
break;
end
end
end
k=k+1;
end
iffzbreak;
end
k=k+1;
end
xx=xk（:
xl）;
yy=f0（n0,xx）;
end
9、内部循环一维调用函数
functionxx=fmin（m1,n0,x,d,e）
x0=0;%初始步长默认为0
[amin,amax]=range1（n0,x,d,x0）;
ifm1==1
xx=gold（n0,x,d,amin,amax,e）;
elseifm1==2
xx=newton（n0,x,d,amax,e）;
elseifm1==3
xx=wolfe（n0,x,d）;
end
end
10、一维搜索确定区间函数
function[amin,amax]=range1（n0,x,d,x0）
%UNTITLED5此处显示有关此函数的摘要
%此处显示详细说明
h=1;
a1=x0;y1=f_1（n0,x,d,a1）;
a2=a1+h;y2=f_1（n0,x,d,a2）;
ify2>y1
h=-h;a3=a1;y3=y1;
a1=a2;
a2=a3;y2=y3;
end
a3=a2+h;y3=f_1（n0,x,d,a3）;
whiley3h=h*2;
a1=a2;
a2=a3;y2=y3;
a3=a2+h;y3=f_1（n0,x,d,a3）;
end
amin=min（a1,a3）;
amax=max（a1,a3）;
end
11、黄金一维法
functionxx=gold（no,x,d,amin,amax,e）
%UNTITLED6此处显示有关此函数的摘要
%此处显示详细说明
a1=amax-0.618*（amax-amin）;
y1=f_1（no,x,d,a1）;
a2=amin+0.618*（amax-amin）;
y2=f_1（no,x,d,a2）;
whileabs（amax-amin）>=e
ify1>=y2
amin=a1;
a1=a2;
y1=y2;
a2=amin+0.618*（amax-amin）;
y2=f_1（no,x,d,a2）;
else
amax=a2;
a2=a1;
y2=y1;
a1=amax-0.618*（amax-amin）;
y1=f_1（no,x,d,a1）;
end
end
xx=（amax+amin）/2;
end
12、牛顿一维法
functionxx=newton（n0,x,d,amax,e）
%UNTITLED9此处显示有关此函数的摘要
%此处显示详细说明
symss;
z=x+s*d;
ifn0==1
a=z
（1）;
b=z
（2）;
f=a^2+2*b^2-2*a*b-4*a;
elseifn0==2
a=z
（1）;
b=z
（2）;
f=a^4-2*a*a*b-2*a*b+a*a+2*b*b+4.5*a-4*b+4;
elseifn0==3
a=z
（1）;
b=z
（2）;
c=z（3）;
f=a*a+b*b+c*c;
end
x0=amax;
while
（1）
ifsubs（diff（diff（f,s）,s）,x0）==0
break;
end
x0=x0-double（
升级会员 