版塊導(dǎo)航: 正在加載中...

登錄注冊

應(yīng)《網(wǎng)絡(luò)安全法》要求，自2017年10月1日起，未進(jìn)行實(shí)名認(rèn)證將不得使用互聯(lián)網(wǎng)跟帖服務(wù)。為保障您的帳號(hào)能夠正常使用，請盡快對帳號(hào)進(jìn)行手機(jī)號(hào)驗(yàn)證，感謝您的理解與支持！

24小時(shí)熱門版塊排行榜

北京石油化工學(xué)院2026年研究生招生接收調(diào)劑公告

返回列表

apolloking

金蟲 (小有名氣)

應(yīng)助: 0 (幼兒園)
金幣: 1083.1
散金: 35
紅花: 2
帖子: 156
在線: 60.5小時(shí)
蟲號(hào): 75702
注冊: 2005-06-20
專業(yè): 能源化工

[交流] 【求助】使用Matlab預(yù)估動(dòng)力學(xué)方程問題已有4人參與

大家好，問題如下：
使用固定床積分反應(yīng)器來測量甲烷催化燃燒的本征動(dòng)力學(xué)實(shí)驗(yàn)數(shù)據(jù)。
整個(gè)計(jì)算過程見圖片附件

本征動(dòng)力學(xué)方程為ra=k0*exp（-E/(R*T)）*Ca
ra:反應(yīng)速率 mol/(kgcat*s)
k0:指前因子 m3/(kgcat*s)
E:活化能 kJ/mol
Ca：甲烷摩爾濃度 mol/m3

待定參數(shù)為k0和E

文獻(xiàn)使用了Matlab的半隱式龍格庫塔法和最小二乘法去約束優(yōu)化
最后得出的k0為3.77*10E8
E為25kJ/mol

在積分反應(yīng)器上有;
-ra=dxa/d(w/Fa0)=-k0*exp（-E/(R*T)）*Ca
xa為轉(zhuǎn)化率 w為催化劑的質(zhì)量 Fa0為流速
又因?yàn)镃a=C0*(1-xa)
所以上式轉(zhuǎn)化為
dxa/d(w/Fa0)=-k0*exp（-E/(R*T)）*C0（1-xa）
因?yàn)閷?shí)驗(yàn)所用的甲烷體積濃度為1%，上式可轉(zhuǎn)化為
1%*22.4L/mol

dxa/d(w/Fa0)=-k0*exp（-E/(R*T)）*（1-xa）/2.24
（我的理解是w/Fa0是空速的倒數(shù),但文獻(xiàn)中給的空速值單位是ml/(g*h)，我不知道需不需要在進(jìn)行單位轉(zhuǎn)化，還是直接取倒數(shù)就行了）

文獻(xiàn)中說，“采用Matlab程序?qū)Ρ菊鲃?dòng)力學(xué)實(shí)驗(yàn)數(shù)據(jù)進(jìn)行搜索，選優(yōu)”求的k0和E
按照這個(gè)思路，是不是針對每一個(gè)溫度，認(rèn)為溫度值已知，空速和轉(zhuǎn)化率為變量求的一組k0和E。然后將所有溫度得到的k0和E值取平均值？
我自己寫了一個(gè)550度的計(jì)算程序，但是得不到正確的結(jié)果，希望大家能幫忙看一下。

function jifen
%積分法動(dòng)力學(xué)參數(shù)估算
close all,
clear,
clc,
%實(shí)驗(yàn)數(shù)據(jù)
t=[1/5000 1/10000 1/15000 1/20000 1/25000];
ca=[0.989578 0.972126 0.956323 0.897261 0.870735];
%最小二乘法估計(jì)動(dòng)力學(xué)參數(shù)
x0=[0 1];
x=lsqnonlin(@myfun,x0,[],[],optimset('Display','iter'),ca,t);
k0=x(1);
E=x(2);
%給出估算結(jié)果
fprintf('動(dòng)力學(xué)參數(shù)為：\n'),
fprintf('k=%f\ta=%f\n',k0,E),
%繪制反應(yīng)速率擬合曲線
[T,CA]=ode45(@reaction_rate,[0,1/5000],ca(1),[],k0,a);
plot(t,ca,'*k'),
hold on,
plot(T,CA,'-k'),
xlabel('t/min')
ylabel('C_A/mol/L')
legend('實(shí)驗(yàn)數(shù)據(jù)','模擬曲線')
%實(shí)驗(yàn)導(dǎo)數(shù)與計(jì)算導(dǎo)數(shù)之差
function F=myfun(x,ca,tspan)
k0=x(1);
E=x(2);
[t,cac]=ode45(@reaction_rate,tspan,ca(1),[],k0,E);
F=cac-ca';
%反應(yīng)速率方程
function r=reaction_rate(t,ca,k0,E)
R=8.314;
T=823.15;
r=k0*exp(-E/R/T)*(1-ca)/2.24;

