24小时热门版块排行榜    

北京石油化工学院2026年研究生招生接收调剂公告
查看: 2680  |  回复: 5

2011207156

金虫 (小有名气)

[求助] matlab求模型参数有错误,求指导 已有2人参与

以下的程序是我照葫芦画瓢改的,不知道哪里错了。。。。。。本人刚开始接触matlab,求大神指导啊!

function KineticsEst1_Int  % 用积分法分析得到kh,b,和kd
clear all; clc
global nm
t=[31 32 33 34 35 36 37 38 39 40 41 42 43 44 45];
nm=[20.41220 22.73032 24.55624 25.88679 26.69019 27.16334 27.45866 27.64919 27.57615 27.70635 27.77303 27.86195 27,91911 27.95721 27.98897]';

% 非线性拟合
beta0=[0.0001 2 0.0001]; tspan=[31 32 33 34 35 36 37 38 39 40 41 42 43 44 45]; n0=20.4122;
[beta,resnorm,resid,exitflag,output,lambda,jacobian] = ...
    lsqnonlin(@OptObjFunc,beta0,[],[],[],tspan,n0,nm)
ci=nlparci(beta,resid,jacobian)

% 拟合效果图(实验与拟合的比较)
[t4plot n4plot]=ode45(@KineticsEqs,[tspan(1) tspan(end)],n0,[],beta);
plot(t,nm,'bo',t4plot,n4plot,'k-')
legend('Exp','Model'), xlabel('时间t,s'), ylabel('累积量n,mol')

%残差关于拟合值的残差图
[t nc]=ode45(@KineticsEqs,tspan,nO,[],beta);
figure, plot(nc,resid,'*'), xlabel('累积值(mol)'), ylabel('残差R(mol)'), refline(0,0)

%参数辨识结果
fprintf('\n\nEstimated Parameters:\n'), fprintf('\tkh=%.4f±%.4f\n',beta(1),ci(1,2)-beta(1))
fprintf('\tb = %.2f ± %.2f\n',beta(2),ci(2,2)-beta(2))
fprintf('\tkd = %.2f ± %.2f\n',beta(3),ci(3,2)-beta(3))
% ------------------------------------------------------------------
function f = OptObjFunc(beta,tspan,n0,nm)
[t nc] = ode45(@KineticsEqs,tspan,n0,[],beta);
f = nc - nm;
% ------------------------------------------------------------------
function dndt = KineticsEqs(t,n,beta)
dndt = 0.001414/(1/(4.5727*10^4*beta(1)*(1-0.02119*n)^beta(2))+1/beta(3)/45727);    % kh= beta(1), b= beta(2),kd=beta(3)


报错:
??? Error using ==> minus
Matrix dimensions must agree.

Error in ==> KineticsEst1_Int4>OptObjFunc at 29
f = nc - nm;

Error in ==> lsqnonlin at 203
            initVals.F = feval(funfcn{3},xCurrent,varargin{:});

Error in ==> KineticsEst1_Int4 at 9
[beta,resnorm,resid,exitflag,output,lambda,jacobian] = ...

Caused by:
    Failure in initial user-supplied objective function
    evaluation. LSQNONLIN cannot continue.
回复此楼

» 猜你喜欢

» 本主题相关价值贴推荐,对您同样有帮助:

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

spiderone

木虫 (小有名气)

【答案】应助回帖

感谢参与,应助指数 +1
矩阵维数不对,nc和nm的个数应该不对,可以仔细看看。
2楼2014-06-30 16:47:09
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

月只蓝

主管区长 (职业作家)

【答案】应助回帖

★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★
感谢参与,应助指数 +1
2011207156: 金币+20, ★★★★★最佳答案, 非常谢谢哦~ 2014-07-01 16:01:45
CODE:
function KineticsEst1_Int_modified_by_Yuezhilan
clear all;clc
format long
tspan=[31 32 33 34 35 36 37 38 39 40 41 42 43 44 45]-31;
yexp=[22.73032 24.55624 25.88679 26.69019 27.16334 27.45866 27.64919 27.57615 27.70635 27.77303 27.86195 27.91911 27.95721 27.98897]';

