| 查看: 685 | 回复: 0 | ||
[求助]
求1stOpt程序代跑，求解微分方程组系数，可有偿，捉急捉急
|
|
需要拟合k1-k8的参数值 以下是我编写的matlabm文件: function QQQ_4_13_NJ_cs_2 format long clear all clc tspan = [0 0.25 1 3 9 22 35 91 135 ]; x0 = [19750 0 0 0 5250]; k0 = [0.147;0.19;0.03 ;0.02;0.01;0.000001184;3.7;13]; lb = [0;0;0;0;0;0;0;0]; ub = []; data=... [ 0.25 21300 100.225875 46.63575 16.168 3536.970375 1 22050 40.657785 14.123025 8.6688 2874.46475 3 21750 101.43456 48.06252 19.9383 3052.32539 9 19500 279.2715 190.106 23.3415 4908.7372 22 20200 390.3344 403.5 38.5308 3759.5435 35 18550 245.784 262.899 11.2602 5756.7424 91 17950 296.3025 363.12 20.31 6015.16875 135 17500 279.4176 467.418 11.8272 6891.531067 ]; yexp = data(:,2:6); [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('\tk(1) = %.9f ± %.9f\n',k(1),ci(1,2)-k(1)) fprintf('\tk(2) = %.9f ± %.9f\n',k(2),ci(2,2)-k(2)) fprintf('\tk(3) = %.9f ± %.9f\n',k(3),ci(3,2)-k(3)) fprintf('\tk(4) = %.9f ± %.9f\n',k(4),ci(4,2)-k(4)) fprintf('\tk(5) = %.9f ± %.9f\n',k(5),ci(5,2)-k(5)) fprintf('\tk(6) = %.9f ± %.9f\n',k(6),ci(6,2)-k(6)) fprintf('\tk(7) = %.9f ± %.9f\n',k(7),ci(7,2)-k(7)) fprintf('\tk(8) = %.9f ± %.9f\n',k(8),ci(8,2)-k(8)) 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); Xsim2=Xsim(:,2); Xsim3=Xsim(:,3); Xsim4=Xsim(:,4); Xsim4=Xsim(:,5); ysim(:,1) = Xsim1(2:end); ysim(:,2) = Xsim2(2:end); ysim(:,3) = Xsim3(2:end); ysim(:,4) = Xsim4(2:end); ysim(:,5) = Xsim4(2:end); size(ysim(:,1)); size(ysim(:,2)); size(ysim(:,3)); size(ysim(:,4)); size(ysim(:,5)); size(yexp(:,1)); size(yexp(:,2)); size(yexp(:,3)); size(yexp(:,4)); size(yexp(:,5)); f1=ysim(:,1)-yexp(:,1); f2=ysim(:,2)-yexp(:,2); f3=ysim(:,3)-yexp(:,3); f4=ysim(:,4)-yexp(:,4); f5=ysim(:,5)-yexp(:,5); f = [f1;f2;f3;f4;f5]; function dCdt = KineticsEqs(t,C,k) % ODE模型方程 dCAdt =-(k(3)+k(5)+k(7))*C(1)+k(1)*C(2)+k(2)*C(3)+k(2)*C(4)+k(8)*C(5); dCBdt = k(3)*C(1)-(k(1)+k(6))*C(2)+k(4)*C(3)+k(4)*C(4); dCCdt = k(5)*C(1)+k(6)*C(2)-(k(2)+k(4))*C(3); dCDdt = k(5)*C(1)+k(6)*C(2)-(k(2)+k(4))*C(4); dCEdt = k(7)*C(1)-k(8)*C(5); dCdt = [dCAdt;dCBdt;dCCdt;dCDdt;dCEdt]; |
回复此楼» 猜你喜欢
中国科学院东莞材料科学与技术研究所-2026年博士招生-吴昊研究员-磁学与自旋电子学
已经有0人回复
-
《电磁学》教材推荐
已经有1人回复
-
物理学I论文润色/翻译怎么收费?
已经有62人回复
-
【急招】合肥工大核聚变材料计算方向2026级工程博士生
已经有4人回复
-
大豆异黄酮分离
已经有0人回复
-
湖南大学材料学院急招2026年博士生，临时增加一名博士联培指标
已经有10人回复
-
天津理工大学晶体材料全国重点实验室刘红军教授课题组招收博士生1-2名
已经有1人回复
-
中国科学院物理研究所谌志国研究员团队招收2027年博士研究生
已经有3人回复
-
2026年中德博士后交流项目 - 新型量子和磁性材料：材料制备表征和中子散射研究
已经有12人回复
10