24小时热门版块排行榜    

查看: 538  |  回复: 2

laojin

木虫 (小有名气)

[求助] 用matlab求动力学参数代码有错误 已有1人参与

初学者,所以编了一个最简单的方程,但是有一点问题。
反应方程:-dc/dt=k*c^a

function KA
format long
clear all
clc
tspan = [0  2  4  6 8 10];
x0 = 0.846812;
k0 = [0.02  0];   
lb = [0  -1];
ub = [1  1];

data=...
    [

0        0.846812
2        0.802778
4        0.757897
6        0.702854
8        0.653739
10        0.609704

];
yexp = data(:,2);

[k,resnorm,residual,exitflag,output,lambda,jacobian] = ...
    lsqnonlin(@ObjFunc,k0,lb,ub,[],tspan,x0,yexp);      
ci = nlparci(k,residual,jacobian);
fprintf('\n\n使用函数lsqnonlin()估计得到的参数值为:\n')
fprintf('\tk1 = %.9f ± %.9f\n',k(1),ci(1,2)-k(1))
fprintf('\tk2 = %.9f ± %.9f\n',k(2),ci(2,2)-k(2))
fprintf('  The sum of the squares is: %.9e\n\n',resnorm)

function f = ObjFunc(k,tspan,x0,yexp)           % 目标函数
[t Xsim] = ode45(@KineticsEqs,tspan,x0,[],k);
Xsim1=Xsim(:,1);
ysim(:,1) = Xsim1(1:end);

size(ysim(:,1));
size(yexp(:,1));

f = (ysim(:,1)-yexp(:,1));


function dCdt = KineticsEqs(t,C,k)              % ODE模型方程
dCAdt = -k(1)*(C^k(2));
dCdt = dCAdt;

运行之后,得到结果
使用函数lsqnonlin()估计得到的参数值为:
        k1 = 0.021862398 ± 0.005117281
        k2 = -0.282935990 ± 0.797237745
  The sum of the squares is: 4.866526201e-05
K1还好,但是K2显然有问题,我是按照论坛里的Code改的,不知道出了什么问题,请各位帮忙看看。
我自己编了一个零级反应的数据,但是也是K2出现问题。
不知道是不是数据太少的缘故。
烦请大神修正,不胜感激。
回复此楼

» 猜你喜欢

» 本主题相关价值贴推荐,对您同样有帮助:

已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

dingd

铁杆木虫 (职业作家)

【答案】应助回帖

★ ★ ★ ★ ★ ★ ★ ★ ★ ★
感谢参与,应助指数 +1
laojin: 金币+10, ★★★很有帮助 2015-01-23 11:44:11
应该没错吧,1stOpt求解:
CODE:
Variable t,c;
ODEFunction c'=-k*c^a;
Data;
0        0.846812
2        0.802778
4        0.757897
6        0.702854
8        0.653739
10        0.609704

均方差(RMSE):0.00311978140170499
残差平方和(SSE):4.86651799721219E-5
相关系数(R): 0.999116092954796
相关系数之平方(R^2): 0.998232967201257
确定系数(DC): 0.997978305044057
F统计(F-Statistic): 1483.90339065772

参数                  最佳估算
--------------------        -------------
k        0.0218659845981955
a        -0.282434977032699

====== 结果输出 ======

文件: 数据文件-1
No        t        目标 c        计算 c
1        2        0.802778        0.800616656747102
2        4        0.757897        0.753655380090554
3        6        0.702854        0.705852003645687
4        8        0.653739        0.657115339412452
5        10        0.609704        0.607334455281825
2楼2015-01-23 10:34:55
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

laojin

木虫 (小有名气)

引用回帖:
2楼: Originally posted by dingd at 2015-01-23 10:34:55
应该没错吧,1stOpt求解:

Variable t,c;
ODEFunction c'=-k*c^a;
Data;
0        0.846812
2        0.802778
4        0.757897
6        0.702854
8        0.653739
10        0.609704


均方差 ...

只是觉得K2的偏差太大了,不知道是不是数据量太少的缘故。
3楼2015-01-23 11:44:41
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
相关版块跳转 我要订阅楼主 laojin 的主题更新
最具人气热帖推荐 [查看全部] 作者 回/看 最后发表
[考研] 材料与化工考研调剂 +9 孅華 2026-03-22 9/450 2026-03-25 13:09 by cmz0325
[考研] 求调剂323材料与化工 +4 1124361 2026-03-24 4/200 2026-03-25 11:19 by shulmg
[考研] 300分,材料,求调剂,英一数二 +5 超赞的 2026-03-24 5/250 2026-03-24 21:07 by 星空星月
[考研] 0854AI CV方向招收调剂 +3 章小鱼567 2026-03-23 3/150 2026-03-24 20:25 by 汪!?!
[考研] 085404电子信息284分求调剂 +4 13659058978 2026-03-24 4/200 2026-03-24 12:15 by syl20081243
[考研] 335分 | 材料与化工专硕 | GPA 4.07 | 有科研经历 +4 cccchenso 2026-03-23 4/200 2026-03-23 23:00 by 徐ckkk
[考研] 316求调剂 +7 梁茜雯 2026-03-19 7/350 2026-03-23 16:21 by lingjue
[考研] 263求调剂 +6 yqdszhdap- 2026-03-22 9/450 2026-03-23 12:57 by yqdszhdap-
[考研] 280分求调剂 一志愿085802 +4 PUMPT 2026-03-22 7/350 2026-03-22 22:13 by 星空星月
[考研] 材料学硕301分求调剂 +7 Liyouyumairs 2026-03-21 7/350 2026-03-21 22:31 by peike
[考研] 一志愿东华大学控制学硕320求调剂 +3 Grand777 2026-03-21 3/150 2026-03-21 19:23 by 简之-
[考研] 材料与化工(0856)304求B区调剂 +3 邱gl 2026-03-20 7/350 2026-03-21 19:05 by 15709483992
[考研] 求调剂 +4 要好好无聊 2026-03-21 4/200 2026-03-21 18:57 by 学员8dgXkO
[考研] 0703化学调剂 +4 妮妮ninicgb 2026-03-21 4/200 2026-03-21 18:39 by 学员8dgXkO
[考研] 一志愿重庆大学085700资源与环境总分308求调剂 +7 墨墨漠 2026-03-20 7/350 2026-03-21 16:36 by barlinike
[考研] 279求调剂 +5 红衣隐官 2026-03-21 5/250 2026-03-21 14:59 by lature00
[考研] 考研调剂求学校推荐 +3 伯乐29 2026-03-18 5/250 2026-03-20 22:59 by JourneyLucky
[考研] 求调剂一志愿南京航空航天大学289分 +3 @taotao 2026-03-19 3/150 2026-03-20 21:34 by JourneyLucky
[考研] A区线材料学调剂 +5 周周无极 2026-03-20 5/250 2026-03-20 21:33 by laoshidan
[考研] 一志愿西安交通大学 学硕 354求调剂211或者双一流 +3 我想要读研究生 2026-03-20 3/150 2026-03-20 20:13 by JourneyLucky
信息提示
请填处理意见