• Home
  • Line#
  • Scopes#
  • Navigate#
  • Raw
  • Download
1<html>
2<head>
3<meta http-equiv="Content-Type" content="text/html; charset=UTF-8">
4<title>Using boost::units</title>
5<link rel="stylesheet" href="../../../../../../../doc/src/boostbook.css" type="text/css">
6<meta name="generator" content="DocBook XSL Stylesheets V1.79.1">
7<link rel="home" href="../../index.html" title="Chapter 1. Boost.Numeric.Odeint">
8<link rel="up" href="../tutorial.html" title="Tutorial">
9<link rel="prev" href="ensembles_of_oscillators.html" title="Ensembles of oscillators">
10<link rel="next" href="using_matrices_as_state_types.html" title="Using matrices as state types">
11</head>
12<body bgcolor="white" text="black" link="#0000FF" vlink="#840084" alink="#0000FF">
13<table cellpadding="2" width="100%"><tr>
14<td valign="top"><img alt="Boost C++ Libraries" width="277" height="86" src="../../logo.jpg"></td>
15<td align="center"><a href="../../../../../../../index.html">Home</a></td>
16<td align="center"><a href="../../../../../../../libs/libraries.htm">Libraries</a></td>
17<td align="center"><a href="http://www.boost.org/users/people.html">People</a></td>
18<td align="center"><a href="http://www.boost.org/users/faq.html">FAQ</a></td>
19<td align="center"><a href="../../../../../../../more/index.htm">More</a></td>
20</tr></table>
21<hr>
22<div class="spirit-nav">
23<a accesskey="p" href="ensembles_of_oscillators.html"><img src="../../../../../../../doc/src/images/prev.png" alt="Prev"></a><a accesskey="u" href="../tutorial.html"><img src="../../../../../../../doc/src/images/up.png" alt="Up"></a><a accesskey="h" href="../../index.html"><img src="../../../../../../../doc/src/images/home.png" alt="Home"></a><a accesskey="n" href="using_matrices_as_state_types.html"><img src="../../../../../../../doc/src/images/next.png" alt="Next"></a>
24</div>
25<div class="section">
26<div class="titlepage"><div><div><h3 class="title">
27<a name="boost_numeric_odeint.tutorial.using_boost__units"></a><a class="link" href="using_boost__units.html" title="Using boost::units">Using
28      boost::units</a>
29</h3></div></div></div>
30<p>
31        odeint also works well with <a href="http://www.boost.org/doc/libs/release/libs/units/" target="_top">Boost.Units</a>
32        - a library for compile time unit and dimension analysis. It works by decoding
33        unit information into the types of values. For a one-dimensional unit you
34        can just use the Boost.Unit types as state type, deriv type and time type
35        and hand the <code class="computeroutput"><span class="identifier">vector_space_algebra</span></code>
36        to the stepper definition and everything works just fine:
37      </p>
38<p>
39</p>
40<pre class="programlisting"><span class="keyword">typedef</span> <span class="identifier">units</span><span class="special">::</span><span class="identifier">quantity</span><span class="special">&lt;</span> <span class="identifier">si</span><span class="special">::</span><span class="identifier">time</span> <span class="special">,</span> <span class="keyword">double</span> <span class="special">&gt;</span> <span class="identifier">time_type</span><span class="special">;</span>
41<span class="keyword">typedef</span> <span class="identifier">units</span><span class="special">::</span><span class="identifier">quantity</span><span class="special">&lt;</span> <span class="identifier">si</span><span class="special">::</span><span class="identifier">length</span> <span class="special">,</span> <span class="keyword">double</span> <span class="special">&gt;</span> <span class="identifier">length_type</span><span class="special">;</span>
42<span class="keyword">typedef</span> <span class="identifier">units</span><span class="special">::</span><span class="identifier">quantity</span><span class="special">&lt;</span> <span class="identifier">si</span><span class="special">::</span><span class="identifier">velocity</span> <span class="special">,</span> <span class="keyword">double</span> <span class="special">&gt;</span> <span class="identifier">velocity_type</span><span class="special">;</span>
43
44<span class="keyword">typedef</span> <span class="identifier">runge_kutta4</span><span class="special">&lt;</span> <span class="identifier">length_type</span> <span class="special">,</span> <span class="keyword">double</span> <span class="special">,</span> <span class="identifier">velocity_type</span> <span class="special">,</span> <span class="identifier">time_type</span> <span class="special">,</span>
45                      <span class="identifier">vector_space_algebra</span> <span class="special">&gt;</span> <span class="identifier">stepper_type</span><span class="special">;</span>
46</pre>
47<p>
48      </p>
49<p>
50        If you want to solve more-dimensional problems the individual entries typically
51        have different units. That means that the <code class="computeroutput"><span class="identifier">state_type</span></code>
52        is now possibly heterogeneous, meaning that every entry might have a different
53        type. To solve this problem, compile-time sequences from <a href="http://www.boost.org/doc/libs/release/libs/fusion/" target="_top">Boost.Fusion</a>
54        can be used.
55      </p>
56<p>
57        To illustrate how odeint works with <a href="http://www.boost.org/doc/libs/release/libs/units/" target="_top">Boost.Units</a>
58        we use the harmonic oscillator as primary example. We start with defining
59        all quantities
60      </p>
61<p>
62</p>
63<pre class="programlisting"><span class="preprocessor">#include</span> <span class="special">&lt;</span><span class="identifier">boost</span><span class="special">/</span><span class="identifier">numeric</span><span class="special">/</span><span class="identifier">odeint</span><span class="special">.</span><span class="identifier">hpp</span><span class="special">&gt;</span>
64<span class="preprocessor">#include</span> <span class="special">&lt;</span><span class="identifier">boost</span><span class="special">/</span><span class="identifier">numeric</span><span class="special">/</span><span class="identifier">odeint</span><span class="special">/</span><span class="identifier">algebra</span><span class="special">/</span><span class="identifier">fusion_algebra</span><span class="special">.</span><span class="identifier">hpp</span><span class="special">&gt;</span>
65<span class="preprocessor">#include</span> <span class="special">&lt;</span><span class="identifier">boost</span><span class="special">/</span><span class="identifier">numeric</span><span class="special">/</span><span class="identifier">odeint</span><span class="special">/</span><span class="identifier">algebra</span><span class="special">/</span><span class="identifier">fusion_algebra_dispatcher</span><span class="special">.</span><span class="identifier">hpp</span><span class="special">&gt;</span>
66
67<span class="preprocessor">#include</span> <span class="special">&lt;</span><span class="identifier">boost</span><span class="special">/</span><span class="identifier">units</span><span class="special">/</span><span class="identifier">systems</span><span class="special">/</span><span class="identifier">si</span><span class="special">/</span><span class="identifier">length</span><span class="special">.</span><span class="identifier">hpp</span><span class="special">&gt;</span>
68<span class="preprocessor">#include</span> <span class="special">&lt;</span><span class="identifier">boost</span><span class="special">/</span><span class="identifier">units</span><span class="special">/</span><span class="identifier">systems</span><span class="special">/</span><span class="identifier">si</span><span class="special">/</span><span class="identifier">time</span><span class="special">.</span><span class="identifier">hpp</span><span class="special">&gt;</span>
69<span class="preprocessor">#include</span> <span class="special">&lt;</span><span class="identifier">boost</span><span class="special">/</span><span class="identifier">units</span><span class="special">/</span><span class="identifier">systems</span><span class="special">/</span><span class="identifier">si</span><span class="special">/</span><span class="identifier">velocity</span><span class="special">.</span><span class="identifier">hpp</span><span class="special">&gt;</span>
70<span class="preprocessor">#include</span> <span class="special">&lt;</span><span class="identifier">boost</span><span class="special">/</span><span class="identifier">units</span><span class="special">/</span><span class="identifier">systems</span><span class="special">/</span><span class="identifier">si</span><span class="special">/</span><span class="identifier">acceleration</span><span class="special">.</span><span class="identifier">hpp</span><span class="special">&gt;</span>
71<span class="preprocessor">#include</span> <span class="special">&lt;</span><span class="identifier">boost</span><span class="special">/</span><span class="identifier">units</span><span class="special">/</span><span class="identifier">systems</span><span class="special">/</span><span class="identifier">si</span><span class="special">/</span><span class="identifier">io</span><span class="special">.</span><span class="identifier">hpp</span><span class="special">&gt;</span>
72
73<span class="preprocessor">#include</span> <span class="special">&lt;</span><span class="identifier">boost</span><span class="special">/</span><span class="identifier">fusion</span><span class="special">/</span><span class="identifier">container</span><span class="special">.</span><span class="identifier">hpp</span><span class="special">&gt;</span>
74
75<span class="keyword">using</span> <span class="keyword">namespace</span> <span class="identifier">std</span><span class="special">;</span>
76<span class="keyword">using</span> <span class="keyword">namespace</span> <span class="identifier">boost</span><span class="special">::</span><span class="identifier">numeric</span><span class="special">::</span><span class="identifier">odeint</span><span class="special">;</span>
77<span class="keyword">namespace</span> <span class="identifier">fusion</span> <span class="special">=</span> <span class="identifier">boost</span><span class="special">::</span><span class="identifier">fusion</span><span class="special">;</span>
78<span class="keyword">namespace</span> <span class="identifier">units</span> <span class="special">=</span> <span class="identifier">boost</span><span class="special">::</span><span class="identifier">units</span><span class="special">;</span>
79<span class="keyword">namespace</span> <span class="identifier">si</span> <span class="special">=</span> <span class="identifier">boost</span><span class="special">::</span><span class="identifier">units</span><span class="special">::</span><span class="identifier">si</span><span class="special">;</span>
80
81<span class="keyword">typedef</span> <span class="identifier">units</span><span class="special">::</span><span class="identifier">quantity</span><span class="special">&lt;</span> <span class="identifier">si</span><span class="special">::</span><span class="identifier">time</span> <span class="special">,</span> <span class="keyword">double</span> <span class="special">&gt;</span> <span class="identifier">time_type</span><span class="special">;</span>
82<span class="keyword">typedef</span> <span class="identifier">units</span><span class="special">::</span><span class="identifier">quantity</span><span class="special">&lt;</span> <span class="identifier">si</span><span class="special">::</span><span class="identifier">length</span> <span class="special">,</span> <span class="keyword">double</span> <span class="special">&gt;</span> <span class="identifier">length_type</span><span class="special">;</span>
83<span class="keyword">typedef</span> <span class="identifier">units</span><span class="special">::</span><span class="identifier">quantity</span><span class="special">&lt;</span> <span class="identifier">si</span><span class="special">::</span><span class="identifier">velocity</span> <span class="special">,</span> <span class="keyword">double</span> <span class="special">&gt;</span> <span class="identifier">velocity_type</span><span class="special">;</span>
84<span class="keyword">typedef</span> <span class="identifier">units</span><span class="special">::</span><span class="identifier">quantity</span><span class="special">&lt;</span> <span class="identifier">si</span><span class="special">::</span><span class="identifier">acceleration</span> <span class="special">,</span> <span class="keyword">double</span> <span class="special">&gt;</span> <span class="identifier">acceleration_type</span><span class="special">;</span>
85<span class="keyword">typedef</span> <span class="identifier">units</span><span class="special">::</span><span class="identifier">quantity</span><span class="special">&lt;</span> <span class="identifier">si</span><span class="special">::</span><span class="identifier">frequency</span> <span class="special">,</span> <span class="keyword">double</span> <span class="special">&gt;</span> <span class="identifier">frequency_type</span><span class="special">;</span>
86
87<span class="keyword">typedef</span> <span class="identifier">fusion</span><span class="special">::</span><span class="identifier">vector</span><span class="special">&lt;</span> <span class="identifier">length_type</span> <span class="special">,</span> <span class="identifier">velocity_type</span> <span class="special">&gt;</span> <span class="identifier">state_type</span><span class="special">;</span>
88<span class="keyword">typedef</span> <span class="identifier">fusion</span><span class="special">::</span><span class="identifier">vector</span><span class="special">&lt;</span> <span class="identifier">velocity_type</span> <span class="special">,</span> <span class="identifier">acceleration_type</span> <span class="special">&gt;</span> <span class="identifier">deriv_type</span><span class="special">;</span>
89</pre>
90<p>
91      </p>
92<p>
93        Note, that the <code class="computeroutput"><span class="identifier">state_type</span></code>
94        and the <code class="computeroutput"><span class="identifier">deriv_type</span></code> are now
95        a compile-time fusion sequences. <code class="computeroutput"><span class="identifier">deriv_type</span></code>
96        represents <span class="emphasis"><em>x'</em></span> and is now different from the state type
97        as it has different unit definitions. Next, we define the ordinary differential
98        equation which is completely equivalent to the example in <a class="link" href="harmonic_oscillator.html" title="Harmonic oscillator">Harmonic
99        Oscillator</a>:
100      </p>
101<p>
102</p>
103<pre class="programlisting"><span class="keyword">struct</span> <span class="identifier">oscillator</span>
104<span class="special">{</span>
105    <span class="identifier">frequency_type</span> <span class="identifier">m_omega</span><span class="special">;</span>
106
107    <span class="identifier">oscillator</span><span class="special">(</span> <span class="keyword">const</span> <span class="identifier">frequency_type</span> <span class="special">&amp;</span><span class="identifier">omega</span> <span class="special">=</span> <span class="number">1.0</span> <span class="special">*</span> <span class="identifier">si</span><span class="special">::</span><span class="identifier">hertz</span> <span class="special">)</span> <span class="special">:</span> <span class="identifier">m_omega</span><span class="special">(</span> <span class="identifier">omega</span> <span class="special">)</span> <span class="special">{</span> <span class="special">}</span>
108
109    <span class="keyword">void</span> <span class="keyword">operator</span><span class="special">()(</span> <span class="keyword">const</span> <span class="identifier">state_type</span> <span class="special">&amp;</span><span class="identifier">x</span> <span class="special">,</span> <span class="identifier">deriv_type</span> <span class="special">&amp;</span><span class="identifier">dxdt</span> <span class="special">,</span> <span class="identifier">time_type</span> <span class="identifier">t</span> <span class="special">)</span> <span class="keyword">const</span>
110    <span class="special">{</span>
111        <span class="identifier">fusion</span><span class="special">::</span><span class="identifier">at_c</span><span class="special">&lt;</span> <span class="number">0</span> <span class="special">&gt;(</span> <span class="identifier">dxdt</span> <span class="special">)</span> <span class="special">=</span> <span class="identifier">fusion</span><span class="special">::</span><span class="identifier">at_c</span><span class="special">&lt;</span> <span class="number">1</span> <span class="special">&gt;(</span> <span class="identifier">x</span> <span class="special">);</span>
112        <span class="identifier">fusion</span><span class="special">::</span><span class="identifier">at_c</span><span class="special">&lt;</span> <span class="number">1</span> <span class="special">&gt;(</span> <span class="identifier">dxdt</span> <span class="special">)</span> <span class="special">=</span> <span class="special">-</span> <span class="identifier">m_omega</span> <span class="special">*</span> <span class="identifier">m_omega</span> <span class="special">*</span> <span class="identifier">fusion</span><span class="special">::</span><span class="identifier">at_c</span><span class="special">&lt;</span> <span class="number">0</span> <span class="special">&gt;(</span> <span class="identifier">x</span> <span class="special">);</span>
113    <span class="special">}</span>
114<span class="special">};</span>
115</pre>
116<p>
117      </p>
118<p>
119        Next, we instantiate an appropriate stepper. We must explicitly parametrize
120        the stepper with the <code class="computeroutput"><span class="identifier">state_type</span></code>,
121        <code class="computeroutput"><span class="identifier">deriv_type</span></code>, <code class="computeroutput"><span class="identifier">time_type</span></code>.
122      </p>
123<p>
124</p>
125<pre class="programlisting"><span class="keyword">typedef</span> <span class="identifier">runge_kutta_dopri5</span><span class="special">&lt;</span> <span class="identifier">state_type</span> <span class="special">,</span> <span class="keyword">double</span> <span class="special">,</span> <span class="identifier">deriv_type</span> <span class="special">,</span> <span class="identifier">time_type</span> <span class="special">&gt;</span> <span class="identifier">stepper_type</span><span class="special">;</span>
126
127<span class="identifier">state_type</span> <span class="identifier">x</span><span class="special">(</span> <span class="number">1.0</span> <span class="special">*</span> <span class="identifier">si</span><span class="special">::</span><span class="identifier">meter</span> <span class="special">,</span> <span class="number">0.0</span> <span class="special">*</span> <span class="identifier">si</span><span class="special">::</span><span class="identifier">meter_per_second</span> <span class="special">);</span>
128
129<span class="identifier">integrate_const</span><span class="special">(</span> <span class="identifier">make_dense_output</span><span class="special">(</span> <span class="number">1.0e-6</span> <span class="special">,</span> <span class="number">1.0e-6</span> <span class="special">,</span> <span class="identifier">stepper_type</span><span class="special">()</span> <span class="special">)</span> <span class="special">,</span> <span class="identifier">oscillator</span><span class="special">(</span> <span class="number">2.0</span> <span class="special">*</span> <span class="identifier">si</span><span class="special">::</span><span class="identifier">hertz</span> <span class="special">)</span> <span class="special">,</span>
130                 <span class="identifier">x</span> <span class="special">,</span> <span class="number">0.0</span> <span class="special">*</span> <span class="identifier">si</span><span class="special">::</span><span class="identifier">second</span> <span class="special">,</span> <span class="number">100.0</span> <span class="special">*</span> <span class="identifier">si</span><span class="special">::</span><span class="identifier">second</span> <span class="special">,</span> <span class="number">0.1</span> <span class="special">*</span> <span class="identifier">si</span><span class="special">::</span><span class="identifier">second</span> <span class="special">,</span> <span class="identifier">streaming_observer</span><span class="special">(</span> <span class="identifier">cout</span> <span class="special">)</span> <span class="special">);</span>
131</pre>
132<p>
133      </p>
134<div class="note"><table border="0" summary="Note">
135<tr>
136<td rowspan="2" align="center" valign="top" width="25"><img alt="[Note]" src="../../../../../../../doc/src/images/note.png"></td>
137<th align="left">Note</th>
138</tr>
139<tr><td align="left" valign="top"><p>
140          When using compile-time sequences, the iteration over vector elements is
141          done by the <code class="computeroutput"><span class="identifier">fusion_algebra</span></code>,
142          which is automatically chosen by odeint. For more on the state types /
143          algebras see chapter <a class="link" href="../odeint_in_detail/state_types__algebras_and_operations.html" title="State types, algebras and operations">Adapt
144          your own state types</a>.
145        </p></td></tr>
146</table></div>
147<p>
148        It is quite easy but the compilation time might take very long. Furthermore,
149        the observer is defined a bit different
150      </p>
151<p>
152</p>
153<pre class="programlisting"><span class="keyword">struct</span> <span class="identifier">streaming_observer</span>
154<span class="special">{</span>
155    <span class="identifier">std</span><span class="special">::</span><span class="identifier">ostream</span><span class="special">&amp;</span> <span class="identifier">m_out</span><span class="special">;</span>
156
157    <span class="identifier">streaming_observer</span><span class="special">(</span> <span class="identifier">std</span><span class="special">::</span><span class="identifier">ostream</span> <span class="special">&amp;</span><span class="identifier">out</span> <span class="special">)</span> <span class="special">:</span> <span class="identifier">m_out</span><span class="special">(</span> <span class="identifier">out</span> <span class="special">)</span> <span class="special">{</span> <span class="special">}</span>
158
159    <span class="keyword">struct</span> <span class="identifier">write_element</span>
160    <span class="special">{</span>
161        <span class="identifier">std</span><span class="special">::</span><span class="identifier">ostream</span> <span class="special">&amp;</span><span class="identifier">m_out</span><span class="special">;</span>
162        <span class="identifier">write_element</span><span class="special">(</span> <span class="identifier">std</span><span class="special">::</span><span class="identifier">ostream</span> <span class="special">&amp;</span><span class="identifier">out</span> <span class="special">)</span> <span class="special">:</span> <span class="identifier">m_out</span><span class="special">(</span> <span class="identifier">out</span> <span class="special">)</span> <span class="special">{</span> <span class="special">};</span>
163
164        <span class="keyword">template</span><span class="special">&lt;</span> <span class="keyword">class</span> <span class="identifier">T</span> <span class="special">&gt;</span>
165        <span class="keyword">void</span> <span class="keyword">operator</span><span class="special">()(</span> <span class="keyword">const</span> <span class="identifier">T</span> <span class="special">&amp;</span><span class="identifier">t</span> <span class="special">)</span> <span class="keyword">const</span>
166        <span class="special">{</span>
167            <span class="identifier">m_out</span> <span class="special">&lt;&lt;</span> <span class="string">"\t"</span> <span class="special">&lt;&lt;</span> <span class="identifier">t</span><span class="special">;</span>
168        <span class="special">}</span>
169    <span class="special">};</span>
170
171    <span class="keyword">template</span><span class="special">&lt;</span> <span class="keyword">class</span> <span class="identifier">State</span> <span class="special">,</span> <span class="keyword">class</span> <span class="identifier">Time</span> <span class="special">&gt;</span>
172    <span class="keyword">void</span> <span class="keyword">operator</span><span class="special">()(</span> <span class="keyword">const</span> <span class="identifier">State</span> <span class="special">&amp;</span><span class="identifier">x</span> <span class="special">,</span> <span class="keyword">const</span> <span class="identifier">Time</span> <span class="special">&amp;</span><span class="identifier">t</span> <span class="special">)</span> <span class="keyword">const</span>
173    <span class="special">{</span>
174        <span class="identifier">m_out</span> <span class="special">&lt;&lt;</span> <span class="identifier">t</span><span class="special">;</span>
175        <span class="identifier">fusion</span><span class="special">::</span><span class="identifier">for_each</span><span class="special">(</span> <span class="identifier">x</span> <span class="special">,</span> <span class="identifier">write_element</span><span class="special">(</span> <span class="identifier">m_out</span> <span class="special">)</span> <span class="special">);</span>
176        <span class="identifier">m_out</span> <span class="special">&lt;&lt;</span> <span class="string">"\n"</span><span class="special">;</span>
177    <span class="special">}</span>
178<span class="special">};</span>
179</pre>
180<p>
181      </p>
182<div class="caution"><table border="0" summary="Caution">
183<tr>
184<td rowspan="2" align="center" valign="top" width="25"><img alt="[Caution]" src="../../../../../../../doc/src/images/caution.png"></td>
185<th align="left">Caution</th>
186</tr>
187<tr><td align="left" valign="top"><p>
188          Using <a href="http://www.boost.org/doc/libs/release/libs/units/" target="_top">Boost.Units</a>
189          works nicely but compilation can be very time and memory consuming. For
190          example the unit test for the usage of <a href="http://www.boost.org/doc/libs/release/libs/units/" target="_top">Boost.Units</a>
191          in odeint take up to 4 GB of memory at compilation.
192        </p></td></tr>
193</table></div>
194<p>
195        The full cpp file for this example can be found here <a href="https://github.com/headmyshoulder/odeint-v2/blob/master/examples/harmonic_oscillator_units.cpp" target="_top">harmonic_oscillator_units.cpp</a>.
196      </p>
197</div>
198<table xmlns:rev="http://www.cs.rpi.edu/~gregod/boost/tools/doc/revision" width="100%"><tr>
199<td align="left"></td>
200<td align="right"><div class="copyright-footer">Copyright © 2009-2015 Karsten Ahnert and Mario Mulansky<p>
201        Distributed under the Boost Software License, Version 1.0. (See accompanying
202        file LICENSE_1_0.txt or copy at <a href="http://www.boost.org/LICENSE_1_0.txt" target="_top">http://www.boost.org/LICENSE_1_0.txt</a>)
203      </p>
204</div></td>
205</tr></table>
206<hr>
207<div class="spirit-nav">
208<a accesskey="p" href="ensembles_of_oscillators.html"><img src="../../../../../../../doc/src/images/prev.png" alt="Prev"></a><a accesskey="u" href="../tutorial.html"><img src="../../../../../../../doc/src/images/up.png" alt="Up"></a><a accesskey="h" href="../../index.html"><img src="../../../../../../../doc/src/images/home.png" alt="Home"></a><a accesskey="n" href="using_matrices_as_state_types.html"><img src="../../../../../../../doc/src/images/next.png" alt="Next"></a>
209</div>
210</body>
211</html>
212