beta0=[100 10 0.1 ];
y0=20.4122;
lb=[-1 -1 -1]*1e6;
ub=[1e9 1e9 1e2];   
yy=[y0 yexp'];
k0=beta0;
[k,resnorm,residual,exitflag,output,lambda,jacobian] = ...
    lsqnonlin(@ObjFunc,k0,lb,ub,[],tspan,y0,yexp);
fprintf('\n\n使用函数lsqnonlin()估计得到的参数值为:\n')
fprintf('\t待拟合参数 k1 = %.6f\n',k(1))
fprintf('\t待拟合参数 k2 = %.6f\n',k(2))
fprintf('\t待拟合参数 k3 = %.6f\n',k(3))
fprintf('  残差平方和= %.6f\n\n',resnorm)
ts=0:1:max(tspan);

[ts ys]=ode45(@KineticsEqs,ts,y0,[],k);
[ttt XXsim] = ode45(@KineticsEqs,tspan,y0,[],k);
y=XXsim(2:end);
xexp=yexp;
R2=1-sum((xexp-y).^2)./sum((xexp-mean(y)).^2);
fprintf('\n\t决定系数R-Square = %.6f',R2);
figure(1)
plot(ts,ys,'b',tspan,yy,'or'),legend('计算值','实验值','Location','best');
yr=y-yexp;
figure(2)
plot(tspan(2:end),yr,'r*',[-1 15],[0 0]),axis([-1 15 -0.5 0.5]);
figure(3)
plot(yexp,y,'ro',[21 29],[21 29],'b-');

%---------------------------------------------------------
function f = ObjFunc(k,tspan,y0,yexp)           % 目标函数
[t Xsim] = ode45(@KineticsEqs,tspan,y0,[],k) ;
ysim = Xsim(2:end);
size(ysim);
size(yexp);
f=ysim-yexp;
%----------------------------------------------------------

function dydt = KineticsEqs(t,y,k)
beta(1)=k(1);
beta(2)=k(2);
beta(3)=k(3);
dydt = 0.001414/(1/(4.5727*10^4*beta(1)*(1-0.02119*y)^beta(2))+1/beta(3)/45727);

计算结果显示:
使用函数lsqnonlin()估计得到的参数值为:
        待拟合参数 k1 = 73143442.357330
        待拟合参数 k2 = 27.497659
        待拟合参数 k3 = 0.035271
  残差平方和= 0.241151

决定系数R-Square = 0.992227
matlab求模型参数有错误,求指导
附图1.png

MATLAB、MS小问题、普通问题请发帖求助!时间精力有限,恕不接受无偿私信求助。
3楼2014-06-30 20:05:14
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

2011207156

金虫 (小有名气)

引用回帖:
3楼: Originally posted by 月只蓝 at 2014-06-30 20:05:14
function KineticsEst1_Int_modified_by_Yuezhilan
clear all;clc
format long
tspan=-31;
yexp=';

beta0=;
y0=20.4122;
lb=*1e6;
ub=;   
yy=;
k0=beta0;
= ...
    lsqnonlin(@ObjFunc,k0,lb,u ...

真心感谢大神~
4楼2014-07-01 16:02:35
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

你的CEO

铁虫 (小有名气)

引用回帖:
3楼: Originally posted by 月只蓝 at 2014-06-30 20:05:14
function KineticsEst1_Int_modified_by_Yuezhilan
clear all;clc
format long
tspan=-31;
yexp=';

beta0=;
y0=20.4122;
lb=*1e6;
ub=;   
yy=;
k0=beta0;
= ...
    lsqnonlin(@ObjFunc,k0,lb,u ...

大神你好,能否请教下你关于动力学模型参数拟合的问题
未来为你而来!
5楼2015-09-22 13:50:07
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

Irene_JY

新虫 (初入文坛)

引用回帖:
3楼: Originally posted by 月只蓝 at 2014-06-30 20:05:14
function KineticsEst1_Int_modified_by_Yuezhilan
clear all;clc
format long
tspan=-31;
yexp=';

beta0=;
y0=20.4122;
lb=*1e6;
ub=;   
yy=;
k0=beta0;
= ...
    lsqnonlin(@ObjFunc,k0,lb,u ...

你好,请问你这个怎么生成的,为什么我用你的代码生成不了?求助
6楼2017-01-20 17:25:45
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
相关版块跳转 我要订阅楼主 2011207156 的主题更新
最具人气热帖推荐 [查看全部] 作者 回/看 最后发表
[考研] 一志愿 江南大学 085602 化工专硕 338分求调剂 +13 路痴小琪 2026-04-05 13/650 2026-04-06 06:19 by houyaoxu
[考研] 材料专硕(0856) 339分求调剂 +9 哈哈哈鹅哈哈哈 2026-04-05 9/450 2026-04-05 22:24 by dongzh2009
[考研] 材料调剂 +7 一样YWY 2026-04-05 8/400 2026-04-05 21:20 by 学员8dgXkO
[考研] 求调剂,一志愿厦门大学,生物与医药,总分272,本科211 +3 Electron1cc 2026-04-01 4/200 2026-04-05 20:24 by lys0704
[考研] 326求调剂 +3 顾若浮生 2026-04-05 3/150 2026-04-05 18:32 by 蓝云思雨
[考研] 工科08专硕机械275求调剂 +3 AaAa7420 2026-04-02 3/150 2026-04-05 13:26 by jp9609
[考研] 0854求调剂 +4 assdll 2026-04-04 4/200 2026-04-05 09:44 by zhq0425
[考研] 材料调剂 +9 革微桂 2026-04-04 9/450 2026-04-05 08:27 by 544594351
[考研] 材料与化工306分找调剂 +12 沧海轻舟e 2026-04-03 13/650 2026-04-04 23:45 by lqwchd
[考研] 278求调剂 +14 范婷娜 2026-04-04 15/750 2026-04-04 22:15 by lqwchd
[考研] 0835学硕299求调剂 08大类可接受 +5 useryy 2026-04-03 5/250 2026-04-04 20:07 by 蓝云思雨
[考研] 085602化工求调剂(331分) +9 111@127 2026-03-30 9/450 2026-04-02 20:00 by dick_runner
[考研] 22408 266求调剂 +3 masss11222 2026-04-02 3/150 2026-04-02 18:11 by 笔落锦州
[考研] 求调剂 +7 Aniyaio 2026-04-02 7/350 2026-04-02 16:42 by zzsw+
[考研] 材料专业求调剂 +10 月月鸟木 2026-04-01 10/500 2026-04-02 12:57 by wxiongid
[考研] 298求B区调剂 +4 zzz,,r 2026-04-02 5/250 2026-04-02 12:17 by 土木硕士招生
[考研] 【求调剂】新能源材料本科,一志愿211,初试321 +6 求调剂学校, 2026-04-02 6/300 2026-04-02 09:41 by 晴空210210
[考研] 求调剂,一志愿南京师范大学计算机专硕,初试373,六级通过, +3 计算机追梦人 2026-04-01 3/150 2026-04-02 07:57 by fxue1114
[考研] 求调剂0703 +5 周嘉尧 2026-03-31 8/400 2026-04-01 20:32 by ltltkkk
[考研] 省双一流重点一本大学招收调剂 +4 wwwwffffff 2026-03-31 7/350 2026-04-01 15:23 by wwwwffffff
信息提示
请填处理意见