ceph/src/boost/libs/numeric/odeint/doc/getting_started.qbk

   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 Getting started]
  16
  17 [section Overview]
  18
  19 odeint is a library for solving  initial value problems (IVP) of ordinary
  20 differential equations. Mathematically, these problems are formulated as
  21 follows:
  22
  23 ['x'(t) = f(x,t)], ['x(0) = x0].
  24
  25 ['x] and ['f] can be vectors and the solution is some function ['x(t)] fulfilling both equations above. In the following we will refer to ['x'(t)] also `dxdt` which is also our notation for the derivative in the source code.
  26
  27 Ordinary differential equations occur nearly everywhere in natural sciences. For example, the whole Newtonian mechanics are described by second order differential equations. Be sure, you will find them in every discipline. They also occur if partial differential equations (PDEs) are discretized. Then, a system of coupled ordinary differential occurs, sometimes also referred as lattices ODEs.
  28
  29 Numerical approximations for the solution ['x(t)] are calculated iteratively. The easiest algorithm is the Euler scheme, where starting at ['x(0)] one finds ['x(dt) = x(0) + dt f(x(0),0)]. Now one can use ['x(dt)] and obtain ['x(2dt)] in a similar way and so on. The Euler method is of order 1, that means the error at each step is ['~ dt[super 2]]. This is, of course, not very satisfying, which is why the Euler method is rarely used for real life problems and serves just as illustrative example.
  30
  31 The main focus of odeint is to provide numerical methods implemented in a way where the algorithm is completely independent on the data structure used to represent the state /x/.
  32 In doing so, odeint is applicable for a broad variety of situations and it can be used with many other libraries. Besides the usual case where the state is defined as a `std::vector` or a `boost::array`, we provide native support for the following libraries:
  33
  34 * __ublas
  35 * __thrust, making odeint naturally running on CUDA devices
  36 * gsl_vector for compatibility with the many numerical function in the GSL
  37 * __boost_range
  38 * __boost_fusion (the state type can be a fusion vector)
  39 * __boost_units
  40 * __intel_mkl for maximum performance
  41 * __vexcl for OpenCL
  42 * __boost_graph (still experimentally)
  43
  44 In odeint, the following algorithms are implemented:
  45
  46 [include stepper_table.qbk]
  47
  48 [endsect]
  49
  50
  51
  52 [section Usage, Compilation, Headers]
  53
  54 odeint is a header-only library, no linking against pre-compiled code is required. It can be included by
  55
  56 ``
  57 #include <boost/numeric/odeint.hpp>
  58 ``
  59
  60 which includes all headers of the library. All functions and classes from odeint live in the namespace
  61
  62 ``
  63  using namespace boost::numeric::odeint;
  64 ``
  65
  66 It is also possible to include only parts of the library. This is the recommended way since it saves a lot of compilation time.
  67
  68 * `#include <boost/numeric/odeint/stepper/XYZ.hpp>` - the include path for all steppers, XYZ is a placeholder for a stepper.
  69 * `#include <boost/numeric/odeint/algebra/XYZ.hpp>` - all algebras.
  70 * `#include <boost/numeric/odeint/util/XYZ.hpp>` - the utility functions like `is_resizeable`, `same_size`, or `resize`.
  71 * `#include <boost/numeric/odeint/integrate/XYZ.hpp>` - the integrate routines.
  72 * `#include <boost/numeric/odeint/iterator/XYZ.hpp>` - the range and iterator functions.
  73 * `#include <boost/numeric/odeint/external/XYZ.hpp>` - any binders to external libraries.
  74
  75
  76
  77 [endsect]
  78
  79 [section Short Example]
  80
  81 Imaging, you want to numerically integrate a harmonic oscillator with
  82 friction. The equations of motion are given by ['x'' = -x + __gamma x'].
  83 Odeint only deals with first order ODEs that have no higher derivatives than
  84 x' involved. However, any higher order ODE can be transformed to a system of
  85 first order ODEs by introducing the new variables ['q=x] and ['p=x'] such that ['w=(q,p)]. To apply numerical integration one first has to design the right hand side of the equation ['w' = f(w) = (p,-q+__gamma p)]:
  86
  87 [import ../examples/harmonic_oscillator.cpp]
  88 [rhs_function]
  89
  90 Here we chose `vector<double>` as the state type, but others are also
  91 possible, for example `boost::array<double,2>`. odeint is designed in such a
  92 way that you can easily use your own state types. Next, the ODE is defined
  93 which is in this case a simple function calculating ['f(x)]. The parameter
  94 signature of this function is crucial: the integration methods will always
  95 call them in the form `f(x, dxdt, t)` (there are exceptions for some special routines). So, even if there is no explicit time dependence, one has to define `t` as a function parameter.
  96
  97 Now, we have to define the initial state from which the integration should start:
  98
  99 [state_initialization]
 100
 101 For the integration itself we'll use the `integrate` function, which is a  convenient way to get quick results. It is based on the error-controlled `runge_kutta54_cash_karp` stepper (5th order) and uses adaptive step-size.
 102
 103 [integration]
 104
 105 The integrate function expects as parameters the rhs of the ode as defined
 106 above, the initial state `x`, the start-and end-time of the integration as
 107 well as the initial time step=size. Note, that `integrate` uses an adaptive step-size during
 108 the integration steps so the time points will not be equally spaced. The
 109 integration returns the number of steps that were applied and updates x which
 110 is set to the approximate solution of the ODE at the end of integration.
 111
 112 It is also possible to represent the ode system as a class. The
 113 rhs must then be implemented as a functor - a class with an overloaded function call operator:
 114
 115 [rhs_class]
 116
 117 which can be used via
 118
 119 [integration_class]
 120
 121 In order to observe the solution
 122 during the integration steps all you have to do is
 123 to provide a reasonable observer. An example is
 124
 125 [integrate_observer]
 126
 127 which stores the intermediate steps in a container.
 128 Note, the argument structure of the ()-operator: odeint calls the observer
 129 exactly in this way, providing the current state and time. Now, you only have
 130 to pass this container to the integration function:
 131
 132 [integrate_observ]
 133
 134 That is all. You can use functional libraries like __boost_lambda or __boost_phoenix to ease the creation of observer functions.
 135
 136 The full cpp file for this example can be found here: [github_link examples/harmonic_oscillator.cpp harmonic_oscillator.cpp]
 137
 138
 139
 140 [endsect]
 141
 142 [endsect]