最速下降法MATLAB程序实现.docx
《最速下降法MATLAB程序实现.docx》由会员分享，可在线阅读，更多相关《最速下降法MATLAB程序实现.docx（3页珍藏版）》请在冰点文库上搜索。
function[drag,ddrag]=fxy（a0,s0）
%要优化的目标函数
%返回值：
drag，函数在（a0，s0）点处的值；ddrag，函数在（a0，s0，drag）处的梯度（行向量）
%参数：
（a0,s0）
functiond=f（a,s）
w=0.12754*a^1.5+（0.12754^2*a^3+4*（4940/s+45.42））^0.5;
w=w^2*s/4;
cl=2*w/（1506.75*s）;
cd=0.03062702/s+0.00962623*a^0.1*s^-0.1+cl^2*a^-1/（pi*0.96）;
d=w*cd/cl;
end
ddrag=zeros（1,2）;
drag=f（a0,s0）;
h=1e-6;
ddrag
（1）=（f（a0+h,s0）-drag）/h;%函数对a的偏导数
ddrag
（2）=（f（a0,s0+h）-drag）/h;%函数对s的偏导数
end
%最速下降法
%-------------------------------------------------
%-------------------------------------------------
clear;
%绘制等高线
x=5:
0.1:
30;%展弦比A
y=5:
0.1:
30;%机翼面积S
[X,Y]=meshgrid（x,y）;
dim=size（X）;
fori=1:
dim
（1）
forj=1:
dim
（2）
[Z（i,j）,null]=fxy（X（i,j）,Y（i,j））;
end
end
%mesh（X,Y,Z）;
holdon;
[c,h]=contourf（X,Y,Z,10）;clabel（c,h）;colormapcool
xlabel（'A'）;ylabel（'S'）;
%-------------------------------------------------
%-------------------------------------------------
eg=1e-6;%精度要求
ea=1e-6;
er=1e-2;
num=0;%迭代次数
num1=0;%确定最优步长迭代的次数
x0=15;y0=30;%起始点
[z0,dz0]=fxy（x0,y0）;
whilenorm（dz0）>eg
z00=z0;
x00=x0;
y00=y0;
h=norm（dz0）/1e4;%初始步长
dz0=dz0/norm（dz0）;%沿梯度方向的单位向量
x1=x0-dz0
（1）*h;
y1=y0-dz0
（2）*h;
[z1,dz1]=fxy（x1,y1）;
step0=0;%[step0,step1]，步长极小值区间
step1=step0+h;
step=step0;
ifz1whilez1step=step0;
step0=step1;
z0=z1;
h=2*h;
step1=step1+h;
x1=x0-dz0
（1）*step1;
y1=y0-dz0
（2）*step1;
[z1,dz1]=fxy（x1,y1）;
end
step1=step0;
step0=step;
end
while（step1-step0）>eg%黄金分割法求最优步长
x2=x0-dz0
（1）*（step0+0.618*（step1-step0））;
y2=y0-dz0
（2）*（step0+0.618*（step1-step0））;
[z2,dz2]=fxy（x2,y2）;
x1=x0-dz0
（1）*（step0+0.382*（step1-step0））;
y1=y0-dz0
（2）*（step0+0.382*（step1-step0））;
[z1,dz1]=fxy（x1,y1）;
ifz1>z2
step0=step0+0.382*（step1-step0）;
elseifz1step1=step0+0.618*（step1-step0）;
else
step0=step0+0.382*（step1-step0）;
step1=step0+0.618*（step1-step0）;
end
end
num1=num1+1;
end
step_opt=（step0+step1）/2;%最优步长
x0=x0-dz0
（1）*step_opt;
y0=y0-dz0
（2）*step_opt;
[z0,dz0]=fxy（x0,y0）;
%line（[x00,x0],[y00,y0],[z00,z0],'Color','r','LineWidth',2）;
line（[x00,x0],[y00,y0],'Color','r','LineWidth',2）;
plot（x0,y0,'r+'）;
num=num+1;
ifabs（z0-z00）<=ea+er*abs（z00）
break;
end
end
holdoff;
disp（['Opt_Piont=（',num2str（x0）,',',num2str（y0）,'）']）;
disp（['Opt_Value=',num2str（z0）]）;
disp（['迭代次数：
',num2str（num）]）;
disp（['线性搜索迭代次数：
',num2str（num1）]）;
升级会员 