Symaptic: os/ossrv/ossrv_pub/boost_apis/boost/numeric/ublas/blas.hpp@bde4ae8d615e

     1 //

     2 //  Copyright (c) 2000-2002

     3 //  Joerg Walter, Mathias Koch

     4 //

     5 //  Permission to use, copy, modify, distribute and sell this software

     6 //  and its documentation for any purpose is hereby granted without fee,

     7 //  provided that the above copyright notice appear in all copies and

     8 //  that both that copyright notice and this permission notice appear

     9 //  in supporting documentation.  The authors make no representations

    10 //  about the suitability of this software for any purpose.

    11 //  It is provided "as is" without express or implied warranty.

    12 //

    13 //  The authors gratefully acknowledge the support of

    14 //  GeNeSys mbH & Co. KG in producing this work.

    15 //

    17 #ifndef _BOOST_UBLAS_BLAS_

    18 #define _BOOST_UBLAS_BLAS_

    20 #include <boost/numeric/ublas/traits.hpp>

    22 namespace boost { namespace numeric { namespace ublas {

    24     namespace blas_1 {

    26           /** \namespace boost::numeric::ublas::blas_1

    27                   \brief wrapper functions for level 1 blas

    28           */

    31           /** \brief 1-Norm: \f$\sum_i |x_i|\f$

    32                   \ingroup blas1

    33            */

    34         template<class V>

    35         typename type_traits<typename V::value_type>::real_type

    36         asum (const V &v) {

    37             return norm_1 (v);

    38         }

    39           /** \brief 2-Norm: \f$\sum_i |x_i|^2\f$

    40                   \ingroup blas1

    41            */

    42         template<class V>

    43         typename type_traits<typename V::value_type>::real_type

    44         nrm2 (const V &v) {

    45             return norm_2 (v);

    46         }

    47           /** \brief element with larges absolute value: \f$\max_i |x_i|\f$

    48                   \ingroup blas1

    49           */

    50         template<class V>

    51         typename type_traits<typename V::value_type>::real_type

    52         amax (const V &v) {

    53             return norm_inf (v);

    54         }

    56           /** \brief inner product of vectors \a v1 and \a v2

    57                   \ingroup blas1

    58           */

    59         template<class V1, class V2>

    60         typename promote_traits<typename V1::value_type, typename V2::value_type>::promote_type

    61         dot (const V1 &v1, const V2 &v2) {

    62             return inner_prod (v1, v2);

    63         }

    65           /** \brief copy vector \a v2 to \a v1

    66                   \ingroup blas1

    67           */

    68         template<class V1, class V2>

    69         V1 &

    70         copy (V1 &v1, const V2 &v2) {

    71             return v1.assign (v2);

    72         }

    74           /** \brief swap vectors \a v1 and \a v2

    75                   \ingroup blas1

    76           */

    77         template<class V1, class V2>

    78         void swap (V1 &v1, V2 &v2) {

    79             v1.swap (v2);

    80         }

    82           /** \brief scale vector \a v with scalar \a t

    83                   \ingroup blas1

    84           */

    85         template<class V, class T>

    86         V &

    87         scal (V &v, const T &t) {

    88             return v *= t;

    89         }

    91           /** \brief compute \a v1 = \a v1 + \a t * \a v2

    92                   \ingroup blas1

    93           */

    94         template<class V1, class T, class V2>

    95         V1 &

    96         axpy (V1 &v1, const T &t, const V2 &v2) {

    97             return v1.plus_assign (t * v2);

    98         }

   100           /** \brief apply plane rotation

   101                   \ingroup blas1

   102           */

   103         template<class T1, class V1, class T2, class V2>

   104         void

   105         rot (const T1 &t1, V1 &v1, const T2 &t2, V2 &v2) {

   106             typedef typename promote_traits<typename V1::value_type, typename V2::value_type>::promote_type promote_type;

   107             vector<promote_type> vt (t1 * v1 + t2 * v2);

   108             v2.assign (- t2 * v1 + t1 * v2);

   109             v1.assign (vt);

   110         }

   112     }

   114     namespace blas_2 {

   116           /** \namespace boost::numeric::ublas::blas_2

   117                   \brief wrapper functions for level 2 blas

   118           */

   120           /** \brief multiply vector \a v with triangular matrix \a m

   121                   \ingroup blas2

   122                   \todo: check that matrix is really triangular

   123           */

   124         template<class V, class M>

   125         V &

   126         tmv (V &v, const M &m) {

   127             return v = prod (m, v);

   128         }

   130           /** \brief solve \a m \a x = \a v in place, \a m is triangular matrix

   131                   \ingroup blas2

   132           */

   133         template<class V, class M, class C>

   134         V &

   135         tsv (V &v, const M &m, C) {

   136             return v = solve (m, v, C ());

   137         }

   139           /** \brief compute \a v1 = \a t1 * \a v1 + \a t2 * (\a m * \a v2)

   140                   \ingroup blas2

   141           */

   142         template<class V1, class T1, class T2, class M, class V2>

   143         V1 &

   144         gmv (V1 &v1, const T1 &t1, const T2 &t2, const M &m, const V2 &v2) {

   145             return v1 = t1 * v1 + t2 * prod (m, v2);

   146         }

   148           /** \brief rank 1 update: \a m = \a m + \a t * (\a v1 * \a v2<sup>T</sup>)

   149                   \ingroup blas2

   150           */

   151         template<class M, class T, class V1, class V2>

   152         M &

   153         gr (M &m, const T &t, const V1 &v1, const V2 &v2) {

   154 #ifndef BOOST_UBLAS_SIMPLE_ET_DEBUG

   155             return m += t * outer_prod (v1, v2);

   156 #else

   157             return m = m + t * outer_prod (v1, v2);

   158 #endif

   159         }

   161           /** \brief symmetric rank 1 update: \a m = \a m + \a t * (\a v * \a v<sup>T</sup>)

   162                   \ingroup blas2

   163           */

   164         template<class M, class T, class V>

   165         M &

   166         sr (M &m, const T &t, const V &v) {

   167 #ifndef BOOST_UBLAS_SIMPLE_ET_DEBUG

   168             return m += t * outer_prod (v, v);

   169 #else

   170             return m = m + t * outer_prod (v, v);

   171 #endif

   172         }

   173           /** \brief hermitian rank 1 update: \a m = \a m + \a t * (\a v * \a v<sup>H</sup>)

   174                   \ingroup blas2

   175           */

   176         template<class M, class T, class V>

   177         M &

   178         hr (M &m, const T &t, const V &v) {

   179 #ifndef BOOST_UBLAS_SIMPLE_ET_DEBUG

   180             return m += t * outer_prod (v, conj (v));

   181 #else

   182             return m = m + t * outer_prod (v, conj (v));

   183 #endif

   184         }

   186           /** \brief symmetric rank 2 update: \a m = \a m + \a t *

   187                   (\a v1 * \a v2<sup>T</sup> + \a v2 * \a v1<sup>T</sup>)

   188                   \ingroup blas2

   189           */

   190         template<class M, class T, class V1, class V2>

   191         M &

   192         sr2 (M &m, const T &t, const V1 &v1, const V2 &v2) {

   193 #ifndef BOOST_UBLAS_SIMPLE_ET_DEBUG

   194             return m += t * (outer_prod (v1, v2) + outer_prod (v2, v1));

   195 #else

   196             return m = m + t * (outer_prod (v1, v2) + outer_prod (v2, v1));

   197 #endif

   198         }

   199           /** \brief hermitian rank 2 update: \a m = \a m +

   200                   \a t * (\a v1 * \a v2<sup>H</sup>)

   201                   + \a v2 * (\a t * \a v1)<sup>H</sup>)

   202                   \ingroup blas2

   203           */

   204         template<class M, class T, class V1, class V2>

   205         M &

   206         hr2 (M &m, const T &t, const V1 &v1, const V2 &v2) {

   207 #ifndef BOOST_UBLAS_SIMPLE_ET_DEBUG

   208             return m += t * outer_prod (v1, conj (v2)) + type_traits<T>::conj (t) * outer_prod (v2, conj (v1));

   209 #else

   210             return m = m + t * outer_prod (v1, conj (v2)) + type_traits<T>::conj (t) * outer_prod (v2, conj (v1));

   211 #endif

   212         }

   214     }

   216     namespace blas_3 {

   218           /** \namespace boost::numeric::ublas::blas_3

   219                   \brief wrapper functions for level 3 blas

   220           */

   222           /** \brief triangular matrix multiplication

   223                   \ingroup blas3

   224           */

   225         template<class M1, class T, class M2, class M3>

   226         M1 &

   227         tmm (M1 &m1, const T &t, const M2 &m2, const M3 &m3) {

   228             return m1 = t * prod (m2, m3);

   229         }

   231           /** \brief triangular solve \a m2 * \a x = \a t * \a m1 in place,

   232                   \a m2 is a triangular matrix

   233                   \ingroup blas3

   234           */

   235         template<class M1, class T, class M2, class C>

   236         M1 &

   237         tsm (M1 &m1, const T &t, const M2 &m2, C) {

   238             return m1 = solve (m2, t * m1, C ());

   239         }

   241           /** \brief general matrix multiplication

   242                   \ingroup blas3

   243           */

   244         template<class M1, class T1, class T2, class M2, class M3>

   245         M1 &

   246         gmm (M1 &m1, const T1 &t1, const T2 &t2, const M2 &m2, const M3 &m3) {

   247             return m1 = t1 * m1 + t2 * prod (m2, m3);

   248         }

   250           /** \brief symmetric rank k update: \a m1 = \a t * \a m1 +

   251                   \a t2 * (\a m2 * \a m2<sup>T</sup>)

   252                   \ingroup blas3

   253                   \todo use opb_prod()

   254           */

   255         template<class M1, class T1, class T2, class M2>

   256         M1 &

   257         srk (M1 &m1, const T1 &t1, const T2 &t2, const M2 &m2) {

   258             return m1 = t1 * m1 + t2 * prod (m2, trans (m2));

   259         }

   260           /** \brief hermitian rank k update: \a m1 = \a t * \a m1 +

   261                   \a t2 * (\a m2 * \a m2<sup>H</sup>)

   262                   \ingroup blas3

   263                   \todo use opb_prod()

   264           */

   265         template<class M1, class T1, class T2, class M2>

   266         M1 &

   267         hrk (M1 &m1, const T1 &t1, const T2 &t2, const M2 &m2) {

   268             return m1 = t1 * m1 + t2 * prod (m2, herm (m2));

   269         }

   271           /** \brief generalized symmetric rank k update:

   272                   \a m1 = \a t1 * \a m1 + \a t2 * (\a m2 * \a m3<sup>T</sup>)

   273                   + \a t2 * (\a m3 * \a m2<sup>T</sup>)

   274                   \ingroup blas3

   275                   \todo use opb_prod()

   276           */

   277         template<class M1, class T1, class T2, class M2, class M3>

   278         M1 &

   279         sr2k (M1 &m1, const T1 &t1, const T2 &t2, const M2 &m2, const M3 &m3) {

   280             return m1 = t1 * m1 + t2 * (prod (m2, trans (m3)) + prod (m3, trans (m2)));

   281         }

   282           /** \brief generalized hermitian rank k update:

   283                   \a m1 = \a t1 * \a m1 + \a t2 * (\a m2 * \a m3<sup>H</sup>)

   284                   + (\a m3 * (\a t2 * \a m2)<sup>H</sup>)

   285                   \ingroup blas3

   286                   \todo use opb_prod()

   287           */

   288         template<class M1, class T1, class T2, class M2, class M3>

   289         M1 &

   290         hr2k (M1 &m1, const T1 &t1, const T2 &t2, const M2 &m2, const M3 &m3) {

   291             return m1 = t1 * m1 + t2 * prod (m2, herm (m3)) + type_traits<T2>::conj (t2) * prod (m3, herm (m2));

   292         }

   294     }

   296 }}}

   298 #endif

author	sl@SLION-WIN7.fritz.box
	Fri, 15 Jun 2012 03:10:57 +0200
changeset 0	bde4ae8d615e
permissions	-rw-r--r--