| 5 | 1/1 | 返回列表 |
| 查看: 549 | 回復(fù): 3 | ||
| 當(dāng)前只顯示滿足指定條件的回帖,點(diǎn)擊這里查看本話題的所有回帖 | ||
你的CEO鐵蟲 (小有名氣)
|
[求助]
有關(guān)反應(yīng)動(dòng)力學(xué)推導(dǎo)的matlab程序
|
|
|
請(qǐng)蟲友大神幫小弟看看為啥誤差這么大,小弟在此謝過了 *****程序行數(shù)較多帖子不能完全 顯示,只能添加附件了,煩請(qǐng)大神幫看下,小弟在此謝過了!********* |

鐵蟲 (小有名氣)

鐵蟲 (小有名氣)
|
function KineticsEst1_int % 動(dòng)力學(xué)參數(shù)辨識(shí): 用積分法進(jìn)行反應(yīng)速率分析得到速率常數(shù)k和反應(yīng)級(jí)數(shù)n % Analysis of kinetic rate data by using the integral method % Reaction of the type -- rate = kCA^order % order - reaction order % rate -- reaction rate vector % YA -- yield vector for reactant A % T -- vector of reaction time % N -- number of data points % k- reacion rate constant clear all clc format short global Y_exp Y_sim tspan = [0,21.065]; % t=w/F(MEOH.0),W:催化劑的質(zhì)量,F(xiàn)(MEOH.0):甲醇的初始進(jìn)料流率 k0 = [0 0 0 0]; lb = [0 0 0 0]; ub = [+inf +inf +inf +inf ]; Y0 = [0 0 0 0 ]; Y_exp =[0.2245 0.0086 0.0013 0.0038; 0.2503 0.0082 0.0009 0.0018; 0.267 0.0098 0.0013 0.003; 0.2814 0.0094 0.0015 0.0012; ]; %同一空速四個(gè)不同組成配比下的反應(yīng)管出口生成物的收率Y % 使用函數(shù)lsqnonlin()進(jìn)行參數(shù)估計(jì) [k,resnorm,residual,exitflag,output,lambda,jacobian] = lsqnonlin(@ObjFunc,k0,lb,ub,optimset('TolFun',1.0000e-20),tspan,Y0,Y_exp); ci = nlparci(k,residual,jacobian) %[k,resnorm,residual,exitflag,output,lambda,jacobian] = lsqnonlin(@ObjFunc4LNL,k0,lb,ub,optimset('TolFun',1.0000e-6),Y0,Y_exp); %ci = nlparci(k,residual,jacobian) fprintf('\n\n使用函數(shù)lsqnonlin()估計(jì)得到的參數(shù)值為:\n') fprintf('\tk1 = %.4f ± %.4f\n',k(1),ci(1,2)-k(1)) fprintf('\tk2 = %.4f ± %.4f\n',k(2),ci(2,2)-k(2)) fprintf('\tk3 = %.4f ± %.4f\n',k(3),ci(3,2)-k(3)) fprintf('\tk4 = %.4f ± %.4f\n',k(4),ci(4,2)-k(4)) fprintf('\tThe sum of the squares is: %.1e\n\n',resnorm) %擬合效果圖(實(shí)驗(yàn)與擬合的比較) %[t4plot Y4plot] = ode45(@KineticsEqs1,[tspan(1) tspan(end)],Y0,[],k0) %figure %plot(tspan,Y_exp(:,1),'bo',t4plot,Y4plot(:,1),'r--'); %legend('Exp','Model') %xlabel('空時(shí)t=w/F_A_0, h') %ylabel('收率Y_D_M_M') %title('擬合效果圖') %figure %plot(tspan,Y_exp(:,2),'bo',t4plot,Y4plot(:,2),'r--'); %legend('Exp','Model') %xlabel('空時(shí)t=w/F_A_0, h') %ylabel('收率Y_M_F') %title('擬合效果圖') %figure %plot(tspan,Y_exp(:,3),'bo',t4plot,Y4plot(:,3),'r--'); %legend('Exp','Model') %xlabel('空時(shí)t=w/F_A_0, h') %ylabel('收率Y_F_A') %title('擬合效果圖') %figure %plot(tspan,Y_exp(:,4),'bo',t4plot,Y4plot(:,4),'r--'); %legend('Exp','Model') %xlabel('空時(shí)t=w/F_A_0, h') %ylabel('收率Y_D_M_E') %title('擬合效果圖') function f = ObjFunc(k,tspan,Y0,Y_exp) % 目標(biāo)函數(shù) [t,Y_sim1] = ode45(@KineticsEqs1,tspan,Y0,[],k) f1 = 1*(Y_sim1(end,1)-Y_exp(1,1)); f2= 10*(Y_sim1(end,2)-Y_exp(1,2)); f3= 10*(Y_sim1(end,3)-Y_exp(1,3)); f4= 10*(Y_sim1(end,4)-Y_exp(1,4)); [t,Y_sim2] = ode45(@KineticsEqs2,tspan,Y0,[],k) f5 = 1*(Y_sim2(end,1)-Y_exp(2,1)); f6 =10*(Y_sim2(end,2)-Y_exp(2,2)); f7 =10*(Y_sim2(end,3)-Y_exp(2,3)); f8 =10*(Y_sim2(end,4)-Y_exp(2,4)); [t,Y_sim3] = ode45(@KineticsEqs3,tspan,Y0,[],k) f9 = 1*(Y_sim3(end,1)-Y_exp(3,1)); f10 =10*(Y_sim3(end,2)-Y_exp(3,2)); f11 =10*(Y_sim3(end,3)-Y_exp(3,3)); f12 =10*(Y_sim3(end,4)-Y_exp(3,4)); [t,Y_sim4] = ode45(@KineticsEqs4,tspan,Y0,[],k) f13 = 1*(Y_sim4(end,1)-Y_exp(4,1)); f14 =10*(Y_sim4(end,2)-Y_exp(4,2)); f15 =10*(Y_sim4(end,3)-Y_exp(4,3)); f16 =10*(Y_sim4(end,4)-Y_exp(4,4)); f=[f1 f2 f3 f4; f5 f6 f7 f8; f9 f10 f11 f12; f13 f14 f15 f16; ] % ------------------------------------------------------------------ function dYdt = KineticsEqs1(t,Y,k) %p=[0.9 0.5 0.7 0.8 0.6]; %m=[0.6 0.4 0.5 0.3 0.25]; m:氧醇比 %n=[5.4 2 2 2.1 2]; % n:氮醇比 yMEOH1=(1+0.6)/(0.6+2.5*Y(1)+2*Y(2)+1.5*Y(3)+2*Y(4)+5.4); yO21=(0.6-0.5*Y(1)-Y(2)-0.5*Y(3))/(0.6+2.5*Y(1)+2*Y(2)+1.5*Y(3)+2*Y(4)+5.4); a=Y(1) dYdt=... [k(1)*(0.9^2)*yMEOH1^0.5*yO21^1.5; k(2)*(0.9^2)*yMEOH1^1*yO21^1; k(3)*(0.9^2)*yMEOH1^1*yO21^1; k(4)*(0.9^2)*yMEOH1^1.0*yO21^1; ]; function dYdt = KineticsEqs2(t,Y,k) %p=[0.9 0.5 0.7 0.8 0.6]; %m=[0.6 0.4 0.5 0.3 0.25]; m:氧醇比 %n=[5.4 2 2 2.1 2]; % n:氮醇比 yMEOH2=(1+0.4)/(0.6+2.5*Y(1)+2*Y(2)+1.5*Y(3)+2*Y(4)+2); yO22=(0.4-0.5*Y(1)-Y(2)-0.5*Y(3))/(0.4+2.5*Y(1)+2*Y(2)+1.5*Y(3)+2*Y(4)+2); dYdt=... [k(1)*(0.5^2)*yMEOH2^0.5*yO22^1.5; k(2)*(0.5^2)*yMEOH2^1*yO22^1; k(3)*(0.5^2)*yMEOH2^1*yO22^1; k(4)*(0.5^2)*yMEOH2^1.0*yO22^1; ]; function dYdt = KineticsEqs3(t,Y,k) %p=[0.9 0.5 0.7 0.8 0.6]; %m=[0.6 0.4 0.5 0.3 0.25]; m:氧醇比 %n=[5.4 2 2 2.1 2]; % n:氮醇比 yMEOH3=(1+0.5)/(0.6+2.5*Y(1)+2*Y(2)+1.5*Y(3)+2*Y(4)+2); yO23=(0.5-0.5*Y(1)-Y(2)-0.5*Y(3))/(0.5+2.5*Y(1)+2*Y(2)+1.5*Y(3)+2*Y(4)+2); dYdt=... [k(1)*(0.7^2)*yMEOH3^0.5*yO23^1.5; k(2)*(0.7^2)*yMEOH3^1*yO23^1; k(3)*(0.7^2)*yMEOH3^1*yO23^1; k(4)*(0.7^2)*yMEOH3^1.0*yO23^1; ]; function dYdt = KineticsEqs4(t,Y,k) %p=[0.9 0.5 0.7 0.8 0.6]; %m=[0.6 0.4 0.5 0.3 0.25]; m:氧醇比 %n=[5.4 2 2 2.1 2]; % n:氮醇比 yMEOH4=(1+0.8)/(0.6+2.5*Y(1)+2*Y(2)+1.5*Y(3)+2*Y(4)+2.1); yO24=(0.3-0.5*Y(1)-Y(2)-0.5*Y(3))/(0.3+2.5*Y(1)+2*Y(2)+1.5*Y(3)+2*Y(4)+2.1); dYdt=... [k(1)*(0.8^2)*yMEOH4^0.5*yO24^1.5; k(2)*(0.8^2)*yMEOH4^1*yO24^1; k(3)*(0.8^2)*yMEOH4^1*yO24^1; k(4)*(0.8^2)*yMEOH4^1.0*yO24^1; ]; |

