24小时热门版块排行榜    

查看: 237  |  回复: 0
当前主题已经存档。

wlf815795

银虫 (初入文坛)

[交流] 【金币求助】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 ]
回复此楼
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
相关版块跳转 我要订阅楼主 wlf815795 的主题更新
普通表情 高级回复 (可上传附件)
最具人气热帖推荐 [查看全部] 作者 回/看 最后发表
[公派出国] 售SCI文章,我:8O.5.5.1O.54,科目全,可十急 +3 EZ67DG6CvcKL 2026-08-31 3/150 2026-08-31 22:22 by o3EMxHUAQDtS
[基金申请] 基金未中,这种答复是模板吗? +6 zhaosm1982 2026-08-27 7/350 2026-08-31 21:18 by qdxxmc
[基金申请] 面上意见出来了 +12 黄鸟于飞Chao 2026-08-29 23/1150 2026-08-31 18:57 by 黄鸟于飞Chao
[基金申请] 能否申诉? +7 echo8914667 2026-08-30 8/400 2026-08-31 17:00 by yihongxu
[基金申请] 基金不中,共勉 +12 eulota 2026-08-26 12/600 2026-08-31 08:42 by ZJTJZ
[基金申请] 为什么到现在没收到通知? +5 tannykie 2026-08-29 5/250 2026-08-30 21:05 by purplejack
[基金申请] 我就是申请一个面上项目而已,这评审意见是按照杰青的条件评的吧? +6 gouxfjh 2026-08-28 11/550 2026-08-30 07:57 by gouxfjh
[基金申请] 国自然面上复盘~欢迎讨论 (金币+15) +15 晴天加油 2026-08-26 16/800 2026-08-29 18:28 by symmetry
[基金申请] 2026年叶企孙基金 +4 bud_bud 2026-08-27 7/350 2026-08-29 07:23 by foolishmani
[基金申请] 为什么资助数各大高校都创新高,自己申请怎么就这么难 +12 Kittylucky 2026-08-27 13/650 2026-08-29 00:04 by superceng
[基金申请] 系统查不到 +11 董八千 2026-08-26 11/550 2026-08-28 18:06 by Leogzhya
[基金申请] 基金系统什么内容也没有 30+4 winsaint 2026-08-27 9/450 2026-08-28 11:06 by maolC
[基金申请] 看板上这么多中的,有点像50人群里49个人都是骗子的那种感觉…… +5 a089 2026-08-26 6/300 2026-08-27 14:05 by jonewore
[基金申请] 为什么 国际(地区)合作与交流项目 没有放榜? 10+3 majunge000 2026-08-26 11/550 2026-08-27 08:42 by 北京莱茵编辑
[基金申请] 我不理解! +15 Edward_pc 2026-08-26 23/1150 2026-08-26 20:34 by zzuzxg
[基金申请] 2026年8月25日国自然放榜前突然收到列入评审专家邮件,有关系吗? +25 木水思豆 2026-08-25 28/1400 2026-08-26 14:53 by draco1987
[基金申请] 出来了 +9 trojank 2026-08-26 9/450 2026-08-26 14:25 by 宝贝虫子
[基金申请] 国合里面能看到了 +7 一怀馨秋 2026-08-26 7/350 2026-08-26 11:23 by zhaosm1982
[基金申请] 国际合作可查了,中了面上 (EPI+1)(金币+50) +18 Ldrop2023 2026-08-26 18/900 2026-08-26 11:15 by cmrandy
[基金申请] 国合可查了 +3 paperzjh 2026-08-26 3/150 2026-08-26 10:41 by LemmonTr
信息提示
请填处理意见