libs/numeric/odeint/doc/html/boost_numeric_odeint/tutorial/using_opencl_via_vexcl.html

   1 <html>
   2 <head>
   3 <meta http-equiv="Content-Type" content="text/html; charset=US-ASCII">
   4 <title>Using OpenCL via VexCL</title>
   5 <link rel="stylesheet" href="../../../../../../../doc/src/boostbook.css" type="text/css">
   6 <meta name="generator" content="DocBook XSL Stylesheets V1.78.1">
   7 <link rel="home" href="../../index.html" title="Chapter&#160;1.&#160;Boost.Numeric.Odeint">
   8 <link rel="up" href="../tutorial.html" title="Tutorial">
   9 <link rel="prev" href="using_cuda__or_openmp__tbb_______via_thrust.html" title="Using CUDA (or OpenMP, TBB, ...) via Thrust">
  10 <link rel="next" href="parallel_computation_with_openmp_and_mpi.html" title="Parallel computation with OpenMP and MPI">
  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="using_cuda__or_openmp__tbb_______via_thrust.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="parallel_computation_with_openmp_and_mpi.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_opencl_via_vexcl"></a><a class="link" href="using_opencl_via_vexcl.html" title="Using OpenCL via VexCL">Using
  28       OpenCL via VexCL</a>
  29 </h3></div></div></div>
  30 <p>
  31         In the previous section the usage of odeint in combination with <a href="http://code.google.com/p/thrust/" target="_top">Thrust</a>
  32         was shown. In this section we show how one can use OpenCL with odeint. The
  33         point of odeint is not to implement its own low-level data structures and
  34         algorithms, but to use high level libraries doing this task. Here, we will
  35         use the <a href="https://github.com/ddemidov/vexcl" target="_top">VexCL</a> framework
  36         to use OpenCL. <a href="https://github.com/ddemidov/vexcl" target="_top">VexCL</a>
  37         is a nice library for general computations and it uses heavily expression
  38         templates. With the help of <a href="https://github.com/ddemidov/vexcl" target="_top">VexCL</a>
  39         it is possible to write very compact and expressive application.
  40       </p>
  41 <div class="note"><table border="0" summary="Note">
  42 <tr>
  43 <td rowspan="2" align="center" valign="top" width="25"><img alt="[Note]" src="../../../../../../../doc/src/images/note.png"></td>
  44 <th align="left">Note</th>
  45 </tr>
  46 <tr><td align="left" valign="top"><p>
  47           vexcl needs C++11 features! So you have to compile with C++11 support enabled.
  48         </p></td></tr>
  49 </table></div>
  50 <p>
  51         To use <a href="https://github.com/ddemidov/vexcl" target="_top">VexCL</a> one needs
  52         to include one additional header which includes the data-types and algorithms
  53         from vexcl and the adaption to odeint. Adaption to odeint means here only
  54         to adapt the resizing functionality of <a href="https://github.com/ddemidov/vexcl" target="_top">VexCL</a>
  55         to odeint.
  56       </p>
  57 <p>
  58 </p>
  59 <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">external</span><span class="special">/</span><span class="identifier">vexcl</span><span class="special">/</span><span class="identifier">vexcl</span><span class="special">.</span><span class="identifier">hpp</span><span class="special">&gt;</span>
  60 </pre>
  61 <p>
  62       </p>
  63 <p>
  64         To demonstrate the use of <a href="https://github.com/ddemidov/vexcl" target="_top">VexCL</a>
  65         we integrate an ensemble of Lorenz system. The example is very similar to
  66         the parameter study of the Lorenz system in the previous section except that
  67         we do not compute the Lyapunov exponents. Again, we vary the parameter R
  68         of the Lorenz system an solve a whole ensemble of Lorenz systems in parallel
  69         (each with a different parameter R). First, we define the state type and
  70         a vector type
  71       </p>
  72 <p>
  73 </p>
  74 <pre class="programlisting"><span class="keyword">typedef</span> <span class="identifier">vex</span><span class="special">::</span><span class="identifier">vector</span><span class="special">&lt;</span> <span class="keyword">double</span> <span class="special">&gt;</span>    <span class="identifier">vector_type</span><span class="special">;</span>
  75 <span class="keyword">typedef</span> <span class="identifier">vex</span><span class="special">::</span><span class="identifier">multivector</span><span class="special">&lt;</span> <span class="keyword">double</span><span class="special">,</span> <span class="number">3</span> <span class="special">&gt;</span> <span class="identifier">state_type</span><span class="special">;</span>
  76 </pre>
  77 <p>
  78       </p>
  79 <p>
  80         The <code class="computeroutput"><span class="identifier">vector_type</span></code> is used to
  81         represent the parameter R. The <code class="computeroutput"><span class="identifier">state_type</span></code>
  82         is a multi-vector of three sub vectors and is used to represent. The first
  83         component of this multi-vector represent all <code class="computeroutput"><span class="identifier">x</span></code>
  84         components of the Lorenz system, while the second all <code class="computeroutput"><span class="identifier">y</span></code>
  85         components and the third all <code class="computeroutput"><span class="identifier">z</span></code>
  86         components. The components of this vector can be obtained via
  87       </p>
  88 <p>
  89 </p>
  90 <pre class="programlisting"><span class="keyword">auto</span> <span class="special">&amp;</span><span class="identifier">x</span> <span class="special">=</span> <span class="identifier">X</span><span class="special">(</span><span class="number">0</span><span class="special">);</span>
  91 <span class="keyword">auto</span> <span class="special">&amp;</span><span class="identifier">y</span> <span class="special">=</span> <span class="identifier">X</span><span class="special">(</span><span class="number">1</span><span class="special">);</span>
  92 <span class="keyword">auto</span> <span class="special">&amp;</span><span class="identifier">z</span> <span class="special">=</span> <span class="identifier">X</span><span class="special">(</span><span class="number">2</span><span class="special">);</span>
  93 </pre>
  94 <p>
  95       </p>
  96 <p>
  97         As already mentioned <a href="https://github.com/ddemidov/vexcl" target="_top">VexCL</a>
  98         supports expression templates and we will use them to implement the system
  99         function for the Lorenz ensemble:
 100       </p>
 101 <p>
 102 </p>
 103 <pre class="programlisting"><span class="keyword">const</span> <span class="keyword">double</span> <span class="identifier">sigma</span> <span class="special">=</span> <span class="number">10.0</span><span class="special">;</span>
 104 <span class="keyword">const</span> <span class="keyword">double</span> <span class="identifier">b</span> <span class="special">=</span> <span class="number">8.0</span> <span class="special">/</span> <span class="number">3.0</span><span class="special">;</span>
 105
 106 <span class="keyword">struct</span> <span class="identifier">sys_func</span>
 107 <span class="special">{</span>
 108     <span class="keyword">const</span> <span class="identifier">vector_type</span> <span class="special">&amp;</span><span class="identifier">R</span><span class="special">;</span>
 109
 110     <span class="identifier">sys_func</span><span class="special">(</span> <span class="keyword">const</span> <span class="identifier">vector_type</span> <span class="special">&amp;</span><span class="identifier">_R</span> <span class="special">)</span> <span class="special">:</span> <span class="identifier">R</span><span class="special">(</span> <span class="identifier">_R</span> <span class="special">)</span> <span class="special">{</span> <span class="special">}</span>
 111
 112     <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">state_type</span> <span class="special">&amp;</span><span class="identifier">dxdt</span> <span class="special">,</span> <span class="keyword">double</span> <span class="identifier">t</span> <span class="special">)</span> <span class="keyword">const</span>
 113     <span class="special">{</span>
 114         <span class="identifier">dxdt</span><span class="special">(</span><span class="number">0</span><span class="special">)</span> <span class="special">=</span> <span class="special">-</span><span class="identifier">sigma</span> <span class="special">*</span> <span class="special">(</span> <span class="identifier">x</span><span class="special">(</span><span class="number">0</span><span class="special">)</span> <span class="special">-</span> <span class="identifier">x</span><span class="special">(</span><span class="number">1</span><span class="special">)</span> <span class="special">);</span>
 115         <span class="identifier">dxdt</span><span class="special">(</span><span class="number">1</span><span class="special">)</span> <span class="special">=</span> <span class="identifier">R</span> <span class="special">*</span> <span class="identifier">x</span><span class="special">(</span><span class="number">0</span><span class="special">)</span> <span class="special">-</span> <span class="identifier">x</span><span class="special">(</span><span class="number">1</span><span class="special">)</span> <span class="special">-</span> <span class="identifier">x</span><span class="special">(</span><span class="number">0</span><span class="special">)</span> <span class="special">*</span> <span class="identifier">x</span><span class="special">(</span><span class="number">2</span><span class="special">);</span>
 116         <span class="identifier">dxdt</span><span class="special">(</span><span class="number">2</span><span class="special">)</span> <span class="special">=</span> <span class="special">-</span> <span class="identifier">b</span> <span class="special">*</span> <span class="identifier">x</span><span class="special">(</span><span class="number">2</span><span class="special">)</span> <span class="special">+</span> <span class="identifier">x</span><span class="special">(</span><span class="number">0</span><span class="special">)</span> <span class="special">*</span> <span class="identifier">x</span><span class="special">(</span><span class="number">1</span><span class="special">);</span>
 117     <span class="special">}</span>
 118 <span class="special">};</span>
 119 </pre>
 120 <p>
 121       </p>
 122 <p>
 123         It's very easy, isn't it? These three little lines do all the computations
 124         for you. There is no need to write your own OpenCL kernels. <a href="https://github.com/ddemidov/vexcl" target="_top">VexCL</a>
 125         does everything for you. Next we have to write the main application. We initialize
 126         the vector of parameters (R) and the initial state. Note that <a href="https://github.com/ddemidov/vexcl" target="_top">VexCL</a>
 127         requires the <code class="computeroutput"><span class="identifier">vector_space_algebra</span></code>,
 128         but that is automatically deduced and configured by odeint internally, so
 129         we only have to specify the <code class="computeroutput"><span class="identifier">state_type</span></code>
 130         when instantiating the stepper and we are done:
 131       </p>
 132 <p>
 133 </p>
 134 <pre class="programlisting"><span class="comment">// setup the opencl context</span>
 135 <span class="identifier">vex</span><span class="special">::</span><span class="identifier">Context</span> <span class="identifier">ctx</span><span class="special">(</span> <span class="identifier">vex</span><span class="special">::</span><span class="identifier">Filter</span><span class="special">::</span><span class="identifier">Type</span><span class="special">(</span><span class="identifier">CL_DEVICE_TYPE_GPU</span><span class="special">)</span> <span class="special">);</span>
 136 <span class="identifier">std</span><span class="special">::</span><span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="identifier">ctx</span> <span class="special">&lt;&lt;</span> <span class="identifier">std</span><span class="special">::</span><span class="identifier">endl</span><span class="special">;</span>
 137
 138 <span class="comment">// set up number of system, time step and integration time</span>
 139 <span class="keyword">const</span> <span class="identifier">size_t</span> <span class="identifier">n</span> <span class="special">=</span> <span class="number">1024</span> <span class="special">*</span> <span class="number">1024</span><span class="special">;</span>
 140 <span class="keyword">const</span> <span class="keyword">double</span> <span class="identifier">dt</span> <span class="special">=</span> <span class="number">0.01</span><span class="special">;</span>
 141 <span class="keyword">const</span> <span class="keyword">double</span> <span class="identifier">t_max</span> <span class="special">=</span> <span class="number">1000.0</span><span class="special">;</span>
 142
 143 <span class="comment">// initialize R</span>
 144 <span class="keyword">double</span> <span class="identifier">Rmin</span> <span class="special">=</span> <span class="number">0.1</span> <span class="special">,</span> <span class="identifier">Rmax</span> <span class="special">=</span> <span class="number">50.0</span> <span class="special">,</span> <span class="identifier">dR</span> <span class="special">=</span> <span class="special">(</span> <span class="identifier">Rmax</span> <span class="special">-</span> <span class="identifier">Rmin</span> <span class="special">)</span> <span class="special">/</span> <span class="keyword">double</span><span class="special">(</span> <span class="identifier">n</span> <span class="special">-</span> <span class="number">1</span> <span class="special">);</span>
 145 <span class="identifier">std</span><span class="special">::</span><span class="identifier">vector</span><span class="special">&lt;</span><span class="keyword">double</span><span class="special">&gt;</span> <span class="identifier">x</span><span class="special">(</span> <span class="identifier">n</span> <span class="special">*</span> <span class="number">3</span> <span class="special">)</span> <span class="special">,</span> <span class="identifier">r</span><span class="special">(</span> <span class="identifier">n</span> <span class="special">);</span>
 146 <span class="keyword">for</span><span class="special">(</span> <span class="identifier">size_t</span> <span class="identifier">i</span><span class="special">=</span><span class="number">0</span> <span class="special">;</span> <span class="identifier">i</span><span class="special">&lt;</span><span class="identifier">n</span> <span class="special">;</span> <span class="special">++</span><span class="identifier">i</span> <span class="special">)</span> <span class="identifier">r</span><span class="special">[</span><span class="identifier">i</span><span class="special">]</span> <span class="special">=</span> <span class="identifier">Rmin</span> <span class="special">+</span> <span class="identifier">dR</span> <span class="special">*</span> <span class="keyword">double</span><span class="special">(</span> <span class="identifier">i</span> <span class="special">);</span>
 147 <span class="identifier">vector_type</span> <span class="identifier">R</span><span class="special">(</span> <span class="identifier">ctx</span><span class="special">.</span><span class="identifier">queue</span><span class="special">()</span> <span class="special">,</span> <span class="identifier">r</span> <span class="special">);</span>
 148
 149 <span class="comment">// initialize the state of the lorenz ensemble</span>
 150 <span class="identifier">state_type</span> <span class="identifier">X</span><span class="special">(</span><span class="identifier">ctx</span><span class="special">.</span><span class="identifier">queue</span><span class="special">(),</span> <span class="identifier">n</span><span class="special">);</span>
 151 <span class="identifier">X</span><span class="special">(</span><span class="number">0</span><span class="special">)</span> <span class="special">=</span> <span class="number">10.0</span><span class="special">;</span>
 152 <span class="identifier">X</span><span class="special">(</span><span class="number">1</span><span class="special">)</span> <span class="special">=</span> <span class="number">10.0</span><span class="special">;</span>
 153 <span class="identifier">X</span><span class="special">(</span><span class="number">2</span><span class="special">)</span> <span class="special">=</span> <span class="number">10.0</span><span class="special">;</span>
 154
 155 <span class="comment">// create a stepper</span>
 156 <span class="identifier">runge_kutta4</span><span class="special">&lt;</span> <span class="identifier">state_type</span> <span class="special">&gt;</span> <span class="identifier">stepper</span><span class="special">;</span>
 157
 158 <span class="comment">// solve the system</span>
 159 <span class="identifier">integrate_const</span><span class="special">(</span> <span class="identifier">stepper</span> <span class="special">,</span> <span class="identifier">sys_func</span><span class="special">(</span> <span class="identifier">R</span> <span class="special">)</span> <span class="special">,</span> <span class="identifier">X</span> <span class="special">,</span> <span class="number">0.0</span> <span class="special">,</span> <span class="identifier">t_max</span> <span class="special">,</span> <span class="identifier">dt</span> <span class="special">);</span>
 160 </pre>
 161 <p>
 162       </p>
 163 </div>
 164 <table xmlns:rev="http://www.cs.rpi.edu/~gregod/boost/tools/doc/revision" width="100%"><tr>
 165 <td align="left"></td>
 166 <td align="right"><div class="copyright-footer">Copyright &#169; 2009-2012 Karsten
 167       Ahnert and Mario Mulansky<p>
 168         Distributed under the Boost Software License, Version 1.0. (See accompanying
 169         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>)
 170       </p>
 171 </div></td>
 172 </tr></table>
 173 <hr>
 174 <div class="spirit-nav">
 175 <a accesskey="p" href="using_cuda__or_openmp__tbb_______via_thrust.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="parallel_computation_with_openmp_and_mpi.html"><img src="../../../../../../../doc/src/images/next.png" alt="Next"></a>
 176 </div>
 177 </body>
 178 </html>