| 最具人氣熱帖推薦 [查看全部] | 作者 | 回/看 | 最后發(fā)表 | |
|---|---|---|---|---|
|
[考研] 307求調(diào)劑 +17 | 超級(jí)伊昂大王 2026-03-24 | 18/900 |
|
|---|---|---|---|---|
|
[考研] 288資源與環(huán)境專碩求調(diào)劑,不限專業(yè),有學(xué)上就行 +11 | lllllos 2026-03-30 | 11/550 |
|
|
[考研] 285求調(diào)劑 +6 | AZMK 2026-03-29 | 9/450 |
|
|
[考研] 求調(diào)劑,一志愿 南京航空航天大學(xué) ,080500材料科學(xué)與工程學(xué)碩,總分289分 +8 | @taotao 2026-03-29 | 8/400 |
|
|
[考研] 279求調(diào)劑 +12 | j的立方 2026-03-29 | 12/600 |
|
|
[考研] 346求調(diào)劑 一志愿070303有機(jī)化學(xué) +7 | 蘿卜燉青菜 2026-03-28 | 8/400 |
|
|
[考研] 材料專碩 085600求調(diào)劑 +7 | BBQ233 2026-03-30 | 7/350 |
|
|
[考研] 329求調(diào)劑,一志愿西北工業(yè)大學(xué),材料工程(085601) +5 | 小小機(jī)靈蟲 2026-03-29 | 11/550 |
|
|
[考研] 071010 323 分求調(diào)劑 +3 | Baekzhy 2026-03-27 | 3/150 |
|
|
[考研] 283求調(diào)劑(080500) +14 | A child 2026-03-27 | 14/700 |
|
|
[考研] 282求調(diào)劑 +4 | wcq131415 2026-03-24 | 4/200 |
|
|
[考研] 一志愿中南大學(xué)化學(xué)0703總分337求調(diào)劑 +6 | niko- 2026-03-27 | 6/300 |
|
|
[考研] 294分080500材料科學(xué)與工程求調(diào)劑 +8 | 柳溪邊 2026-03-26 | 8/400 |
|
|
[考研] 2026年華南師范大學(xué)歡迎化學(xué),化工,生物,生醫(yī)工等專業(yè)優(yōu)秀學(xué)子加入! +3 | llss0711 2026-03-28 | 6/300 |
|
|
[考研] 266分,求材料冶金能源化工等調(diào)劑 +7 | 哇呼哼呼哼 2026-03-27 | 9/450 |
|
|
[考研] 265求調(diào)劑11408 +3 | 劉小鹿lu 2026-03-27 | 3/150 |
|
|
[考研] 266分求材料化工冶金礦業(yè)等專業(yè)的調(diào)劑 +4 | 哇呼哼呼哼 2026-03-26 | 4/200 |
|
|
[考研] 機(jī)械學(xué)碩310分,數(shù)一英一,一志愿211本科雙非找調(diào)劑信息 +3 | @357 2026-03-25 | 3/150 |
|
|
[考研] 281求調(diào)劑 +6 | Koxui 2026-03-24 | 7/350 |
|
|
[考研] 一志愿天津大學(xué)339材料與化工求調(diào)劑 +3 | 江往賣魚 2026-03-26 | 3/150 |
|