Problen with optimization of two parameters, and with differential equations

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

Answer 1

[Torsten]

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.

[Jorge Armando Vazquez]

But I need to find the values of k0 and Ea, is there no other way?

[Torsten]

No. They cannot be estimated independently if T does not change during your integration.

[Jorge Armando Vazquez]

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

Answer 2

[Torsten]

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