|
马上注册,结交更多好友,享用更多功能,让你轻松玩转社区。
您需要 登录 才可以下载或查看,没有账号?我要加入
x
这是我编的一个BFGS程序,但是运行时好像是死循环,用Ctrl+C中断程序后出错信息为:
Error in ==> num2cell at 33
c{i} = a(i);
Error in ==> sym.subs at 68
VaR = num2cell(sym([ '[' vars ']' ]));
Error in ==> cwdBFGS at 29
fz0=subs(f,{x1,x2},{x0(1,1),x0(2,1)});%计算函数值
Error in ==> cwdBFGS at 120
p0=-g0;
Error in ==> cwdBFGS at 115
if fz>=fz0 %进入算法第五步
小弟虚心求教,请高手帮忙指点一下。
%BFGS算法求解
%求梯度
syms x1 x2 t
f=21.5+x1*sin(4.0*pi*x1)+x2*sin(20.0*pi*x2); %目标函数,
fx=diff(f,x1); %对x求梯度
fy=diff(f,x2); %对y求梯度
fxy=[fx;fy];
%初始化赋值
dim=2; %维数
e=0.01; % H的终止准则限
%算法第一步
x0=[1;5]; %给初值
fz0=subs(f,{x1,x2},{x0(1,1),x0(2,1)});%计算函数值
g0=subs(fxy,{x1,x2},{x0(1,1),x0(2,1)});%计算梯度值
%算法第二步
k=0; %迭代次数
H0=[1,0;0,1]; %迭代矩阵H初值为单位矩阵
p0=-g0;
%算法第三步,一维搜索
m=1; %标志量
iter=0;
while (m)
iter=iter+1;
x=x0+t*p0;
ft=subs(f,{x1,x2},{x(1,1),x(2,1)});
%调用进退法算法,确定范围
d0=0;%初始点
h=1;%初始步长
s=2; %加倍系数
y0=subs(ft,t,d0);
n=0;
p=1;%标志量
while (p)
d1=d0+h;
y1=subs(ft,t,d1);
if y1<y0
h=s*h;
d=d0;
d0=d1;
y0=y1;
n=n+1;
elseif n==0;
h=-h;
d=d1;
else
p=0;
break;
end
end
a=min(d,d1);%得到搜索区间
b=max(d,d1);
%0.618算法确定t的值
r=0.618;
h1=a+(1-r)*(b-a);%两个边界点
u1=a+r*(b-a);
y1=subs(ft,t,h1);%边界点的值
y2=subs(ft,t,u1);
while (abs(h1-u1)>0.01)
if y1<y2
a=u1;
u1=h1;
y2=y1;
h1=a+(1-r)*(b-a);
y1=subs(ft,t,h1);
else
a=h1;
h1=u1;
y1=y2;
u1=a+r*(b-a);
y2=subs(ft,t,u1);
end
end
t=(h1+u1)/2; %得到最优值
%一维搜索后使新的迭代点不超过x的取值范围,限幅处理
x=x0+t*p0;
x=double(x);
if x(1,1)>12.1
x(1,1)=12.1;
end
if x(1,1)<-3.0
x(1,1)=-3.0;
end
if x(2,1)>5.8
x(2,1)=5.8;
end
if x(2,1)<4.1
x(2,1)=4.1;
end
%计算f(k+1)和g(k+1)的值
fz=subs(f,{x1,x2},{x(1,1),x(2,1)});%计算新的值
g=subs(fxy,{x1,x2},{x(1,1),x(2,1)});% 计算新的梯度
%算法第四步
while(norm(g(:,1),2)>=e) %判断H是否满足终止条件,计算2-范数
if fz>=fz0 %进入算法第五步
x0=x0; %迭代后函数值没有下降,重新迭代
fz0=fz0;
g0=g0;
H0=[1,0;0,1]; %迭代矩阵H初值为单位矩阵
p0=-g0;
m=1;
k=0;
continue; %转算法第三步
elseif k==dim %进入算法第六步
x0=x;
fz0=fz;
g0=g;
H0=[1,0;0,1]; %迭代矩阵H初值为单位矩阵
p0=-g0;
m=1;
k=0;
continue; %转算法第三步
else %进入算法第七步
y=g-g0;
s=x-x0;
H=H0+[(1+y'*H0*y/(s'*y))*s*s'-H0*y*s'-s*y'*H0]/(s'*y);
p=-H*g;
H0=H;
p0=p;
fz0=fz;
g0=g;
x0=x;
k=k+1;
continue; %转算法第三步
end
end
m=0;
end
xpso=x %输出x,f的值,BFGS算法结束
fpso=fz
iter |
|