Parameter estimation and fitting kinetic model

Can you help me with a code to estimate the parameters in this kinetic model (dX_A)/dt=Ae^((-E_a)/RT) C_Ao^(n-1) (1-X_A )^n. I want to fit my experimental data and determine the values of Ae, E_a, and n.

Answer 1

[praguna manvi]

Hi @BONGO,

The kinematic system in non-linear, hence “lsqonlin” function can be used in estimating “Ae”, “E_a” and “n” parameters. Here is a link of the documentation of “lsqonlin” :

https://www.mathworks.com/help/optim/ug/lsqcurvefit.html

Below is a small example at solving the problem using an “ode45” solver to estimate "X_A" from the kinetic model:

% Define the kinetic model function
function dxdt = kineticModel(t, x, params, C_Ao, R, T)
    Ae = params(1);
    Ea = params(2);
    n = params(3);
    dxdt = Ae * exp(-Ea / (R * T)) * C_Ao^(n-1) * (1 - x)^n;
end

% Experimental data
time_exp = [0.1, 0.2, 0.3, 0.4, 0.5]; % Example time data
X_A_exp = [0.05, 0.15, 0.25, 0.35, 0.45]; % Example conversion data

% Constants
C_Ao = 1.0; % Initial concentration (example value)
R = 8.314; % Universal gas constant (J/(mol*K))
T = 298; % Temperature in Kelvin (example value)

% Initial guesses for parameters [Ae, Ea, n]
initial_guesses = [1e6, 50000, 1];

% Define the objective function for lsqnonlin
objectiveFunction = @(params) arrayfun(@(t) solveODE(params, t, C_Ao, R, T), time_exp) - X_A_exp;

% Perform the curve fitting
options = optimoptions('lsqnonlin', 'Display', 'iter');
[params_fit, resnorm] = lsqnonlin(objectiveFunction, initial_guesses, [], [], options);

% Extract the fitted parameters
Ae_fit = params_fit(1);
Ea_fit = params_fit(2);
n_fit = params_fit(3);

fprintf('Fitted parameters:\n');
fprintf('Ae = %.4e\n', Ae_fit);
fprintf('Ea = %.4f J/mol\n', Ea_fit);
fprintf('n = %.4f\n', n_fit);

% Helper function to solve the ODE
function X_A = solveODE(params, t, C_Ao, R, T)
    % Initial condition
    x0 = 0; % Assuming conversion starts from 0
    
    % Solve the ODE
    [~, X_A] = ode45(@(t, x) kineticModel(t, x, params, C_Ao, R, T), [0 t], x0);
    
    % Return the conversion at time t
    X_A = X_A(end);
end

                                            Norm of      First-order 
 Iteration  Func-count      Resnorm            step       optimality
     0          4          0.410868                         3.29e-07

Initial point is a local minimum.

Optimization completed because the size of the gradient at the initial point 
is less than the value of the optimality tolerance.

Fitted parameters:

Ae = 1.0000e+06

Ea = 50000.0000 J/mol

n = 1.0000