[ Last edited by apolloking on 2009-7-18 at 18:11 ]

回復(fù)此樓

1樓 2009-07-18 17:24:51

已閱回復(fù)此樓關(guān)注TA 給TA發(fā)消息送TA紅花 TA的回帖

snipher950

木蟲 (正式寫手)

應(yīng)助: 2 (幼兒園)
金幣: 3024.6
紅花: 3
帖子: 797
在線: 22.5小時(shí)
蟲號(hào): 522854
注冊: 2008-03-12
性別: MM
專業(yè): 化工系統(tǒng)工程

私下解答。

回復(fù)此樓

2樓2009-07-19 14:48:56

已閱回復(fù)此樓關(guān)注TA 給TA發(fā)消息送TA紅花 TA的回帖

apolloking

金蟲 (小有名氣)

應(yīng)助: 0 (幼兒園)
金幣: 1083.1
散金: 35
紅花: 2
帖子: 156
在線: 60.5小時(shí)
蟲號(hào): 75702
注冊: 2005-06-20
專業(yè): 能源化工

function jifen
%積分法動(dòng)力學(xué)參數(shù)估算
close all,
clear,
clc,
%實(shí)驗(yàn)數(shù)據(jù)
t=[350 375 400 425 450 475 500 525 550 350 375 400 425 450 475 500 525 550 350 375 400 425 450 475 500 525 550 350 375 400 425 450 475 500 525 550 350 375 400 425 450 475 500 525 550];
w=[22.4*3600/5000 22.4*3600/5000 22.4*3600/5000 22.4*3600/5000 22.4*3600/5000 22.4*3600/5000 22.4*3600/5000 22.4*3600/5000 22.4*3600/5000 22.4*3600/10000 22.4*3600/10000 22.4*3600/10000 22.4*3600/10000 22.4*3600/10000 22.4*3600/10000 22.4*3600/10000 22.4*3600/10000 22.4*3600/10000 22.4*3600/15000 22.4*3600/15000 22.4*3600/15000 22.4*3600/15000 22.4*3600/15000 22.4*3600/15000 22.4*3600/15000 22.4*3600/15000 22.4*3600/15000 22.4*3600/20000 22.4*3600/20000 22.4*3600/20000 22.4*3600/20000 22.4*3600/20000 22.4*3600/20000 22.4*3600/20000 22.4*3600/20000 22.4*3600/20000 22.4*3600/25000 22.4*3600/25000 22.4*3600/25000 22.4*3600/25000 22.4*3600/25000 22.4*3600/25000 22.4*3600/25000 22.4*3600/25000 22.4*3600/25000];
xa=[0.224067 0.308978 0.418433 0.583276 0.712133 0.834063 0.927528 0.971279 0.989578 0.167333 0.26823 0.359565 0.492762 0.605102 0.706325 0.870419 0.936043 0.972126 0.161258 0.229717 0.305041 0.445914 0.562024 0.702027 0.847299 0.901723 0.956323 0.703751 0.413094 0.227758 0.297666 0.456384 0.592814 0.743624 0.825191 0.897261 0.089297 0.149173 0.213542 0.304452 0.433972 0.563372 0.703383 0.801953 0.870735];
%最小二乘法估計(jì)動(dòng)力學(xué)參數(shù)
x0=[3E8 20000];
x=lsqnonlin(@myfun,x0,[],[],optimset('Display','iter'),t,xa,w);
k0=x(1);
e=x(2);
%給出估算結(jié)果
fprintf('動(dòng)力學(xué)參數(shù)為：\n'),
fprintf('k=%f\ta=%f\n',k0,e),
%繪制反應(yīng)速率擬合曲線
[W,XA]=ode45(@reaction_rate,[0,22.4*3600/5000],xa(1),[],k0,e);
plot(w,xa,'*r'),
hold on,
plot(W,XA,'-r'),
xlabel('t/min')
ylabel('C_A/mol/L')
legend('實(shí)驗(yàn)數(shù)據(jù)','模擬曲線')
%實(shí)驗(yàn)導(dǎo)數(shù)與計(jì)算導(dǎo)數(shù)之差
function F=myfun(x,xa,tspan)
k0=x(1);
e=x(2);
[w,xac]=ode45(@reaction_rate,tspan,xa(1),[],k0,e);
F=xac-xa';
%反應(yīng)速率方程
function r=reaction_rate(t,w,xa,k0,e)
r=-k0*exp(-e/8.314/(t+273.15))(1-xa)/2.24;

