cblas2.hpp

Go to the documentation of this file.
00001 /*
00002  * 
00003  * Copyright (c) Kresimir Fresl 2002, 2003
00004  *
00005  * Distributed under the Boost Software License, Version 1.0.
00006  * (See accompanying file LICENSE_1_0.txt or copy at
00007  * http://www.boost.org/LICENSE_1_0.txt)
00008  *
00009  * Author acknowledges the support of the Faculty of Civil Engineering, 
00010  * University of Zagreb, Croatia.
00011  *
00012  */
00013 
00014 #ifndef BOOST_NUMERIC_BINDINGS_CBLAS_LEVEL_2_HPP
00015 #define BOOST_NUMERIC_BINDINGS_CBLAS_LEVEL_2_HPP
00016 
00017 #include <cassert>
00018 
00019 #include <boost/numeric/bindings/traits/traits.hpp>
00020 #include <boost/numeric/bindings/traits/type_traits.hpp>
00021 #include <boost/numeric/bindings/atlas/cblas2_overloads.hpp>
00022 #include <boost/numeric/bindings/atlas/cblas_enum.hpp>
00023 
00024 #ifndef BOOST_NUMERIC_BINDINGS_NO_STRUCTURE_CHECK
00025 #  include <boost/static_assert.hpp>
00026 #  include <boost/type_traits/same_traits.hpp>
00027 #endif
00028 
00029 
00030 namespace boost { namespace numeric { namespace bindings { 
00031 
00032   namespace atlas {
00033 
00034     // y <- alpha * op (A) * x + beta * y
00035     // op (A) == A || A^T || A^H
00036     template <typename T, typename Matr, typename VctX, typename VctY>
00037     inline 
00038     void gemv (CBLAS_TRANSPOSE const TransA, 
00039                T const& alpha, Matr const& a, VctX const& x, 
00040                T const& beta, VctY& y
00041                )
00042     {
00043 #ifndef BOOST_NUMERIC_BINDINGS_NO_STRUCTURE_CHECK 
00044       BOOST_STATIC_ASSERT((boost::is_same<
00045         typename traits::matrix_traits<Matr>::matrix_structure, 
00046         traits::general_t
00047       >::value)); 
00048 #endif 
00049 
00050       int const m = traits::matrix_size1 (a);
00051       int const n = traits::matrix_size2 (a);
00052       assert (traits::vector_size (x) >= (TransA == CblasNoTrans ? n : m)); 
00053       assert (traits::vector_size (y) >= (TransA == CblasNoTrans ? m : n)); 
00054       // .. what about AtlasConj? 
00055 
00056       CBLAS_ORDER const stor_ord
00057         = enum_cast<CBLAS_ORDER const>
00058         (storage_order<
00059 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
00060            typename traits::matrix_traits<Matr>::ordering_type
00061 #else
00062            typename Matr::orientation_category 
00063 #endif 
00064          >::value); 
00065 
00066       detail::gemv (stor_ord, TransA, m, n, alpha, 
00067 #ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
00068                     traits::matrix_storage (a), 
00069 #else
00070                     traits::matrix_storage_const (a), 
00071 #endif
00072                     traits::leading_dimension (a),
00073 #ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
00074                     traits::vector_storage (x), 
00075 #else
00076                     traits::vector_storage_const (x), 
00077 #endif
00078                     traits::vector_stride (x),
00079                     beta, 
00080                     traits::vector_storage (y), 
00081                     traits::vector_stride (y));
00082     }
00083 
00084     // y <- alpha * A * x + beta * y 
00085     template <typename T, typename Matr, typename VctX, typename VctY>
00086     inline 
00087     void gemv (T const& alpha, Matr const& a, VctX const& x, 
00088                T const& beta, VctY& y) {
00089       gemv (CblasNoTrans, alpha, a, x, beta, y); 
00090     }
00091 
00092     // y <- A * x
00093     template <typename Matr, typename VctX, typename VctY>
00094     inline 
00095     void gemv (Matr const& a, VctX const& x, VctY& y) {
00096 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
00097       typedef typename traits::matrix_traits<Matr>::value_type val_t; 
00098 #else
00099       typedef typename Matr::value_type val_t; 
00100 #endif 
00101       gemv (CblasNoTrans, (val_t) 1, a, x, (val_t) 0, y);
00102     }
00103 
00104 
00105     // y <- alpha * A * x + beta * y
00106     // A real symmetric matrix (T == float | double)
00107     // [from 'dsymv.f':]
00108     /* [...] with UPLO = 'U' or 'u', the  leading n by n upper 
00109      * triangular part of the array A must contain the upper triangular
00110      * part of the symmetric matrix and the strictly lower triangular 
00111      * part of A is not referenced.
00112      * [...] with UPLO = 'L' or 'l', the leading n by n lower
00113      * triangular part of the array A must contain the lower triangular  
00114      * part of the symmetric matrix and the strictly  upper
00115      * triangular part of A is not referenced.
00116      */ 
00117     template <typename T, typename SymmMatr, typename VctX, typename VctY>
00118     inline 
00119     void symv (CBLAS_UPLO const uplo, T const& alpha, SymmMatr const& a, 
00120                VctX const& x, T const& beta, VctY& y)
00121     {
00122 #ifndef BOOST_NUMERIC_BINDINGS_NO_STRUCTURE_CHECK 
00123       BOOST_STATIC_ASSERT((boost::is_same<
00124         typename traits::matrix_traits<SymmMatr>::matrix_structure, 
00125         traits::general_t
00126       >::value)); 
00127 #endif 
00128 
00129       int const n = traits::matrix_size1 (a);
00130       assert (n == traits::matrix_size2 (a)); 
00131       assert (traits::vector_size (x) >= n); 
00132       assert (traits::vector_size (y) >= n); 
00133 
00134       CBLAS_ORDER const stor_ord
00135         = enum_cast<CBLAS_ORDER const>
00136         (storage_order<
00137 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
00138            typename traits::matrix_traits<SymmMatr>::ordering_type
00139 #else
00140            typename SymmMatr::orientation_category 
00141 #endif
00142          >::value); 
00143 
00144       detail::symv (stor_ord, uplo, n, alpha, 
00145 #ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
00146                     traits::matrix_storage (a), 
00147 #else
00148                     traits::matrix_storage_const (a), 
00149 #endif
00150                     traits::leading_dimension (a),
00151 #ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
00152                     traits::vector_storage (x), 
00153 #else
00154                     traits::vector_storage_const (x), 
00155 #endif
00156                     traits::vector_stride (x),
00157                     beta, 
00158                     traits::vector_storage (y), 
00159                     traits::vector_stride (y));
00160     }
00161 
00162     template <typename T, typename SymmMatr, typename VctX, typename VctY>
00163     inline 
00164     void symv (T const& alpha, SymmMatr const& a, VctX const& x, 
00165                T const& beta, VctY& y)
00166     {
00167 #ifndef BOOST_NUMERIC_BINDINGS_NO_STRUCTURE_CHECK 
00168       BOOST_STATIC_ASSERT((boost::is_same<
00169         typename traits::matrix_traits<SymmMatr>::matrix_structure, 
00170         traits::symmetric_t
00171       >::value)); 
00172 #endif 
00173 
00174       int const n = traits::matrix_size1 (a);
00175       assert (n == traits::matrix_size2 (a)); 
00176       assert (traits::vector_size (x) >= n); 
00177       assert (traits::vector_size (y) >= n); 
00178 
00179       CBLAS_ORDER const stor_ord
00180         = enum_cast<CBLAS_ORDER const>
00181         (storage_order<
00182 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
00183            typename traits::matrix_traits<SymmMatr>::ordering_type
00184 #else
00185            typename SymmMatr::orientation_category 
00186 #endif
00187          >::value); 
00188 
00189       CBLAS_UPLO const uplo
00190         = enum_cast<CBLAS_UPLO const>
00191         (uplo_triang<
00192 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
00193            typename traits::matrix_traits<SymmMatr>::uplo_type
00194 #else
00195            typename SymmMatr::packed_category 
00196 #endif 
00197          >::value); 
00198 
00199       detail::symv (stor_ord, uplo, n, alpha, 
00200 #ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
00201                     traits::matrix_storage (a), 
00202 #else
00203                     traits::matrix_storage_const (a), 
00204 #endif
00205                     traits::leading_dimension (a),
00206 #ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
00207                     traits::vector_storage (x), 
00208 #else
00209                     traits::vector_storage_const (x), 
00210 #endif
00211                     traits::vector_stride (x),
00212                     beta, 
00213                     traits::vector_storage (y), 
00214                     traits::vector_stride (y));
00215     }
00216 
00217     // y <- A * x
00218     template <typename SymmMatr, typename VctX, typename VctY>
00219     inline 
00220     void symv (SymmMatr const& a, VctX const& x, VctY& y) {
00221 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
00222       typedef typename traits::matrix_traits<SymmMatr>::value_type val_t; 
00223 #else
00224       typedef typename SymmMatr::value_type val_t; 
00225 #endif
00226       symv ((val_t) 1, a, x, (val_t) 0, y);
00227     }
00228 
00229 
00230     // y <- alpha * A * x + beta * y
00231     // A real symmetric matrix (T == float | double) in packed form
00232     // [from 'dspmv.f' (description assumes column major order):]
00233     /* A is an array of DIMENSION ( ( n*( n + 1 ) )/2 ).
00234      * Before entry with UPLO = 'U' or 'u', the array A must contain
00235      * the upper triangular  part of the symmetric matrix packed 
00236      * sequentially, column by column, so that A( 1 ) contains a(1,1),
00237      * A( 2 ) and A( 3 ) contain a(1,2) and a(2,2) respectively, and
00238      * so on.
00239      * Before entry with UPLO = 'L' or 'l', the array A must contain
00240      * the lower triangular  part of the symmetric matrix packed 
00241      * sequentially, column by column, so that A( 1 ) contains a(1,1),
00242      * A( 2 ) and A( 3 ) contain a(2,1) and a(3,1) respectively, and
00243      * so on. 
00244      */ 
00245     template <typename T, typename SymmMatr, typename VctX, typename VctY>
00246     inline 
00247     void spmv (T const& alpha, SymmMatr const& a, VctX const& x, 
00248                T const& beta, VctY& y)
00249     {
00250 #ifndef BOOST_NUMERIC_BINDINGS_NO_STRUCTURE_CHECK 
00251       BOOST_STATIC_ASSERT((boost::is_same<
00252         typename traits::matrix_traits<SymmMatr>::matrix_structure, 
00253         traits::symmetric_packed_t
00254       >::value)); 
00255 #endif 
00256 
00257       int const n = traits::matrix_size1 (a);
00258       assert (n == traits::matrix_size2 (a)); 
00259       assert (traits::vector_size (x) >= n); 
00260       assert (traits::vector_size (y) >= n); 
00261 
00262       CBLAS_ORDER const stor_ord
00263         = enum_cast<CBLAS_ORDER const>
00264         (storage_order<
00265 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
00266            typename traits::matrix_traits<SymmMatr>::ordering_type
00267 #else
00268            typename SymmMatr::orientation_category 
00269 #endif
00270          >::value); 
00271 
00272       CBLAS_UPLO const uplo
00273         = enum_cast<CBLAS_UPLO const>
00274         (uplo_triang<
00275 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
00276            typename traits::matrix_traits<SymmMatr>::uplo_type
00277 #else
00278            typename SymmMatr::packed_category 
00279 #endif 
00280          >::value); 
00281 
00282       detail::spmv (stor_ord, uplo, n, alpha, 
00283 #ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
00284                     traits::matrix_storage (a), 
00285                     traits::vector_storage (x), 
00286 #else
00287                     traits::matrix_storage_const (a), 
00288                     traits::vector_storage_const (x), 
00289 #endif
00290                     traits::vector_stride (x),
00291                     beta, 
00292                     traits::vector_storage (y), 
00293                     traits::vector_stride (y));
00294     }
00295 
00296     // y <- A * x
00297     template <typename SymmMatr, typename VctX, typename VctY>
00298     inline 
00299     void spmv (SymmMatr const& a, VctX const& x, VctY& y) {
00300 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
00301       typedef typename traits::matrix_traits<SymmMatr>::value_type val_t; 
00302 #else
00303       typedef typename SymmMatr::value_type val_t; 
00304 #endif
00305       spmv ((val_t) 1, a, x, (val_t) 0, y);
00306     }
00307 
00308 
00309     // y <- alpha * A * x + beta * y
00310     // A complex hermitian matrix 
00311     // (T == std::complex<float> | std::complex<double>)
00312     template <typename T, typename HermMatr, typename VctX, typename VctY>
00313     inline 
00314     void hemv (CBLAS_UPLO const uplo, T const& alpha, HermMatr const& a, 
00315                VctX const& x, T const& beta, VctY& y)
00316     {
00317 #ifndef BOOST_NUMERIC_BINDINGS_NO_STRUCTURE_CHECK 
00318       BOOST_STATIC_ASSERT((boost::is_same<
00319         typename traits::matrix_traits<HermMatr>::matrix_structure, 
00320         traits::general_t
00321       >::value)); 
00322 #endif 
00323 
00324       int const n = traits::matrix_size1 (a);
00325       assert (n == traits::matrix_size2 (a)); 
00326       assert (traits::vector_size (x) >= n); 
00327       assert (traits::vector_size (y) >= n); 
00328 
00329       CBLAS_ORDER const stor_ord
00330         = enum_cast<CBLAS_ORDER const>
00331         (storage_order<
00332 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
00333            typename traits::matrix_traits<HermMatr>::ordering_type
00334 #else
00335            typename HermMatr::orientation_category 
00336 #endif
00337          >::value); 
00338 
00339       detail::hemv (stor_ord, uplo, n, alpha, 
00340 #ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
00341                     traits::matrix_storage (a), 
00342 #else
00343                     traits::matrix_storage_const (a), 
00344 #endif
00345                     traits::leading_dimension (a),
00346 #ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
00347                     traits::vector_storage (x), 
00348 #else
00349                     traits::vector_storage_const (x), 
00350 #endif
00351                     traits::vector_stride (x),
00352                     beta, 
00353                     traits::vector_storage (y), 
00354                     traits::vector_stride (y));
00355     }
00356 
00357     template <typename T, typename HermMatr, typename VctX, typename VctY>
00358     inline 
00359     void hemv (T const& alpha, HermMatr const& a, VctX const& x, 
00360                T const& beta, VctY& y)
00361     {
00362 #ifndef BOOST_NUMERIC_BINDINGS_NO_STRUCTURE_CHECK 
00363       BOOST_STATIC_ASSERT((boost::is_same<
00364         typename traits::matrix_traits<HermMatr>::matrix_structure, 
00365         traits::hermitian_t
00366       >::value)); 
00367 #endif 
00368 
00369       int const n = traits::matrix_size1 (a);
00370       assert (n == traits::matrix_size2 (a)); 
00371       assert (traits::vector_size (x) >= n); 
00372       assert (traits::vector_size (y) >= n); 
00373 
00374       CBLAS_ORDER const stor_ord
00375         = enum_cast<CBLAS_ORDER const>
00376         (storage_order<
00377 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
00378            typename traits::matrix_traits<HermMatr>::ordering_type
00379 #else
00380            typename HermMatr::orientation_category 
00381 #endif
00382          >::value); 
00383 
00384       CBLAS_UPLO const uplo
00385         = enum_cast<CBLAS_UPLO const>
00386         (uplo_triang<
00387 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
00388            typename traits::matrix_traits<HermMatr>::uplo_type
00389 #else
00390            typename HermMatr::packed_category 
00391 #endif 
00392          >::value); 
00393 
00394       detail::hemv (stor_ord, uplo, n, alpha, 
00395 #ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
00396                     traits::matrix_storage (a), 
00397 #else
00398                     traits::matrix_storage_const (a), 
00399 #endif
00400                     traits::leading_dimension (a),
00401 #ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
00402                     traits::vector_storage (x), 
00403 #else
00404                     traits::vector_storage_const (x), 
00405 #endif
00406                     traits::vector_stride (x),
00407                     beta, 
00408                     traits::vector_storage (y), 
00409                     traits::vector_stride (y));
00410     }
00411 
00412     // y <- A * x
00413     template <typename HermMatr, typename VctX, typename VctY>
00414     inline 
00415     void hemv (HermMatr const& a, VctX const& x, VctY& y) {
00416 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
00417       typedef typename traits::matrix_traits<HermMatr>::value_type val_t; 
00418 #else
00419       typedef typename HermMatr::value_type val_t; 
00420 #endif
00421       hemv ((val_t) 1, a, x, (val_t) 0, y);
00422     }
00423 
00424 
00425     // y <- alpha * A * x + beta * y
00426     // A complex hermitian matrix in packed form 
00427     // (T == std::complex<float> | std::complex<double>)
00428     template <typename T, typename HermMatr, typename VctX, typename VctY>
00429     inline 
00430     void hpmv (T const& alpha, HermMatr const& a, VctX const& x, 
00431                T const& beta, VctY& y)
00432     {
00433 #ifndef BOOST_NUMERIC_BINDINGS_NO_STRUCTURE_CHECK 
00434       BOOST_STATIC_ASSERT((boost::is_same<
00435         typename traits::matrix_traits<HermMatr>::matrix_structure, 
00436         traits::hermitian_packed_t
00437       >::value)); 
00438 #endif 
00439 
00440       int const n = traits::matrix_size1 (a);
00441       assert (n == traits::matrix_size2 (a)); 
00442       assert (traits::vector_size (x) >= n); 
00443       assert (traits::vector_size (y) >= n); 
00444 
00445       CBLAS_ORDER const stor_ord
00446         = enum_cast<CBLAS_ORDER const>
00447         (storage_order<
00448 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
00449            typename traits::matrix_traits<HermMatr>::ordering_type
00450 #else
00451            typename HermMatr::orientation_category 
00452 #endif
00453          >::value); 
00454 
00455       CBLAS_UPLO const uplo
00456         = enum_cast<CBLAS_UPLO const>
00457         (uplo_triang<
00458 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
00459            typename traits::matrix_traits<HermMatr>::uplo_type
00460 #else
00461            typename HermMatr::packed_category 
00462 #endif 
00463          >::value); 
00464 
00465       detail::hpmv (stor_ord, uplo, n, alpha, 
00466 #ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
00467                     traits::matrix_storage (a), 
00468                     traits::vector_storage (x), 
00469 #else
00470                     traits::matrix_storage_const (a), 
00471                     traits::vector_storage_const (x), 
00472 #endif
00473                     traits::vector_stride (x),
00474                     beta, 
00475                     traits::vector_storage (y), 
00476                     traits::vector_stride (y));
00477     }
00478 
00479     // y <- A * x
00480     template <typename HermMatr, typename VctX, typename VctY>
00481     inline 
00482     void hpmv (HermMatr const& a, VctX const& x, VctY& y) {
00483 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
00484       typedef typename traits::matrix_traits<HermMatr>::value_type val_t; 
00485 #else
00486       typedef typename HermMatr::value_type val_t; 
00487 #endif
00488       hpmv ((val_t) 1, a, x, (val_t) 0, y);
00489     }
00490 
00491 
00492     // A <- alpha * x * y^T + A 
00493     // .. real & complex types
00494     template <typename T, typename Matr, typename VctX, typename VctY>
00495     inline 
00496     void ger (T const& alpha, VctX const& x, VctY const& y, Matr& a) {
00497 #ifndef BOOST_NUMERIC_BINDINGS_NO_STRUCTURE_CHECK 
00498       BOOST_STATIC_ASSERT((boost::is_same<
00499         typename traits::matrix_traits<Matr>::matrix_structure, 
00500         traits::general_t
00501       >::value)); 
00502 #endif 
00503 
00504       int const m = traits::matrix_size1 (a);
00505       int const n = traits::matrix_size2 (a);
00506       assert (traits::vector_size (x) >= m); 
00507       assert (traits::vector_size (y) >= n); 
00508 
00509       CBLAS_ORDER const stor_ord
00510         = enum_cast<CBLAS_ORDER const>
00511         (storage_order<
00512 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
00513            typename traits::matrix_traits<Matr>::ordering_type
00514 #else
00515            typename Matr::orientation_category 
00516 #endif 
00517          >::value); 
00518 
00519       detail::ger (stor_ord, m, n, alpha, 
00520 #ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
00521                    traits::vector_storage (x), 
00522 #else
00523                    traits::vector_storage_const (x), 
00524 #endif
00525                    traits::vector_stride (x),
00526 #ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
00527                    traits::vector_storage (y), 
00528 #else
00529                    traits::vector_storage_const (y), 
00530 #endif
00531                    traits::vector_stride (y),
00532                    traits::matrix_storage (a), 
00533                    traits::leading_dimension (a)); 
00534     }
00535 
00536     // A <- x * y^T + A 
00537     template <typename Matr, typename VctX, typename VctY>
00538     inline 
00539     void ger (VctX const& x, VctY const& y, Matr& a) {
00540 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
00541       typedef typename traits::matrix_traits<Matr>::value_type val_t; 
00542 #else
00543       typedef typename Matr::value_type val_t; 
00544 #endif 
00545       ger ((val_t) 1, x, y, a); 
00546     }
00547 
00548     // A <- alpha * x * y^T + A 
00549     // .. complex types only
00550     template <typename T, typename Matr, typename VctX, typename VctY>
00551     inline 
00552     void geru (T const& alpha, VctX const& x, VctY const& y, Matr& a) {
00553 #ifndef BOOST_NUMERIC_BINDINGS_NO_STRUCTURE_CHECK 
00554       BOOST_STATIC_ASSERT((boost::is_same<
00555         typename traits::matrix_traits<Matr>::matrix_structure, 
00556         traits::general_t
00557       >::value)); 
00558 #endif 
00559 
00560       int const m = traits::matrix_size1 (a);
00561       int const n = traits::matrix_size2 (a);
00562       assert (traits::vector_size (x) >= m); 
00563       assert (traits::vector_size (y) >= n); 
00564 
00565       CBLAS_ORDER const stor_ord
00566         = enum_cast<CBLAS_ORDER const>
00567         (storage_order<
00568 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
00569            typename traits::matrix_traits<Matr>::ordering_type
00570 #else
00571            typename Matr::orientation_category 
00572 #endif 
00573          >::value); 
00574 
00575       detail::geru (stor_ord, m, n, alpha, 
00576 #ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
00577                     traits::vector_storage (x), 
00578 #else
00579                     traits::vector_storage_const (x), 
00580 #endif
00581                     traits::vector_stride (x),
00582 #ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
00583                     traits::vector_storage (y), 
00584 #else
00585                     traits::vector_storage_const (y), 
00586 #endif
00587                     traits::vector_stride (y),
00588                     traits::matrix_storage (a), 
00589                     traits::leading_dimension (a)); 
00590     }
00591 
00592     // A <- x * y^T + A 
00593     template <typename Matr, typename VctX, typename VctY>
00594     inline 
00595     void geru (VctX const& x, VctY const& y, Matr& a) {
00596 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
00597       typedef typename traits::matrix_traits<Matr>::value_type val_t; 
00598 #else
00599       typedef typename Matr::value_type val_t; 
00600 #endif 
00601       geru ((val_t) 1, x, y, a); 
00602     }
00603 
00604     // A <- alpha * x * y^H + A 
00605     // .. complex types only
00606     template <typename T, typename Matr, typename VctX, typename VctY>
00607     inline 
00608     void gerc (T const& alpha, VctX const& x, VctY const& y, Matr& a) {
00609 #ifndef BOOST_NUMERIC_BINDINGS_NO_STRUCTURE_CHECK 
00610       BOOST_STATIC_ASSERT((boost::is_same<
00611         typename traits::matrix_traits<Matr>::matrix_structure, 
00612         traits::general_t
00613       >::value)); 
00614 #endif 
00615 
00616       int const m = traits::matrix_size1 (a);
00617       int const n = traits::matrix_size2 (a);
00618       assert (traits::vector_size (x) >= m); 
00619       assert (traits::vector_size (y) >= n); 
00620 
00621       CBLAS_ORDER const stor_ord
00622         = enum_cast<CBLAS_ORDER const>
00623         (storage_order<
00624 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
00625            typename traits::matrix_traits<Matr>::ordering_type
00626 #else
00627            typename Matr::orientation_category 
00628 #endif 
00629          >::value); 
00630 
00631       detail::gerc (stor_ord, m, n, alpha, 
00632 #ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
00633                     traits::vector_storage (x), 
00634 #else
00635                     traits::vector_storage_const (x), 
00636 #endif
00637                     traits::vector_stride (x),
00638 #ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
00639                     traits::vector_storage (y), 
00640 #else
00641                     traits::vector_storage_const (y), 
00642 #endif
00643                     traits::vector_stride (y),
00644                     traits::matrix_storage (a), 
00645                     traits::leading_dimension (a)); 
00646     }
00647 
00648     // A <- x * y^H + A 
00649     template <typename Matr, typename VctX, typename VctY>
00650     inline 
00651     void gerc (VctX const& x, VctY const& y, Matr& a) {
00652 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
00653       typedef typename traits::matrix_traits<Matr>::value_type val_t; 
00654 #else
00655       typedef typename Matr::value_type val_t; 
00656 #endif 
00657       gerc ((val_t) 1, x, y, a); 
00658     }
00659 
00660 
00661     // A <- alpha * x * x^T + A 
00662     // A real symmetric (see leading comments for 'symv()') 
00663     template <typename T, typename SymmM, typename VctX>
00664     inline 
00665     void syr (CBLAS_UPLO const uplo, T const& alpha, VctX const& x, SymmM& a) {
00666 #ifndef BOOST_NUMERIC_BINDINGS_NO_STRUCTURE_CHECK 
00667       BOOST_STATIC_ASSERT((boost::is_same<
00668         typename traits::matrix_traits<SymmM>::matrix_structure, 
00669         traits::general_t
00670       >::value)); 
00671 #endif 
00672 
00673       int const n = traits::matrix_size1 (a);
00674       assert (n == traits::matrix_size2 (a)); 
00675       assert (traits::vector_size (x) >= n); 
00676 
00677       CBLAS_ORDER const stor_ord
00678         = enum_cast<CBLAS_ORDER const>
00679         (storage_order<
00680 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
00681            typename traits::matrix_traits<SymmM>::ordering_type
00682 #else
00683            typename SymmM::orientation_category 
00684 #endif
00685          >::value); 
00686 
00687       detail::syr (stor_ord, uplo, n, alpha, 
00688 #ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
00689                    traits::vector_storage (x), 
00690 #else
00691                    traits::vector_storage_const (x), 
00692 #endif
00693                    traits::vector_stride (x),
00694                    traits::matrix_storage (a), 
00695                    traits::leading_dimension (a)); 
00696     }
00697 
00698     template <typename T, typename SymmM, typename VctX>
00699     inline 
00700     void syr (T const& alpha, VctX const& x, SymmM& a) {
00701 #ifndef BOOST_NUMERIC_BINDINGS_NO_STRUCTURE_CHECK 
00702       BOOST_STATIC_ASSERT((boost::is_same<
00703         typename traits::matrix_traits<SymmM>::matrix_structure, 
00704         traits::symmetric_t
00705       >::value)); 
00706 #endif 
00707 
00708       int const n = traits::matrix_size1 (a);
00709       assert (n == traits::matrix_size2 (a)); 
00710       assert (traits::vector_size (x) >= n); 
00711 
00712       CBLAS_ORDER const stor_ord
00713         = enum_cast<CBLAS_ORDER const>
00714         (storage_order<
00715 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
00716            typename traits::matrix_traits<SymmM>::ordering_type
00717 #else
00718            typename SymmM::orientation_category 
00719 #endif
00720          >::value); 
00721 
00722       CBLAS_UPLO const uplo
00723         = enum_cast<CBLAS_UPLO const>
00724         (uplo_triang<
00725 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
00726            typename traits::matrix_traits<SymmM>::uplo_type
00727 #else
00728            typename SymmM::packed_category 
00729 #endif 
00730          >::value); 
00731 
00732       detail::syr (stor_ord, uplo, n, alpha, 
00733 #ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
00734                    traits::vector_storage (x), 
00735 #else
00736                    traits::vector_storage_const (x), 
00737 #endif
00738                    traits::vector_stride (x),
00739                    traits::matrix_storage (a), 
00740                    traits::leading_dimension (a)); 
00741     }
00742 
00743     // A <- x * x^T + A 
00744     template <typename SymmM, typename VctX>
00745     inline 
00746     void syr (VctX const& x, SymmM& a) { 
00747 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
00748       typedef typename traits::matrix_traits<SymmM>::value_type val_t; 
00749 #else
00750       typedef typename SymmM::value_type val_t; 
00751 #endif 
00752       syr ((val_t) 1, x, a); 
00753     }
00754 
00755 
00756     // A <- alpha * x * x^T + A 
00757     // A real symmetric in packed form (see leading comments for 'spmv()') 
00758     template <typename T, typename SymmM, typename VctX>
00759     inline 
00760     void spr (T const& alpha, VctX const& x, SymmM& a) {
00761 #ifndef BOOST_NUMERIC_BINDINGS_NO_STRUCTURE_CHECK 
00762       BOOST_STATIC_ASSERT((boost::is_same<
00763         typename traits::matrix_traits<SymmM>::matrix_structure, 
00764         traits::symmetric_packed_t
00765       >::value)); 
00766 #endif 
00767 
00768       int const n = traits::matrix_size1 (a);
00769       assert (n == traits::matrix_size2 (a)); 
00770       assert (traits::vector_size (x) >= n); 
00771 
00772       CBLAS_ORDER const stor_ord
00773         = enum_cast<CBLAS_ORDER const>
00774         (storage_order<
00775 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
00776            typename traits::matrix_traits<SymmM>::ordering_type
00777 #else
00778            typename SymmM::orientation_category 
00779 #endif
00780          >::value); 
00781 
00782       CBLAS_UPLO const uplo
00783         = enum_cast<CBLAS_UPLO const>
00784         (uplo_triang<
00785 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
00786            typename traits::matrix_traits<SymmM>::uplo_type
00787 #else
00788            typename SymmM::packed_category 
00789 #endif 
00790          >::value); 
00791 
00792       detail::spr (stor_ord, uplo, n, alpha, 
00793 #ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
00794                    traits::vector_storage (x), 
00795 #else
00796                    traits::vector_storage_const (x), 
00797 #endif
00798                    traits::vector_stride (x),
00799                    traits::matrix_storage (a));
00800     }
00801 
00802     // A <- x * x^T + A 
00803     template <typename SymmM, typename VctX>
00804     inline 
00805     void spr (VctX const& x, SymmM& a) { 
00806 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
00807       typedef typename traits::matrix_traits<SymmM>::value_type val_t; 
00808 #else
00809       typedef typename SymmM::value_type val_t; 
00810 #endif 
00811       spr ((val_t) 1, x, a); 
00812     }
00813 
00814 
00815     // A <- alpha * x * y^T + alpha * y * x^T + A 
00816     // A real symmetric (see leading comments for 'symv()') 
00817     template <typename T, typename SymmM, typename VctX, typename VctY>
00818     inline 
00819     void syr2 (CBLAS_UPLO const uplo, T const& alpha, 
00820                VctX const& x, VctY const& y, SymmM& a) 
00821     {
00822 #ifndef BOOST_NUMERIC_BINDINGS_NO_STRUCTURE_CHECK 
00823       BOOST_STATIC_ASSERT((boost::is_same<
00824         typename traits::matrix_traits<SymmM>::matrix_structure, 
00825         traits::general_t
00826       >::value)); 
00827 #endif 
00828 
00829       int const n = traits::matrix_size1 (a);
00830       assert (n == traits::matrix_size2 (a)); 
00831       assert (traits::vector_size (x) >= n); 
00832       assert (traits::vector_size (y) >= n); 
00833 
00834       CBLAS_ORDER const stor_ord
00835         = enum_cast<CBLAS_ORDER const>
00836         (storage_order<
00837 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
00838            typename traits::matrix_traits<SymmM>::ordering_type
00839 #else
00840            typename SymmM::orientation_category 
00841 #endif
00842          >::value); 
00843 
00844       detail::syr2 (stor_ord, uplo, n, alpha, 
00845 #ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
00846                     traits::vector_storage (x), 
00847 #else
00848                     traits::vector_storage_const (x), 
00849 #endif
00850                     traits::vector_stride (x),
00851 #ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
00852                     traits::vector_storage (y), 
00853 #else
00854                     traits::vector_storage_const (y), 
00855 #endif
00856                     traits::vector_stride (y),
00857                     traits::matrix_storage (a), 
00858                     traits::leading_dimension (a)); 
00859     }
00860 
00861     template <typename T, typename SymmM, typename VctX, typename VctY>
00862     inline 
00863     void syr2 (T const& alpha, VctX const& x, VctY const& y, SymmM& a) {
00864 #ifndef BOOST_NUMERIC_BINDINGS_NO_STRUCTURE_CHECK 
00865       BOOST_STATIC_ASSERT((boost::is_same<
00866         typename traits::matrix_traits<SymmM>::matrix_structure, 
00867         traits::symmetric_t
00868       >::value)); 
00869 #endif 
00870 
00871       int const n = traits::matrix_size1 (a);
00872       assert (n == traits::matrix_size2 (a)); 
00873       assert (traits::vector_size (x) >= n); 
00874       assert (traits::vector_size (y) >= n); 
00875 
00876       CBLAS_ORDER const stor_ord
00877         = enum_cast<CBLAS_ORDER const>
00878         (storage_order<
00879 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
00880            typename traits::matrix_traits<SymmM>::ordering_type
00881 #else
00882            typename SymmM::orientation_category 
00883 #endif
00884          >::value); 
00885 
00886       CBLAS_UPLO const uplo
00887         = enum_cast<CBLAS_UPLO const>
00888         (uplo_triang<
00889 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
00890            typename traits::matrix_traits<SymmM>::uplo_type
00891 #else
00892            typename SymmM::packed_category 
00893 #endif 
00894          >::value); 
00895 
00896       detail::syr2 (stor_ord, uplo, n, alpha, 
00897 #ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
00898                     traits::vector_storage (x), 
00899 #else
00900                     traits::vector_storage_const (x), 
00901 #endif
00902                     traits::vector_stride (x),
00903 #ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
00904                     traits::vector_storage (y), 
00905 #else
00906                     traits::vector_storage_const (y), 
00907 #endif
00908                     traits::vector_stride (y),
00909                     traits::matrix_storage (a), 
00910                     traits::leading_dimension (a)); 
00911     }
00912 
00913     // A <- x * y^T + y * x^T + A 
00914     template <typename SymmM, typename VctX, typename VctY>
00915     inline 
00916     void syr2 (VctX const& x, VctY const& y, SymmM& a) { 
00917 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
00918       typedef typename traits::matrix_traits<SymmM>::value_type val_t; 
00919 #else
00920       typedef typename SymmM::value_type val_t; 
00921 #endif 
00922       syr2 ((val_t) 1, x, y, a); 
00923     }
00924 
00925 
00926     // A <- alpha * x * y^T + alpha * y * x^T + A 
00927     // A real symmetric in packed form (see leading comments for 'spmv()') 
00928     template <typename T, typename SymmM, typename VctX, typename VctY>
00929     inline 
00930     void spr2 (T const& alpha, VctX const& x, VctY const& y, SymmM& a) {
00931 #ifndef BOOST_NUMERIC_BINDINGS_NO_STRUCTURE_CHECK 
00932       BOOST_STATIC_ASSERT((boost::is_same<
00933         typename traits::matrix_traits<SymmM>::matrix_structure, 
00934         traits::symmetric_packed_t
00935       >::value)); 
00936 #endif 
00937 
00938       int const n = traits::matrix_size1 (a);
00939       assert (n == traits::matrix_size2 (a)); 
00940       assert (traits::vector_size (x) >= n); 
00941       assert (traits::vector_size (y) >= n); 
00942 
00943       CBLAS_ORDER const stor_ord
00944         = enum_cast<CBLAS_ORDER const>
00945         (storage_order<
00946 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
00947            typename traits::matrix_traits<SymmM>::ordering_type
00948 #else
00949            typename SymmM::orientation_category 
00950 #endif
00951          >::value); 
00952 
00953       CBLAS_UPLO const uplo
00954         = enum_cast<CBLAS_UPLO const>
00955         (uplo_triang<
00956 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
00957            typename traits::matrix_traits<SymmM>::uplo_type
00958 #else
00959            typename SymmM::packed_category 
00960 #endif 
00961          >::value); 
00962 
00963       detail::spr2 (stor_ord, uplo, n, alpha, 
00964 #ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
00965                     traits::vector_storage (x), 
00966 #else
00967                     traits::vector_storage_const (x), 
00968 #endif
00969                     traits::vector_stride (x),
00970 #ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
00971                     traits::vector_storage (y), 
00972 #else
00973                     traits::vector_storage_const (y), 
00974 #endif
00975                     traits::vector_stride (y),
00976                     traits::matrix_storage (a));
00977     }
00978 
00979     // A <- x * y^T + y * x^T + A 
00980     template <typename SymmM, typename VctX, typename VctY>
00981     inline 
00982     void spr2 (VctX const& x, VctY const& y, SymmM& a) { 
00983 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
00984       typedef typename traits::matrix_traits<SymmM>::value_type val_t; 
00985 #else
00986       typedef typename SymmM::value_type val_t; 
00987 #endif 
00988       spr2 ((val_t) 1, x, y, a); 
00989     }
00990 
00991 
00992     // A <- alpha * x * x^H + A 
00993     // A hermitian 
00994     template <typename T, typename HermM, typename VctX>
00995     inline 
00996     void her (CBLAS_UPLO const uplo, T const& alpha, VctX const& x, HermM& a) {
00997 #ifndef BOOST_NUMERIC_BINDINGS_NO_STRUCTURE_CHECK 
00998       BOOST_STATIC_ASSERT((boost::is_same<
00999         typename traits::matrix_traits<HermM>::matrix_structure, 
01000         traits::general_t
01001       >::value)); 
01002 #endif 
01003 
01004       int const n = traits::matrix_size1 (a);
01005       assert (n == traits::matrix_size2 (a)); 
01006       assert (traits::vector_size (x) >= n); 
01007 
01008       CBLAS_ORDER const stor_ord
01009         = enum_cast<CBLAS_ORDER const>
01010         (storage_order<
01011 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
01012            typename traits::matrix_traits<HermM>::ordering_type
01013 #else
01014            typename HermM::orientation_category 
01015 #endif
01016          >::value); 
01017 
01018       detail::her (stor_ord, uplo, n, alpha, 
01019 #ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
01020                    traits::vector_storage (x), 
01021 #else
01022                    traits::vector_storage_const (x), 
01023 #endif
01024                    traits::vector_stride (x),
01025                    traits::matrix_storage (a), 
01026                    traits::leading_dimension (a)); 
01027     }
01028 
01029     template <typename T, typename HermM, typename VctX>
01030     inline 
01031     void her (T const& alpha, VctX const& x, HermM& a) {
01032 #ifndef BOOST_NUMERIC_BINDINGS_NO_STRUCTURE_CHECK 
01033       BOOST_STATIC_ASSERT((boost::is_same<
01034         typename traits::matrix_traits<HermM>::matrix_structure, 
01035         traits::hermitian_t
01036       >::value)); 
01037 #endif 
01038 
01039       int const n = traits::matrix_size1 (a);
01040       assert (n == traits::matrix_size2 (a)); 
01041       assert (traits::vector_size (x) >= n); 
01042 
01043       CBLAS_ORDER const stor_ord
01044         = enum_cast<CBLAS_ORDER const>
01045         (storage_order<
01046 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
01047            typename traits::matrix_traits<HermM>::ordering_type
01048 #else
01049            typename HermM::orientation_category 
01050 #endif
01051          >::value); 
01052 
01053       CBLAS_UPLO const uplo
01054         = enum_cast<CBLAS_UPLO const>
01055         (uplo_triang<
01056 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
01057            typename traits::matrix_traits<HermM>::uplo_type
01058 #else
01059            typename HermM::packed_category 
01060 #endif 
01061          >::value); 
01062 
01063       detail::her (stor_ord, uplo, n, alpha, 
01064 #ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
01065                    traits::vector_storage (x), 
01066 #else
01067                    traits::vector_storage_const (x), 
01068 #endif
01069                    traits::vector_stride (x),
01070                    traits::matrix_storage (a), 
01071                    traits::leading_dimension (a)); 
01072     }
01073 
01074     // A <- x * x^H + A 
01075     template <typename HermM, typename VctX>
01076     inline 
01077     void her (VctX const& x, HermM& a) { 
01078 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
01079       typedef typename traits::matrix_traits<HermM>::value_type val_t; 
01080 #else
01081       typedef typename HermM::value_type val_t; 
01082 #endif 
01083       typedef typename traits::type_traits<val_t>::real_type real_t; 
01084       her ((real_t) 1, x, a); 
01085     }
01086 
01087 
01088     // A <- alpha * x * x^H + A 
01089     // A hermitian in packed form 
01090     template <typename T, typename HermM, typename VctX>
01091     inline 
01092     void hpr (T const& alpha, VctX const& x, HermM& a) {
01093 #ifndef BOOST_NUMERIC_BINDINGS_NO_STRUCTURE_CHECK 
01094       BOOST_STATIC_ASSERT((boost::is_same<
01095         typename traits::matrix_traits<HermM>::matrix_structure, 
01096         traits::hermitian_packed_t
01097       >::value)); 
01098 #endif 
01099 
01100       int const n = traits::matrix_size1 (a);
01101       assert (n == traits::matrix_size2 (a)); 
01102       assert (traits::vector_size (x) >= n); 
01103 
01104       CBLAS_ORDER const stor_ord
01105         = enum_cast<CBLAS_ORDER const>
01106         (storage_order<
01107 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
01108            typename traits::matrix_traits<HermM>::ordering_type
01109 #else
01110            typename HermM::orientation_category 
01111 #endif
01112          >::value); 
01113 
01114       CBLAS_UPLO const uplo
01115         = enum_cast<CBLAS_UPLO const>
01116         (uplo_triang<
01117 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
01118            typename traits::matrix_traits<HermM>::uplo_type
01119 #else
01120            typename HermM::packed_category 
01121 #endif 
01122          >::value); 
01123 
01124       detail::hpr (stor_ord, uplo, n, alpha, 
01125 #ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
01126                    traits::vector_storage (x), 
01127 #else
01128                    traits::vector_storage_const (x), 
01129 #endif
01130                    traits::vector_stride (x),
01131                    traits::matrix_storage (a));
01132     }
01133 
01134     // A <- x * x^H + A 
01135     template <typename HermM, typename VctX>
01136     inline 
01137     void hpr (VctX const& x, HermM& a) { 
01138 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
01139       typedef typename traits::matrix_traits<HermM>::value_type val_t; 
01140 #else
01141       typedef typename HermM::value_type val_t; 
01142 #endif 
01143       typedef typename traits::type_traits<val_t>::real_type real_t; 
01144       hpr ((real_t) 1, x, a); 
01145     }
01146 
01147 
01148     // A <- alpha * x * y^H + y * (alpha * x)^H + A 
01149     // A hermitian
01150     template <typename T, typename HermM, typename VctX, typename VctY>
01151     inline 
01152     void her2 (CBLAS_UPLO const uplo, T const& alpha, 
01153                VctX const& x, VctY const& y, HermM& a) 
01154     {
01155 #ifndef BOOST_NUMERIC_BINDINGS_NO_STRUCTURE_CHECK 
01156       BOOST_STATIC_ASSERT((boost::is_same<
01157         typename traits::matrix_traits<HermM>::matrix_structure, 
01158         traits::general_t
01159       >::value)); 
01160 #endif 
01161 
01162       int const n = traits::matrix_size1 (a);
01163       assert (n == traits::matrix_size2 (a)); 
01164       assert (traits::vector_size (x) >= n); 
01165       assert (traits::vector_size (y) >= n); 
01166 
01167       CBLAS_ORDER const stor_ord
01168         = enum_cast<CBLAS_ORDER const>
01169         (storage_order<
01170 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
01171            typename traits::matrix_traits<HermM>::ordering_type
01172 #else
01173            typename HermM::orientation_category 
01174 #endif
01175          >::value); 
01176 
01177       detail::her2 (stor_ord, uplo, n, alpha, 
01178 #ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
01179                     traits::vector_storage (x), 
01180 #else
01181                     traits::vector_storage_const (x), 
01182 #endif
01183                     traits::vector_stride (x),
01184 #ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
01185                     traits::vector_storage (y), 
01186 #else
01187                     traits::vector_storage_const (y), 
01188 #endif
01189                     traits::vector_stride (y),
01190                     traits::matrix_storage (a), 
01191                     traits::leading_dimension (a)); 
01192     }
01193 
01194     template <typename T, typename HermM, typename VctX, typename VctY>
01195     inline 
01196     void her2 (T const& alpha, VctX const& x, VctY const& y, HermM& a) {
01197 #ifndef BOOST_NUMERIC_BINDINGS_NO_STRUCTURE_CHECK 
01198       BOOST_STATIC_ASSERT((boost::is_same<
01199         typename traits::matrix_traits<HermM>::matrix_structure, 
01200         traits::hermitian_t
01201       >::value)); 
01202 #endif 
01203 
01204       int const n = traits::matrix_size1 (a);
01205       assert (n == traits::matrix_size2 (a)); 
01206       assert (traits::vector_size (x) >= n); 
01207       assert (traits::vector_size (y) >= n); 
01208 
01209       CBLAS_ORDER const stor_ord
01210         = enum_cast<CBLAS_ORDER const>
01211         (storage_order<
01212 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
01213            typename traits::matrix_traits<HermM>::ordering_type
01214 #else
01215            typename HermM::orientation_category 
01216 #endif
01217          >::value); 
01218 
01219       CBLAS_UPLO const uplo
01220         = enum_cast<CBLAS_UPLO const>
01221         (uplo_triang<
01222 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
01223            typename traits::matrix_traits<HermM>::uplo_type
01224 #else
01225            typename HermM::packed_category 
01226 #endif 
01227          >::value); 
01228 
01229       detail::her2 (stor_ord, uplo, n, alpha, 
01230 #ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
01231                     traits::vector_storage (x), 
01232 #else
01233                     traits::vector_storage_const (x), 
01234 #endif
01235                     traits::vector_stride (x),
01236 #ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
01237                     traits::vector_storage (y), 
01238 #else
01239                     traits::vector_storage_const (y), 
01240 #endif
01241                     traits::vector_stride (y),
01242                     traits::matrix_storage (a), 
01243                     traits::leading_dimension (a)); 
01244     }
01245 
01246     // A <- x * y^H + y * x^H + A 
01247     template <typename HermM, typename VctX, typename VctY>
01248     inline 
01249     void her2 (VctX const& x, VctY const& y, HermM& a) { 
01250 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
01251       typedef typename traits::matrix_traits<HermM>::value_type val_t; 
01252 #else
01253       typedef typename HermM::value_type val_t; 
01254 #endif 
01255       her2 ((val_t) 1, x, y, a); 
01256     }
01257 
01258 
01259     // A <- alpha * x * y^H + y * (alpha * x)^H + A 
01260     // A hermitian in packed form 
01261     template <typename T, typename HermM, typename VctX, typename VctY>
01262     inline 
01263     void hpr2 (T const& alpha, VctX const& x, VctY const& y, HermM& a) {
01264 #ifndef BOOST_NUMERIC_BINDINGS_NO_STRUCTURE_CHECK 
01265       BOOST_STATIC_ASSERT((boost::is_same<
01266         typename traits::matrix_traits<HermM>::matrix_structure, 
01267         traits::hermitian_packed_t
01268       >::value)); 
01269 #endif 
01270 
01271       int const n = traits::matrix_size1 (a);
01272       assert (n == traits::matrix_size2 (a)); 
01273       assert (traits::vector_size (x) >= n); 
01274       assert (traits::vector_size (y) >= n); 
01275 
01276       CBLAS_ORDER const stor_ord
01277         = enum_cast<CBLAS_ORDER const>
01278         (storage_order<
01279 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
01280            typename traits::matrix_traits<HermM>::ordering_type
01281 #else
01282            typename HermM::orientation_category 
01283 #endif
01284          >::value); 
01285 
01286       CBLAS_UPLO const uplo
01287         = enum_cast<CBLAS_UPLO const>
01288         (uplo_triang<
01289 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
01290            typename traits::matrix_traits<HermM>::uplo_type
01291 #else
01292            typename HermM::packed_category 
01293 #endif 
01294          >::value); 
01295 
01296       detail::hpr2 (stor_ord, uplo, n, alpha, 
01297 #ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
01298                     traits::vector_storage (x), 
01299 #else
01300                     traits::vector_storage_const (x), 
01301 #endif
01302                     traits::vector_stride (x),
01303 #ifndef BOOST_NO_FUNCTION_TEMPLATE_ORDERING
01304                     traits::vector_storage (y), 
01305 #else
01306                     traits::vector_storage_const (y), 
01307 #endif
01308                     traits::vector_stride (y),
01309                     traits::matrix_storage (a));
01310     }
01311 
01312     // A <- x * y^H + y * x^H + A 
01313     template <typename HermM, typename VctX, typename VctY>
01314     inline 
01315     void hpr2 (VctX const& x, VctY const& y, HermM& a) { 
01316 #ifndef BOOST_NUMERIC_BINDINGS_POOR_MANS_TRAITS
01317       typedef typename traits::matrix_traits<HermM>::value_type val_t; 
01318 #else
01319       typedef typename HermM::value_type val_t; 
01320 #endif 
01321       hpr2 ((val_t) 1, x, y, a); 
01322     }
01323 
01324 
01325   } // namespace atlas
01326 
01327 }}} 
01328 
01329 #endif // BOOST_NUMERIC_BINDINGS_CBLAS_LEVEL_2_HPP