24小时热门版块排行榜

返回列表

zhaoshazhu

新虫 (小有名气)

应助: 0 (幼儿园)
金币: 125.7
散金: 10
帖子: 51
在线: 38小时
虫号: 1812564
注册: 2012-05-12
专业: 环境化工

[求助] Matlab求助

谁帮我看一下下面的程序哪里有问题，我运行不出来了。
function DMMxin
clear all
clc
k0 = [1 1  1  1 1 1 1 1 1 1];
lb = [0 0 0 0 0 0 0 0 0 0];
ub = [inf  inf  inf  inf inf inf inf inf inf inf];

P0 =[0.015481 4.644213 0.340306 0  0  0;
   0.018316 4.579047 0.402637 0  0  0;
   0.022423 4.484655 0.492922 0  0  0;
   0.040656 4.065616 0.893728 0  0  0;
   0.015481 4.644213 0.340306 0  0  0;
   0.018316 4.579047 0.402637 0  0  0;
   0.022423 4.484655 0.492922 0  0  0;
   0.040656 4.065616 0.893728 0  0  0;
   0.015481 4.644213 0.340306 0  0  0;
   0.018316 4.579047 0.402637 0  0  0;
   0.022423 4.484655 0.492922 0  0  0
];  % 初始分压,MPa

Pi=[0.008122 4.641505 0.345489 0.002907 0.001773 0.000203;
0.008845 4.576146 0.409071 0.003090 0.002578 0.000270;
0.011420 4.480752 0.500745 0.003683 0.003051 0.000350;
0.019618 4.062645 0.907326 0.005871 0.003936 0.000604;
0.006237 4.639595 0.347298 0.002693 0.003667 0.000511;
0.007462 4.574544 0.410503 0.003129 0.003859 0.000502;
0.008862 4.478527 0.503080 0.003801 0.005014 0.000716;
0.014593 4.064102 0.909456 0.005734 0.005297 0.000819;
0.004422 4.640240 0.347865 0.002427 0.004377 0.000669;
0.005254 4.574451 0.411560 0.002858 0.005010 0.000867;
0.006444 4.480676 0.503265 0.003476 0.005285 0.000856;
];
% 经过Wc/F0后，各物质分压,MPa

% 使用函数lsqnonlin()进行参数估计
opt=optimset('Algorithm','levenberg-marquardt');
[k,resnorm,residual,exitflag,output,lambda,jacobian] = lsqnonlin(@ObjFunc,k0,[],[],opt,P0,Pi);
ci = nlparci(k,residual,jacobian);

fprintf('\n\n使用函数lsqnonlin()估计得到的参数值为:\n')
fprintf('\tk1 = %.4f ± %.4f\n',k(1),ci(1,2)-k(1))
fprintf('\tk2 = %.4f ± %.4f\n',k(2),ci(2,2)-k(2))
fprintf('\tk3 = %.4f ± %.4f\n',k(3),ci(3,2)-k(3))
fprintf('\tk4 = %.4f ± %.4f\n',k(4),ci(4,2)-k(4))
fprintf('\tk5 = %.4f ± %.4f\n',k(5),ci(5,2)-k(5))
fprintf('\tk6 = %.4f ± %.4f\n',k(6),ci(6,2)-k(6))
fprintf('\tk7 = %.4f ± %.4f\n',k(7),ci(7,2)-k(7))
fprintf('\tk8 = %.4f ± %.4f\n',k(8),ci(8,2)-k(8))
fprintf('\tk9 = %.4f ± %.4f\n',k(9),ci(9,2)-k(9))
fprintf('\tk10 = %.4f ± %.4f\n',k(10),ci(10,2)-k(10))
fprintf('\t残差平方和 = %.4f\n',resnorm)
fprintf('\texitflag = %.4f\n',exitflag)
fprintf('\tresidual = %.4f\n',residual)

% ------------------------------------------------------------------
function f = ObjFunc(k,P0,Pi)          % 目标函数
[m,n] = size(P0);
Pcal = zeros(m,n);
tspan = [0  264;
      0  264;
      0  264;
      0  264;
      0  330;
      0  330;
      0  330;
      0  330;
      0  440;
      0  440;
      0  440
];       % 即Wc/F0,g.h/mol
for i = 1:m
[t PP] = ode45(@Euqations,tspan(i,

,P0(i,

,[],k);
Pcal(i,

= PP(end,

;
end
f= Pcal-Pi;

% ------------------------------------------------------------------
function dPdt = Euqations(t, P, k) % here t = Wc / F0
denom = 1+k(4)*P(1)+k(5)*P(3)+k(6)*P(4)+k(7)*P(5)+k(8)*P(6); % k(4) = KDMM,k(5) = KME ,k(6)=KHPM,k(7)=KPDO,k(8)=KNPA,k(9)=Kp1,k(10)=Kp2
theA =k(3)*k(4)*P(1)*P(2)*(1-P(4)*P(3)/k(9)*P(1)*P(2)^2) / denom^2;
theB =k(6)*k(4)* P(4)*P(2)*(1-P(5)*P(3)/k(10)*P(4)*P(2)^2)/ denom^2;
theC =k(7)*P(5)*P(2)/denom^2;
r1 = k(1)*theA;
r2 = k(2)*theB;
r3 = k(3)*theC;

dPDMMdt = -r1;
dPHdt = -2*r1-2*r2;
dPMEdt = r1+r2;
dPHPMdt = r1-r2;
dPPDOdt = r2-r3;
dPNPAdt = r3;

dPdt = [dPDMMdt;dPHdt;dPMEdt;dPHPMdt;dPPDOdt;dPNPAdt];

回复此楼

» 猜你喜欢

【复旦大学】二维材料方向招收2026年博士研究生1名已经有0人回复
北京纳米能源与系统研究所王中林院士/曹南颖研究员课题组2026级硕/博/博后招生已经有10人回复
物理学I论文润色/翻译怎么收费? 已经有87人回复
荷兰Utrecht University超快太赫兹光谱王海教授课题招收2026 CSC博士生已经有23人回复
反铁磁体中的磁性切换：两种不同的机制已成功可视化已经有0人回复
求标准粉末衍射卡号 ICDD 01-076-1802 已经有0人回复
新西兰Robinson研究所招收全奖PhD 已经有0人回复
石墨烯转移--二氧化硅衬底石墨烯已经有0人回复

1楼 2014-11-22 14:42:42

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

相关版块跳转我要订阅楼主 zhaoshazhu 的主题更新

返回列表