| 查看: 1083 | 回复: 5 | ||
| 当前只显示满足指定条件的回帖，点击这里查看本话题的所有回帖 | ||
[求助]
matlab lsqnonlin函数进行非线性动力学参数拟合，求高手指点.已有1人参与
|
||
|
function kinetics1 clear all clc a0=[10^11 81 0.3]; lb=[0 0 0];ub=[+inf +inf +inf]; c0=[0.0525564 100]; c16exp=... [0.0525564 100 0.0403533 73.3083 0.0267938 50.3759 0.015037 42.1053 0.0082474 24.812 0.00415801 13.9098 0.00278195 11.6541 0.00185881 10.1504 0.00183435 9.77444 0.0015879 9.77444]; c14exp=... [0.0525564 100 0.0493758 89.8496 0.0441658 81.203 0.0378303 72.9323 0.0357841 67.6692 0.0337337 61.6541 0.0276196 54.8872 0.0239832 47.7444 0.018766 40.6015 0.0160359 35.7143 0.0153455 34.5865]; c12exp=... [0.0526 100 0.0516 98.1203 0.0507 94.7368 0.0484 93.2331 0.0464 89.8496 0.0439 80.4511 0.0429 80.8271 0.0407 77.8195 0.0393 76.3158 0.0381 75.5639 0.0365 72.1805]; [a,resnorm,residual,exitflag,output,lambda,jacobian]=lsqnonlin(@objlm,a0,lb,ub,[],c0,c16exp,c14exp,c12exp); ci=nlparci(a,residual,jacobian); fprintf('\n\n使用函数lsqnonlin()估计得到的参数值为：\n'),a function f=objlm(a,c0,c16exp,c14exp,c12exp) %160度 tspan1=[0 0.1005 0.2143 0.3360 0.4365 0.5661 0.6667 0.7593 0.8677 0.9603]; [t,c16]=ode45(@kinetic16,tspan1,c0,[],a); f1=c16(:,1)-c16exp(:,1); f2=c16(:,2)-c16exp(:,2); %140度 tspan2=[0 0.1005 0.1984 0.2857 0.3571 0.4471 0.5529 0.6746 0.8042 0.9074 0.9683]; [t,c14]=ode45(@kinetic14,tspan2,c0,[],a); f3=c14(:,1)-c14exp(:,1); f4=c14(:,2)-c14exp(:,2); %120度 tspan3=[0 0.1005 0.2011 0.2937 0.4021 0.5000 0.5979 0.7090 0.8069 0.8862 0.9894]; [t,c12]=ode45(@kinetic12,tspan3,c0,[],a); f5=c12(:,1)-c12exp(:,1); f6=c12(:,2)-c12exp(:,2); f=[f1;f2;f3;f4;f5;f6]; %构造目标函数 function dc16dt=kinetic16(t,c16,a) %160度 T1=433.15;R=0.008314; dc16dt=... [-(a(1).*exp(-a(2)./(R.*T1)).*((c16(2)./100).^(a(3))).*c16(1).*0.943) -(7.*a(1).*exp(-a(2)./(R.*T1)).*((c16(2)./100).^(a(3))).*c16(1).*0.943)]; function dc14dt=kinetic14(t,c14,a) %140度 T2=413.15;R=0.008314; dc14dt=... [-(a(1).*exp(-a(2)./(R.*T2)).*((c14(2)./100).^(a(3))).*c14(1).*0.926) -(7.*a(1).*exp(-a(2)./(R.*T2)).*((c14(2)./100).^(a(3))).*c14(1).*0.926)]; function dc12dt=kinetic12(t,c12,a) %120度 T3=393.15;R=0.008314; dc12dt=... [-(a(1).*exp(-a(2)./(R.*T3)).*((c12(2)./100).^(a(3))).*c12(1).*0.907) -(7.*a(1).*exp(-a(2)./(R.*T3)).*((c12(2)./100).^(a(3))).*c12(1).*0.907)]; |
回复此楼» 猜你喜欢
博士读完未来一定会好吗
已经有14人回复
-
心脉受损
已经有4人回复
-
Springer期刊投稿求助
已经有4人回复
-
读博
已经有3人回复
-
小论文投稿
已经有3人回复
-
Bioresource Technology期刊，第一次返修的时候被退回好几次了
已经有9人回复
-
到新单位后，换了新的研究方向，没有团队，持续积累2区以上论文，能申请到面上吗
已经有8人回复
-
申请2026年博士
已经有6人回复
-
请问哪里可以有青B申请的本子可以借鉴一下。
已经有5人回复
6楼2019-12-22 18:21:34