Difference between revisions of "Software tutorial/Integration of ODEs"
m |
m |
||
| Line 1: | Line 1: | ||
In our course we have learned several ways of integrating a single equation \(\displaystyle \frac{dy(t)}{dt} = f(t, y) \) with a given initial condition \(y(t=0)=y_0\); and we have shown that a system of \(n\) first order ODE's can be integrated: \[ \displaystyle \frac{d{\bf y}(t)}{dt} = {\bf f}(t, {\bf y}) \]given a set of \(n\) initial conditions \({\bf y}(t=0) = {\bf y}_0\), i.e. \( {\bf y}_0\) is an \(n \times 1\) vector. The [[Ordinary_differential_equations | course notes covered]] Euler's method, Heun's method and Runge-Kutta methods. | In our course we have learned several ways of integrating a single equation \(\displaystyle \frac{dy(t)}{dt} = f(t, y) \) with a given initial condition \(y(t=0)=y_0\); and we have shown that a system of \(n\) first order ODE's can be integrated: \[ \displaystyle \frac{d{\bf y}(t)}{dt} = {\bf f}(t, {\bf y}) \]given a set of \(n\) initial conditions \({\bf y}(t=0) = {\bf y}_0\), i.e. \( {\bf y}_0\) is an \(n \times 1\) vector. The [[Ordinary_differential_equations | course notes covered]] Euler's method, Heun's method and Runge-Kutta methods. | ||
However, we only coded Euler's method (because it was simple!), but not the others. These other methods have been (reliably) coded in software packages and sophisticated error correction tools built into their algorithms. You should always use these toolbox functions to | However, we only coded Euler's method (because it was simple!), but not the others. These other methods have been (reliably) coded in software packages and sophisticated error correction tools built into their algorithms. You should always use these toolbox functions to integrate your differential equation models. In this section we will show how to use MATLAB and Python's built-in functions to integrate: | ||
* a single differential equation | * a single differential equation | ||
* a system of differential and algebraic equations. | * a system of differential and algebraic equations. | ||
| Line 17: | Line 17: | ||
== Example problem == | == Example problem == | ||
The example we will work with is a common modelling reaction: a liquid-based stirred tank reactor, with (initially) constant physical properties, a second order chemical reaction where species A is converted to B according to \({\sf A} \rightarrow {\sf B} \), with reaction rate \(r = k C_{\sf A}^2\). | The example we will work with is a common modelling reaction: a liquid-based stirred tank reactor, with (initially) constant physical properties, a second order chemical reaction where species A is converted to B according to \({\sf A} \rightarrow {\sf B} \), with reaction rate \(r = k C_{\sf A}^2\). One can find the time-dependent mass balance for this system to be: | ||
\[ \frac{dC_{\sf A}}{dt} = \frac{F(t)}{V} \bigg( C_{A,in} - C_{\sf A} \bigg) - k C_{\sf A}^2 \] | |||
{| class="wikitable" | {| class="wikitable" | ||
|- | |- | ||
Revision as of 14:51, 18 November 2010
In our course we have learned several ways of integrating a single equation \(\displaystyle \frac{dy(t)}{dt} = f(t, y) \) with a given initial condition \(y(t=0)=y_0\); and we have shown that a system of \(n\) first order ODE's can be integrated: \[ \displaystyle \frac{d{\bf y}(t)}{dt} = {\bf f}(t, {\bf y}) \]given a set of \(n\) initial conditions \({\bf y}(t=0) = {\bf y}_0\), i.e. \( {\bf y}_0\) is an \(n \times 1\) vector. The course notes covered Euler's method, Heun's method and Runge-Kutta methods.
However, we only coded Euler's method (because it was simple!), but not the others. These other methods have been (reliably) coded in software packages and sophisticated error correction tools built into their algorithms. You should always use these toolbox functions to integrate your differential equation models. In this section we will show how to use MATLAB and Python's built-in functions to integrate:
- a single differential equation
- a system of differential and algebraic equations.
We will only look at initial value problems (IVPs) in this tutorial.
MATLAB
MATLAB's has several ODE solvers available. You can read up more about the implementation details in this technical document.
Python
Like MATLAB, several integrators are available in Python. The integrator I will use in this tutorial is one of the most recent additions to SciPy - the VODE integrator developed at Lawrence Livermore National Laboratories in 1988. It is a good general-purpose solver for both stiff and non-stiff systems.
Example problem
The example we will work with is a common modelling reaction: a liquid-based stirred tank reactor, with (initially) constant physical properties, a second order chemical reaction where species A is converted to B according to \({\sf A} \rightarrow {\sf B} \), with reaction rate \(r = k C_{\sf A}^2\). One can find the time-dependent mass balance for this system to be:
\[ \frac{dC_{\sf A}}{dt} = \frac{F(t)}{V} \bigg( C_{A,in} - C_{\sf A} \bigg) - k C_{\sf A}^2 \]
| MATLAB | Python |
|---|---|