24小时热门版块排行榜

返回列表

当前只显示满足指定条件的回帖，点击这里查看本话题的所有回帖

heisekaishui

银虫 (小有名气)

应助: 0 (幼儿园)
金币: 215
帖子: 52
在线: 59.2小时
虫号: 2445260

[交流] Matlab同时拟合活化能，反应级数，速率常数

各位，您好!
本人研究的方向是催化化学，最近在评估催化剂效果时想引入动力学的内容，对其中的一些地方有一些疑问。
受月之蓝版主的这个帖子http://muchong.com/bbs/viewthread.php?tid=6425538&target=self&page=1以及其它前辈的启发，对于反应级数以及反应速率常数的拟合有一定的了解。
但是如果我的动力学模型是这种形式：-dC/dt=k0*exp(-E/R*T)*C^m*P^n（只有一个反应物和氢气），其中P与T分别为压力与时间，就是说是一个加氢反应，试验中可以得到不同压力、不同温度下浓度随时间的变化值，应该如何拟合k0,E,m,n。
我的一个想法是先拟合不同压力下的k1（k1=k0*exp(-E/R*T)*P^n)（例如以上述帖子的代码方法），再根据k1=k2*P^n（k2=k0*exp(-E/R*T)）拟合出k2和n，此后拟合出不同温度下的k2，根据k2=k0*exp(-E/R*T)）拟合出E和k0
但是我觉得这种方法较为复杂，一个matlab里有好几个拟合程序，不知道是否有更简洁的CODE。

PS:不知道版里有没有相似的帖子有类似的CODE可以让我学习一下（我没搜到啊

）。谢谢

回复此楼

» 猜你喜欢

散金币祈福已经有186人回复
求助，如何提取ELK的rt-TDDFT在某一时刻的自旋密度分布和电子密度分布已经有1人回复
物理学I论文润色/翻译怎么收费? 已经有179人回复
澳门大学IAPME李海峰教授常年招收博士研究生已经有0人回复
澳门大学IAPME李海峰教授常年招收博士研究生已经有0人回复
澳门大学IAPME李海峰教授常年招收博士研究生已经有0人回复
深圳大学应用技术学院招聘凝聚态物理博士后已经有114人回复
深圳大学2027级光电信息工程/物理学研究生推免已经有70人回复
求助中医院临床中医中标的青年标书学习，谢谢已经有4人回复
浙大物院表面物理课题组长期招收普博及直博学生已经有257人回复
【征稿】Springer Nature热电专题 Discover Materials / Discover Chemistry 已经有39人回复

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

动力学反应级数已经有4人回复
求助用matlab三次样条法算一个反应的反应速率已经有13人回复
求助用matlab拟合动力学方程已经有9人回复
速率方程的问题求教！已经有4人回复
有关反应速率常数的计算已经有9人回复
SCR反应动力学方程求解的问题已经有4人回复
催化剂增大了速率常数但没有改变活化能-来自一篇经典文献，如何理解？已经有3人回复
动力学参数拟合已经有26人回复
matlab 拟合反应动力学参数结果很差。大家帮忙看一下已经有14人回复
根据活化能怎么判断一级反应还是零级反应已经有12人回复
焓变活化能键离解能能垒反应速率的巩固求教已经有14人回复
帮我看看聚合反应三个速率常数怎么算？已经有9人回复
【求助】拟合动力学方程求助已经有13人回复
【求助】催化反应动力学matlab计算各基元反应的速率常数时，该如何避免较小量被忽略？已经有3人回复
【求助】反应速率常数怎么算已经有5人回复
【求助】关于催化反应中的活化能已经有24人回复
【求助】使用Matlab拟合反应动力学方程问题已经有7人回复
【求助】使用Matlab预估动力学方程问题已经有13人回复

1楼 2014-08-25 17:30:13

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

whqs8426212

铜虫 (正式写手)

应助: 5 (幼儿园)
金币: 13.2
帖子: 782
在线: 685.7小时
虫号: 2572504

★
小木虫: 金币+0.5, 给个红包，谢谢回帖

引用回帖:

16楼: Originally posted by heisekaishui at 2014-08-27 11:00:36
非常感谢，对，数据格式是这样的。
主要是数据还没有做出来，如果可以的话可以让我看一下程序吗，我想看一下大概的思路。...

function k1k2k3k4
format long
clear all
clc
tspan = [ 0.6 1.2  2  3 ];
C0 = [2 20 200]; %初值只有一个C0
k0 = [0  0 0 0]; %k0 E M n 参数的初值
lb = [0  0 0 0]; %上届
%ub = [100  10000 1000 1000 ]; %下界

dataTP1=...
[
%t       C
XX    XX
XXX    XX
XXX    X
XX    X    ];    %实验数据1
dataTP2=...
[
%t       C
XX    XX
XXX    XX
XXX    X
XX    X    ]; %实验数据2
dataTP3=...
[
%t       C
XX    XX
XXX    XX
XXX    X
XX    X    ]; %实验数据3

[k,resnorm,residual,exitflag,output,lambda,jacobian] = ...
lsqnonlin(@ObjFunc,k0,lb,[],options,tspan,C0,dataTP1,dataTP2,dataTP3);
ci = nlparci(k,residual,jacobian);
fprintf('\n\n使用函数lsqnonlin()估计得到的参数值为:\n')
fprintf('\tk0 = %.9f \n',k(1))
fprintf('\tE = %.9f \n',k(2))
fprintf('\tm = %.9f \n',K(3))
fprintf('\tn = %.9f \n',k(4))

fprintf('  The sum of the squares is: %.9e\n\n',resnorm)

%-----------------------------------------------------
function f = ObjFunc(k,tspan,C0,dataTP1,dataTP2,dataTP3)          % 目标函数
T1=?;P1=?;
[t XsimTP1] = ode23s(@KineticsEqs,tspan,C0(1),[],k,T1,P1);
T2=? ;P2=?;
[t XsimTP2] = ode23s(@KineticsEqs,tspan,C0(2),[],k,T2,P2);
T3=?;P3=?;
[t XsimTP3] = ode23s(@KineticsEqs,tspan,C0(3),[],k,T3,P3);
f = [(XsimTP1(:,1)-dataTP1(:,2)) (XsimTP2(:,1)-dataTP2(:,2)) (XsimTP3(:,1)-dataTP3(:,2))];

%----------------------------------------------------------
function dCdt = KineticsEqs(t,C,k,T,P)             % ODE模型方程
R=8;
dC=-k(1)*exp(-k(2)/R*T)*(C^k(3))*(P^k(4));       %k(1)=k0,k(2)=E,k(3)=m;K(4)=n
dCdt = [dC];