錯(cuò)誤提示
??? Error: File: D:\My Documents\Untitled3.m Line: 34 Column: 7
()-indexing must appear last in an index expression.

第34行是
r=-k0*exp(-e/8.314/(t+273.15))(1-xa)/2.24;

[ Last edited by apolloking on 2009-7-20 at 11:53 ]

贊一下

回復(fù)此樓

3樓2009-07-20 11:51:15

已閱回復(fù)此樓關(guān)注TA 給TA發(fā)消息送TA紅花 TA的回帖

change0618

鐵桿木蟲 (著名寫手)

方丈大師

應(yīng)助: 44 (小學(xué)生)
金幣: 17724.5
紅花: 17
帖子: 2413
在線: 546.7小時(shí)
蟲號(hào): 496517
注冊: 2008-01-19
專業(yè): 化學(xué)反應(yīng)工程

★
小木蟲(金幣+0.5):給個(gè)紅包，謝謝回帖交流

r=-k0*exp(-e/8.314/(t+273.15))*(1-xa)/2.24;

贊一下(1人)

回復(fù)此樓

4樓2009-07-20 13:01:01

已閱回復(fù)此樓關(guān)注TA 給TA發(fā)消息送TA紅花 TA的回帖

apolloking

金蟲 (小有名氣)

應(yīng)助: 0 (幼兒園)
金幣: 1083.1
散金: 35
紅花: 2
帖子: 156
在線: 60.5小時(shí)
蟲號(hào): 75702
注冊: 2005-06-20
專業(yè): 能源化工

改了這個(gè)也不對，運(yùn)行還是有錯(cuò)誤

贊一下

回復(fù)此樓

5樓2009-07-20 19:13:06

已閱回復(fù)此樓關(guān)注TA 給TA發(fā)消息送TA紅花 TA的回帖

hitzhang

木蟲 (正式寫手)

仿真EPI: 1
應(yīng)助: 0 (幼兒園)
貴賓: 2.15
金幣: 1376.7
散金: 969
紅花: 8
帖子: 862
在線: 226.4小時(shí)
蟲號(hào): 390575
注冊: 2007-06-02
性別: GG
專業(yè): 無機(jī)非金屬類電介質(zhì)與電解

★
小木蟲(金幣+0.5):給個(gè)紅包，謝謝回帖交流

這句x=lsqnonlin(@myfun,x0,[],[],optimset('Display','iter'),t,xa,w);，t, xa, w應(yīng)該是myfun的額外輸入?yún)?shù)，而你的myfun函數(shù)有x, xa, tspan，不懂作者是什么意思

贊一下(1人)

回復(fù)此樓

6樓2009-07-20 22:03:01

已閱回復(fù)此樓關(guān)注TA 給TA發(fā)消息送TA紅花 TA的回帖

apolloking

金蟲 (小有名氣)

應(yīng)助: 0 (幼兒園)
金幣: 1083.1
散金: 35
紅花: 2
帖子: 156
在線: 60.5小時(shí)
蟲號(hào): 75702
注冊: 2005-06-20
專業(yè): 能源化工

實(shí)話說這個(gè)程序是我從另外一個(gè)程序改過來的，我也不是很清楚這個(gè) tspan具體的含義

贊一下

回復(fù)此樓

7樓2009-07-20 22:05:23

已閱回復(fù)此樓關(guān)注TA 給TA發(fā)消息送TA紅花 TA的回帖

apolloking

金蟲 (小有名氣)

應(yīng)助: 0 (幼兒園)
金幣: 1083.1
散金: 35
紅花: 2
帖子: 156
在線: 60.5小時(shí)
蟲號(hào): 75702
注冊: 2005-06-20
專業(yè): 能源化工

snipher950 網(wǎng)友說先固定一個(gè)溫度，然后輸入空速和轉(zhuǎn)化率，但是之后我就不知道怎么辦了。我看相關(guān)的文獻(xiàn)好像都是一次性就做出來一個(gè)k0和E的

