← Back to branch summary

~ubuntu-branches/ubuntu/maverick/freecad/maverick

« back to all changes in this revision

Viewing changes to src/3rdParty/boost/numeric/bindings/atlas/cblas1.hpp

  • Committer: Bazaar Package Importer
  • Author(s): Teemu Ikonen
  • Date: 2009-07-16 18:37:41 UTC
  • Revision ID: james.westby@ubuntu.com-20090716183741-oww9kcxqrk991i1n
Tags: upstream-0.8.2237
Import upstream version 0.8.2237

Show diffs side-by-side

added added

removed removed

Lines of Context:
src/3rdParty/boost/numeric/bindings/atlas/cblas1.hpp
 
1
/*
 
2
 * 
 
3
 * Copyright (c) Kresimir Fresl 2002 
 
4
 *
 
5
 * Distributed under the Boost Software License, Version 1.0.
 
6
 * (See accompanying file LICENSE_1_0.txt or copy at
 
7
 * http://www.boost.org/LICENSE_1_0.txt)
 
8
 *
 
9
 * Author acknowledges the support of the Faculty of Civil Engineering, 
 
10
 * University of Zagreb, Croatia.
 
11
 *
 
12
 */
 
13
 
 
14
#ifndef BOOST_NUMERIC_BINDINGS_CBLAS_LEVEL_1_HPP
 
15
#define BOOST_NUMERIC_BINDINGS_CBLAS_LEVEL_1_HPP
 
16
 
 
17
#include <cassert>
 
18
 
 
19
#include <boost/numeric/bindings/traits/type_traits.hpp>
 
20
#include <boost/numeric/bindings/traits/vector_traits.hpp>
 
21
#include <boost/numeric/bindings/atlas/cblas1_overloads.hpp>
 
22
 
 
23
#ifndef BOOST_NUMERIC_BINDINGS_NO_TYPE_CHECK
 
24
#  include <boost/type_traits/same_traits.hpp>
 
25
#  include <boost/static_assert.hpp>
 
26
#endif 
 
