24小时热门版块排行榜    

查看: 2356  |  回复: 5

lxyy

银虫 (小有名气)

[求助] matlab拟合模型参数

CA,CB,CC,CW,CR每个时间点t下都有具体数值,要怎么用matlab编程啊
回复此楼

» 本帖附件资源列表

  • 欢迎监督和反馈:小木虫仅提供交流平台,不对该内容负责。
    本内容由用户自主发布,如果其内容涉及到知识产权问题,其责任在于用户本人,如对版权有异议,请联系邮箱:xiaomuchong@tal.com
  • 附件 1 : 新建 Microsoft Word 文档.doc
  • 2012-03-23 11:05:51, 20.5 K

» 猜你喜欢

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

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

dingd

铁杆木虫 (职业作家)

【答案】应助回帖

★ ★ ★ ★ ★ ★ ★ ★ ★ ★
感谢参与,应助指数 +1
lxyy: 金币+10 2012-03-23 19:15:43
这种微分方程拟合不是Matlab的强项,用1stOpt试试!
2楼2012-03-23 15:34:22
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

change0618

铁杆木虫 (著名写手)

方丈大师

【答案】应助回帖

★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ...
感谢参与,应助指数 +1
dbb627: 金币+2, 谢谢应助! 2012-03-23 17:16:55
lxyy: 金币+60, 如果可以,剩下的金币会再给你,呵呵 2012-03-23 19:17:04
function IsothermTR
clear all
clc
global keq1 keq2 t0 c0

keq1 = ?;         % 平衡常数k1
keq2 = ?;   % 平衡常数k2

t0 = [0.5, 1.0, 1.5, 2.0, 2.5, 3.0];
c0 = [?, ? ?]        % A C W 组成的初始浓度
ci=[0.184         0.2443        0.1348;
    0.4457        0.3214        0.1733;
    0.7023        0.4183        0.2713;
    0.7758        0.4789        0.1056;
    0.3365        0.7035        0.1951;
    0.7946        0.8864        0.2235];
% cA cC cR 对应时间数值

k0=[2 2 3];   % k1+, k2, k3+ 初始向值
lb = [0 0 0];
ub = [1000  1000 500];  % 上下限

% 使用函数lsqnonlin()进行参数估计
[k,resnorm,residual,exitflag] = lsqnonlin(@ObjFunc,k0,lb,ub,[],ci);      

k1plus=k(1);
k1minus=keq1/k1plus;
k2=k(2);
k3plus=k(3);
k3minus=keq1/k3plus;

% ------------------------------------------------------------------
function f = ObjFunc(k,ci)    % 目标函数
global t0  c0

[t,c_cal] = ode45(@Euqations,t0,c0,[],k);
f= c_cal-ci;

% ------------------------------------------------------------------
function dcdt = Euqations(t,c,k)
global keq1 keq2

k1plus=k(1);
k1minus=keq1/k1plus;
k2=k(2);
k3plus=k(3);
k3minus=keq1/k3plus;

cA=c(1);
cC=c(2);
cR=c(3);
cB=?;                 % cB物料衡算表达
cW=?;                % cW物料衡算表达

dcAdt = k1plus*cA*cB-k1minus*cC+k2*cA*cW;
dccdt = -k1plus*cA*cB+k1minus*cC+k3plus*cB*cR-k3minus*cC*cW;
dcRdt = k2*cA*cW-k3plus*cB*cR+k3minus*cC*cW;
dcdt = [dcAdt; dcCdt;dcRdt];
3楼2012-03-23 17:13:43
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

lxyy

银虫 (小有名气)

引用回帖:
2楼: Originally posted by dingd at 2012-03-23 15:34:22:
这种微分方程拟合不是Matlab的强项,用1stOpt试试!

可是我不怎么会
4楼2012-03-23 19:15:18
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

lxyy

银虫 (小有名气)

引用回帖:
3楼: Originally posted by change0618 at 2012-03-23 17:13:43:
function IsothermTR
clear all
clc
global keq1 keq2 t0 c0

keq1 = ?;         % 平衡常数k1
keq2 = ?;   % 平衡常数k2

t0 = ;
c0 =         % A C W 组成的初始浓度
ci=;
% cA cC cR 对应时间数值

