24小时热门版块排行榜    

查看: 703  |  回复: 2

llgcty

金虫 (正式写手)

[求助] 在线等,请版主看看如何运算不了。总是报错

小弟实在是一窍不通对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];
回复此楼

» 本帖@通知

» 猜你喜欢

已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

llgcty

金虫 (正式写手)

快点来吧,帮忙一下!!!
2楼2015-11-29 21:52:01
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

llgcty

金虫 (正式写手)

各位大神帮帮忙
3楼2015-11-30 11:21:13
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
相关版块跳转 我要订阅楼主 llgcty 的主题更新
最具人气热帖推荐 [查看全部] 作者 回/看 最后发表
[基金申请] 投票:  有多少人是今天查系统知道结果的? +16 爱看书的可乐 2026-08-26 18/900 2026-08-29 20:07 by winsaint
[基金申请] 国自然面上复盘~欢迎讨论 (金币+15) +15 晴天加油 2026-08-26 16/800 2026-08-29 18:28 by symmetry
[论文投稿] 售SCI一区T0P文章,我:8O.55.1.O.5.4,科目齐全,可+急 +3 gy1nBQXYQJqL 2026-08-29 4/200 2026-08-29 18:15 by 4FFAWE8HcgUD
[博后之家] 售SCI文章,我:8O5.5.1.O.54,科目齐全,可+急 +3 gy1nBQXYQJqL 2026-08-29 4/200 2026-08-29 18:04 by 4FFAWE8HcgUD
[博后之家] 售SCI文章,我:8O.5.5.1O.54,科目全,可十急 +6 ASdOkHsho7FD 2026-08-28 7/350 2026-08-29 17:31 by 4FFAWE8HcgUD
[考博] 售SCI一区T0P文章,我:8.O.55.1.O54,科目全,可伽急 +5 ASdOkHsho7FD 2026-08-28 8/400 2026-08-29 17:22 by 4FFAWE8HcgUD
[论文投稿] 售SCI一区文章,我:8.O.551.O.5.4,科目全,可伽急 +3 ASdOkHsho7FD 2026-08-28 4/200 2026-08-29 14:13 by jCd0dEvKHShX
[公派出国] 售SCI一区T0P文章,我:8.O.55.1.O.5.4,科目全,可+急 +3 ASdOkHsho7FD 2026-08-28 4/200 2026-08-29 14:02 by jCd0dEvKHShX
[找工作] 售SCI文章,我:8O5.5.1.O.54,科目齐全,可+急 +3 ASdOkHsho7FD 2026-08-28 6/300 2026-08-29 11:51 by jCd0dEvKHShX
[考研] 售一区SCI文章T0P,我:8O.551.O54,科目全,可十急 +4 ASdOkHsho7FD 2026-08-28 6/300 2026-08-29 11:40 by jCd0dEvKHShX
[找工作] 售SCI文章,我:8O5.5.1.O.54,科目齐全,可+急 +4 ASdOkHsho7FD 2026-08-28 6/300 2026-08-29 11:29 by jCd0dEvKHShX
[教师之家] 售SCI-T0P文章,我:8O.5.5.1.O.54,科目齐全,可+急 +4 ASdOkHsho7FD 2026-08-28 7/350 2026-08-29 09:26 by G6APbkg8SA6w
[基金申请] 系统查不到 +11 董八千 2026-08-26 11/550 2026-08-28 18:06 by Leogzhya
[基金申请] 基金系统什么内容也没有 30+4 winsaint 2026-08-27 9/450 2026-08-28 11:06 by maolC
[基金申请] 申请删除本帖 +6 lyz123lyz 2026-08-27 7/350 2026-08-27 17:31 by 宁静致远sy
[基金申请] 怎么看青基中了没有啊 +5 叶九微 2026-08-26 5/250 2026-08-27 10:35 by l_zh2008
[文学芳草园] 梦想 +3 myrtle 2026-08-26 3/150 2026-08-27 10:01 by angelyueyi
[基金申请] 为什么国自然不能直接公布 +4 bjdxyxy 2026-08-26 4/200 2026-08-26 13:12 by qingmu1201
[基金申请] 牛来!米来!面来! +8 beefly 2026-08-26 8/400 2026-08-26 08:37 by xuzhipiao
[基金申请] 在坚冰还盖着北海的时候,我看到了怒放的梅花。 (金币+10) +6 ziyangfang 2026-08-25 9/450 2026-08-25 20:26 by huagongfeihu
信息提示
请填处理意见