How to find fitting optimized parameters to fit a system of non-linear ODEs to experiment.

18 Ansichten (letzte 30 Tage)
Hi
I have a set of ODEs (attached), I have been able to solve them using ode45, however, my issue now is my experimental results don't match the integrated values of the equations. So, I am looking to fit only the solution for epsilon with it's experimental results to find the best parameters A, B, (A0/alpha), k0, Q, and QG. Attached is my code based on an answer from another thread, but it just runs continuously but I couln't figure out what the problem is. Could it be that the there are too many parameters to fit? Any help is greatly appreciated. Thank you.
  7 Kommentare
Continuum
Continuum am 27 Jan. 2024
If you observe the data, there is a period where the temp stayed almost constant for some time.
Torsten
Torsten am 27 Jan. 2024
Bearbeitet: Torsten am 27 Jan. 2024
But it's one long experiment with a dynamic development in time. So your measurements have a history. It's usually necessary that you start at t = 0 with a fixed temperature which is kept constant over time until the experiment has finished.
But we are in a MATLAB forum here ...

Melden Sie sich an, um zu kommentieren.

Akzeptierte Antwort

Star Strider
Star Strider am 27 Jan. 2024
You will need to post the ‘modified_data.xlsx’ file as well to run this here.
My edit of your code —
Model_Exp_Fit
function Model_Exp_Fit
% 2016 12 03
% NOTES:
%
% 1. The ‘theta’ (parameter) argument has to be first in your
% ‘kinetics’ funciton,
% 2. You need to return ALL the values from ‘DifEq’ since you are fitting
% all the values
dataTable1=readtable('modified_data.xlsx','Sheet','Sheet1', 'VariableNamigRule','preserve');
dataTable1 = dataTable1(:,4:end);
tempFunc = griddedInterpolant(dataTable1.time1_min_,dataTable1.T1_K_);
% Temp = tempFunc(t);
% tT = [t Temp];
tT = dataTable1{:,[1 2]}; % Time And Temperature
t = tT(:,1);
x = dataTable1{:,3}; % Epsilon
% return
function X = evolution(theta,tT)
T = tT(:,2);
tv = tT(:,1);
x0=[20;0.4363;0.0000001];
[T,Xv] = ode15s(@DifEq,t,x0);
k0 = theta(1);
QG = theta(2);
A0_alpha = theta(3);
Q = theta(4);
A = theta(5);
B = theta(6);
function dX = DifEq(t,x)
%thet_c = 0.508;
rho_c = 0.948;
R = 8.314;
Temp = interp1(tv, T, t); % Interpolates 'T' For Each Value Of 't'
% dataTable1=readtable('modified_data.xlsx','Sheet','Sheet1')
% tempFunc = griddedInterpolant(dataTable1.time1_min_,dataTable1.T1_K_);
% Temp = tempFunc(t);
% tT = [t Temp];
dxdt = zeros(3,1);
phi1 = (theta(1)/3)*(1-rho_c)^1.5;
phi2 = 9/(2*theta(3));
phi3 = 3/(2*theta(3));
dxdt(1) = (exp(-theta(2)./(R*Temp))).*(phi1./x(1).^2).*1./((1-rho_c+x(2)).^1.5);
dxdt(2) = (-phi2./(x(1).*Temp)).*(1./exp(theta(4)./(R*Temp))).*x(2).^theta(6).*(((1-x(2)).^3)./(1-x(2)).^theta(5));
dxdt(3) = (-phi3./(x(1).*Temp)).*(1./exp(theta(4)./(R*Temp))).*x(2).^theta(6).*(((1-x(2)).^2)./(1-x(2)).^theta(5));
dX = dxdt;
end
X=Xv(:,3);
end
% return
% eps_exp = load("exp_epsilon.txt");
% t = eps_exp(:,1);
% t = [0;t];
% % x = eps_exp(:,2)/100;
% x = [0; x];
theta0 = [(29.65e-5)/(1/60);164.8e3;2.03*((1/60)/(10^6));217.2e3;11.35;0.49];
[theta,Rsdnrm,Rsd,ExFlg,OptmInfo,Lmda,Jmat]=lsqcurvefit(@evolution,theta0,tT,x);
fprintf(1,'\tRate Constants:\n')
for k1 = 1:length(theta)
fprintf(1, '\t\tTheta(%d) = %8.5f\n', k1, theta(k1))
end
% tv = linspace(min(t), max(t));
Xfit = evolution(theta, tT);
figure(1)
plot(t, x, 'p')
hold on
hlp = plot(t, Xfit);
hold off
grid
xlabel('Time')
ylabel('Concentration')
% legend(hlp, 'G(t)', 'theta(t)', 'epsilon(t)', 'Location','N')
end
I left in the sections that I commented-out so that you can understand my edits. I will give this a shot with some standardized code that I use to run genetic algorithm parameter estimations to see if I can improve on these results. If I get good results with it, I will post the results and the code here. For the time being, this code runs, and you can experiment with it.
.
  2 Kommentare
Continuum
Continuum am 27 Jan. 2024
Thanks very much, it works now, I will play around with it a bit and see if I can get some relatively good results. Thanks very much for your input.

Melden Sie sich an, um zu kommentieren.

Weitere Antworten (1)

William Rose
William Rose am 27 Jan. 2024
  2 Kommentare
William Rose
William Rose am 27 Jan. 2024
@Buhari Ibrahim, You are fitting six parameters. That is not too many. I have fitted ODE models with more parameters. You say you have attached code, but none is attached.
Continuum
Continuum am 27 Jan. 2024
Thanks for the pointers, I reattached the code, sorry about that.

Melden Sie sich an, um zu kommentieren.

Kategorien

Mehr zu Least Squares finden Sie in Help Center und File Exchange

Community Treasure Hunt

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

Start Hunting!

Translated by