An important class of ordinary differential equations are so called stiff
system which are characterized by two or more time scales of different order.
Examples of such systems are found in chemical systems where reaction rates
of individual sub-reaction might differ over large ranges, for example:

The state type has to be a ublas::vector
and the matrix type must by a ublas::matrix
since the stiff integrator only accepts these types. However, you might want
use non-stiff integrators on this system, too - we will do so later for demonstration.
Therefore we want to use the same function also with other state_types, realized
by templatizing the operator():

Now you can use stiff_system
in combination with std::vector or boost::array.
In the example the explicit time derivative of f(x,t)
is introduced separately in the Jacobian. If df / dt = 0
simply fill dfdt with zeros.

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:

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.