MATLAB Answers

0

How to calculate the rate using a given data set?

Asked by Danny Helwegen on 10 Sep 2019 at 15:59
Latest activity Commented on by Star Strider
on 17 Sep 2019 at 0:29
Hey guys,
I need to estimate the reaction rate from given data in excel. To answer this I set up the differential equations and appointed the given data (for the sake of convenience I briefly reduced the data set here). My first thought was that the most simple way to estimate the rate was to use ODE45, but since it doesn't works I assume I made a mistake somewhere.
This was the code I made:
function [dxdt] = data(t)
A = [0.98 0.97 0.95 0.92 0.90]';
B = [0.99 0.97 0.95 0.93 0.90]';
C = [0.01 0.03 0.07 0.11 0.16]';
A_rad = [0.01 0.02 0.02 0.03 0.03]';
B_rad = [0.01 0.02 0.02 0.03 0.03]';
dAdt = -1 .* A - 100 .* B_rad .* A + 1 .* (A_rad).^2;
dBdt = -100 .* A_rad .* B + 1000 .* B_rad .* C;
dCdt = (20000 ((A).^0.5) .* B) ./(100 + 1000 .*(C./A)); %This is the rate I want to estimate
%dBr_radical = 0; Since equal to zero it can be left out
%dH_radical = 0;
[dxdt] = [dAdt; dBdt; dCdt]
And next to call it:
init = [1 1 0 0 0];
tspan = [0 100];
[t,c] = ode45(@data, tspan, init, []);
plot(t,c)
Can someone help me improving my code?
Thanks in advance,
Kind regards,
Danny

  1 Comment

darova
on 12 Sep 2019 at 12:37
You don't need to use ode45 since you already have table data
This is a result you want (as you wrote above):
dCdt = (20000 ((A).^0.5) .* B) ./(100 + 1000 .*(C./A)); %This is the rate I want to estimate

Sign in to comment.

1 Answer

Answer by Star Strider
on 10 Sep 2019 at 16:31
 Accepted Answer

You are not coding your differential equations and data correctly. See for example: Parameter Estimation for a System of Differential Equations.

  13 Comments

Okay, so I have these differentials and I have to estimate the outcome for A, B, C, A_rad and B_rad and these outcomes have to be compared to the dataset, so they have to be approximatly similar.
Given are the differentials:
dAdt = A^2 - A - 100 * A * B_rad
dBdt = 1000 * B_rad * C - 100 * A_rad * B
dCdt = 100 * A_rad * B + 100 * B_rad * A - 1000 * B_rad * C
dA_raddt = dB_raddt = 0
A maximum time of 10, the begin values of A, B, C, A_rad and B_rad (=[1 1 0 0 0]) and the data set for comparison.
I have found the anser, it was easier than I thought:
This was my final code:
function [dxdt] = Code(t, v)
A = v(1);
B = v(2);
C = v(3);
D = v(4); %B_rad
E = v(5); %H_rad
dAdt = - A - 100 * E * A + D^2 ;
dBdt = - 100 * D * B + 100 * E * C;
dCdt = 100 * D * B + 100 * E * A - 1000 * E * C);
dDdt = 2 * A - 100 * D * B + 100 * E * A + E * C - 2000 * D^2 ;
dEdt = 100 * D * B - 100 * E * A - E * C;
dxdt = [dAdt; dBdt; dCdt; dDdt; dEdt;];
And for calling:
init = [1 1 0 0 0 ];
tspan = [0:0.01:10];
[t,v] = ode45(@AJA, tspan, init, []);
plot(t,v)
Thanks for your help star strider.
As always, my pleasure!

Sign in to comment.