| 查看: 432 | 回复: 0 | |||
ws89金虫 (著名写手)
悠然书生
|
[求助]
参数拟合程序，不知道是不是陷入了死循环，运行两天还在跑。求帮助
|
|
参数拟合程序，不知道是不是陷入了死循环，运行两天还在跑。师兄的程序，没有得到最后跟师兄一致的结果。求帮助。大神给看看是怎磨回事。 clc clear allT2=[-20 -10 -10 0 0 0 0 10 10 20]'; allf=1.0e+003 *[2.1304 0.9212 0.9212 0.4141 0.4141 0.4141 0.4141 0.193 0.193 0.093]'; allw2=[1 1 2 0.1 0.5 1 5 1 2 1]'; allG2=[23.46 8.062 14 1 2.97 4.18 13.483 2.13 3.128 1.08]'*10^6; allG1=[17 5.8 10 2.5 3.3 3.8 9.7 3.0 3.4 2.7]'*10^6; for i=1:10 f2(i,1)=allf(i); w2(i,1)=allw2(i); G1(i,1)=allG1(i); G2(i,1)=allG2(i); end x1=[f2 w2]; GG2=@(p,x)(p(1).*p(2).^2./p(4).^2+2.*p(1).*p(2)./p(4).*x1(:,1).^p(5).*x1(:,2).^p(5).*cos(pi/2.*p(5))+p(1).*x1(:,1).^(2.*p(5)).*x1(:,2).^(2.*p(5))+p(3).*p(2).^2./p(4)^2.*x1(:,1).^p(7).*x1(:,2).^p(7).*cos... (pi/2.*p(7))+2.*p(3).*p(2)./p(4).*x1(:,1).^(p(5)+p(7)).*x1(:,2).^(p(5)+p(7)).*cos(pi/2.*p(5)).*cos(pi/2.*p(7))+p(3).*x1(:,1).^(2.*p(5)+p(7)).*x1(:,2).^(2.*p(5)+p(7)).*cos(pi/2.*p(7))+p(2).^2.... ./p(4).*cos(pi/2.*p(6))+p(2).*x1(:,1).^p(5).*x1(:,2).^p(5).*cos(pi/2.*(p(5)+p(6)))); GG1=@(p,x)(p(3).*p(2).^2./p(4).^2.*x1(:,1).^p(7).*x1(:,2).^p(7).*sin(pi/2.*p(7))+2.*p(3).*p(2)./p(4).*x1(:,1).^(p(5)+p(7)).*x1(:,2).^(p(5)+p(7)).*sin(pi/2.*p(5)).*sin(pi/2.*p(7))+p(3).*x1(:,1).^... (2.*p(5)+p(7)).*x1(:,2).^(2*p(5)+p(7)).*sin(pi/2.*p(7))+p(2).^2./p(4).*sin(pi/2.*p(6))+p(2).*x1(:,1).^p(5).*x1(:,2).^p(5).*sin(pi/2.*(p(5)+p(6)))); p0=[97.06669468363,272.111066749,60321.3626849,164.3298821730,0.517346991136036,0.20053385863749,0.26009044292186]; p=p0; GG22=p(1).*p(2).^2./p(4).^2+2.*p(1).*p(2)./p(4).*x1(:,1).^p(5).*x1(:,2).^p(5).*cos(pi/2.*p(5))+p(1).*x1(:,1).^(2.*p(5)).*x1(:,2).^(2.*p(5))+p(3).*p(2).^2./p(4)^2.*x1(:,1).^p(7).*x1(:,2).^p(7).*cos... (pi/2.*p(7))+2.*p(3).*p(2)./p(4).*x1(:,1).^(p(5)+p(7)).*x1(:,2).^(p(5)+p(7)).*cos(pi/2.*p(5)).*cos(pi/2.*p(7))+p(3).*x1(:,1).^(2.*p(5)+p(7)).*x1(:,2).^(2.*p(5)+p(7)).*cos(pi/2.*p(7))+p(2).^2.... ./p(4).*cos(pi/2.*p(6))+p(2).*x1(:,1).^p(5).*x1(:,2).^p(5).*cos(pi/2.*(p(5)+p(6))); %p0=p; GG11=p(3).*p(2).^2./p(4).^2.*x1(:,1).^p(7).*x1(:,2).^p(7).*sin(pi/2.*p(7))+2.*p(3).*p(2)./p(4).*x1(:,1).^(p(5)+p(7)).*x1(:,2).^(p(5)+p(7)).*sin(pi/2.*p(5)).*sin(pi/2.*p(7))+p(3).*x1(:,1).^... (2.*p(5)+p(7)).*x1(:,2).^(2*p(5)+p(7)).*sin(pi/2.*p(7))+p(2).^2./p(4).*sin(pi/2.*p(6))+p(2).*x1(:,1).^p(5).*x1(:,2).^p(5).*sin(pi/2.*(p(5)+p(6))); a1=0.5;a2=1-a1; while(1) if (norm(GG11-allG1,inf)<1.0e6&&norm(GG22-allG2,inf)<1.0e6) break end options = optimset('MaxFunEvals',4000,'MaxIter',40000); p=fminsearch(@(p)(a1*(sum(((G1-GG1(p,x1))./G1).^2)+a2*sum(((G2-GG2(p,x1))./G2).^2))),p0); GG22=p(1).*p(2).^2./p(4).^2+2.*p(1).*p(2)./p(4).*x1(:,1).^p(5).*x1(:,2).^p(5).*cos(pi/2.*p(5))+p(1).*x1(:,1).^(2.*p(5)).*x1(:,2).^(2.*p(5))+p(3).*p(2).^2./p(4)^2.*x1(:,1).^p(7).*x1(:,2).^p(7).*cos... (pi/2.*p(7))+2.*p(3).*p(2)./p(4).*x1(:,1).^(p(5)+p(7)).*x1(:,2).^(p(5)+p(7)).*cos(pi/2.*p(5)).*cos(pi/2.*p(7))+p(3).*x1(:,1).^(2.*p(5)+p(7)).*x1(:,2).^(2.*p(5)+p(7)).*cos(pi/2.*p(7))+p(2).^2.... ./p(4).*cos(pi/2.*p(6))+p(2).*x1(:,1).^p(5).*x1(:,2).^p(5).*cos(pi/2.*(p(5)+p(6))); GG11=p(3).*p(2).^2./p(4).^2.*x1(:,1).^p(7).*x1(:,2).^p(7).*sin(pi/2.*p(7))+2.*p(3).*p(2)./p(4).*x1(:,1).^(p(5)+p(7)).*x1(:,2).^(p(5)+p(7)).*sin(pi/2.*p(5)).*sin(pi/2.*p(7))+p(3).*x1(:,1).^... (2.*p(5)+p(7)).*x1(:,2).^(2*p(5)+p(7)).*sin(pi/2.*p(7))+p(2).^2./p(4).*sin(pi/2.*p(6))+p(2).*x1(:,1).^p(5).*x1(:,2).^p(5).*sin(pi/2.*(p(5)+p(6))); p0=p; end |
回复此楼» 猜你喜欢
听说今天filecode变了
已经有48人回复
-
FileCode能看出啥？
已经有17人回复
-
奇怪，两个人的filecode固定段从头到尾一模一样
已经有6人回复
-
filecode
已经有9人回复
-
应该是93bebmhtak前后十一个字符比较关键
已经有23人回复
-
我的国基提前知道中了，可是同事的操作让我实在接受不了，怎么会有这样的人
已经有13人回复
-
小木虫上这么多卖论文的，真有人买论文么？感觉没必要啊
已经有6人回复
-
为什么网上很多人说本周 12号出结果
已经有5人回复
-
帮忙看看fileCode
已经有11人回复
-
基金中了
已经有15人回复
