| 查看: 228 | 回复: 0 | |||
| 当前主题已经存档。 | |||
[交流]
【金币求助】matlab求解CO变换反应动力学方程
|
|||
|
我要求的是CO变换反应动力学方程，用的是幂函数形式的模型 采用龙格库塔法求解微分方程，用单纯性法优化参数 我编的程序如下： function DEF clear all;clc data yco=yco./(1.+yr); tspan= [0 15 30]; yco1=[0.005 0.098 0.033 0.12 0.085 0.013 0.068 0.094 0.005 0.046 0.017 0.05 0.037] beta0=[0.9 0.3 -0.6 0.1 1.00 1.0000] lb=-10.0; ub=10.0; [beta,resnorm,resid,exitflag,output,lambda] = lsqnonlin(@OptObjFunc,beta0,lb,ub,[],tspan,yco,yco1) %..................... %[beta,resnorm]=lsqnonlin(@wfun,beta0,yco1) %[beta,resnorm] =lsqnonlin(@wfun,beta0,[],[],[],yco,yco1) function f=OptObjFunc(beta,tspan,yco,yco1) [t yco3]=ode23s(@feng,tspan,yco,[],beta); data yhho=yr./(1.+yr); yco4=yco3(3)./(1-(yhho-yco+yco3(3)));%yco3---用龙格库塔法算出的CO湿基含量 yco4--是计算出的CO干基含量 f=(yco1-yco4).^(2);%yco1--是实验测得的CO干基含量 yco4--是计算出的CO干基含量 %.................. function dydx=feng(t,y,beta) data ycoo=ycoo./(1.+yr); ycoo=ycoo+yco; yhho=yr./(1.+yr); yhho=yhho-yco; yhh=yhh./(1.+yr); yhh=yhh+yco; yv=yv.*(1.+yr); kp=exp(5025.163./yt-0.0936.*log(yt)+1.4555/1000*yt-2.4887/10000000*yt-5.2894); %b=(ycoo-y').*(yhh-y')./kp./(yhho+y'); dydx=(22.4*beta(5)*100./yv).*exp(-beta(6)*10000/8.314./yt).*(y.^beta(1))'.*((yhho'+y).^beta(2))'.*((ycoo'-y).^beta(3))'.*((yhh'-y).^beta(4))'.*(1.0-((ycoo-y').*(yhh-y')./kp./(yhho+y')./y')); dydx=dydx'; 最后算得的参数是复数，盼望大家指点 [ Last edited by wlf815795 on 2009-5-5 at 08:59 ] |
回复此楼» 猜你喜欢
filecode
已经有3人回复
-
我的国基提前知道中了，可是同事的操作让我实在接受不了，怎么会有这样的人
已经有10人回复
-
静等基金结果
已经有13人回复
-
关于Filecode分析方法
已经有3人回复
-
国自然结果
已经有8人回复
-
FileCode能看出啥？
已经有8人回复
-
应该是93bebmhtak前后十一个字符比较关键
已经有22人回复
-
奇怪，两个人的filecode固定段从头到尾一模一样
已经有4人回复
-
8月时间戳变的，举个手。玩一下，释放压力
已经有13人回复
-
好奇怪的filecode
已经有5人回复
100