| 查看: 630 | 回复: 2 | ||
[求助]
在线等，请版主看看如何运算不了。总是报错
|
|
小弟实在是一窍不通对matlab，不要怪小弟提的问题比较傻，一个是报错中的问题，这个[k,fval,flag] = fmincon(@ObjFunc7Fmincon一句的@ObjFunc7Fmincon是怎么设定的？我看有的人写@ObjFunc4Fmincon. 然后麻烦大家看看这个我自己模拟的代码，是按照版主 @dbb627的代码写的，到底是哪里不对麻烦改下谢谢了。 function parafit % % r1 = k(1)*C(1); % r2 = k(2)*C(1); % r3 = k(3)*C(2); % r4 = k(4)*C(3); % r5 = k(5)*C(2); % r6 = k(6)*C(3); % % dCAdt = - r1 - r2; % dCBdt = r1 + r4 - r3 - r5; % dCCdt = r2 + r3 - r4 – r6 ; % dCDdt = r5 + r6; clear all clc % t/min CA-1 CB-2 CC-3 CD-4 / mol/L Kinetics=[0 100 0 0 0 10 0.2441 0.1390 0.28.31 0.0171 20 0.432 0.1629 0.2370 0.0192 30 0.756 0.1992 0.2509 0.0349 40 0.522 0.2091 0.2737 0.0531 50 0.426 0.1967 0.3231 0.0634 60 0.397 0.1803 0.3536 0.0822 70 0.431 0.1732 0.3330 0.0945 80 0.531 0.1645 0.3310 0.0987]; k0 = [0.0000000005 0.0000000005 0.0000000005 0.00000000005 0.00000005 0.00000005]; % 参数初值 lb = [0 0 0 0 0 0]; % 参数下限 ub = [1 1 1 1 1 1]; % 参数上限 x0 = [100 0 0 0]; yexp = Kinetics; % yexp: 实验数据[x1 x4 x5 x6] warning off % 使用函数fmincon()进行参数估计 [k,fval,flag] = fmincon(@ObjFunc7Fmincon,k0,[],[],[],[],lb,ub,[],[],x0,yexp); fprintf('\n使用函数fmincon()估计得到的参数值为:\n') fprintf('\tk1 = %.11f\n',k(1)) fprintf('\tk2 = %.11f\n',k(2)) fprintf('\tk3 = %.11f\n',k(3)) fprintf('\tk4 = %.11f\n',k(4)) fprintf('\tk5 = %.11f\n',k(5)) fprintf('\tk6 = %.11f\n',k(6)) fprintf(' The sum of the squares is: %.1e\n\n',fval) k_fm= k; warning off % 使用函数lsqnonlin()进行参数估计 [k,resnorm,residual,exitflag,output,lambda,jacobian] = ... lsqnonlin(@ObjFunc7LNL,k0,lb,ub,[],x0,yexp); ci = nlparci(k,residual,jacobian); fprintf('\n\n使用函数lsqnonlin()估计得到的参数值为:\n') fprintf('\tk1 = %.11f\n',k(1)) fprintf('\tk2 = %.11f\n',k(2)) fprintf('\tk3 = %.11f\n',k(3)) fprintf('\tk4 = %.11f\n',k(4)) fprintf('\tk5 = %.11f\n',k(5)) fprintf('\tk6 = %.11f\n',k(6)) fprintf(' The sum of the squares is: %.1e\n\n',resnorm) k_ls = k; output warning off % 以函数fmincon()估计得到的结果为初值，使用函数lsqnonlin()进行参数估计 k0 = k_fm; [k,resnorm,residual,exitflag,output,lambda,jacobian] = ... lsqnonlin(@ObjFunc7LNL,k0,lb,ub,[],x0,yexp); ci = nlparci(k,residual,jacobian); fprintf('\n\n以fmincon()的结果为初值，使用函数lsqnonlin()估计得到的参数值为:\n') fprintf('\tk1 = %.11f\n',k(1)) fprintf('\tk2 = %.11f\n',k(2)) fprintf('\tk3 = %.11f\n',k(3)) fprintf('\tk4 = %.11f\n',k(4)) fprintf('\tk5 = %.11f\n',k(5)) fprintf('\tk6 = %.11f\n',k(6)) fprintf(' The sum of the squares is: %.1e\n\n',resnorm) k_fmls = k; output tspan = [0 : 10 : 80]; [t x] = ode45(@KineticEqs,tspan,x0,[],k_fmls); figure; plot(t,x(:,1),t,yexp(:,2),'*');legend('ca-pr','ca-real') figure;plot(t,x(:,2:5)); hold on plot(t,yexp(:,3:6),'o');legend('cb-pr','cc-pr','cd-pr','ce-pr','cb-real','cc-real','cd-real','ce-real') % ------------------------------------------------------------------ function f = ObjFunc7Fmincon(k,x0,yexp) tspan = [0 : 10 : 80]; [t x] = ode45(@KineticEqs,tspan,x0,[],k); y(:,2) = x(:,1); y(:,3:6) = x(:,2:5); f = sum((y(:,2)-yexp(:,2)).^2) + sum((y(:,3)-yexp(:,3)).^2) ... + sum((y(:,4)-yexp(:,4)).^2) + sum((y(:,5)-yexp(:,5)).^2) ... + sum((y(:,6)-yexp(:,6)).^2) ; % ------------------------------------------------------------------ function f = ObjFunc7LNL(k,x0,yexp) tspan = [0.0 : 10 : 80]; [t x] = ode45(@KineticEqs,tspan,x0,[],k); y(:,2) = x(:,1); y(:,3:6) = x(:,2:5); f1 = y(:,2) - yexp(:,2); f2 = y(:,3) - yexp(:,3); f3 = y(:,4) - yexp(:,4); f4 = y(:,5) - yexp(:,5); f5 = y(:,6) - yexp(:,6); f = [f1; f2; f3; f4; f5]; % ------------------------------------------------------------------ function dxdt = KineticEqs(t,x,k) r1 = k(1)*x(1); r2 = k(2)*x(1); r3 = k(3)*x(2); r4 = k(4)*x(3); r5 = k(5)*x(2); r6 = k(6)*x(3); dCAdt = - r1 - r2; dCBdt = r2 + r4 - r3 - r5; dCCdt = r2 + r3 - r4 - r6; dCDdt = r5 + r6; dxdt = [dCAdt; dCBdt; dCCdt; dCDdt]; |
回复此楼» 本帖@通知
» 猜你喜欢
光学工程学硕调剂信息
已经有32人回复
-
欢迎加入课题组
已经有0人回复
-
物理学I论文润色/翻译怎么收费?
已经有80人回复
-
散金币，求好运，祝面上顺利！
已经有79人回复
-
华东师范大学芯片设计徐珑真实水平咋样
已经有1人回复
-
调剂
已经有0人回复
-
【新加坡】纳米电子器件项目组有“联合培养博士生”名额
已经有0人回复
-
中德博士后2026项目5月30日截止- 国内博士后入站、德国两年全职
已经有0人回复
-
面向全球申请人的中德博士后2026项目5月30日截止- 国内博后入站德国两年全职
已经有0人回复
-
【新加坡】纳米电子器件项目组有“联合培养博士生”名额
已经有2人回复
-
重庆交大材料学院26年硕士招生复试实施细则已出。欢迎加入光子学微结构与器件课题组。
已经有0人回复
2楼2015-11-29 21:52:01
3楼2015-11-30 11:21:13
