24小时热门版块排行榜

返回列表

llgcty

金虫 (正式写手)

应助: 5 (幼儿园)
金币: 985.4
散金: 5
红花: 1
帖子: 337
在线: 143.9小时
虫号: 1073112
注册: 2010-08-11
专业: 能源化工

[求助] 在线等，请版主看看如何运算不了。总是报错

小弟实在是一窍不通对matlab，不要怪小弟提的问题比较傻，一个是报错中的问题，这个[k,fval,flag] = fmincon(@ObjFunc7Fmincon一句的@ObjFunc7Fmincon是怎么设定的？我看有的人写@ObjFunc4Fmincon.
然后麻烦大家看看这个我自己模拟的代码，是按照版主 @dbb627的代码写的，到底是哪里不对麻烦改下谢谢了。
function parafit
%
% r1 = k(1)*C(1);
% r2 = k(2)*C(1);
% r3 = k(3)*C(2);
% r4 = k(4)*C(3);
% r5 = k(5)*C(2);
% r6 = k(6)*C(3);
%
% dCAdt = - r1 - r2;
% dCBdt = r1 + r4 - r3 - r5;
% dCCdt = r2 + r3 - r4 – r6 ;
% dCDdt = r5 + r6;
clear all
clc
%       t/min CA-1    CB-2       CC-3 CD-4 / mol/L
  Kinetics=[0    100 0       0       0
      10       0.2441  0.1390 0.28.31 0.0171
      20       0.432  0.1629 0.2370 0.0192
      30       0.756 0.1992  0.2509 0.0349
      40       0.522 0.2091 0.2737 0.0531
      50       0.426 0.1967 0.3231 0.0634
      60       0.397 0.1803 0.3536 0.0822
      70       0.431 0.1732 0.3330 0.0945
      80       0.531 0.1645 0.3310  0.0987];
k0 = [0.0000000005  0.0000000005  0.0000000005  0.00000000005  0.00000005  0.00000005];       % 参数初值
lb = [0  0  0  0  0  0];                % 参数下限
ub = [1 1  1  1  1  1]; % 参数上限
x0 = [100  0  0  0];
yexp = Kinetics;                % yexp: 实验数据[x1       x4       x5       x6]
warning off
% 使用函数fmincon()进行参数估计
[k,fval,flag] = fmincon(@ObjFunc7Fmincon,k0,[],[],[],[],lb,ub,[],[],x0,yexp);
fprintf('\n使用函数fmincon()估计得到的参数值为:\n')
fprintf('\tk1 = %.11f\n',k(1))
fprintf('\tk2 = %.11f\n',k(2))
fprintf('\tk3 = %.11f\n',k(3))
fprintf('\tk4 = %.11f\n',k(4))
fprintf('\tk5 = %.11f\n',k(5))
fprintf('\tk6 = %.11f\n',k(6))
fprintf('  The sum of the squares is: %.1e\n\n',fval)
k_fm= k;
warning off
% 使用函数lsqnonlin()进行参数估计
[k,resnorm,residual,exitflag,output,lambda,jacobian] = ...
lsqnonlin(@ObjFunc7LNL,k0,lb,ub,[],x0,yexp);
ci = nlparci(k,residual,jacobian);
fprintf('\n\n使用函数lsqnonlin()估计得到的参数值为:\n')
fprintf('\tk1 = %.11f\n',k(1))
fprintf('\tk2 = %.11f\n',k(2))
fprintf('\tk3 = %.11f\n',k(3))
fprintf('\tk4 = %.11f\n',k(4))
fprintf('\tk5 = %.11f\n',k(5))
fprintf('\tk6 = %.11f\n',k(6))
fprintf('  The sum of the squares is: %.1e\n\n',resnorm)
k_ls = k;
output
warning off
% 以函数fmincon()估计得到的结果为初值，使用函数lsqnonlin()进行参数估计
k0 = k_fm;
[k,resnorm,residual,exitflag,output,lambda,jacobian] = ...
lsqnonlin(@ObjFunc7LNL,k0,lb,ub,[],x0,yexp);
ci = nlparci(k,residual,jacobian);
fprintf('\n\n以fmincon()的结果为初值，使用函数lsqnonlin()估计得到的参数值为:\n')
fprintf('\tk1 = %.11f\n',k(1))
fprintf('\tk2 = %.11f\n',k(2))
fprintf('\tk3 = %.11f\n',k(3))
fprintf('\tk4 = %.11f\n',k(4))
fprintf('\tk5 = %.11f\n',k(5))
fprintf('\tk6 = %.11f\n',k(6))
fprintf('  The sum of the squares is: %.1e\n\n',resnorm)
k_fmls = k;
output
tspan = [0 : 10 : 80];
[t x] = ode45(@KineticEqs,tspan,x0,[],k_fmls);
figure;
plot(t,x(:,1),t,yexp(:,2),'*');legend('ca-pr','ca-real')
figure;plot(t,x(:,2:5));
hold on
plot(t,yexp(:,3:6),'o');legend('cb-pr','cc-pr','cd-pr','ce-pr','cb-real','cc-real','cd-real','ce-real')
% ------------------------------------------------------------------
function f = ObjFunc7Fmincon(k,x0,yexp)
tspan = [0 : 10 : 80];
[t x] = ode45(@KineticEqs,tspan,x0,[],k);
y(:,2) = x(:,1);
y(:,3:6) = x(:,2:5);
f =  sum((y(:,2)-yexp(:,2)).^2) + sum((y(:,3)-yexp(:,3)).^2) ...
+ sum((y(:,4)-yexp(:,4)).^2) + sum((y(:,5)-yexp(:,5)).^2) ...
+ sum((y(:,6)-yexp(:,6)).^2) ;

% ------------------------------------------------------------------
function f = ObjFunc7LNL(k,x0,yexp)
tspan = [0.0 : 10 : 80];
[t x] = ode45(@KineticEqs,tspan,x0,[],k);
y(:,2) = x(:,1);
y(:,3:6) = x(:,2:5);
f1 = y(:,2) - yexp(:,2);
f2 = y(:,3) - yexp(:,3);
f3 = y(:,4) - yexp(:,4);
f4 = y(:,5) - yexp(:,5);
f5 = y(:,6) - yexp(:,6);
f = [f1; f2; f3; f4; f5];

% ------------------------------------------------------------------
function dxdt = KineticEqs(t,x,k)
r1 = k(1)*x(1);
r2 = k(2)*x(1);
r3 = k(3)*x(2);
r4 = k(4)*x(3);
r5 = k(5)*x(2);
r6 = k(6)*x(3);

dCAdt = - r1 - r2;
dCBdt = r2 + r4 - r3 - r5;
dCCdt = r2 + r3 - r4 - r6;
dCDdt = r5 + r6;

dxdt = [dCAdt; dCBdt; dCCdt; dCDdt];

回复此楼