贊一下

回復(fù)此樓

8樓2009-07-20 22:07:53

已閱回復(fù)此樓關(guān)注TA 給TA發(fā)消息送TA紅花 TA的回帖

hitzhang

木蟲 (正式寫手)

仿真EPI: 1
應(yīng)助: 0 (幼兒園)
貴賓: 2.15
金幣: 1376.7
散金: 969
紅花: 8
帖子: 862
在線: 226.4小時(shí)
蟲號(hào): 390575
注冊: 2007-06-02
性別: GG
專業(yè): 無機(jī)非金屬類電介質(zhì)與電解

是這樣，能把你程序里的函數(shù)調(diào)用關(guān)系以及各個(gè)函數(shù)的作用描述一下嗎？

贊一下

回復(fù)此樓

9樓2009-07-20 22:21:52

已閱回復(fù)此樓關(guān)注TA 給TA發(fā)消息送TA紅花 TA的回帖

apolloking

金蟲 (小有名氣)

應(yīng)助: 0 (幼兒園)
金幣: 1083.1
散金: 35
紅花: 2
帖子: 156
在線: 60.5小時(shí)
蟲號(hào): 75702
注冊: 2005-06-20
專業(yè): 能源化工

function jifen
%積分法動(dòng)力學(xué)參數(shù)估算
close all,
clear,
clc,
%實(shí)驗(yàn)數(shù)據(jù)
t=[350 375 400 425 450 475 500 525 550 350 375 400 425 450 475 500 525 550 350 375 400 425 450 475 500 525 550 350 375 400 425 450 475 500 525 550 350 375 400 425 450 475 500 525 550];
w=[22.4*3600/5000 22.4*3600/5000 22.4*3600/5000 22.4*3600/5000 22.4*3600/5000 22.4*3600/5000 22.4*3600/5000 22.4*3600/5000 22.4*3600/5000 22.4*3600/10000 22.4*3600/10000 22.4*3600/10000 22.4*3600/10000 22.4*3600/10000 22.4*3600/10000 22.4*3600/10000 22.4*3600/10000 22.4*3600/10000 22.4*3600/15000 22.4*3600/15000 22.4*3600/15000 22.4*3600/15000 22.4*3600/15000 22.4*3600/15000 22.4*3600/15000 22.4*3600/15000 22.4*3600/15000 22.4*3600/20000 22.4*3600/20000 22.4*3600/20000 22.4*3600/20000 22.4*3600/20000 22.4*3600/20000 22.4*3600/20000 22.4*3600/20000 22.4*3600/20000 22.4*3600/25000 22.4*3600/25000 22.4*3600/25000 22.4*3600/25000 22.4*3600/25000 22.4*3600/25000 22.4*3600/25000 22.4*3600/25000 22.4*3600/25000];
xa=[0.224067 0.308978 0.418433 0.583276 0.712133 0.834063 0.927528 0.971279 0.989578 0.167333 0.26823 0.359565 0.492762 0.605102 0.706325 0.870419 0.936043 0.972126 0.161258 0.229717 0.305041 0.445914 0.562024 0.702027 0.847299 0.901723 0.956323 0.703751 0.413094 0.227758 0.297666 0.456384 0.592814 0.743624 0.825191 0.897261 0.089297 0.149173 0.213542 0.304452 0.433972 0.563372 0.703383 0.801953 0.870735];

%以上是我要輸入的參數(shù)，一共是45組，每組參數(shù)分別為溫度t 空速的倒數(shù)w 轉(zhuǎn)化率xa
%最小二乘法估計(jì)動(dòng)力學(xué)參數(shù)
x0=[3E8 20000];
x=lsqnonlin(@myfun,x0,[],[],optimset('Display','iter'),t,xa,w);

%這個(gè)是非線性最小二乘法來調(diào)用myfun這個(gè)函數(shù)來判斷我計(jì)算出來的轉(zhuǎn)化率和實(shí)際測得的轉(zhuǎn)化率之差
k0=x(1);
e=x(2);
%給出估算結(jié)果
fprintf('動(dòng)力學(xué)參數(shù)為：\n'),
fprintf('k=%f\ta=%f\n',k0,e),
%實(shí)驗(yàn)導(dǎo)數(shù)與計(jì)算導(dǎo)數(shù)之差
function F=myfun(x,xa,tspan)
k0=x(1);
e=x(2);
[w,xac]=ode45(@reaction_rate,tspan,xa(1),[],k0,e);
F=xac-xa';

