(%i1) | kill(all)$ |
ode : part6,intro: ricatti
[email protected],05/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); |
(%i7) |
eq4:subst(rhs(eq1),CA,eq3); eq5: lhs(eq3)=rhs(eq4); |
(%i9) |
eq6:contrib_ode(eq5,CB,t)[1]; method; |
function :build in maxima = 'expintegral_e(r,z)'
info :
https://en.wikipedia.org/wiki/Exponential_integral
http://maxima-online.org/help/index/expintegral_e
http://functions.wolfram.com/GammaBetaErf/ExpIntegralE/
http://encyclopedia2.thefreedictionary.com/Exponential+integral
(%i10) | eq7:expintegral_e(2,-(k2*(CA0*k1*t+1))/(CA0*k1)); |
find %c,by t=0 and CB=0 in eq6
(%i12) |
eq8:ev(eq6,t=0,CB=0); eq9:solve(eq8,%c)[1]; |
solution : n=1,m=2 for eq5
(%i13) | eq10:ev(eq6,eq9); |