How to calculate the rate using a given data set?

19 次查看(过去 30 天)
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 个评论
darova
darova 2019-9-12
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

请先登录,再进行评论。

采纳的回答

Star Strider
Star Strider 2019-9-10
You are not coding your differential equations and data correctly. See for example: Parameter Estimation for a System of Differential Equations.
  13 个评论
Danny Helwegen
Danny Helwegen 2019-9-16
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.

请先登录,再进行评论。

更多回答(0 个)

类别

Help CenterFile Exchange 中查找有关 Polymers & Soft Matter 的更多信息

Community Treasure Hunt

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

Start Hunting!

Translated by