24小时热门版块排行榜

返回列表

1135725495

铁杆木虫 (著名写手)

应助: 1 (幼儿园)
金币: 6868.7
散金: 18
红花: 1
帖子: 1492
在线: 209.6小时
虫号: 2756771
注册: 2013-10-27
性别: MM
专业: 有机化工

[求助] 运行matlab出现的问题

Caused by:
Failure in initial user-supplied objective function evaluation. FMINCON cannot continue.
谁知道这是怎么回事啊？
function parafit
%
% r1 =(k(1)*c^2-0.5*k(2)*(9.404-c))/(1+k(3)*c+0.5*k(4)*(9.404-c)+2*k(5)*c^2/(k(4)*(9.404-c)));
%
% dCAdt = - r1;
clear all
clc
%       t/min CA    / mol/L
  Kinetics=[0 7.9556
3 5.0759
6 4.3477
9 3.9301
12 3.6489
15 3.4439
18 3.2875
21 3.1645
31 2.8989
41 2.7540
51 2.6694
61 2.6182];
k0 = [0.1  0.1  0.01  0.01  0.01];       % 参数初值
lb = [0  0  0  0 0];                % 参数下限
ub = [100 100  100  100  100]; % 参数上限
x0 = [7.9556];
yexp = Kinetics;
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('  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('  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('  The sum of the squares is: %.1e\n\n',resnorm)
k_fmls = k;
output
tspan = [0 3 6 9 12 15 18 21 31 41 51 61];
[t x] = ode45(@KineticEqs,tspan,x0,[],k_fmls);
figure;
plot(t,x,'o-',t,yexp(:,2),'*');legend('ca-pr','ca-real')
x

% ------------------------------------------------------------------
function f = ObjFunc7Fmincon(k,x0,yexp)
tspan = [0 3 6 9 12 15 18 21 31 41 51 61];
[t x] = ode45(@KineticEqs,tspan,x0,[],k);
f =  sum((x-yexp(:,2)).^2) ;
% ------------------------------------------------------------------
function f = ObjFunc7LNL(k,x0,yexp)
tspan = [0 3 6 9 12 15 18 21 31 41 51 61];
[t x] = ode45(@KineticEqs,tspan,x0,[],k);
f1 = x - yexp(:,2);
f = [f1];
% ------------------------------------------------------------------
function dxdt = KineticEqs(t,x,k)
dxdt =-(k(1)*x^2-0.5*k(2)*(9.404-x))/(1+k(3)*c+0.5*k(4)*(9.404-x)+2*k(5)*x^2/(k(4)*(9.404-x)));

回复此楼