Differential model won't calculate with powers <1

1 次查看（过去 30 天）
Jort Puiman
Jort Puiman 2021-3-8
评论： Jort Puiman 2021-3-10
Hi everyone,
I am working on a model which describes multiple reactions over time, using simple power law models (R = k*[C]^n).
ode1 = diff(Ca) == -k1*Ca;
cond1 = Ca(0) == 5;
ode2 = diff(Cb) == k1*Ca - k2*Cb;
cond2 = Cb(0) == 0;
ode3 = diff(Cc) == k2*Cb - k3*Cc;
cond3 = Cc(0) == 0;
ode4 = diff(Cd) == k3*Cc - k4*Cd;
cond4 = Cd(0) == 0;
ode5 = diff(Ce) == k4*Cd;
cond5 = Ce(0) == 0;
t works just like I want to with n = 1, however, our data suggests that n < 1. I tried adding powers to my concentrations, but then, Matlab has a hard time calculating it, and it never finishes.
I want to calculate the concentrations of all components over time. All constants (k1, k2, k3, k4) and reaction orders will be added later. However, the model should be able to calculate the concentrations when the order is lower than 1.
Do you have any suggestions on how to tackle this problem? Do I have to use other functions?
  2 个评论
Jan
Jan 2021-3-8
You forgot to mention how you try to calculate what. These are the equations only.
Jort Puiman
Jort Puiman 2021-3-8
Thank you for your quick remark, I edited my question. Hope it is clear right now!

请先登录，再进行评论。

请先登录，再回答此问题。

采纳的回答

Bjorn Gustavsson
Bjorn Gustavsson 2021-3-8
It is not given that a general nonlinear system of ODEs have analytical solutions. If you cannot get the symbolic toolbox to find one for you, the best (possibly?) advice is to turn to numerical solutions. For that have a look at ode45 and its siblings, and the large number of ODE-examples. This should be very straightforward to implement for numerical integration. The ODE-function might look something like this:
function dCdt = your_power_reactions_ODE(t,C,k,g)
Ca = C(1);
Cb = C(2);
Cc = C(3);
Cd = C(4);
Ce = C(5);
dCadt = -k(1)*Ca^g(1);
dCbdt = k(1)*Ca^g(1) - k(2)*Cb^g(2);
dCcdt = k(2)*Cb^g(2) - k(3)*Cc^g(3);
dCddt = k(3)*Cc^g(3) - k(4)*Cd^g(4);
dCedt = k(4)*Cd^g(4);
% etc...
end
That should be pretty standard to integrate. It might be a numerical stiff set of ODEs - I've not spent any time to judge that. You then simply set the initial conditions and the time-period and then run:
C0 = [5 0 0 0 0];
t_span = [0 17];
k_all = [1 log(2),2^pi,12];
g_all = 1./[1:4];
[t_out,C_out] = ode45(@(t,C) your_power_reactions_ODE(t,C,k_all,g_all),...
t_span,C0);
If it is too stiff, try ode15s, ode113, or ode23s...
HTH
  3 个评论
显示 1更早的评论
Bjorn Gustavsson
Bjorn Gustavsson 2021-3-8
But to my inderstanding these are all "decay-type reactions" - like for the concentrations of elements in a radioactive decay-chain?
Jort Puiman
Jort Puiman 2021-3-10
Yes you could compare it to a radioactive decay-chain. We start with a set concentration of Ca, and it oxidises to Cb, then Cc, etc... I however want to get as much of Cd, and as few of Ce. That's the general idea.

请先登录，再进行评论。

更多回答（0 个）

类别

在 Help Center 和 File Exchange 中查找有关 Particle & Nuclear Physics 的更多信息

标签

Community Treasure Hunt

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

Start Hunting!

Translated by