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可以让我学习一下（我没搜到啊

）。谢谢

回复此楼

» 猜你喜欢

反铁磁体中的磁性切换：两种不同的机制已成功可视化已经有0人回复
求标准粉末衍射卡号 ICDD 01-076-1802 已经有0人回复
物理学I论文润色/翻译怎么收费? 已经有116人回复
新西兰Robinson研究所招收全奖PhD 已经有0人回复
石墨烯转移--二氧化硅衬底石墨烯已经有0人回复
笼目材料中量子自旋液体基态的证据已经有0人回复
数学教学论硕士可以读数学物理博士吗？已经有0人回复
德国亥姆霍兹Hereon中心汉堡分部招镁合金腐蚀裂变SCC课题方向2026公派博士生已经有4人回复
澳门大学应用物理及材料工程研究院潘晖教授课题组诚招博士后已经有11人回复

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

动力学反应级数已经有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];