k0=;   % k ...

我试试看,谢谢
5楼2012-03-23 19:16:10
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

lxyy

银虫 (小有名气)

引用回帖:
3楼: Originally posted by change0618 at 2012-03-23 17:13:43:
function IsothermTR
clear all
clc
global keq1 keq2 t0 c0

keq1 = ?;         % 平衡常数k1
keq2 = ?;   % 平衡常数k2

t0 = ;
c0 =         % A C W 组成的初始浓度
ci=;
% cA cC cR 对应时间数值

k0=;   % k ...

%30℃
function IsothermTR
clear all
clc
global keq1 keq3 t0 c0

keq1 = 3.94;         % ƽºâ³£Êýk1
keq3 = 1.08;   % ƽºâ³£Êýk3

t0 = [5,10, 15, 30, 45, 60,90,120,160,200,240,300,360];
c0 = [5.939601,5.203973,0,0.022484,0.078956]        % A B C R W ×é³ÉµÄ³õʼŨ¶È
ci=[3.742508811 4.213566482 0.990406938 0.252498405 0.090415512;
3.045716838 3.387013322 1.816960097 0.215546741 0.081072022;
2.504754523 2.901991595 2.301981824 0.193427102 0.099462604;
1.65098702  2.036686653 3.167286766 0.170236582 0.095144044;
1.271633833 1.655460403 3.548513016 0.124730781 0.10334349;
1.129887499 1.630815394 3.573158025 0.129798695 0.097764543;
1.01239208  1.460748979 3.74322444  0.144649391 0.073614958;
1.007417603 1.429180292 3.774793127 0.123477124 0.094249307;
0.990005723 1.499216939 3.704756481 0.135524067 0.076072022;
0.970073726 1.493966669 3.71000675  0.130791851 0.076387812;
0.96903024  1.429705429 3.77426799  0.119510017 0.101202216;
0.980873385 1.457988124 3.745985296 0.124189274 0.097022161;
0.98477262  1.453402965 3.750570455 0.1200936   0.0762608];

% cA cC cR ¶ÔӦʱ¼äÊýÖµ

k0=[1 1 1];   % k1+, k2, k3+ ³õʼÏòÖµ
lb = [0 0 0];
ub = [1000  1000 500];  % ÉÏÏÂÏÞ

% ʹÓú¯Êýlsqnonlin()½øÐвÎÊý¹À¼Æ
[k,resnorm,residual,exitflag] = lsqnonlin(@ObjFunc,k0,lb,ub,[],ci);      

k1plus=k(1);
k1minus=keq1/k1plus;
k2=k(2);
k3plus=k(3);
k3minus=keq3/k3plus;

% ------------------------------------------------------------------
function f = ObjFunc(k,ci)    % Ä¿±êº¯Êý
global t0  c0

[t,c_cal] = ode45(@Euqations,t0,c0,[],k);
f= c_cal-ci;

% ------------------------------------------------------------------
function dcdt = Euqations(t,c,k)
global keq1 keq3

k1plus=k(1);
k1minus=keq1/k1plus;
k2=k(2);
k3plus=k(3);
k3minus=keq3/k3plus;

cA=c(1);
cB=c(2);
cC=c(3);
cR=c(4);
cW=c(5);

dcAdt = k1plus*cA*cB-k1minus*cC+k2*cA*cW;
dcCdt = -k1plus*cA*cB+k1minus*cC+k3plus*cB*cR-k3minus*cC*cW;
dcRdt = k2*cA*cW-k3plus*cB*cR+k3minus*cC*cW;
dcdt = [dcAdt; dcCdt;dcRdt];

