- Submitted by
- SciPy Central on 22 July 2011
- Download
- File
- Update history
- Revision 2 of 2: previous
- Updated by
- SciPy Central on 06 August 2011
- Tags

We extend the example from http://scpyce.org/12/ here to integrate 2 coupled ODEs and include several algebraic equations.

Consider a chemical reactor with a second order reaction \({\sf A} \rightarrow {\sf B}\). Our aim is calculate the concentration, \(C_{\sf A}(t)\), and tank temperature, \(T(t)\), as functions of time using an ODE integrator.

We know the reaction rate is \(r = k C_{\sf A}^2\), and is an algebraic function of temperature: \(k = 0.15 e^{- E_a/(RT)}\) L/(mol.min). Furthermore, the reaction is known to release heat, with \(\Delta H_r = -590\) J/mol. The time-dependent mass and energy balance for this system can be derived:

Notice that these are coupled ODEs, as the first ODE is a function of \(T(t)\), while the second ODE is a function of \(C_{\sf A}(t)\).

We also know that:

- \(C_{{\sf A},\text{in}} = 2.5\) mol/L,
- \(T_\text{in} = 288\) K
- a constant volume of \(V = 100\) L
- constant inlet flow of \(F(t) = 20.1\) L/min, though we could easily make it a function of time and adjust the function
`tank`

to use the time variable,`t`

- molar heat capacity is a weak function of the system temperature: \(C_p(T) = 4.184 - 0.002(T-273)\) J/(kg.K),
- liquid phase density is \(\rho=1.05\) kg/L.
- \(E_a = 5000\) J/mol
- \(R = 8.314\) J/(mol.K),

We need two initial conditions, one for each of the ODEs:

- \(C_{\sf A}(t=0) = 0.5\) mol/L
- \(T(t=0) = 295\) K

In the code we will integrate the ODE from \(t_\text{start} = 0.0\) minutes up to \(t_\text{final} = 45.0\) minutes

and plot the time-varying trajectories of temperature and concentration in the tank using `matplotlib`

.

Read the SciPy documentation for `ode`

.