马上注册,结交更多好友,享用更多功能,让你轻松玩转社区。
您需要 登录 才可以下载或查看,没有账号?我要加入
x
唉,以前用Fortran编写过最小二乘法的程序,也没有感觉什么;后来用maple,出了点问题;现在用matlab还是麻烦。所以,特地向各位高人求助。
在书中,学到的最小二乘法的编程例子(函数自己适当变化过,权当练习了)如下:
function y=c8f3(a,x)
y=(a(1)*3+a(3))*3*x+(a(2)*2-a(1))*x.^2.*exp(-a(3)*x)+a(4);
x=0.1:0.1:1;
y=[2.3201,2.6470,2.9707,3.2885,3.6008,3.9090,4.2147,4.5191,4.8232,5.1275];
a=lsqcurvefit('c8f3',[1;2;2;3],x,y);a'
运行没有问题,可是下面我自己的程序,就不能运行呢?
function y=c8f3(a,x)
y=53802.86163*1/x+(-.3072209409*a(1)-.5888875307e-3*a(3)-.1176917018e-7*a(6)-4886632.223-.4703299774e-6*a(5)-.1647355068e-1*a(2)-.1754269266e-4*a(4))/x^2+(.2355550123e-2*a(4)+.6139942431e-4*a(5)+.1505055928e-5*a(6)+.8236775338e-1*a(3)+2.457767527*a(2)+196117179.3+51.56522571*a(1))/x^3+(-3.072209409*a(3)-2936843549.-103.1304514*a(2)-.2061106358e-2*a(5)-.4911953945e-4*a(6)-2884.973156*a(1)-.8236775338e-1*a(4))/x^4+53802.86163/x^5*a(1)+53802.86163/x^6*a(2)+53802.86163/x^7*a(3)+53802.86163/x^8*a(4)+53802.86163/x^9*a(5)+53802.86163/x^10*a(6);
x=[1708.4,807.03,502.72,346.79,248.93,177.67,113.22,55.948037,36.285,33.483,32.01,31,30.229,29.606,29.082,28.63,28.233,27.879];
y=[30,60,90,120,150,180,210,220.64,240,270,300,330,360,390,420,450,480,510];
a=lsqcurvefit('c8f3',[-.94e7; .218265e9; -.16440790e10;-732184771.; -82569287.4; -6088025.940],x,y);a'
计算机提示为:
???Error using ==> lsqncommon at 101
LSQCURVEFITcannot continue because user supplied objective function failed with thefollowing error:
Errorusing ==> mrdivide
Matrixdimensions must agree.
Errorin ==> lsqcurvefit at 182
[x,Resnorm,FVAL,EXITFLAG,OUTPUT,LAMBDA,JACOB]= ...
Errorin ==> mai at 6
a=lsqcurvefit('c8f3',[-.94e7;.218265e9; -.16440790e10; -732184771.; -82569287.4; -6088025.940],x,y);a'
这是为什么呢,谁能帮帮我呀,不胜感激!!!
另外,还有一个问题,如果上面成功的例子,函数变成下面形式会如何呢?
function y=c8f3(a,x)
y=(a(1)*3+a(3) +2*tm-3*tn)*3*x+(a(2)*2-a(1) -2*tm+3*tn)*x.^2.*exp(-a(3)*x)+a(4);
也就是说,除了待拟合的参数a以外,函数里面含有已知的参数tm和tn(例如分别等于8和9),这样的怎么写程序呢。
[ 本帖最后由 ChaChing 于 2009-10-12 11:18 编辑 ] |