]>
Commit | Line | Data |
---|---|---|
7c673cae FG |
1 | [/============================================================================ |
2 | Boost.odeint | |
3 | ||
4 | Copyright 2011-2013 Karsten Ahnert | |
5 | Copyright 2011-2012 Mario Mulansky | |
6 | Copyright 2012 Sylwester Arabas | |
7 | ||
8 | Use, modification and distribution is subject to the Boost Software License, | |
9 | Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at | |
10 | http://www.boost.org/LICENSE_1_0.txt) | |
11 | =============================================================================/] | |
12 | ||
13 | ||
14 | ||
15 | [section Stiff systems] | |
16 | ||
17 | [import ../examples/stiff_system.cpp] | |
18 | ||
19 | An important class of ordinary differential equations are so called stiff | |
20 | system which are characterized by two or more time scales of different | |
21 | order. Examples of such systems are found in chemical systems where reaction | |
22 | rates of individual sub-reaction might differ over large ranges, for example: | |
23 | ||
24 | ['d S[subl 1] / dt = - 101 S[subl 2] - 100 S[subl 1]] | |
25 | ||
26 | ['d S[subl 2] / dt = S[subl 1]] | |
27 | ||
28 | ||
29 | In order to efficiently solve stiff systems numerically the Jacobian | |
30 | ||
31 | ['J = d f[subl i] / d x[subl j]] | |
32 | ||
33 | is needed. Here is the definition of the above example | |
34 | ||
35 | [stiff_system_definition] | |
36 | ||
37 | The state type has to be a `ublas::vector` and the matrix type must by a | |
38 | `ublas::matrix` since the stiff integrator only accepts these types. | |
39 | However, you might want use non-stiff integrators on this system, too - we will | |
40 | do so later for demonstration. Therefore we want to use the same function also | |
41 | with other state_types, realized by templatizing the `operator()`: | |
42 | ||
43 | [stiff_system_alternative_definition] | |
44 | ||
45 | Now you can use `stiff_system` in combination with `std::vector` or | |
46 | `boost::array`. In the example the explicit time derivative of ['f(x,t)] is | |
47 | introduced separately in the Jacobian. If ['df / dt = 0] simply fill `dfdt` with zeros. | |
48 | ||
49 | A well know solver for stiff systems is the Rosenbrock method. It has a step size control and dense output facilities and can be used like all the other steppers: | |
50 | ||
51 | [integrate_stiff_system] | |
52 | ||
53 | During the integration 71 steps have been done. Comparing to a classical Runge-Kutta solver this is a very good result. For example the Dormand-Prince 5 method with step size control and dense output yields 1531 steps. | |
54 | ||
55 | [integrate_stiff_system_alternative] | |
56 | ||
57 | Note, that we have used __boost_phoenix, a great functional programming library, to create and compose the observer. | |
58 | ||
59 | The full example can be found here: [github_link examples/stiff_system.cpp stiff_system.cpp] | |
60 | ||
61 | ||
62 | [endsect] |