how to solve a non linear ode equation ?

Hello,
i have the following ODE that i want to solve
where k_r and kaq are constant values.
i am really stuck on what approach should i take to implement the equation in matalb, i have tried to linearize the equation using odetovectorfield function to solve it using ode45 or another solver but it won't work since the equation is not a quasi-linear differential equation as explained in the function description.
I would be grateful for any guidance or help
thanks in advance !

답변 (1개)

Bjorn Gustavsson
Bjorn Gustavsson 2021년 9월 14일

0 개 추천

You have a quadratic equation for dCAdt, solve that one as a standard quadratic equation, you will end up with two solutions for dCAdt that will be possible to solve. If you're somewhat lucky you might end up with solutions that are simple enough to solve analyticalls. Depending on your initial condition you might have to discard one of the solutions and select one to proceed with.
HTH

댓글 수: 4

Thanks for your fast response Bjorn !
You are right and that is the first thing i have tried. However, the obtained results are far from what i would expect as an answer. That's why, i was wondering if there is another appraoch to solve it.
I have attached a pdf file where i demonstrate how to obtain two expressions for dCA/dt and the following code to solve them but i don't understand why my CA variable is not changing with time.
clear all
%
Kr=8.6876E-10; % reaction transfer coefficient (m4.mol-1.s-1)
C0=116.9; % initial uranium concentration (mol.m-3)
kaq=[1.24E-04 1.43E-04 1.57E-04 1.63E-04 1.69E-04 1.73E-04]; % aqueous phase mass transfer coefficient for difference residence time (m/s)
ts_aq=[0.033135652 0.019419467 0.013737619 0.011995912 0.010648619 0.009601691]; % aqueous phase residence time (s)
S=1.25E04; % A/V (m-1)
a1=(S/2*Kr).*kaq.^2;
a2=(2*Kr)./kaq;
a3=2*a2;
Caq_out1=zeros(1,length(kaq));
Caq_out2=zeros(1,length(kaq));
for i=1:length(kaq)
tspan=[0 ts_aq(i)];
dCAdt=@(t,CA) a1(i)*(a2(i)*CA-1+sqrt(1-a3(i)*CA));
[t,CA]=ode45(dCAdt,tspan,C0);
Caq_out1(i)=CA(length(CA));
end
for i=1:length(kaq)
tspan=[0 ts_aq(i)];
dCAdt=@(t,CA) a1(i)*(a2(i)*CA-1-sqrt(1-a3(i)*CA));
[t,CA]=ode45(dCAdt,tspan,C0);
Caq_out2(i)=CA(length(CA));
end
Efficiency1=Caq_out1/C0;
Efficiency2=Caq_out2/C0;
Are you sure about the units of your coefficients and CA? When try:
dCAdt(tspan(1),C0)
For the first case (i=1) it returns:
-1.1219e-19
Which is a very small change from an initial condition of 116.9.
My guess is that you're off by something like an elementary charge, q_e, somewhere (not that I have personal experience with exactly that goof or anything...).
hamza karim
hamza karim 2021년 9월 16일
편집: hamza karim 2021년 9월 16일
thanks again for your reply Bjorn,
Regarding the units, yes i am sure. I even have done dimensional analysis to verify but it is correct. I will look further in the problem
thanks again
Regarding units again: you might have forgotten to correct reaction-rates from cm³/s to m³/s or some such easy to forget step. That would give you something like a 6 orders of magnitude change in reaction-rates, but would be hard to capture with dimensional analysis...

댓글을 달려면 로그인하십시오.

카테고리

도움말 센터File Exchange에서 Programming에 대해 자세히 알아보기

태그

질문:

2021년 9월 14일

댓글:

2021년 9월 16일

Community Treasure Hunt

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

Start Hunting!

Translated by