这个是我根据你的程序把自己的数据放进去了,可是运行不了,会出错,麻烦您再帮忙看看哪里有问题
6楼2012-03-26 11:40:05
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
相关版块跳转 我要订阅楼主 lxyy 的主题更新
最具人气热帖推荐 [查看全部] 作者 回/看 最后发表
[基金申请] 师弟论文见刊大半年才想起来申请专利,还能抢救一下吗? +5 13108017953 2026-09-21 5/250 2026-09-27 16:46 by Leogzhya
[教师之家] 某top大学教授说“能够在市场中兑现的能力才是真能力”无比同意! +3 zju2000 2026-09-26 3/150 2026-09-27 16:26 by lala1008
[硕博家园] 上海工程技术大学 激光智能制造课题组 2027级博士研究生招生 +5 oxidpl 2026-09-20 7/350 2026-09-27 10:43 by Equinoxhua
[找工作] 售SCI一区T0P文章,我:8O.55.1.O.54,科目全,可伽急 +3 ViqFlLSxHDGt 2026-09-26 3/150 2026-09-27 10:26 by 4IW0sJvtEqX8
[博后之家] 售SCI一区T0P文章,我:8O.55.1.O.5.4,科目齐全,可+急 +3 ViqFlLSxHDGt 2026-09-26 3/150 2026-09-27 10:21 by 4IW0sJvtEqX8
[公派出国] 售SCI一区T0P文章,我:8.O.55.1.O54,科目全,可伽急 +3 ViqFlLSxHDGt 2026-09-26 3/150 2026-09-27 10:15 by 4IW0sJvtEqX8
[考博] 售SCI文章,我:8O5.5.1.O.54,科目齐全,可+急 +3 ViqFlLSxHDGt 2026-09-26 3/150 2026-09-27 10:03 by 4IW0sJvtEqX8
[公派出国] 售SCI文章,我:8O5.5.1.O.54,科目齐全,可+急 +3 GxESXJmptOtQ 2026-09-25 4/200 2026-09-27 06:44 by XgCC8uTcwILl
[考博] 售SCI一区T0P文章,我:8.O.55.1.O.5.4,科目全,可+急 +3 GxESXJmptOtQ 2026-09-25 4/200 2026-09-27 06:42 by XgCC8uTcwILl
[论文投稿] 售SCI一区T0P文章,我:8O.55.1.O.5.4,科目齐全,可+急 +3 GxESXJmptOtQ 2026-09-25 4/200 2026-09-27 06:39 by XgCC8uTcwILl
[找工作] 售SCI一区T0P文章,我:8.O.55.1.O54,科目全,可伽急 +3 GxESXJmptOtQ 2026-09-25 3/150 2026-09-27 06:32 by XgCC8uTcwILl
[公派出国] 售SCI一区T0P文章,我:8.O55.1.O.54,科目全,可十急 +3 vQZoDrrm7mUF 2026-09-25 3/150 2026-09-27 06:05 by XgCC8uTcwILl
[考研] 售一区SCI文章T0P,我:8O.551.O54,科目全,可十急 +3 vQZoDrrm7mUF 2026-09-25 3/150 2026-09-27 05:56 by XgCC8uTcwILl
[论文投稿] 售SCI一区文章,我:8.O.55.1.O.54,科目齐全,可伽急 +4 vQZoDrrm7mUF 2026-09-25 4/200 2026-09-27 05:56 by XgCC8uTcwILl
[找工作] 售SCI一区文章,我:8.O.551.O.5.4,科目全,可伽急 +5 vQZoDrrm7mUF 2026-09-25 5/250 2026-09-27 05:53 by XgCC8uTcwILl
[博后之家] 售SCI一区T0P文章,我:8.O55.1.O.54,科目全,可十急 +4 vQZoDrrm7mUF 2026-09-25 4/200 2026-09-27 05:50 by XgCC8uTcwILl
[公派出国] 售SCI文章,我:8O5.5.1.O.54,科目齐全,可+急 +4 vQZoDrrm7mUF 2026-09-25 4/200 2026-09-27 05:44 by XgCC8uTcwILl
[论文投稿] 售SCI一区T0P文章,我:8O.55.1.O.54,科目全,可伽急 +3 vQZoDrrm7mUF 2026-09-25 6/300 2026-09-27 05:42 by XgCC8uTcwILl
[公派出国] 售SCI-T0P文章,我:8O.5.5.1.O.54,科目齐全,可+急 +4 PGTSwSC3F6nQ 2026-09-24 7/350 2026-09-25 11:18 by yhq9807
[教师之家] 课题组招2027级博士 上海工程技术大学 激光智能制造方向 +4 oxidpl 2026-09-20 6/300 2026-09-24 20:37 by oxidpl
信息提示
请填处理意见