Chris@16: /*
Chris@16:  [auto_generated]
Chris@16:  boost/numeric/odeint/util/ublas_wrapper.hpp
Chris@16: 
Chris@16:  [begin_description]
Chris@16:  Resizing for ublas::vector and ublas::matrix
Chris@16:  [end_description]
Chris@16: 
Chris@101:  Copyright 2011-2013 Mario Mulansky
Chris@101:  Copyright 2011-2013 Karsten Ahnert
Chris@16: 
Chris@16:  Distributed under the Boost Software License, Version 1.0.
Chris@16:  (See accompanying file LICENSE_1_0.txt or
Chris@16:  copy at http://www.boost.org/LICENSE_1_0.txt)
Chris@16:  */
Chris@16: 
Chris@16: 
Chris@16: #ifndef BOOST_NUMERIC_ODEINT_UTIL_UBLAS_WRAPPER_HPP_INCLUDED
Chris@16: #define BOOST_NUMERIC_ODEINT_UTIL_UBLAS_WRAPPER_HPP_INCLUDED
Chris@16: 
Chris@16: 
Chris@16: #include <boost/type_traits/integral_constant.hpp>
Chris@16: #include <boost/numeric/ublas/vector.hpp>
Chris@16: #include <boost/numeric/ublas/matrix.hpp>
Chris@16: #include <boost/numeric/ublas/lu.hpp>
Chris@101: #include <boost/numeric/ublas/vector_expression.hpp>
Chris@101: #include <boost/numeric/ublas/matrix_expression.hpp>
Chris@101: 
Chris@101: #include <boost/numeric/odeint/algebra/vector_space_algebra.hpp>
Chris@101: #include <boost/numeric/odeint/algebra/default_operations.hpp>
Chris@16: 
Chris@16: #include <boost/numeric/odeint/util/is_resizeable.hpp>
Chris@16: #include <boost/numeric/odeint/util/state_wrapper.hpp>
Chris@16: 
Chris@101: 
Chris@101: /* extend ublas by a few operations */
Chris@101: 
Chris@101: /* map norm_inf onto reduce( v , default_operations::maximum ) */
Chris@101: namespace boost { namespace numeric { namespace odeint {
Chris@101: 
Chris@101:     template< typename T , typename A >
Chris@101:     struct vector_space_norm_inf< boost::numeric::ublas::vector<T,A> >
Chris@101:     {
Chris@101:         typedef T result_type;
Chris@101: 
Chris@101:         result_type operator()( const boost::numeric::ublas::vector<T,A> &x ) const
Chris@101:         {
Chris@101:             return boost::numeric::ublas::norm_inf( x );
Chris@101:         }
Chris@101:     };
Chris@101: 
Chris@101: 
Chris@101:     template< class T , class L , class A >
Chris@101:     struct vector_space_norm_inf< boost::numeric::ublas::matrix<T,L,A> >
Chris@101:     {
Chris@101:         typedef T result_type;
Chris@101: 
Chris@101:         result_type operator()( const boost::numeric::ublas::matrix<T,L,A> &x ) const
Chris@101:         {
Chris@101:             return boost::numeric::ublas::norm_inf( x );
Chris@101:         }
Chris@101:     };
Chris@101: } } }
Chris@101: 
Chris@101: /* additional operations:
Chris@101:  * abs( v )
Chris@101:  * v / w
Chris@101:  * a + v
Chris@101:  */
Chris@101: namespace boost { namespace numeric { namespace ublas {
Chris@101: 
Chris@101: 
Chris@101:     // elementwise abs - calculates absolute values of the elements
Chris@101:     template<class T>
Chris@101:     struct scalar_abs: public scalar_unary_functor<T> {
Chris@101:         typedef typename scalar_unary_functor<T>::value_type value_type;
Chris@101:         typedef typename scalar_unary_functor<T>::argument_type argument_type;
Chris@101:         typedef typename scalar_unary_functor<T>::result_type result_type;
Chris@101: 
Chris@101:         static BOOST_UBLAS_INLINE
Chris@101:         result_type apply (argument_type t) {
Chris@101:             using std::abs;
Chris@101:             return abs (t);
Chris@101:         }
Chris@101:     };
Chris@101: 
Chris@101: 
Chris@101:     // (abs v) [i] = abs (v [i])
Chris@101:     template<class E>
Chris@101:     BOOST_UBLAS_INLINE
Chris@101:     typename vector_unary_traits<E, scalar_abs<typename E::value_type> >::result_type
Chris@101:     abs (const vector_expression<E> &e) {
Chris@101:         typedef typename vector_unary_traits<E, scalar_abs<typename E::value_type> >::expression_type expression_type;
Chris@101:         return expression_type (e ());
Chris@101:     }
Chris@101: 
Chris@101:     // (abs m) [i] = abs (m [i])
Chris@101:     template<class E>
Chris@101:     BOOST_UBLAS_INLINE
Chris@101:     typename matrix_unary1_traits<E, scalar_abs<typename E::value_type> >::result_type
Chris@101:     abs (const matrix_expression<E> &e) {
Chris@101:         typedef typename matrix_unary1_traits<E, scalar_abs<typename E::value_type> >::expression_type expression_type;
Chris@101:         return expression_type (e ());
Chris@101:     }
Chris@101: 
Chris@101: 
Chris@101:     // elementwise division (v1 / v2) [i] = v1 [i] / v2 [i]
Chris@101:     template<class E1, class E2>
Chris@101:     BOOST_UBLAS_INLINE
Chris@101:     typename vector_binary_traits<E1, E2, scalar_divides<typename E1::value_type,
Chris@101:                                                          typename E2::value_type> >::result_type
Chris@101:     operator / (const vector_expression<E1> &e1,
Chris@101:                 const vector_expression<E2> &e2) {
Chris@101:         typedef typename vector_binary_traits<E1, E2, scalar_divides<typename E1::value_type,
Chris@101:                                                                      typename E2::value_type> >::expression_type expression_type;
Chris@101:         return expression_type (e1 (), e2 ());
Chris@101:     }
Chris@101: 
Chris@101: 
Chris@101:     // elementwise division (m1 / m2) [i] = m1 [i] / m2 [i]
Chris@101:     template<class E1, class E2>
Chris@101:     BOOST_UBLAS_INLINE
Chris@101:     typename matrix_binary_traits<E1, E2, scalar_divides<typename E1::value_type,
Chris@101:     typename E2::value_type> >::result_type
Chris@101:     operator / (const matrix_expression<E1> &e1,
Chris@101:                 const matrix_expression<E2> &e2) {
Chris@101:         typedef typename matrix_binary_traits<E1, E2, scalar_divides<typename E1::value_type,
Chris@101:                                                                      typename E2::value_type> >::expression_type expression_type;
Chris@101:         return expression_type (e1 (), e2 ());
Chris@101:     }
Chris@101: 
Chris@101:     // addition with scalar
Chris@101:     // (t + v) [i] = t + v [i]
Chris@101:     template<class T1, class E2>
Chris@101:     BOOST_UBLAS_INLINE
Chris@101:     typename enable_if< is_convertible<T1, typename E2::value_type >,
Chris@101:     typename vector_binary_scalar1_traits<const T1, E2, scalar_plus<T1, typename E2::value_type> >::result_type
Chris@101:     >::type
Chris@101:     operator + (const T1 &e1,
Chris@101:                 const vector_expression<E2> &e2) {
Chris@101:         typedef typename vector_binary_scalar1_traits<const T1, E2, scalar_plus<T1, typename E2::value_type> >::expression_type expression_type;
Chris@101:         return expression_type (e1, e2 ());
Chris@101:     }
Chris@101: 
Chris@101:     // addition with scalar
Chris@101:     // (t + m) [i] = t + m [i]
Chris@101:     template<class T1, class E2>
Chris@101:     BOOST_UBLAS_INLINE
Chris@101:     typename enable_if< is_convertible<T1, typename E2::value_type >,
Chris@101:     typename matrix_binary_scalar1_traits<const T1, E2, scalar_plus<T1, typename E2::value_type> >::result_type
Chris@101:     >::type
Chris@101:     operator + (const T1 &e1,
Chris@101:                 const matrix_expression<E2> &e2) {
Chris@101:         typedef typename matrix_binary_scalar1_traits<const T1, E2, scalar_plus<T1, typename E2::value_type> >::expression_type expression_type;
Chris@101:         return expression_type (e1, e2 ());
Chris@101:     }
Chris@101: 
Chris@101: } } }
Chris@101: 
Chris@101: 
Chris@101: 
Chris@101: 
Chris@101: /* add resize functionality */
Chris@16: namespace boost {
Chris@16: namespace numeric {
Chris@16: namespace odeint {
Chris@16: 
Chris@16: /*
Chris@16:  * resizeable specialization for boost::numeric::ublas::vector
Chris@16:  */
Chris@16: template< class T , class A >
Chris@16: struct is_resizeable< boost::numeric::ublas::vector< T , A > >
Chris@16: {
Chris@16:     typedef boost::true_type type;
Chris@16:     const static bool value = type::value;
Chris@16: };
Chris@16: 
Chris@16: 
Chris@16: /*
Chris@16:  * resizeable specialization for boost::numeric::ublas::matrix
Chris@16:  */
Chris@16: template< class T , class L , class A >
Chris@16: struct is_resizeable< boost::numeric::ublas::matrix< T , L , A > >
Chris@16: {
Chris@16:     typedef boost::true_type type;
Chris@16:     const static bool value = type::value;
Chris@16: };
Chris@16: 
Chris@16: 
Chris@16: /*
Chris@16:  * resizeable specialization for boost::numeric::ublas::permutation_matrix
Chris@16:  */
Chris@16: template< class T , class A >
Chris@16: struct is_resizeable< boost::numeric::ublas::permutation_matrix< T , A > >
Chris@16: {
Chris@16:     typedef boost::true_type type;
Chris@16:     const static bool value = type::value;
Chris@16: };
Chris@16: 
Chris@16: 
Chris@16: // specialization for ublas::matrix
Chris@16: // same size and resize specialization for matrix-matrix resizing
Chris@16: template< class T , class L , class A , class T2 , class L2 , class A2 >
Chris@16: struct same_size_impl< boost::numeric::ublas::matrix< T , L , A > , boost::numeric::ublas::matrix< T2 , L2 , A2 > >
Chris@16: {
Chris@16:     static bool same_size( const boost::numeric::ublas::matrix< T , L , A > &m1 ,
Chris@16:                            const boost::numeric::ublas::matrix< T2 , L2 , A2 > &m2 )
Chris@16:     {
Chris@16:         return ( ( m1.size1() == m2.size1() ) && ( m1.size2() == m2.size2() ) );
Chris@16:     }
Chris@16: };
Chris@16: 
Chris@16: template< class T , class L , class A , class T2 , class L2 , class A2 >
Chris@16: struct resize_impl< boost::numeric::ublas::matrix< T , L , A > , boost::numeric::ublas::matrix< T2 , L2 , A2 > >
Chris@16: {
Chris@16:     static void resize( boost::numeric::ublas::matrix< T , L , A > &m1 ,
Chris@16:                         const boost::numeric::ublas::matrix< T2 , L2 , A2 > &m2 )
Chris@16:     {
Chris@16:         m1.resize( m2.size1() , m2.size2() );
Chris@16:     }
Chris@16: };
Chris@16: 
Chris@16: 
Chris@16: 
Chris@16: // same size and resize specialization for matrix-vector resizing
Chris@16: template< class T , class L , class A , class T_V , class A_V >
Chris@16: struct same_size_impl< boost::numeric::ublas::matrix< T , L , A > , boost::numeric::ublas::vector< T_V , A_V > >
Chris@16: {
Chris@16:     static bool same_size( const boost::numeric::ublas::matrix< T , L , A > &m ,
Chris@16:                            const boost::numeric::ublas::vector< T_V , A_V > &v )
Chris@16:     {
Chris@16:         return ( ( m.size1() == v.size() ) && ( m.size2() == v.size() ) );
Chris@16:     }
Chris@16: };
Chris@16: 
Chris@16: template< class T , class L , class A , class T_V , class A_V >
Chris@16: struct resize_impl< boost::numeric::ublas::matrix< T , L , A > , boost::numeric::ublas::vector< T_V , A_V > >
Chris@16: {
Chris@16:     static void resize( boost::numeric::ublas::matrix< T , L , A > &m ,
Chris@16:                         const boost::numeric::ublas::vector< T_V , A_V > &v )
Chris@16:     {
Chris@16:         m.resize( v.size() , v.size() );
Chris@16:     }
Chris@16: };
Chris@16: 
Chris@16: 
Chris@16: 
Chris@16: // specialization for ublas::permutation_matrix
Chris@16: // same size and resize specialization for matrix-vector resizing
Chris@16: template< class T , class A , class T_V , class A_V >
Chris@101: struct same_size_impl< boost::numeric::ublas::permutation_matrix< T , A > ,
Chris@101:                        boost::numeric::ublas::vector< T_V , A_V > >
Chris@16: {
Chris@16:     static bool same_size( const boost::numeric::ublas::permutation_matrix< T , A > &m ,
Chris@16:                            const boost::numeric::ublas::vector< T_V , A_V > &v )
Chris@16:     {
Chris@16:         return ( m.size() == v.size() ); // && ( m.size2() == v.size() ) );
Chris@16:     }
Chris@16: };
Chris@16: 
Chris@16: template< class T , class A , class T_V , class A_V >
Chris@101: struct resize_impl< boost::numeric::ublas::vector< T_V , A_V > ,
Chris@101:                     boost::numeric::ublas::permutation_matrix< T , A > >
Chris@16: {
Chris@16:     static void resize( const boost::numeric::ublas::vector< T_V , A_V > &v,
Chris@16:                         boost::numeric::ublas::permutation_matrix< T , A > &m )
Chris@16:     {
Chris@16:         m.resize( v.size() , v.size() );
Chris@16:     }
Chris@16: };
Chris@16: 
Chris@16: 
Chris@16: 
Chris@16: 
Chris@16: 
Chris@16: 
Chris@16: 
Chris@16: template< class T , class A >
Chris@16: struct state_wrapper< boost::numeric::ublas::permutation_matrix< T , A > > // with resizing
Chris@16: {
Chris@16:     typedef boost::numeric::ublas::permutation_matrix< T , A > state_type;
Chris@16:     typedef state_wrapper< state_type > state_wrapper_type;
Chris@16: 
Chris@16:     state_type m_v;
Chris@16: 
Chris@16:     state_wrapper() : m_v( 1 ) // permutation matrix constructor requires a size, choose 1 as default
Chris@16:     { }
Chris@16: 
Chris@16: };
Chris@16: 
Chris@16: 
Chris@16: 
Chris@16: 
Chris@16: } } }
Chris@16: 
Chris@16: 
Chris@16: #endif // BOOST_NUMERIC_ODEINT_UTIL_UBLAS_WRAPPER_HPP_INCLUDED