| 查看: 1486 | 回复: 10 | |||
njhx505木虫 (正式写手)
|
[交流]
【求助】30金币求一程序 已有4人参与
|
 请各位高手帮个忙，CA,CB是不同时间测得浓度值，求k1,k2  数据给出，希望高手可以帮忙，谢谢 [ Last edited by njhx505 on 2010-6-4 at 19:42 ] |
回复此楼» 猜你喜欢
估计是周四
已经有3人回复
-
基金啊基金
已经有3人回复
-
科研孤儿太难了
已经有15人回复
-
让我中一个面上吧！
已经有10人回复
-
言之凿凿，明天官网会通知放榜？
已经有17人回复
-
时间戳又变了
已经有17人回复
-
今天系统多次维护，明天很可能放榜！
已经有11人回复
-
时间戳变了，能看出什么问题？
已经有23人回复
-
只有每年这种时候来逛逛小木虫
已经有22人回复
-
放榜前的不淡定
已经有12人回复
» 本主题相关价值贴推荐，对您同样有帮助:
- 30金币求助：最长路径问题，要求遍历每个节点，总路径最长
已经有7人回复
- 求高人指点用matlab求解非线性方程组，解决了追加100金币；
已经有11人回复
- 【求助】50金币求大家帮忙改一段程序
已经有13人回复
- 【求助】VS2005数值计算程序调试【急求！！！！加金币】
已经有9人回复
- 【求助】100金币求助高人一个程序解方程
已经有15人回复
lijinfeng042
木虫 (小有名气)
Matlab
- 仿真EPI: 2
- 应助: 1 (幼儿园)
- 金币: 2156.1
- 散金: 115
- 帖子: 291
- 在线: 31.5小时
- 虫号: 1019062
- 注册: 2010-05-15
- 性别: GG
- 专业: 通信理论与系统
2楼2010-06-04 17:03:55
change0618
铁杆木虫 (著名写手)
方丈大师
- 应助: 44 (小学生)
- 金币: 17724.5
- 红花: 17
- 帖子: 2413
- 在线: 546.7小时
- 虫号: 496517
- 注册: 2008-01-19
- 专业: 化学反应工程
3楼2010-06-04 17:31:59
change0618
铁杆木虫 (著名写手)
方丈大师
- 应助: 44 (小学生)
- 金币: 17724.5
- 红花: 17
- 帖子: 2413
- 在线: 546.7小时
- 虫号: 496517
- 注册: 2008-01-19
- 专业: 化学反应工程
数据有的恶劣 4楼2010-06-04 20:28:04
njhx505
木虫 (正式写手)
- 应助: 0 (幼儿园)
- 金币: 2163.9
- 散金: 6
- 帖子: 387
- 在线: 16.4小时
- 虫号: 737675
- 注册: 2009-04-01
- 性别: MM
- 专业: 分离过程
5楼2010-06-04 20:39:35
★ ★
zzuwangshilei(金币+2):辛苦了，多谢应助 2010-06-05 10:52:39
njhx505(金币+3):谢谢 2010-06-06 21:57:51
zzuwangshilei(金币+2):辛苦了，多谢应助 2010-06-05 10:52:39
njhx505(金币+3):谢谢 2010-06-06 21:57:51
|
(2)式*2+(1)式得到 dCA(t)/dt+2*dCB(t)/dt=0; 在任意小区间[t0,t]积分得到：CA(t)+2*CB(t)=const ....(*)式 计算一下你表格中的数据，不满足是常数。至于怎么处理，或者重新采集数据或者像你提到的满足某种最优。 另外： 将(*)式代入到上面其中一个方程比如(1)得到： dCA/dt=-k1*CA^2-k2*CA+const; 这个方程建议还是数值求解比较方便，比如向前欧拉，向后欧拉方法等等，要想得到高精度的就用Runge-Kutta。 |
6楼2010-06-05 03:24:37
njhx505
木虫 (正式写手)
- 应助: 0 (幼儿园)
- 金币: 2163.9
- 散金: 6
- 帖子: 387
- 在线: 16.4小时
- 虫号: 737675
- 注册: 2009-04-01
- 性别: MM
- 专业: 分离过程
7楼2010-06-05 12:50:54
change0618
铁杆木虫 (著名写手)
方丈大师
- 应助: 44 (小学生)
- 金币: 17724.5
- 红花: 17
- 帖子: 2413
- 在线: 546.7小时
- 虫号: 496517
- 注册: 2008-01-19
- 专业: 化学反应工程
8楼2010-06-05 14:58:06
njhx505
木虫 (正式写手)
- 应助: 0 (幼儿园)
- 金币: 2163.9
- 散金: 6
- 帖子: 387
- 在线: 16.4小时
- 虫号: 737675
- 注册: 2009-04-01
- 性别: MM
- 专业: 分离过程
9楼2010-06-05 16:38:43
change0618
铁杆木虫 (著名写手)
方丈大师
- 应助: 44 (小学生)
- 金币: 17724.5
- 红花: 17
- 帖子: 2413
- 在线: 546.7小时
- 虫号: 496517
- 注册: 2008-01-19
- 专业: 化学反应工程
★ ★
kuhailangyu(金币+2):欢迎积极参与 2010-06-05 22:36:53
njhx505(金币+24):谢谢你的程序，真是救了我的急 2010-06-06 21:59:25
kuhailangyu(金币+2):欢迎积极参与 2010-06-05 22:36:53
njhx505(金币+24):谢谢你的程序，真是救了我的急 2010-06-06 21:59:25
|
给你一个程序，但是数据不是你的，自己修改去吧 function Kinetics % 动力学ODE方程模型的参数估计 clear all clc ExpData = ... [ 0 0.1883 0.0100 0.2047 0.0200 0.2181 0.0300 0.2291 0.0400 0.2382 0.0500 0.2459 0.0600 0.2523 0.0700 0.2576 0.0800 0.2622 0.0900 0.2660 0.1000 0.2692 0.1100 0.2719 0.1200 0.2742 0.1300 0.2761 0.1400 0.2777 0.1500 0.2790 0.1600 0.2801 0.1700 0.2811 0.1800 0.2819 0.1900 0.2825 0.2000 0.2830 ]; t = ExpData(:,1); % ExpData第一列为时间 CB = ExpData(:,2); % ExpData第二列为组分Ｂ的浓度 CB0 = CB(1); % t=0时，组分B初始浓度 CA0 = 0.8; % t=0时，组分A初始浓度 CA = CA0-2*(ExpData(:,2)-CB0); % 由微分方程式导出来的CA与CB关系式 k0 = [20 50]; % 估值参数的猜想值 lb = [0 0]; % 设定的估值参数上限 ub = [+inf +inf]; % 设定的估值参数下限 [k,resnorm,residual,exitflag,output,lambda,jacobian] = ... lsqnonlin(@ObjFunc,k0,lb,ub,[],t,[CA0,CB0],[CA,CB]); ci = nlparci(k,residual,jacobian); % 计算置信区间 fprintf('\n\n估计参数值为:\n') fprintf('\tk1 = %.6f\t置信区间：[%.6f %.6f]\n',k(1),ci(1,: )) fprintf('\tk2 = %.6f\t置信区间：[%.6f %.6f]\n',k(2),ci(2,: )) tt = linspace(t(1),t(end),101); [tt C] = ode45(@KineticEqs,tt,[CA0,CB0],[],k); % 由估值k计算时间序列tt下的A，B浓度 figure(1) plot(t,CA,'o',tt,C(:,1),'r-') xlabel('t'); ylabel('C_A') figure(2) plot(t,CB,'o',tt,C(:,2),'r-') xlabel('t'); ylabel('C_B') % ------------------------------------------------------------------ function f = ObjFunc(k,tspan,x0,yexp) [t,y] = ode45(@KineticEqs,tspan,x0,[],k); f = y(: ) - yexp(: ); % ------------------------------------------------------------------ function dxdt = KineticEqs(t,x,k) dxdt =[-k(1)*x(1)^2+2*k(2)*x(2); 0.5*k(1)*x(1)^2-k(2)*x(2)]; [ Last edited by change0618 on 2010-6-5 at 18:23 ] |
10楼2010-06-05 18:15:24
40