\( \DeclareMathOperator{\abs}{abs} \newcommand{\ensuremath}[1]{\mbox{$#1$}} \)
(%i1) kill(all)$

ode : part7,intro: ricatti
[email protected],06/02/2017

Theory Ricatti : ODE (reaction kinetics)

(dy/dx) = p(x)*y^2+q(x)*y+r(x)
Why: Consecutive Reactions and Ricatti
                     k1   k2
reaction kinetics : A → B → C
dCA/dt = - k1*CA^n          (1)
dCB/dt = k1*CA^n-k2*CB^m    (2)
dCC/dt = k2*CC^l            (3)

integer values :  n,m,l
concentration  : CA,CB,CC

(%i2) load('contrib_ode)$
assume(t>0,CA>0,CB>0)$

equation (2) n=2,m=1, by 't=0 → CA=CA0'

(%i5) eq1: CA = CA0/(1+k1*CA0*t);
eq2:'diff(CB,t)=k1*CA^2-k2*CB^1;
eq3:lhs(eq2)-part(eq2,2,2)=rhs(eq2)-part(eq2,2,2);
\[\tag{eq1}\label{eq1}\mathit{CA}=\frac{\mathit{CA0}}{\mathit{CA0}\,\mathit{k1}t+1}\] \[\tag{eq2}\label{eq2}\frac{d}{dt}\mathit{CB}={{\mathit{CA}}^{2}}\,\mathit{k1}-\mathit{CB}\,\mathit{k2}\] \[\tag{eq3}\label{eq3}\mathit{CB}\,\mathit{k2}+\frac{d}{dt}\mathit{CB}={{\mathit{CA}}^{2}}\,\mathit{k1}\]
(%i7) eq4:subst(rhs(eq1),CA,eq3);
eq5: lhs(eq3)=rhs(eq4);
\[\tag{eq4}\label{eq4}\mathit{CB}\,\mathit{k2}+\frac{d}{dt}\mathit{CB}=\frac{{{\mathit{CA0}}^{2}}\,\mathit{k1}}{{{\left( \mathit{CA0}\,\mathit{k1}t+1\right) }^{2}}}\] \[\tag{eq5}\label{eq5}\mathit{CB}\,\mathit{k2}+\frac{d}{dt}\mathit{CB}=\frac{{{\mathit{CA0}}^{2}}\,\mathit{k1}}{{{\left( \mathit{CA0}\,\mathit{k1}t+1\right) }^{2}}}\]
(%i8) depends([y,f],[x]);                                                                        $
\[\tag{\%{}o8}\label{o8} [\operatorname{y}(x),\operatorname{f}(x)]\]
(%i9) eq6:'diff(y,x)+alpha*y=f;
\[\tag{eq6}\label{eq6}\frac{d}{dx}y+\alphay=f\]
(%i10) depends([I],[x]);
\[\tag{\%{}o10}\label{o10} [\operatorname{I}(x)]\]

I is a integration factor.

(%i11) eq7:ev('diff(I*y,x)=I*f,nouns,eq7);
\[\tag{eq7}\label{eq7}I\,\left( \frac{d}{dx}y\right) +\left( \frac{d}{dx}I\right) y=If\]
(%i12) eq8:ev(eq7,I=%e^(alpha*x));
\[\tag{eq8}\label{eq8}{{\%{}e}^{\alphax}}\,\left( \frac{d}{dx}y\right) +\left( \frac{d}{dx}{{\%{}e}^{\alphax}}\right) y=f\,{{\%{}e}^{\alphax}}\]
(%i13) eq9:ev(eq8,nouns);
\[\tag{eq9}\label{eq9}{{\%{}e}^{\alphax}}\,\left( \frac{d}{dx}y\right) +\alpha{{\%{}e}^{\alphax}}y=f\,{{\%{}e}^{\alphax}}\]

make : eq5 = eq9

first part : eq5

x = t
y = CB
alpha = k2

(%i16) eq10:subst(t,x,eq9);
eq11:subst(CB,y,eq10);
eq12:subst(k2,alpha,eq11);
\[\tag{eq10}\label{eq10}{{\%{}e}^{\alphat}}\,\left( \frac{d}{dt}y\right) +\alpha{{\%{}e}^{\alphat}}y=f\,{{\%{}e}^{\alphat}}\] \[\tag{eq11}\label{eq11}\mathit{CB}\alpha{{\%{}e}^{\alphat}}+\left( \frac{d}{dt}\mathit{CB}\right) \,{{\%{}e}^{\alphat}}=f\,{{\%{}e}^{\alphat}}\] \[\tag{eq12}\label{eq12}\mathit{CB}\,\mathit{k2}\,{{\%{}e}^{\mathit{k2}t}}+\left( \frac{d}{dt}\mathit{CB}\right) \,{{\%{}e}^{\mathit{k2}t}}=f\,{{\%{}e}^{\mathit{k2}t}}\]

second part :eq5

f=(CA0^2*k1)/(CA0*k1*t+1)^2

(%i18) eq13:factor(subst(rhs(eq5),f,eq12));
eq14:f=rhs(eq5);
\[\tag{eq13}\label{eq13}\left( \mathit{CB}\,\mathit{k2}+\frac{d}{dt}\mathit{CB}\right) \,{{\%{}e}^{\mathit{k2}t}}=\frac{{{\mathit{CA0}}^{2}}\,\mathit{k1}\,{{\%{}e}^{\mathit{k2}t}}}{{{\left( \mathit{CA0}\,\mathit{k1}t+1\right) }^{2}}}\] \[\tag{eq14}\label{eq14}f=\frac{{{\mathit{CA0}}^{2}}\,\mathit{k1}}{{{\left( \mathit{CA0}\,\mathit{k1}t+1\right) }^{2}}}\]

find integration factor I:eq5

(%i19) dependencies;
\[\tag{\%{}o19}\label{o19} [\operatorname{y}(x),\operatorname{f}(x),\operatorname{I}(x)]\]
(%i20) depends([I,f,CB],[t]);
\[\tag{\%{}o20}\label{o20} [\operatorname{I}(t),\operatorname{f}(t),\operatorname{CB}(t)]\]
(%i21) dependencies;
\[\tag{\%{}o21}\label{o21} [\operatorname{y}(x),\operatorname{I}(t),\operatorname{f}(t),\operatorname{CB}(t)]\]
(%i22) eq15:I=part(eq13,1,2);
\[\tag{eq15}\label{eq15}I={{\%{}e}^{\mathit{k2}t}}\]

seperated,find CB=CB(t) from eq5
solve eq5:
I = I(t),f=f(t),%c = integration constant
formule : CB = 1/I*integration(I*f,t)+%c/I

(%i23) eq16:CB =ev(1/rhs(eq15)*integrate(rhs(eq15)*rhs(eq14),t)+%c/rhs(eq15),nouns);
\[\tag{eq16}\label{eq16}\mathit{CB}=\mathit{\%{}c}\,{{\%{}e}^{-\mathit{k2}t}}-\frac{\operatorname{expintegral\_ e}\left( 2,-\frac{\mathit{k2}\,\left( \mathit{CA0}\,\mathit{k1}t+1\right) }{\mathit{CA0}\,\mathit{k1}}\right) \mathit{CA0}\,{{\%{}e}^{-\mathit{k2}t-\frac{\mathit{k2}}{\mathit{CA0}\,\mathit{k1}}}}}{\mathit{CA0}\,\mathit{k1}t+1}\]

find %c , same part6, t=0 and CB=0 in eq16

(%i24) eq17:ev(eq16,t=0,CB=0);
\[\tag{eq17}\label{eq17}0=\mathit{\%{}c}-\operatorname{expintegral\_ e}\left( 2,-\frac{\mathit{k2}}{\mathit{CA0}\,\mathit{k1}}\right) \mathit{CA0}\,{{\%{}e}^{-\frac{\mathit{k2}}{\mathit{CA0}\,\mathit{k1}}}}\]
(%i25) eq18:solve(eq17,%c)[1];
\[\tag{eq18}\label{eq18}\mathit{\%{}c}=\operatorname{expintegral\_ e}\left( 2,-\frac{\mathit{k2}}{\mathit{CA0}\,\mathit{k1}}\right) \mathit{CA0}\,{{\%{}e}^{-\frac{\mathit{k2}}{\mathit{CA0}\,\mathit{k1}}}}\]
(%i26) eq19:ev(eq16,eq18);
\[\tag{eq19}\label{eq19}\mathit{CB}=\operatorname{expintegral\_ e}\left( 2,-\frac{\mathit{k2}}{\mathit{CA0}\,\mathit{k1}}\right) \mathit{CA0}\,{{\%{}e}^{-\mathit{k2}t-\frac{\mathit{k2}}{\mathit{CA0}\,\mathit{k1}}}}-\frac{\operatorname{expintegral\_ e}\left( 2,-\frac{\mathit{k2}\,\left( \mathit{CA0}\,\mathit{k1}t+1\right) }{\mathit{CA0}\,\mathit{k1}}\right) \mathit{CA0}\,{{\%{}e}^{-\mathit{k2}t-\frac{\mathit{k2}}{\mathit{CA0}\,\mathit{k1}}}}}{\mathit{CA0}\,\mathit{k1}t+1}\]

integration :integration(%e^(k2*t)*f,t),f  see eq14

info : gamma_incomplete,see www

(%i27) eq20:integrate(%e^(alpha*x)/x^2,x);
\[\tag{eq20}\label{eq20}\operatorname{gamma\_ incomplete}\left( -1,-\alphax\right) \alpha\]

make taylor expantion :
go to  ' menu: Calculus → Get series.'

(%i28) eq21:taylor(%e^(alpha*x)/x^2, x, 0, 4);
\[\tag{eq21}\label{eq21}\frac{1}{{{x}^{2}}}+\frac{\alpha}{x}+\frac{{{\alpha}^{2}}}{2}+\frac{{{\alpha}^{3}}x}{6}+\frac{{{\alpha}^{4}}\,{{x}^{2}}}{24}+\frac{{{\alpha}^{5}}\,{{x}^{3}}}{120}+\frac{{{\alpha}^{6}}\,{{x}^{4}}}{720}+\mbox{...}\]
(%i29) eq22:integrate(eq21,x);
\[\tag{eq22}\label{eq22}\alpha\log{(x)}+\frac{{{\alpha}^{6}}\,{{x}^{5}}}{3600}+\frac{{{\alpha}^{5}}\,{{x}^{4}}}{480}+\frac{{{\alpha}^{4}}\,{{x}^{3}}}{72}+\frac{{{\alpha}^{3}}\,{{x}^{2}}}{12}+\frac{{{\alpha}^{2}}x}{2}-\frac{1}{x}\]

Integration(I*f,t),use taylor expantion to do integration.

subst :  tau=CA0*k1*t+1,alpha=k2

(%i30) eq23:tau=part(eq14,2,2,1);
\[\tag{eq23}\label{eq23}\tau=\mathit{CA0}\,\mathit{k1}t+1\]
(%i31) eq24:diff(rhs(eq23),t);
\[\tag{eq24}\label{eq24}\mathit{CA0}\,\mathit{k1}\]
(%i32) eq25:k2*solve(eq23,t)[1];
\[\tag{eq25}\label{eq25}\mathit{k2}t=\frac{\mathit{k2}\,\left( \tau-1\right) }{\mathit{CA0}\,\mathit{k1}}\]

build %e^(k2*t) ⇒ function of tau

(%i33) eq26:%e^(lhs(eq25))=%e^(rhs(eq25));
\[\tag{eq26}\label{eq26}{{\%{}e}^{\mathit{k2}t}}={{\%{}e}^{\frac{\mathit{k2}\,\left( \tau-1\right) }{\mathit{CA0}\,\mathit{k1}}}}\]

build : %e^(k2*t)*f, f=f(tau)

(%i34) eq27:ratsubst(tau,rhs(eq23),eq14);
\[\tag{eq27}\label{eq27}f=\frac{{{\mathit{CA0}}^{2}}\,\mathit{k1}}{{{\tau}^{2}}}\]
(%i35) eq28:v=integrate(rhs(eq25)*rhs(eq27),tau);
\[\tag{eq28}\label{eq28}v=\mathit{CA0}\,\mathit{k2}\,\left( \log{\left( \tau\right) }+\frac{1}{\tau}\right) \]
(%i37) eq29:part(eq28,2,3,1);
eq30:ratsubst(rhs(eq23),lhs(eq23),eq29);
\[\tag{eq29}\label{eq29}\log{\left( \tau\right) }\] \[\tag{eq30}\label{eq30}\log{\left( \mathit{CA0}\,\mathit{k1}t+1\right) }\]
(%i38) taylor(eq30, t, 0,3);
\[\tag{\%{}o38)/T}\label{o38)/T} \mathit{CA0}\,\mathit{k1}t-\frac{{{\mathit{CA0}}^{2}}\,{{\mathit{k1}}^{2}}\,{{t}^{2}}}{2}+\frac{{{\mathit{CA0}}^{3}}\,{{\mathit{k1}}^{3}}\,{{t}^{3}}}{3}+\mbox{...}\]
(%i39) eq31:CB=1/rhs(eq15)*rhs(eq28)+%c/rhs(eq15);
\[\tag{eq31}\label{eq31}\mathit{CB}=\mathit{CA0}\,\mathit{k2}\,{{\%{}e}^{-\mathit{k2}t}}\,\left( \log{\left( \tau\right) }+\frac{1}{\tau}\right) +\mathit{\%{}c}\,{{\%{}e}^{-\mathit{k2}t}}\]

t = 0 ,tau = 1+CA0*k1*t,tau = 1 and CB=0

(%i45) eq32:ratsubst(lhs(eq23),rhs(eq23),eq31);
eq33:solve(eq23,t)[1];
eq34:ev(eq32,tau=1,CB=0,lhs(eq33)=rhs(eq33));
eq35:solve(eq34,%c)[1];
eq36:ev(eq32,eq33,eq35);
eq37:factor(ev(eq36,eq23));
\[\tag{eq32}\label{eq32}\mathit{CB}=\frac{{{\%{}e}^{-\mathit{k2}t}}\,\left( \mathit{CA0}\,\mathit{k2}\tau\log{\left( \tau\right) }+\mathit{\%{}c}\tau+\mathit{CA0}\,\mathit{k2}\right) }{\tau}\] \[\tag{eq33}\label{eq33}t=\frac{\tau-1}{\mathit{CA0}\,\mathit{k1}}\] \[\tag{eq34}\label{eq34}0=\mathit{CA0}\,\mathit{k2}+\mathit{\%{}c}\] \[\tag{eq35}\label{eq35}\mathit{\%{}c}=-\mathit{CA0}\,\mathit{k2}\] \[\tag{eq36}\label{eq36}\mathit{CB}=\frac{{{\%{}e}^{-\frac{\mathit{k2}\,\left( \tau-1\right) }{\mathit{CA0}\,\mathit{k1}}}}\,\left( \mathit{CA0}\,\mathit{k2}\tau\log{\left( \tau\right) }-\mathit{CA0}\,\mathit{k2}\tau+\mathit{CA0}\,\mathit{k2}\right) }{\tau}\] \[\tag{eq37}\label{eq37}\mathit{CB}=\frac{\mathit{CA0}\,\mathit{k2}\,{{\%{}e}^{-\mathit{k2}t}}\,\left( \mathit{CA0}\,\mathit{k1}t\,\log{\left( \mathit{CA0}\,\mathit{k1}t+1\right) }+\log{\left( \mathit{CA0}\,\mathit{k1}t+1\right) }-\mathit{CA0}\,\mathit{k1}t\right) }{\mathit{CA0}\,\mathit{k1}t+1}\]
Created with wxMaxima.