24小时热门版块排行榜    

CyRhmU.jpeg
南方科技大学公共卫生及应急管理学院2025级博士研究生招生报考通知
查看: 3398  |  回复: 15
当前只显示满足指定条件的回帖,点击这里查看本话题的所有回帖

lxyy

银虫 (小有名气)

[求助] matlab拟合方程参数时初值的选择

A+B=C  A+W=R  R+B=C+W
先考虑的三个反应都是可逆反应,有以下几个式子:编写程序如下
%30℃
clear all
clc
global keq1 keq2 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= k1plus/keq1;
k2=k(2);
k3plus=k(3);
k3minus= k3plus/keq3;

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

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

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

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

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];
但是k的初值一直影响结果,程序要怎么改一下初值对结果影响不大(这个程序有问题,)
回复此楼

» 猜你喜欢

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

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

lxyy

银虫 (小有名气)

引用回帖:
2楼: Originally posted by dingd at 2012-03-28 10:37:30:
微分方程dcBdt, dcWdt的表达式呢?

我写在附件里面了

» 本帖附件资源列表

  • 欢迎监督和反馈:小木虫仅提供交流平台,不对该内容负责。
    本内容由用户自主发布,如果其内容涉及到知识产权问题,其责任在于用户本人,如对版权有异议,请联系邮箱:xiaomuchong@tal.com
  • 附件 1 : 太残忍doc.doc
  • 2012-03-28 10:52:17, 30.5 K
3楼2012-03-28 10:52:39
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
查看全部 16 个回答

dingd

铁杆木虫 (职业作家)

【答案】应助回帖

感谢参与,应助指数 +1
微分方程dcBdt, dcWdt的表达式呢?
2楼2012-03-28 10:37:30
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

dingd

铁杆木虫 (职业作家)

【答案】应助回帖


zhangguangping: 金币+1, 谢谢指导! 2012-03-28 14:16:22
最好把你的专业问题用通俗的数学语言来描述,非化工专业的很难理解的;
你的问题应该就是一微分方程拟合问题,仅需描述清楚:

1:微分方程组表达式
2:对应的数据
3:要求的参数
4楼2012-03-28 11:07:18
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

lxyy

银虫 (小有名气)

引用回帖:
4楼: Originally posted by dingd at 2012-03-28 11:07:18:
最好把你的专业问题用通俗的数学语言来描述,非化工专业的很难理解的;
你的问题应该就是一微分方程拟合问题,仅需描述清楚:

1:微分方程组表达式
2:对应的数据
3:要求的参数

不好意思,我的微分方程组就是附件里面标(1)(2)(3)的三个式子,对应的数据为

时间: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];
要求的是ki+,k2,k3+,
另外方程中的ki-=ki+/3.14,k3-=k3+/1.08

» 本帖附件资源列表

  • 欢迎监督和反馈:小木虫仅提供交流平台,不对该内容负责。
    本内容由用户自主发布,如果其内容涉及到知识产权问题,其责任在于用户本人,如对版权有异议,请联系邮箱:xiaomuchong@tal.com
  • 附件 1 : 太残忍doc (1)1.doc
  • 2012-03-28 12:18:21, 26 K
5楼2012-03-28 12:18:31
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
信息提示
请填处理意见