%這個(gè)用來計(jì)算轉(zhuǎn)化率
%反應(yīng)速率方程
function r=reaction_rate(t,w,xa,k0,e)
r=-k0*exp(-e/8.314/(t+273.15))(1-xa)/2.24;

%這個(gè)是我的動(dòng)力學(xué)方程

贊一下

回復(fù)此樓

10樓2009-07-20 22:30:53

已閱回復(fù)此樓關(guān)注TA 給TA發(fā)消息送TA紅花 TA的回帖

相關(guān)版塊跳轉(zhuǎn) 我要訂閱樓主 apolloking 的主題更新

返回列表

普通表情龍兔虎貓高級(jí)回復(fù) (可上傳附件)

最具人氣熱帖推薦 [查看全部]		作者	回/看	最后發(fā)表

[考研] 材料化工340求調(diào)劑 +3	jhx777 2026-03-30	3/150	2026-03-30 17:54 by JourneyLucky
[考研] 307求調(diào)劑 +16	超級(jí)伊昂大王 2026-03-24	17/850	2026-03-30 17:15 by 544594351
[考研] 342求調(diào)劑 +4	加油a李zs 2026-03-26	4/200	2026-03-30 16:39 by 晶體之美
[考研] 324求調(diào)劑 +9	hanamiko 2026-03-26	11/550	2026-03-30 14:27 by JourneyLucky
[考研] 0703本科鄭州大學(xué)求調(diào)劑 +7	nhj_ 2026-03-25	7/350	2026-03-30 12:44 by fangnagu
[考研] 282求調(diào)劑 +4	wcq131415 2026-03-24	4/200	2026-03-30 10:39 by Delta2012
[考研] 0856求調(diào)劑 +8	楒桉 2026-03-28	8/400	2026-03-30 10:00 by wzy-lxz
[考研] 085600，材料與化工321分求調(diào)劑 +10	大饞小子 2026-03-28	10/500	2026-03-29 23:35 by 飛行日記西
[考研] 求調(diào)劑 +4	QiMing7 2026-03-25	5/250	2026-03-29 21:10 by 唐沐兒
[考研] 281求調(diào)劑 +4	亞克西good 2026-03-26	6/300	2026-03-29 21:06 by 啊李999
[考研] 275求調(diào)劑 +15	Micky11223 2026-03-25	20/1000	2026-03-29 20:44 by 唐沐兒
[考研] 調(diào)劑考研 +3	王杰一 2026-03-29	3/150	2026-03-29 08:09 by fmesaito
[考研] 081200-314 +3	LILIQQ 2026-03-27	4/200	2026-03-28 09:41 by 保護(hù)地球你我做?/a>
[考研] 085405 考的11408求各位老師帶走 +3	Qiu學(xué)ing 2026-03-28	3/150	2026-03-28 09:19 by 樂呵呵的追夢人
[有機(jī)交流] 高溫高壓反應(yīng)求助 10+4	chibby 2026-03-25	4/200	2026-03-27 21:08 by BT20230424
[考研] 安徽大學(xué)專碩生物與醫(yī)藥專業(yè)(086000)324分，英語已過四六級(jí)，六級(jí)521，求調(diào)劑 +4	美味可樂雞翅 2026-03-26	4/200	2026-03-27 15:27 by 星空星月
[考研] 0856調(diào)劑 +5	求求讓我有書讀?/a> 2026-03-26	6/300	2026-03-27 15:12 by caszguilin
[考研] 322求調(diào)劑 +4	我真的很想學(xué)習(xí) 2026-03-23	4/200	2026-03-27 13:51 by 楊楊楊紫
[考研] 315調(diào)劑 +4	0860求調(diào)劑 2026-03-26	5/250	2026-03-27 11:23 by wangjy2002
[考研] 生物學(xué)學(xué)碩求調(diào)劑 +7	小羊睡著了? 2026-03-23	10/500	2026-03-25 02:24 by 清風(fēng)拂揚(yáng)。 m

24小時(shí)熱門版塊排行榜

[交流] 【求助】使用Matlab預(yù)估動(dòng)力學(xué)方程問題 已有4人參與

» 收錄本帖的淘帖專輯推薦

» 猜你喜歡

» 本主題相關(guān)價(jià)值貼推薦，對您同樣有幫助:

[交流] 【求助】使用Matlab預(yù)估動(dòng)力學(xué)方程問題已有4人參與