27
 
 
28
namespace boost { namespace numeric { namespace bindings { 
 
29
 
 
30
  namespace atlas {
 
31
 
 
32
    // x_i <- alpha for all i
 
33
    template <typename T, typename Vct> 
 
34
    inline 
 
35
    void set (T const& alpha, Vct& x) {
 
36
      detail::set (traits::vector_size (x), alpha, 
 
37
                   traits::vector_storage (x), traits::vector_stride (x)); 
 
38
    }
 
39
 
 
40
    // y <- x
 
41
    template <typename VctX, typename VctY>
 
42
    inline 
 
43
    void copy (VctX const& x, VctY& y) {
 
44
      assert (traits::vector_size (y) >= traits::vector_size (x));
 
45
      detail::copy (traits::vector_size (x),
 
46
#ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
 
47
                    traits::vector_storage (x), 
 
48
#else
 
49
                    traits::vector_storage_const (x), 
 
50
#endif 
 
51
                    traits::vector_stride (x), 
 
52
                    traits::vector_storage (y), traits::vector_stride (y)); 
 
53
    }
 
54
 
 
55
    // x <-> y
 
56
    template <typename VctX, typename VctY>
 
57
    inline 
 
58
    void swap (VctX& x, VctY& y) {
 
59
      assert (traits::vector_size (y) >= traits::vector_size (x));
 
60
      detail::swap (traits::vector_size (x),
 
61
                    traits::vector_storage (x), traits::vector_stride (x), 
 
62
                    traits::vector_storage (y), traits::vector_stride (y)); 
 
63
    }
 
64
 
 
65
    // x <- alpha * x
 
66
    template <typename T, typename Vct> 
 
67
    inline 
 
68
    void scal (T const& alpha, Vct& x) {
 
69
#ifndef BOOST_NUMERIC_BINDINGS_NO_TYPE_CHECK
 
70
      typedef traits::vector_traits<Vct> vtraits;
 
71
      BOOST_STATIC_ASSERT(
 
72
       (boost::is_same<T, typename vtraits::value_type>::value
 
73
        || 
 
74
        boost::is_same<T, 
 
75
          typename traits::type_traits<typename vtraits::value_type>::real_type
 
76
        >::value
 
77
        ));
 
78
#endif 
 
79
      detail::scal (traits::vector_size (x), alpha, 
 
80
                    traits::vector_storage (x), traits::vector_stride (x)); 
 
81
    }
 
82
 
 
83
    // y <- alpha * x + y
 
84
    template <typename T, typename VctX, typename VctY>
 
85
    inline 
 
86
    void axpy (T const& alpha, VctX const& x, VctY& y) {
 
87
      assert (traits::vector_size (y) >= traits::vector_size (x));
 
88
      detail::axpy (traits::vector_size (x), alpha, 
 
89
#ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
 
90
                    traits::vector_storage (x), 
 
91
#else
 
92
                    traits::vector_storage_const (x), 
 
93
#endif
 
94
                    traits::vector_stride (x), 
 
95
                    traits::vector_storage (y), traits::vector_stride (y)); 
 
96
    }
 
97
 
 
98
    // y <- x + y
 
99
    template <typename VctX, typename VctY>
 
100
    inline 
 
101
    void xpy (VctX const& x, VctY& y) {
 
102
#ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
 
103
      typedef typename traits::vector_traits<VctX>::value_type val_t; 
 
104
#else
 
105
      typedef typename VctX::value_type val_t; 
 
106
#endif
 
107
      axpy ((val_t) 1, x, y); 
 
108
    }
 
109
 
 
110
    // y <- alpha * x + beta * y
 
111
    template <typename T, typename VctX, typename VctY>
 
112
    inline 
 
113
    void axpby (T const& alpha, VctX const& x, 
 
114
                T const& beta, VctY& y) { 
 
115
      assert (traits::vector_size (y) >= traits::vector_size (x));
 
116
      detail::axpby (traits::vector_size (x), alpha, 
 
117
#ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
 
118
                     traits::vector_storage (x), 
 
119
#else
 
120
                     traits::vector_storage_const (x), 
 
121
#endif
 
122
                     traits::vector_stride (x), 
 
123
                     beta, 
 
124
                     traits::vector_storage (y), traits::vector_stride (y)); 
 
125
    }
 
126
 
 
127
    ///////////////////////////////////////////
 
128
 
 
129
    // dot <- x^T * y 
 
130
    // .. real & complex types
 
131
    template <typename VctX, typename VctY>
 
132
    inline 
 
133
#ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
 
134
    typename traits::vector_traits<VctX>::value_type 
 
135
#else
 
136
    typename VctX::value_type 
 
137
#endif
 
138
    dot (VctX const& x, VctY const& y) {
 
139
      assert (traits::vector_size (y) >= traits::vector_size (x));
 
140
      return detail::dot (traits::vector_size (x),
 
141
#ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
 
142
                          traits::vector_storage (x), 
 
143
#else
 
144
                          traits::vector_storage_const (x), 
 
145
#endif
 
146
                          traits::vector_stride (x), 
 
147
#ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
 
148
                          traits::vector_storage (y), 
 
149
#else
 
150
                          traits::vector_storage_const (y), 
 
151
#endif
 
152
                          traits::vector_stride (y)); 
 
153
    }
 
154
 
 
155
    // dot <- x^T * y 
 
156
    // .. float only -- with double accumulation
 
157
    template <typename VctX, typename VctY>
 
158
    inline 
 
159
    double dsdot (VctX const& x, VctY const& y) {
 
160
      assert (traits::vector_size (y) >= traits::vector_size (x));
 
161
      return cblas_dsdot (traits::vector_size (x),
 
162
#ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
 
163
                          traits::vector_storage (x), 
 
164
#else
 
165
                          traits::vector_storage_const (x), 
 
166
#endif
 
167
                          traits::vector_stride (x), 
 
168
#ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
 
169
                          traits::vector_storage (y), 
 
170
#else
 
171
                          traits::vector_storage_const (y), 
 
172
#endif
 
173
                          traits::vector_stride (y)); 
 
174
    }
 
175
 
 
176
    // apdot <- alpha + x^T * y    
 
177
    // .. float only -- computation uses double precision 
 
178
    template <typename VctX, typename VctY>
 
179
    inline 
 
180
    float sdsdot (float const alpha, VctX const& x, VctY const& y) {
 
181
      assert (traits::vector_size (y) >= traits::vector_size (x));
 
182
      return cblas_sdsdot (traits::vector_size (x), alpha, 
 
183
#ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
 
184
                           traits::vector_storage (x), 
 
185
#else
 
186
                           traits::vector_storage_const (x), 
 
187
#endif
 
188
                           traits::vector_stride (x), 
 
189
#ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
 
190
                           traits::vector_storage (y), 
 
191
#else
 
192
                           traits::vector_storage_const (y), 
 
193
#endif
 
194
                           traits::vector_stride (y)); 
 
195
    }
 
196
 
 
197
    // dotu <- x^T * y 
 
198
    // .. complex types only
 
199
    // .. function
 
200
    template <typename VctX, typename VctY>
 
201
    inline 
 
202
#ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
 
203
    typename traits::vector_traits<VctX>::value_type 
 
204
#else
 
205
    typename VctX::value_type 
 
206
#endif
 
207
    dotu (VctX const& x, VctY const& y) {
 
208
      assert (traits::vector_size (y) >= traits::vector_size (x));
 
209
#ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
 
210
      typename traits::vector_traits<VctX>::value_type val;
 
211
#else
 
212
      typename VctX::value_type val; 
 
213
#endif
 
214
      detail::dotu (traits::vector_size (x),
 
215
#ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
 
216
                    traits::vector_storage (x), 
 
217
#else
 
218
                    traits::vector_storage_const (x), 
 
219
#endif
 
220
                    traits::vector_stride (x), 
 
221
#ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
 
222
                    traits::vector_storage (y), 
 
223
#else
 
224
                    traits::vector_storage_const (y), 
 
225
#endif
 
226
                    traits::vector_stride (y), 
 
227
                    &val);
 
228
      return val; 
 
229
    }
 
230
    // .. procedure 
 
231
    template <typename VctX, typename VctY>
 
232
    inline 
 
233
    void dotu (VctX const& x, VctY const& y, 
 
234
#ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
 
235
               typename traits::vector_traits<VctX>::value_type& val
 
236
#else
 
237
               typename VctX::value_type& val
 
238
#endif
 
239
    ) {
 
240
      assert (traits::vector_size (y) >= traits::vector_size (x));
 
241
      detail::dotu (traits::vector_size (x),
 
242
#ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
 
243
                    traits::vector_storage (x), 
 
244
#else
 
245
                    traits::vector_storage_const (x), 
 
246
#endif
 
247
                    traits::vector_stride (x), 
 
248
#ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
 
249
                    traits::vector_storage (y), 
 
250
#else
 
251
                    traits::vector_storage_const (y), 
 
252
#endif
 
253
                    traits::vector_stride (y), 
 
254
                    &val);
 
255
    }
 
256
 
 
257
    // dotc <- x^H * y 
 
258
    // .. complex types only
 
259
    // .. function
 
260
    template <typename VctX, typename VctY>
 
261
    inline 
 
262
#ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
 
263
    typename traits::vector_traits<VctX>::value_type 
 
264
#else
 
265
    typename VctX::value_type 
 
266
#endif
 
267
    dotc (VctX const& x, VctY const& y) {
 
268
      assert (traits::vector_size (y) >= traits::vector_size (x));
 
269
#ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
 
270
      typename traits::vector_traits<VctX>::value_type val;
 
271
#else
 
272
      typename VctX::value_type val; 
 
273
#endif
 
274
      detail::dotc (traits::vector_size (x),
 
275
#ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
 
276
                    traits::vector_storage (x), 
 
277
#else
 
278
                    traits::vector_storage_const (x), 
 
279
#endif
 
280
                    traits::vector_stride (x), 
 
281
#ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
 
282
                    traits::vector_storage (y), 
 
283
#else
 
284
                    traits::vector_storage_const (y), 
 
285
#endif
 
286
                    traits::vector_stride (y),
 
287
                    &val);
 
288
      return val; 
 
289
    }
 
290
    // .. procedure 
 
291
    template <typename VctX, typename VctY>
 
292
    inline 
 
293
    void dotc (VctX const& x, VctY const& y, 
 
294
#ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
 
295
               typename traits::vector_traits<VctX>::value_type& val
 
296
#else
 
297
               typename VctX::value_type& val
 
298
#endif
 
299
    ) {
 
300
      assert (traits::vector_size (y) >= traits::vector_size (x));
 
301
      detail::dotc (traits::vector_size (x),
 
302
#ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
 
303
                    traits::vector_storage (x), 
 
304
#else
 
305
                    traits::vector_storage_const (x), 
 
306
#endif
 
307
                    traits::vector_stride (x), 
 
308
#ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
 
309
                    traits::vector_storage (y), 
 
310
#else
 
311
                    traits::vector_storage_const (y), 
 
312
#endif
 
313
                    traits::vector_stride (y), 
 
314
                    &val);
 
315
    }
 
316
 
 
317
    // nrm2 <- ||x||_2
 
318
    template <typename Vct> 
 
319
    inline 
 
320
    typename traits::type_traits<
 
321
#ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
 
322
      typename traits::vector_traits<Vct>::value_type 
 
323
#else 
 
324
      typename Vct::value_type 
 
325
#endif 
 
326
    >::real_type 
 
327
    nrm2 (Vct const& x) {
 
328
      return detail::nrm2 (traits::vector_size (x),
 
329
#ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
 
330
                           traits::vector_storage (x), 
 
331
#else
 
332
                           traits::vector_storage_const (x), 
 
333
#endif
 
334
                           traits::vector_stride (x)); 
 
335
    }
 
336
 
 
337
    // asum <- ||re (x)|| + ||im (x)||
 
338
    template <typename Vct> 
 
339
    inline 
 
340
    typename traits::type_traits<
 
341
#ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
 
342
      typename traits::vector_traits<Vct>::value_type 
 
343
#else 
 
344
      typename Vct::value_type 
 
345
#endif 
 
346
    >::real_type  
 
347
    asum (Vct const& x) {
 
348
      return detail::asum (traits::vector_size (x),
 
349
#ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
 
350
                           traits::vector_storage (x), 
 
351
#else
 
352
                           traits::vector_storage_const (x), 
 
353
#endif
 
354
                           traits::vector_stride (x)); 
 
355
    }
 
356
 
 
357
    // iamax <- 1st i: max (|re (x_i)| + |im (x_i)|)
 
358
    template <typename Vct> 
 
359
    inline 
 
360
    CBLAS_INDEX iamax (Vct const& x) {
 
361
      return detail::iamax (traits::vector_size (x),
 
362
#ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
 
363
                            traits::vector_storage (x), 
 
364
#else
 
365
                            traits::vector_storage_const (x), 
 
366
#endif
 
367
                            traits::vector_stride (x)); 
 
368
    }
 
369
 
 
370
    // TO DO: plane rotations 
 
371
 
 
372
  } // namespace atlas
 
373
 
 
374
}}} 
 
375
 
 
376
#endif // BOOST_NUMERIC_BINDINGS_CBLAS_LEVEL_1_HPP