清空筛选条件
清空筛选条件

Problen with optimization of two parameters, and with differential equations

1 次查看（过去 30 天）
Jorge Armando Vazquez
编辑： Torsten 2023-4-17
I made this program for a chemistry class, but for some reason it shows that the objective function enters optimal values but there really isn't a good similarity between the real data and the optimized data, I don't know if it's a problem of how I wrote the optimizer or because of the resolution of differentials
the data is:
0 0
60 0.28
120 0.36
180 0.42
240 0.44
300 0.445
360 0.47
420 0.5
clc
clear all
%Valores iniciales de parametros
par0 = [690,30 ];
data = load('acido.txt');
%Entrada de prueba de presion
fun_objetivo = @(par)FunObjetivo(par,data(:,:));
t = data(:,1);
xe= data(:,2);
%Argumentos de entrada
A =[]; b =[];
Aeq =[]; beq =[];
lb =[500 0]; ub =[1000 100];
IntCon =[]; nonlcon =[];
nvar =2;
%Nelder_Mead
options = optimset('MaxIter',5000,'MaxFunEvals',...
3000,'FunValCheck','off','Display','iter');
%We use here the Nelder-Mead method
par_optimos = fminsearchbnd(fun_objetivo,par0,lb,ub, options);
function FO =FunObjetivo(par,data)
%Vector tiempo
t = data(:,1);
%Vector conversion real
xe = data(:,2);
%Constantes
cc=0.365851537;
ca0=7.718472259;
keq=3.6021;
T=323.15;
k0=par(1);
Ea=par(2);
%Simulation
y0=0;
dxadt=@(t,xa) k0*exp(-Ea/(8310*T))*cc*ca0*((1-xa^2)-((xa^2)/keq));
tspan=[0:60:420];
[t,xa]=ode45(dxadt,tspan,y0);
%Valor de la Funcion Objetivo
FO = sum((xe-xa).^2);
end

请先登录，再回答此问题。

回答（2 个）

Torsten
Torsten 2023-4-17
移动：Torsten 2023-4-17
The parameters you want to fit are not independent.
Since T remains constant, the complete expression
k0*exp(-Ea/(8310*T))*cc*ca0
can be merged to one parameter to be fitted.
  3 个评论
显示 1更早的评论
Torsten
Torsten 2023-4-17
编辑：Torsten 2023-4-17
No. They cannot be estimated independently if T does not change during your integration.
Jorge Armando Vazquez
Jorge Armando Vazquez 2023-4-17
Ok I change it, but the problem remains, maybe is other thing?
clc
clear all
%Valores iniciales de parametros
par0 = [0.0008];
data = load('acido.txt');
%Entrada de prueba de presion
fun_objetivo = @(par)FunObjetivo(par,data(:,:));
t = data(:,1);
xe= data(:,2);
%Argumentos de entrada
A =[]; b =[];
Aeq =[]; beq =[];
lb =[0]; ub =[1];
IntCon =[]; nonlcon =[];
nvar =1;
%Nelder_Mead
options = optimset('MaxIter',5000,'MaxFunEvals',...
3000,'FunValCheck','off','Display','iter');
%We use here the Nelder-Mead method
par_optimos = fminsearchbnd(fun_objetivo,par0,lb,ub, options);
function FO =FunObjetivo(par,data)
%Vector tiempo
t = data(:,1);
%Vector conversion real
xe = data(:,2);
%Constantes
cc=0.365851537;
ca0=7.718472259;
keq=3.6021;
Te=323.15;
k1=par(1);
%Simulation
y0=0;
dxadt=@(t,xa) k1*cc*ca0*((1-xa^2)-((xa^2)/keq));
%dxadt=@(t,xa) k0*exp(-Ea/(8310*Te))*cc*ca0*((1-xa^2)-((xa^2)/keq));
tspan=[0:60:420];
[t,xa]=ode45(dxadt,tspan,y0);
%Valor de la Funcion Objetivo
FO = sum((xe-xa).^2);
end

请先登录，再进行评论。

Torsten
Torsten 2023-4-17
编辑：Torsten 2023-4-17
Ran in:
All curves that stem from your model tend to sqrt(keq/(1+keq)). Since your measurement data tend to 0.5, my guess is that you have to choose keq such that 0.5 = sqrt(keq/(1+keq)) to get a proper fit. The constant you try to fit gives the velocity with which this equilibrium value of 0.5 is reached.
syms c t keq xa(t)
eqn = diff(xa,t)==c*((1-xa^2)-xa^2/keq);
cond = xa(0)==0;
sol = dsolve(eqn,cond)
sol = 
keqnum = double(solve(sqrt(keq/(keq+1))==0.5,keq))
keqnum = 0.3333
%keqnum=3.6021;
sol = matlabFunction(subs(sol,keq,keqnum));
data=[0 0
60 0.28
120 0.36
180 0.42
240 0.44
300 0.445
360 0.47
420 0.5];
t = data(:,1);
xe = data(:,2);
c = 0:0.001:0.01;
hold on
plot(t,xe,'o')
for i = 1:numel(c)
xa = sol(c(i),t);
plot(t,xa)
end
hold off

类别

在 Help Center 和 File Exchange 中查找有关 Solver-Based Nonlinear Optimization 的更多信息

标签

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by