| 查看: 4025 | 回复: 15 | |||
[求助]
matlab拟合方程参数时初值的选择
|
|
A+B=C A+W=R R+B=C+W 先考虑的三个反应都是可逆反应，有以下几个式子：编写程序如下 %30℃ clear all clc global keq1 keq2 t0 c0 keq1 = 3.94; % 平衡常数k1 keq3 = 1.08; % 平衡常数k3 t0 = [5,10, 15, 30, 45, 60,90,120,160,200,240,300,360]; c0 = [5.939601,5.203973,0,0.022484,0.078956] % A B C R W 组成的初始浓度 ci=[3.742508811 4.213566482 0.990406938 0.252498405 0.090415512; 3.045716838 3.387013322 1.816960097 0.215546741 0.081072022; 2.504754523 2.901991595 2.301981824 0.193427102 0.099462604; 1.65098702 2.036686653 3.167286766 0.170236582 0.095144044; 1.271633833 1.655460403 3.548513016 0.124730781 0.10334349; 1.129887499 1.630815394 3.573158025 0.129798695 0.097764543; 1.01239208 1.460748979 3.74322444 0.144649391 0.073614958; 1.007417603 1.429180292 3.774793127 0.123477124 0.094249307; 0.990005723 1.499216939 3.704756481 0.135524067 0.076072022; 0.970073726 1.493966669 3.71000675 0.130791851 0.076387812; 0.96903024 1.429705429 3.77426799 0.119510017 0.101202216; 0.980873385 1.457988124 3.745985296 0.124189274 0.097022161; 0.98477262 1.453402965 3.750570455 0.1200936 0.0762608]; % cA cC cR 对应时间数值 k0=[1 1 1]; % k1+, k2, k3+ 初始向值 lb = [0 0 0]; ub = [1000 1000 500]; % 上下限 % 使用函数lsqnonlin()进行参数估计 [k,resnorm,residual,exitflag] = lsqnonlin(@ObjFunc,k0,lb,ub,[],ci); k1plus=k(1); k1minus= k1plus/keq1; k2=k(2); k3plus=k(3); k3minus= k3plus/keq3; % ------------------------------------------------------------------ function f = ObjFunc(k,ci) % 目标函数 global t0 c0 [t,c_cal] = ode45(@Euqations,t0,c0,[],k); f= sum((c_cal-ci)^2); % ------------------------------------------------------------------ function dcdt = Euqations(t,c,k) global keq1 keq2 k1plus=k(1); k1minus= k1plus/keq1; k2=k(2); k3plus=k(3); k3minus= k3plus/keq3; cA=c(:,1); cB=c(:,2); cC=c(:,3); cR=c(:,4); cW=c(:,5); dcAdt = k1plus*cA*cB-k1minus*cC+k2*cA*cW; dcCdt = -k1plus*cA*cB+k1minus*cC+k3plus*cB*cR-k3minus*cC*cW; dcRdt = k2*cA*cW-k3plus*cB*cR+k3minus*cC*cW; dcdt = [dcAdt; dcCdt;dcRdt]; 但是k的初值一直影响结果，程序要怎么改一下初值对结果影响不大（这个程序有问题，） |
回复此楼» 猜你喜欢
如何从铁电相到顺电相。
已经有1人回复
-
各位大佬，求一份最新Wien2k版本作为学习用途
已经有0人回复
-
物理学I论文润色/翻译怎么收费?
已经有166人回复
-
南科大活性流体和软物质课题组诚招2027级博士、硕士研究生和博后
已经有42人回复
-
瑞典林雪平大学博士后招聘｜PEC Water Splitting 方向
已经有129人回复
-
PRB投稿，Acknowledgment sent to author状态十天了
已经有3人回复
-
PVA溶解
已经有2人回复
-
基金申请
已经有44人回复
-
CSC与新西兰维多利亚大学PhD奖学金项目
已经有0人回复
-
新西兰Robinson研究所 招聘CSC公派访问人员
已经有0人回复
-
帮我的英语口语老师找学生
已经有0人回复
» 本主题相关价值贴推荐，对您同样有帮助:
- 在matlab中拟合曲面的表达式
已经有29人回复
- matlab 求指点 动力学方程拟合过程中导数的获取
已经有10人回复
- matlab拟合拟均相动力学参数
已经有6人回复
- matlab-常微分方程参数估计
已经有12人回复
- 求教origin拟合简单的渗流方程（又称逾渗方程），做复合材料的应该会遇到吧
已经有20人回复
- 非线性拟合估计方程参数（matlab，spss，origin都可以）
已经有4人回复
- matlab非线性参数拟合问题
已经有7人回复
- MATLAB求解一个比较复杂的方程
已经有14人回复
- 【求助】用matlab最优化方法进行参数拟合
已经有17人回复
- 【求助】使用MATLAB拟合出按照公式的系数
已经有19人回复
- 如何用MATLAB拟合指数函数
已经有3人回复
- 【求助】最小二乘法编写指数拟合方程
已经有19人回复
- 【求助】使用Matlab拟合反应动力学方程问题
已经有7人回复
dingd
铁杆木虫 (职业作家)
- 计算强帖: 4
- 应助: 1641 (讲师)
- 金币: 15037.3
- 散金: 101
- 红花: 234
- 帖子: 3410
- 在线: 1223.7小时
- 虫号: 291104
- 注册: 2006-10-28
2楼2012-03-28 10:37:30
|
我写在附件里面了 |
» 本帖附件资源列表
-
欢迎监督和反馈：小木虫仅提供交流平台，不对该内容负责。
本内容由用户自主发布，如果其内容涉及到知识产权问题，其责任在于用户本人，如对版权有异议，请联系邮箱：xiaomuchong@tal.com - 附件 1 : 太残忍doc.doc
2012-03-28 10:52:17, 30.5 K
3楼2012-03-28 10:52:39
dingd
铁杆木虫 (职业作家)
- 计算强帖: 4
- 应助: 1641 (讲师)
- 金币: 15037.3
- 散金: 101
- 红花: 234
- 帖子: 3410
- 在线: 1223.7小时
- 虫号: 291104
- 注册: 2006-10-28
4楼2012-03-28 11:07:18
|
不好意思，我的微分方程组就是附件里面标（1）（2）（3）的三个式子，对应的数据为 时间：t0 = [5,10, 15, 30, 45, 60,90,120,160,200,240,300,360]; 不同时刻对应的数据 c0 = [5.939601,5.203973,0,0.022484,0.078956] % A B C R W 组成的初始浓度 ci=[3.742508811 4.213566482 0.990406938 0.252498405 0.090415512; 3.045716838 3.387013322 1.816960097 0.215546741 0.081072022; 2.504754523 2.901991595 2.301981824 0.193427102 0.099462604; 1.65098702 2.036686653 3.167286766 0.170236582 0.095144044; 1.271633833 1.655460403 3.548513016 0.124730781 0.10334349; 1.129887499 1.630815394 3.573158025 0.129798695 0.097764543; 1.01239208 1.460748979 3.74322444 0.144649391 0.073614958; 1.007417603 1.429180292 3.774793127 0.123477124 0.094249307; 0.990005723 1.499216939 3.704756481 0.135524067 0.076072022; 0.970073726 1.493966669 3.71000675 0.130791851 0.076387812; 0.96903024 1.429705429 3.77426799 0.119510017 0.101202216; 0.980873385 1.457988124 3.745985296 0.124189274 0.097022161; 0.98477262 1.453402965 3.750570455 0.1200936 0.0762608]; 要求的是ki+,k2,k3+, 另外方程中的ki-=ki+/3.14,k3-=k3+/1.08 |
» 本帖附件资源列表
-
欢迎监督和反馈：小木虫仅提供交流平台，不对该内容负责。
本内容由用户自主发布，如果其内容涉及到知识产权问题，其责任在于用户本人，如对版权有异议，请联系邮箱：xiaomuchong@tal.com - 附件 1 : 太残忍doc (1)1.doc
2012-03-28 12:18:21, 26 K
5楼2012-03-28 12:18:31
dingd
铁杆木虫 (职业作家)
- 计算强帖: 4
- 应助: 1641 (讲师)
- 金币: 15037.3
- 散金: 101
- 红花: 234
- 帖子: 3410
- 在线: 1223.7小时
- 虫号: 291104
- 注册: 2006-10-28
6楼2012-03-28 14:58:43
|
有的，我写在附件里了，真是太谢谢了 |
» 本帖附件资源列表
-
欢迎监督和反馈：小木虫仅提供交流平台，不对该内容负责。
本内容由用户自主发布，如果其内容涉及到知识产权问题，其责任在于用户本人，如对版权有异议，请联系邮箱：xiaomuchong@tal.com - 附件 1 : 太残忍doc (1)1.doc
2012-03-28 15:05:19, 32.5 K
7楼2012-03-28 15:05:21
dingd
铁杆木虫 (职业作家)
- 计算强帖: 4
- 应助: 1641 (讲师)
- 金币: 15037.3
- 散金: 101
- 红花: 234
- 帖子: 3410
- 在线: 1223.7小时
- 虫号: 291104
- 注册: 2006-10-28
8楼2012-03-28 16:30:05
|
dCB/wdt"应该是“dCB/dt，没错 对不起我太心急了，公式已经检查过，我保证这次没错了,真的 |
» 本帖附件资源列表
-
欢迎监督和反馈：小木虫仅提供交流平台，不对该内容负责。
本内容由用户自主发布，如果其内容涉及到知识产权问题，其责任在于用户本人，如对版权有异议，请联系邮箱：xiaomuchong@tal.com - 附件 1 : 太残忍doc(1)2.doc
2012-03-28 16:38:07, 34.5 K
9楼2012-03-28 16:38:11
dingd
铁杆木虫 (职业作家)
- 计算强帖: 4
- 应助: 1641 (讲师)
- 金币: 15037.3
- 散金: 101
- 红花: 234
- 帖子: 3410
- 在线: 1223.7小时
- 虫号: 291104
- 注册: 2006-10-28
10楼2012-03-28 17:10:04
25